积性函数
若定义域为全体正整数的函数(f(x))满足(forall gcd(n,m)=1,f(nm)=f(n)f(m)),则称函数(f(x))为积性函数。
常见的积性函数有(epsilon(n), varphi(n), mu(n), sigma_{k}(n), id_k(n))等。
若(f(n), g(n))都是积性函数,则(f(x^p), f^p(x), f(x)g(x), (f*g)(n))也是积性函数。
线性筛积性函数
线性筛可以在(O(n))的时间复杂度筛出(f(x)),当且仅当(f(x))满足:可以(O(1))求出(f(1),f(p),f(p^k))的函数值。
首先复习一下线性筛素数的代码:
int pri[N], tot, vis[N];
void sieve(int n)
{
for (int i = 2; i <= n; ++i)
{
if (!vis[i])
{
pri[++tot] = i;
//1
}
for (int j = 1; j <= tot && pri[j] * i <= n; ++j)
{
vis[pri[j] * i] = 1;
//2
if (i % pri[j] == 0) break;
}
}
}
显然我们需要在代码中的(1、2)处求(f(x))的值。
(1)处(i)为素数,(f(x))可以(O(1))求出。
(2)处复杂一点,其实是要求出(f(pri_j*i))的值。
注意到:设(i=prod_{i=1}^{k}p_i^{alpha_i}),则(pri_jleq p_1)。因为(pri_j=p_1)后就会(break)。
若(pri_j<p_1),则(f(pri_j*i)=f(pri_j)f(i))。
若(pri_j=p_1),设(low_i=p_1^{alpha_1})(就是上面唯一分解(i)中的第一项),(f(pri_j*i)=f(frac{i}{low_i})f(low_i*pri_j))
有一个例外,当(low_i=i)时,(f(pri_j*i)=f(p_k)),可以(O(1))求出。
细节见代码
void Sieve(int N)
{
f[1] = ...;//求f(1)
low[1] = 1;
for (int i = 2; i <= N; ++i)
{
if (!vis[i])
{
low[i] = pri[++tot] = i;
f[i] = ...;//求f(p)
}
for (int j = 1; j <= tot && pri[j] * i <= N; ++j)
{
vis[pri[j] * i] = 1;
if (i % pri[j] == 0)
{
low[pri[j] * i] = low[i] * pri[j];
if (low[i] == i)
f[pri[j] * i] = ...;//求f(p^k)(一般由f(p^(k-1))推出)
else
f[pri[j] * i] = f[i / low[i]] * f[low[i] * pri[j]];
break;
}
low[pri[j] * i] = pri[j];
f[pri[j] * i] = f[i] * f[pri[j]];
}
}
}
例题
这是一道推柿子题,不过推柿子的过程不是本文重点,仅讨论最后得到结果的求法。
结论是:
设(f(n)=sum_{d|n}dcdotvarphi(d)),不难证(f(n))为积性函数。
线性筛预处理(f(n))即可,时间复杂度(O(n+T))。
#include<bits/stdc++.h>
using namespace std;
#define int long long
inline int read()
{
int x = 0, f = 1; char ch = getchar();
for (; ch < '0' || ch > '9'; ch = getchar()) if (ch == '-') f = -1;
for (; ch >= '0' && ch <= '9'; ch = getchar()) x = (x << 1) + (x << 3) + ch - '0';
return x * f;
}
const int N = 1e6 + 5;
int p[N], tot, f[N], low[N]; //f_n = sum_{d | n} phi_d * d;
bool vis[N];
void Sieve()
{
f[1] = 1;
low[1] = 1;
for (int i = 2; i < N; ++i)
{
if (!vis[i])
{
low[i] = p[++tot] = i;
f[i] = i * (i - 1) + 1;
}
for (int j = 1; j <= tot && p[j] * i < N; ++j)
{
vis[p[j] * i] = 1;
if (i % p[j] == 0)
{
low[p[j] * i] = low[i] * p[j];
if (low[i] == i)
f[p[j] * i] = f[i] + p[j] * i * i * (p[j] - 1);
else
f[p[j] * i] = f[i / low[i]] * f[low[i] * p[j]];
break;
}
low[p[j] * i] = p[j];
f[p[j] * i] = f[i] * f[p[j]];
}
}
}
signed main()
{
Sieve();
int T = read();
while (T--)
{
int n = read();
printf("%lld
", (f[n] + 1) * n / 2);
}
return 0;
}