目录
数论函数
记号与约定
称函数 \(f\) 为数论函数,当且仅当其定义域为 \(\mathbb Z\),陪域为 \(\mathbb C\)。
称数论函数 \(f\) 为积性函数当且仅当对于任意 \(a,b \in \mathbb N_+,\gcd(a,b) = 1\),有 \(f(ab) = f(a)f(b)\)。
称数论函数 \(f\) 为完全积性函数当且仅当对于任意 \(a,b \in \mathbb N_+\),有 \(f(ab) = f(a)f(b)\)。
按照如下记号记录部分数论函数:
- \(I(n) = 1\)
- \(\epsilon(n) = [n = 1]\)
- \(id(n) = n\)
- \(id_k(n) = n^k\)
- \(d(n) = \sum \limits_{d \mid n} 1\)
- \(\sigma(n) = \sum \limits_{d \mid n} d\)
- \(\mu,\varphi\) 定义照常,不再赘述。
其中 \(I,\epsilon,id,id_k\) 完全积性。
另外,文中部分 \(a/b,\frac{a}{b}\) 型分数根据上下文可能有向下取整(\(\left \lfloor \frac ab \right \rfloor\))的含义,请注意甄别。
线性筛
对于任意积性函数,我们可以线性筛求出。
我们只需要关心积性函数在素数幂 \(p^k\) 处的取值。为了检验一个数 \(i\) 是否是 \(p^k\) 型的,我们可以维护 \(low(i)\) 表示 \(i\) 中最小质因子的部分。例如 \(low(6) = 2,low(49) = 49,low(12)=4\)。当且仅当 \(i = low(i)\) 时 \(i\) 是 \(p^k\) 型的。
具体的 \(k\) 和 \(p\) 也是不难维护的。
下面的模板来自 Alex-Wei。
for(int i = 2; i < N; i++) {
if(!vis[i]) pr[++pcnt] = i, f[i] = ..., low[i] = i; // 单独算 f(p)
for(int j = 1; j <= pcnt && i * pr[j] < N; j++) {
vis[i * pr[j]] = 1;
if(i % pr[j] == 0) { // i 与 p 不互质
low[i * pr[j]] = low[i] * pr[j];
if(i == low[i]) f[i * pr[j]] = ...; // i = p ^ k,单独算 f(p ^ {k + 1})
else f[i * pr[j]] = f[i / low[i]] * f[low[i * pr[j]]];
break;
}
low[i * pr[j]] = pr[j];
f[i * pr[j]] = f[i] * f[pr[j]]; // i 与 p 互质,f(ip) = f(i)f(p)
}
}
两例变量分离
-
\(d(ij)\)
\[d(ij) = \sum \limits_{x \mid i} \sum \limits_{y \mid j} [\gcd(x,y) = 1] \]构造双射容易证明。
-
\(\varphi(ij)\)
\[\varphi(ij) =\frac{\varphi(i)\varphi(j)\gcd(i,j)}{\varphi(\gcd(i,j))} \]证明:
\[\varphi(i)\varphi(j)\gcd(i,j) = \gcd(i,j) \left(i \prod \limits_{p \mid i,p \in \mathbb P} \frac{p-1}{p}\right) \left(j \prod \limits_{q \mid j,q \in \mathbb P} \frac{q-1}{q}\right) \\ = ij\gcd(i,j) \left(\prod \limits_{p \mid ij,p \in \mathbb P} \frac{p-1}{p}\right) \left(\prod \limits_{q \mid \gcd(i,j),q \in \mathbb P} \frac{q-1}{q}\right) \\ = \varphi(ij)\varphi(\gcd(i,j)) \]
关于 \(\mu^2(i)\)
有性质:
\[\mu^2(i) = \sum \limits_{d^2 \mid i} \mu(d) \]当 \(\mu(i) \neq 0\) 时,必然不存在平方因子,因此恰有 \(d = 1\) 造成 \(1\) 的贡献。
当 \(\mu(i) = 0\) 时,设 \(p_1,p_2,\cdots,p_k\) 这些质因子在 \(i\) 中的次数至少为 \(2\)。\(d\) 要在右侧造成贡献,只能是这些质因子构成的集合的某个子集的乘积,根据二项式定理,贡献为 \(0\)。
狄利克雷卷积
狄利克雷卷积
狄利克雷卷积:是一个算子。对于数论函数 \(f,g\),我们有
\[(f*g)(n) = \sum \limits_{d \mid n}f(n)g(\frac nd) \]狄利克雷逆:若 \(f* g = \epsilon\) 则称 \(f,g\) 互为狄利克雷逆。
对于 \(f(1) \neq 0\) 的 \(f\),\(f\) 的狄利克雷逆一定存在。显然这条件已经必要。如下推导得到 \(g\):
\[\epsilon(n) = \sum \limits_{d | n} f(d) g(\frac n d) \\ f(1)g(n) = \epsilon(n) - \sum \limits_{d \mid n,d > 1}f(d)g(\frac{n}{d}) \\ g(n) = \frac{1}{f(1)}\left(\epsilon(n) - \sum \limits_{d \mid n,d > 1}f(d)g(\frac{n}{d})\right) \]基本性质
-
交换律
\[f * g= g* f \] -
结合律
\[(f * g) * h = f* (g* h) \] -
分配律
\[(f\pm g) * h = f * h \pm g * h \] -
消去律
\[f * h = g * h \Leftrightarrow f = g \ \ (h(1) \neq 0) \] -
两个积性函数的狄利克雷卷积是积性函数
-
积性函数的狄利克雷逆是积性函数
几个简单的狄利克雷卷积
-
\(\mu * I = \epsilon\)
这也是为什么我们可以把 \([\gcd(x,y) = 1]\) 化开变成 \(\sum \limits_{d \mid \gcd(x,y)} \mu(d)\)。
-
\(\varphi * I = \text{id}\)
这也是为什么我们可以把 \(\gcd(x,y)\) 化开变成 \(\sum \limits_{d \mid \gcd(x,y)} \varphi(d)\)。
-
\(I * I = d\)
根据 \(d\) 的定义。
-
\(\text{id} * I = \sigma\)
根据 \(\sigma\) 的定义。
狄利克雷卷积的求法 (特别鸣谢: pp_orange)
\(f,g\) 已知,每一项可以 \(O(1)\) 求出。设我们要计算 \(h = f*g\) 的前 \(n\) 项。
当 \(f,g\) 是普通函数时
采用定义式 \(O(n \log n)\) 计算。
Code:[V 我 50 我帮你写]。
当 \(f\) 为积性函数时
考察 \(f = I\) 这一特例。此时等价于狄利克雷前缀和,复杂度为 \(O(n \log \log n)\)。
当 \(f\) 为任意积性函数时,我们可以执行类似的操作。考虑枚举质数 \(p\)。初始我们的卷积中只有 \(h'(n) = f(1)g(n)\) 这一项。每枚举一个质数 \(p\),我们就枚举不含因子 \(p\) 的正整数 \(n\),把 \(f(p)h(n),f(p^2)h(n),f(p^3)h(n),\cdots\ (p \nmid n)\) 这样的项贡献到对应的新的 \(h'(np),h'(np^2),h'(np^3),\cdots\) 上。这样我们就相当于在 \(p\) 这个维度做了一遍前缀和,现在 \(h'(n)\) 里面包含了所有 \(f(x)g(n/x)\) 的项,其中 \(x\) 整除 \(n\),且 \(x\) 的质因子集合中只有我们枚举过的质数。
因为高维前缀和是对的,那么这个操作也应当是对的。此时复杂度为 \(\sum \limits_{1 \leq p^k \leq n,k \ge 1,p \in \mathbb P} \left \lfloor \frac{n}{p^k} \right \rfloor = O(n \log \log n)\)。这是因为对于每个 \(p\),后面的部分都可以看作余项(等比数列求和不超过常数倍的第一项),因此其复杂度和狄利克雷前缀和 / 埃氏筛法相同。
Code:[V 我 50 我帮你写]。
当 \(f,g\) 为积性函数时
时间复杂度为 \(O(n)\)。
我们只需要求出 \(h\) 在素数幂处的取值。对于每个素数直接暴力,小于等于 \(\sqrt n\) 的素数复杂度可以看成 \(O\left(\dfrac {\sqrt n}{\log \sqrt n} \times \log^2 n\right) = O(\sqrt n \log n)\),而对于大于 \(\sqrt n\) 的素数只会有一项。因此,计算素数幂处取值的复杂度可以看作低阶项,不影响用 \(O()\) 计算的复杂度。
Code:[V 我 50 我帮你写]。
杜教筛
前置:\(\sum \limits_{i = 1}^{n} i^k \ \ (k\leq 4)\) 公式
- \(k = 1\):\(n(n+1)/2\)。
- \(k = 2\):\(n(n+1)(2n+1)/6\)。
- \(k = 3\):\(n^2(n+1)^2/4\)。
- \(k = 4\):\(n(n+1)(2n+1)(3n^2+3n-1)/30\)。已经很少用到了。
- 对于 \(k > 4\),使用连续点值拉格朗日插值可以 \(O(k)\) 计算一次。
流程
考虑求某数论函数 \(f\) 的前缀和 \(S(n) = \sum \limits_{i = 1}^n f(i)\)。杜教筛的想法是通过一个比较好的数论函数 \(g\) 和比较好的 \(f*g\) 把 \(f\) 的前缀和化简。这和 \(f,g,f*g\) 是否是积性函数无关。
我们知道
\[\sum \limits_{i = 1}^{n} (f*g)(i) = \sum \limits_{i = 1}^{n} g(i) S(\left \lfloor \frac ni \right \rfloor) \\ g(1)S(n ) = \sum \limits_{i = 1}^{n} (f*g)(i) - \sum \limits_{i = 2}^{n} g(i) S(\left \lfloor \frac ni \right \rfloor) \]如果 \(f*g\) 和 \(g\) 的前缀和比较简单,我们就可以快速求出 \(f\) 的前缀和了。
考虑如下伪代码:
inline int SumF(int n){ // 注意:这里计算的是 g(1)S(n). 但大部分时候 g 是积性函数或其点乘,因此 g(1)=1.
int ans=SumFG(n);
for(int l=2,r;l<=n;l=r+1){
r=n/(n/l);
ans-=(SumG(r)-SumG(L-1))*SumF(n/l);
}
return ans;
}
当使用记忆化的时候,根据整除分块经典结论,我们只会用到不超过 \(2 \sqrt n\) 个值。
那么 \(T(n) = O(\sqrt n) + \sum \limits_{i = 1}^{\sqrt n}\left(O(\sqrt i) +O\left(\sqrt {\dfrac{n}{i}}\right)\right)\)。求和改积分容易证得 \(T(n) = O(n^{3/4})\)。
否则时间复杂度不详(OI-Wiki 对不加以记忆化的时间复杂度证明作出了证伪:不能直接将第二层之后的 \(T\left(\sqrt{\dfrac{n}{ij}}\right)\) 看作余项)。
如果我们预处理部分 \(S(i) (i = 1,2,\cdots,m) \ \ (m \geq \sqrt n)\),设预处理复杂度为 \(T_0(m)\),则总复杂度为
\[T(n) = T_0(m) + \sum \limits_{k \in \{\left\lfloor n/x\right\rfloor\},k > m} T(k) \\ T_0(m) + \sum \limits_{k = 1}^{n/m} O(\sqrt{\frac nk}) \\ = O\left(T_0(m) +\frac{n}{\sqrt m}\right) \]那么一般 \(T_0(m) = O(m)\),此时 \(m = \Theta(n^{2/3})\) 最优。
详见 OI-Wiki。
两例简单的杜教筛
-
\(f = \mu\)
根据 \(\mu * I = \epsilon\),取 \(f = \mu,g = I, f*g = \epsilon\)。
-
\(f = \varphi\)
根据 \(\varphi * I = id\),取 \(f = \varphi,g = I,f * g = id\)。
点乘
定义 \(f,g\) 点乘的结果 \(f\cdot g\) 是一个函数,满足 \((f\cdot g)(n) = f(n)g(n)\)。
对于数论函数 \(A,B,C\),当 \(C\) 是完全积性函数的时候,我们有
\[(A\cdot C)*(B \cdot C) = (A*B) \cdot C \]那么当我们要求前缀和的函数 \(f\) 为两个简单函数点乘的时候,我们可以运用这个等式,从而寻求形式较好的 \(g\) 和 \(f*g\),然后杜教筛解决。
下面是两例比较基本的点乘处理方法。
-
\(f = \mu \cdot id_k\):取 \(g = id_k\)。
\[h = (\mu \cdot id_k )* id_k = (\mu \cdot id_k) * (I * id_k) = (\mu * I) * id_k = \epsilon \] -
\(f = \varphi \cdot id_k\):取 \(g = id_k\)。
\[h = (\varphi \cdot id_k )* id_k = (\varphi \cdot id_k) * (I * id_k) = (\varphi * I) * id_k = id_{k+1} \]
下面是抄的一例点乘:
- \(f = \mu^2 * (\mu \cdot id)\):取 \(g = id\)。\[h = \mu^2 * (\mu \cdot id) * id = \mu^2 * ((\mu * I) \cdot id) = \mu ^2 * \epsilon = \mu ^2 \]而 \(\mu^2\) 的前缀和可以通过公式 \(\mu^2(i) = \sum \limits_{d^2 \mid i} \mu(d)\) 计算:简单地交换求和号,前缀和变为 \(\sum \limits_{d = 1}^{\sqrt n} \mu(d) \left \lfloor \frac{n}{d^2} \right \rfloor\)。直接计算即可做到 \(O(\sqrt n)\),求出 \(\mu(i)\) 的前缀和可以进一步做到 \(O(n^{1/3})\)。但是在低于 \(O(\sqrt n)\) 时间复杂度内求出 \(\mu(i)\) 的前缀和需要再写筛法,其实就不如直接 \(O(\sqrt n)\) 求了。
另一例略显奔放的点乘(出自某场联考的最后一部分):
- \(f = id_k \cdot ((\mu \cdot id_k) * I)\):取 \(g = id_{2k}\)。\[h = (id_k\cdot ((\mu \cdot id_k) * I) ) * (id_k \cdot id_k)= ({\color{red}(\mu \cdot id_k) * id_k}* I) \cdot id_k \\ = (\epsilon * I) \cdot id_k \\ = id_k \]红色部分是第一例点乘中出现的结果,它等于 \(\epsilon\)。
小结:欲求 \(f = A\cdot C\) 形的前缀和,我们希望 \(g\) 形如 \(B \cdot C\),且 \(g\) 的前缀和易求。然后,运用上面的等式将 \(f * g\) 化作 \((A*B) \cdot C\)。适当地选取 \(B\) 可以将 \(A*B\) 化作简单的函数,从而 \(f*g\) 的前缀和易求。
简单运用
例 1
求
\[\sum \limits_{i = 1}^{n} \sum \limits_{j = 1}^{m} \varphi(\gcd(i,j)) \]对某个大质数取模
\(1 \leq n \leq m \leq 10^9\)
经过平凡的推导可以得出原式等于
\[\sum \limits_{T = 1}^{n} \left \lfloor \frac nT \right \rfloor \left \lfloor \frac mT \right \rfloor \sum \limits_{i \mid T} \mu(i) \varphi\left(\frac Ti\right) \]想要整除分块,则需要求出 \(\mu * \varphi\) 的前缀和。
令 \(f = \mu * \varphi,g = I,f*g = \varphi\)。利用杜教筛求出 \(\varphi\) 和 \(\mu * \varphi\) 的前缀和即可。注意到复杂度仍然是 \(O(n^{2/3})\) 的,因为我们可以线性筛预处理 \(f\) 的前 \(O(n^{2/3})\) 项。
例 2 (Luogu P3768 简单的数学题)
求
\[\sum \limits_{i = 1}^{n} \sum \limits_{j = 1}^{n} ij\gcd(i,j) \]对某个大质数取模
\(1 \leq n \leq 10^9\)
经过平凡的推导可以得出原式等于
\[\sum \limits_{d = 1}^{n} \varphi(d) d^2 S\left(\left \lfloor \frac ni\right \rfloor\right)^2, \text{where } S(x) = \frac{x(x+1)}{2} \]我们知道了怎么处理点乘,因此杜教筛也是容易的:\(f = \varphi \cdot id_2,g = id_2,(f*g) = id_3\)。
例 3
求
\[\sum \limits_{i = 1}^{n} \sum \limits_{j =1}^{n} \operatorname{lcm}(i,j)^k \]\(1 \leq n \leq 10^9,1 \leq k \leq 1000\)
这碟醋终于出现了!
平凡推导可得原式等于
\[\sum \limits_{T=1}^n S_k\left(\dfrac nT\right)^2 T^k \sum \limits_{t \mid T} \mu(t) \cdot t^k,\\ \text{where }S_k(n) = \sum \limits_{i = 1}^{n} i^k \]按照上面点乘处的推导,令 \(f(n) = n^k \sum \limits_{d \mid n} \mu(d) \cdot d^k\),请出 \(g(n) = n^{2k}\),\(f * g = id_k\)。
一个粗劣的复杂度上界是 \(O((\sqrt n + n^{2/3})k)\),但事实上整个杜教筛过程中也只会用到 \(O(\sqrt n)\) 个 \(S_k\) 的值,对它们记忆化可以分析出更好的复杂度 \(O(\sqrt n k+n^{2/3})\)。
[NOI2016] 循环之美
简要题面:
\[\sum \limits_{i = 1}^{n} \sum \limits_{j = 1}^m [\gcd(i,j) = 1][\gcd(j,k) = 1] \]\(1 \leq n,m \leq 10^9,1 \leq k \leq 10^3\)
这题处理技巧比较奇怪,我们不太能直接把两个艾弗森括号都拆成 \(\mu\),而是先保留下来 \([\gcd(j,k) = 1]\) 的那部分。在这里,和一个给定常数互质这个条件是比较好的。
\[\sum \limits_{i = 1}^{n} \sum \limits_{j = 1}^m [\gcd(i,j) = 1][\gcd(j,k) = 1] \\ = \sum \limits_{j = 1}^m [\gcd(j,k) = 1] \sum \limits_{i = 1}^{n} [\gcd(i,j) = 1] \\ = \sum \limits_{j = 1}^m [\gcd(j,k) = 1] \sum \limits_{i = 1}^{n} \sum \limits_{d \mid \gcd(i,j)} \mu(d) \\\\ = \sum \limits_{i = 1} ^{\min(n,m)} \mu(d) [\gcd(d,k) = 1] \left \lfloor \frac nd \right \rfloor \sum \limits_{j = 1}^{m/d} [\gcd(j,k) = 1] \]然后要求两个东西:
\[f(n) = \sum \limits_{i = 1}^{n} [\gcd(i,k) = 1] \\ g(n,k) = \sum \limits_{i = 1}^{n} [\gcd(i,k) = 1] \mu(i) \]对于 \(f(n)\),根据 \(\gcd(i,k) = \gcd(i+k,k)\),艾弗森括号的取值是 \(k\) - 循环的,我们只需要预处理 \(f(1) \sim f(k)\) 即可。
对于 \(g(n,k)\) 则比较需要有想象力。
\[g(n,k) = \sum \limits_{i = 1}^{n} \mu(i) \sum \limits_{d \mid \gcd(i,k)} \mu(d) \\ = \sum \limits_{1 \leq d \leq n, d \mid k} \mu(d) \sum \limits_{i=1}^{n/d} \mu(id) \\ = \sum \limits_{1 \leq d \leq n, d \mid k} \mu^2(d) \sum \limits_{i = 1}^{n/d} [\gcd(i,d) = 1]\mu(i) \\ = \sum \limits_{1 \leq d \leq n,d \mid k} \mu^2(d) g(n/d,d) \]递归计算,直到递归到 \(d = 1\) 时用杜教筛求 \(\mu\) 前缀和。
注意到最终要算 \(O(\sqrt n)\) 个 \(g(n,...)\),\(...\) 这一维是 \(k\) 的因子,因此总复杂度为 \(O(d(k)\sqrt n+ n^{2/3} + k)\)。
BZOJ3512: DZY Loves Math IV
求
\[\left(\sum \limits_{i = 1}^{n} \sum \limits_{j = 1}^{m} \varphi(ij) \right)\bmod 998244353 \]\(1 \leq n \leq 10^5,1 \leq m \leq 10^9,n \leq m\)
根据《毒瘤之神的考验》一题快进到
\[\sum \limits_{T =1}^{n} \left(\sum \limits_{i = 1}^{n/T} \varphi(iT) \sum \limits_{j = 1}^{m/T} \varphi(jT) \right)\sum \limits_{d \mid T} \frac{d\mu\left(\frac Td\right)}{\varphi(d)} \]记 \(S(a,b) = \sum \limits_{i = 1}^{a} \varphi(ib),f(n) = \sum \limits_{d \mid n} \frac{d\mu\left(\frac Td\right)}{\varphi(d)}\),则上式可写作
\[\sum \limits_{i = 1}^{n} f(i)S(n/i,i)S(m/i,i) \]在《毒瘤之神的考验》一题中,要用到的全体 \(S(x,y)\) 可 \(O(n \log n)\) 预处理,从而最终复杂度也是 \(O(n \log n)\) 的。
本题需要神秘的观察。我们考虑类似扰动法的策略,
\[\text{let } A(n,m) = \sum \limits_{i = 1}^{n} \sum \limits_{j = 1}^{m} \varphi(ij) = \sum \limits_{i = 1}^{n} f(i)S(n/i,i)S(m/i,i), \text{then} \\ S(n,k) = \sum \limits_{i = 1}^{n} \varphi(ik) \\ = A(n,k) - A(n,k-1) \\ = \sum \limits_{i = 1}^{n} f(i)S(n/i,i)\left(S(k/i,i) - S((k-1)/i,i)\right) \\ = \varphi(k) \sum \limits_{i \mid k} f(i) S(n/i,i) \]注意到 \(k\) 这一维是小于等于给定的 \(n\) 的(和此处的变量 \(n\) 作区分),于是对于 \(S(m/i,i)\) 这样的状态来说,计算它们的复杂度不太大。我们可以采用记忆化 + 预处理小的部分来减小递归的复杂度。
Min_25 筛
综述
考虑求一类积性函数的前缀和,其中 \(f\) 是积性函数,满足:
- 对于任意素数 \(p\),\(f(p)\) 处取值是关于 \(p\) 的低次多项式。
- \(f(p^c)\) 处取值容易求得。
则 Min_25 筛可以在 \(O(\frac{n^{3/4}}{\log n})\) 或 \(O(n^{1 - \epsilon})\) 的复杂度内求出 \(f\) 的前缀和。
流程
Min_25 筛的求解分为两步:
- 对于 \(W(n) = \left\{\left\lfloor \dfrac{n}{x}\right\rfloor \mid 1 \leq x \leq n\right\}\),考虑其中每一个元素 \(m\),求出 \(\sum \limits_{1 \leq i \leq m,i \in \mathbb P} f(i)\),即前缀素数点处取值之和;
- 根据第一步求得的 \(O(\sqrt n)\) 个答案进而求解出整个前缀和。
第一步
因为 \(f(p)\) 的取值是关于 \(p\) 的低次多项式,因此,在这一步中,我们只需要考虑 \(f(p) = p^k\) 的情形如何求解(对于项数更多的情形,分别求解后求和即可)。
我们通常找一个 \(f'\) 来替代 \(f\)。这样的 \(f'\) 具有完全积性,且在素数 \(p\) 处和 \(f\) 相同。这里考虑函数 \(f'(x) = x^k\)。
Min_25 筛的第一步,我们要对于每个 \(m\),求出数组 \(g\)。其中 \(g\) 的定义如下:
\[g(m,j) = \sum \limits_{1 \leq i \leq m,(i \in \mathbb P \operatorname {or} L(i) >P_j)} f'(i) \]其中 \(L(i)\) 表示 \(i\) 的最小质因子,\(P_j\) 表示第 \(j\) 个质数,\(P_0 = 0\)。上式即前缀素数点处或最小质因子大于 \(P_j\) 处 \(f'\) 取值之和(注意是 \(f'\) 而非 \(f\)。我们用 \(f'\) 暂时替代了 \(f\))。
对于 \(1\),其不存在最小质因子,故 \([L(1) > P_j]\) 总是 \(0\)。此处可改作 \(2 \leq i \leq m\),但为了推导方便保留可能的 \(i = 1\) 取值,虽然 \(f'(1)\) 永远不能被统计进贡献。
设 \(q(m)\) 表示最大的 \(j\) 使得 \(P_j^2 \leq m\),那么 \(g(m,q(m))\) 即为 \(\sum \limits_{1 \leq i \leq m,i \in \mathbb P} f(i)\)。
考虑递推求出 \(g(m,j)\)。初始有 \(g(m,0) = \sum \limits_{i = 2}^{m} i^k\)。对于 \(j \geq 1\),我们会在 \(g(m,j-1)\) 的基础上删除若干个点。质数是不会被删除的,因此我们会恰好删除最小质因子是 \(P_j\) 的非质数点。
若 \(P_j^2 > m\),\(g(m,j)\) 不会有变化。
若 \(P_j^2 \leq m\),形如 \(P_jx(P_jx \leq m)\),且 \(P_jx\) 的最小质因子恰好是 \(P_j\) 的点会造成贡献。这可以表示为下面的式子:
\[f'(P_j)\left( g\left(\left \lfloor \dfrac{m}{P_j} \right \rfloor,j-1\right)-g(P_{j-1},j-1)\right) \]来解释一下这个式子:首先,我们可以根据完全积性把 \(f'(P_jx)\) 拆开。现在考虑 \(x\) 的取值,很显然是来源于 \(g\left(\left \lfloor \dfrac{m}{P_j} \right \rfloor,j-1\right)\)。但是既小于等于 \(\left \lfloor \dfrac{m}{P_j} \right \rfloor\) 又小于等于 \(P_{j-1}\) 的质数是不能作为 \(x\) 的,因为我们有 \(P_{j-1} < P_j \leq \sqrt m\),所以 \(\left \lfloor \dfrac{m}{P_j} \right \rfloor \geq P_{j-1}\),只需要考虑 \(g(P_{j-1},j-1)\) 就可以了。
因此,我们得到了下式:
\[g(m,j) = \begin{cases} g(m,j-1) & P_j^2 > m \\ g(m,j-1)-f'(P_j) \times \left( g\left(\left \lfloor \dfrac{m}{P_j} \right \rfloor,j-1\right)-g(P_{j-1},j-1)\right) & P_j^2 \leq m\end{cases} \]复杂度分析 · 一
现在考虑复杂度。根据经典结论 \(\left \lfloor \dfrac{n}{ab} \right \rfloor=\left \lfloor \dfrac{\frac{n}{a}}{b} \right \rfloor\),我们对于 \(m \in W(n)\) 有 \(W(m) \subseteq W(n)\)。同时可以证明 \(\forall 1 \leq i \leq \sqrt n\),\(i \in W(n)\),因此对于任意 \(m \in W(m)\) 我们求 \(g(m,j)\) 的时候,只会遍历 \(g\) 数组中第一维在 \(W(n)\) 中的元素。同时,求出单个状态的复杂度是 \(O(1)\) 的,因此状态数即为复杂度。
根据范围素数近似估计,复杂度可写作
\[T(n) = \sum \limits_{i \leq \sqrt n} O(\pi(\sqrt i)) + O\left(\pi\left (\sqrt \frac ni\right) \right) \\ = \sum \limits_{i \leq \sqrt n} O\left(\pi\left (\sqrt \frac ni\right) \right) \\ = \sum \limits_{i \leq \sqrt n} O\left( \dfrac{\sqrt \frac ni }{\ln \sqrt {\frac ni}}\right) \\ \]注意到 \(i \leq \sqrt n\),因此 \(\ln \sqrt \frac ni\) 只和 \(\ln n\) 差常数倍。直接将这部分提出来,上面的部分积分,得到
\[T(n ) =O\left(\frac{n^{3/4}}{\log n}\right) \]第二步
综述
现在我们已经知道所有 \(\sum \limits_{1 \leq i \leq m,i \in \mathbb P} f(i),\ \text{where}\ m \in W(n)\) 了。将其简记为 \(g(m)\)。考虑使用这些 \(g(m)\) 求出原来的解。
设 \(S(m,j)\) 表示 \(\sum \limits_{1 \leq i \leq m}[L(i) > P_j] f(i)\),即最小质因子大于 \(P_j\) 的 \(f(i)\) 之和(我们现在抛弃 \(f'\) 了,同时也不再限定 \(i\) 是质数;只不过 \(1\) 仍然被抛弃在外)。这样再调用 \(S(n,0)\) 就得到答案了。
方法一
这种方法复杂度为 \(O(n^{1 - \epsilon})\),在 \(n \leq 10^{13}\) 的时候也可以看作 \(O\left(\frac{n^{3/4}}{\log n}\right)\),够用。实际测试略快于第二种。这里着重介绍。
仍然考虑如何求解。分 \(i\) 是质数或合数两种情况讨论:
-
\(i\) 是质数
根据 \(g\) 的定义,这部分的贡献是 \(g(n) - g(P_j)\)。
-
\(i\) 是合数
枚举最小质因子 \(P_k (k > j)\),枚举指数 \(e\),利用积性函数的性质进行变量分离,得到
\[\sum \limits_{k > j} \sum \limits_{e \geq 1,P_{k}^e \leq m} f(P_{k}^e)\left(S\left(\left\lfloor\frac{m}{P_k^e}\right\rfloor,k\right)+ [e>1]\right) \]这里 \([e>1]\) 有些耐人寻味。考虑前面的 \(S\left(\left\lfloor\frac{m}{P_k^e}\right\rfloor,k\right)\) 没有考虑 \(P_{k}^e\) 自己的贡献。如果 \(e = 1\),那么 \(P_k^e\) 就是 \(P_k\),应该归到 \(i\) 是质数的分类里,没有贡献。否则,\(P_{k}^e\) 会造成 \(f(P_{k}^e)\) 的贡献。
则 \(S(n,0) + 1\) 就是答案(我们之前没有考虑 \(f(1)\) 的贡献)。
代码 (For Luogu P4213)
略好于之前的版本。
# include <bits/stdc++.h>
const int N=1000010,INF=0x3f3f3f3f,mod=1e9+7,inv6=166666668,inv2=500000004;
typedef long long ll;
int pr[N],pc;
bool vis[N];
int g1[N],g2[N];
ll w[N],idx1[N],idx2[N];
int wtot;
ll MAXN;
ll n;
inline int read(void){
int res,f=1;
char c;
while((c=getchar())<'0'||c>'9')
if(c=='-')f=-1;
res=c-48;
while((c=getchar())>='0'&&c<='9')
res=res*10+c-48;
return res*f;
}
inline void init(void){
for(int i=2;i<=MAXN;++i){
if(!vis[i]) pr[++pc]=i;
for(int j=1;i*pr[j]<=MAXN&&j<=pc;++j){
vis[i*pr[j]]=true;
if(i%pr[j]==0) break;
}
}
return;
}
inline int sum1(ll x){
x%=mod;
return x*(x+1)%mod*inv2%mod;
}
inline int sum2(ll x){
x%=mod;
return x*(x+1)%mod*(2*x+1)%mod*inv6%mod;
}
inline int idx(ll x){
return (x<=MAXN)?idx1[x]:idx2[n/x];
}
inline int adc(int a,int b){
return (a+b>=mod)?(a+b-mod):(a+b);
}
inline int dec(int a,int b){
return (a<b)?(a-b+mod):(a-b);
}
inline void add(int &a,int b){
a=adc(a,b);
}
inline void del(int &a,int b){
a=dec(a,b);
}
inline int mul(int a,int b){
return 1ll*a*b%mod;
}
ll S(ll x,int i){
if(pr[i]>=x) return 0;
int xid=idx(x),pid=idx(pr[i]);
ll res=dec(dec(g2[xid],g1[xid]),dec(g2[pid],g1[pid]));
for(int j=i+1;j<=pc&&1ll*pr[j]*pr[j]<=x;++j){
for(ll e=1,p=pr[j],_p;p<=x;++e,p*=pr[j]){
_p=p%mod;
res=adc(res,mul(mul(dec(_p,1),_p),S(x/p,j)+(e>1)));
}
}
return res;
}
int main(void){
scanf("%lld",&n);
MAXN=1ll*sqrt(n);
init();
for(ll l=1,r;l<=n;l=r+1){
ll val=n/l;
r=n/val,w[++wtot]=val;
g1[wtot]=dec(sum1(val),1),g2[wtot]=dec(sum2(val),1);
if(val<=MAXN) idx1[val]=wtot;
else idx2[n/val]=wtot;
}
for(int i=1;i<=pc;++i){
for(int j=1;1ll*pr[i]*pr[i]<=w[j]&&j<=wtot;++j){
int mid=idx(w[j]/pr[i]),pid=idx(pr[i-1]);
del(g1[j],mul(pr[i],dec(g1[mid],g1[pid])));
del(g2[j],mul(mul(pr[i],pr[i]),dec(g2[mid],g2[pid])));
}
}
printf("%d",adc(S(n,0),1));
return 0;
}
方法二
再度审视 \(S\) 的定义:\(S(m,j)\) 表示 \(\sum \limits_{1 \leq i \leq m}[L(i) > P_j] f(i)\),即最小质因子大于 \(P_j\) 的 \(f(i)\) 之和。
我们直接套用第一步中的递推方式。
\[S(m,j) = S(m,j+1) + \\ [P_{j+1} \leq \sqrt m]\sum \limits_{P_{j+1}^{e} \leq m} f(P_{j+1}^e)\left( S\left( \left \lfloor\frac{m}{P_{j+1}^e} \right \rfloor,j+1 \right) - g\left (\min\left(\left \lfloor\frac{m}{P_{j+1}^e} \right \rfloor,P_{j+1}\right)\right) + [e>1] \right) \]这里的取 \(\min\) 是不好的。我们注意到当 \(P_{j+1}^{e+1} > m\) 的时候,\(S() - g()\) 没有贡献,因为这意味着 \(P_{j+1} > \dfrac{m}{P_{j+1}^e}\)。于是可以改为
\[S(m,j) = S(m,j+1) + \\ [P_{j+1} \leq \sqrt m]\sum \limits_{P_{j+1}^{e+1} \leq m} f(P_{j+1}^e)\left( S\left( \left \lfloor\frac{m}{P_{j+1}^e} \right \rfloor,j+1 \right) - g\left (P_{j+1}\right) \right) + f(P_{j+1}^{e+1}) \]使用和第一步相同的方法进行递推即可。
至此可以使用与第一步类似的方法分析出 \(O\left(\frac{n^{3/4}}{\log n}\right)\)。
方法二的优势是,我们可以求出每个 \(m \in W(n)\) 的前缀 \(m\) 的答案。这无疑在整除分块需要计算 \(O(\sqrt n)\) 个前缀的函数值之和的时候给我们提供了便利。
代码 (For LibreOJ 6784)
# include <bits/stdc++.h>
const int N=2000010,INF=0x3f3f3f3f,mod=1e9+7,inv6=166666668,inv2=500000004;
typedef long long ll;
int pr[N],pc;
bool vis[N];
int g1[N],g0[N];
ll w[N],idx1[N],idx2[N];
int s[N];
int wtot;
ll MAXN;
ll n;
inline int read(void){
int res,f=1;
char c;
while((c=getchar())<'0'||c>'9')
if(c=='-')f=-1;
res=c-48;
while((c=getchar())>='0'&&c<='9')
res=res*10+c-48;
return res*f;
}
inline void init(void){
for(int i=2;i<=MAXN;++i){
if(!vis[i]) pr[++pc]=i;
for(int j=1;i*pr[j]<=MAXN&&j<=pc;++j){
vis[i*pr[j]]=true;
if(i%pr[j]==0) break;
}
}
return;
}
inline int sum1(ll x){
x%=mod;
return x*(x+1)%mod*inv2%mod;
}
inline int sum0(ll x){
return x%mod;
}
inline int idx(ll x){
return (x<=MAXN)?idx1[x]:idx2[n/x];
}
inline int adc(int a,int b){
return (a+b>=mod)?(a+b-mod):(a+b);
}
inline int dec(int a,int b){
return (a<b)?(a-b+mod):(a-b);
}
inline void add(int &a,int b){
a=adc(a,b);
}
inline void del(int &a,int b){
a=dec(a,b);
}
inline int mul(int a,int b){
return 1ll*a*b%mod;
}
inline int f(int x,int k){
return x^k;
}
std::vector <int> vec;
int main(void){
scanf("%lld",&n);
MAXN=1ll*sqrt(n);
init();
for(ll l=1,r;l<=n;l=r+1){
ll val=n/l;
r=n/val,w[++wtot]=val;
g1[wtot]=dec(sum1(val),1),g0[wtot]=dec(sum0(val),1);
if(val<=MAXN) idx1[val]=wtot;
else idx2[n/val]=wtot;
}
for(int i=1;i<=pc;++i){
for(int j=1;1ll*pr[i]*pr[i]<=w[j]&&j<=wtot;++j){
int mid=idx(w[j]/pr[i]),pid=idx(pr[i-1]);
del(g1[j],mul(pr[i],dec(g1[mid],g1[pid])));
del(g0[j],dec(g0[mid],g0[pid]));
}
}
for(int i=1;i<=wtot;++i) s[i]=g1[i]=adc(dec(g1[i],g0[i]),(w[i]>=2)*2);
for(int i=pc;i;--i){
for(int j=1;1ll*pr[i]*pr[i]<=w[j]&&j<=wtot;++j){
for(ll e=1,p=pr[i];p<=w[j]/pr[i];++e,p*=pr[i]){
add(s[j],mul(f(pr[i],e),dec(s[idx(w[j]/p)],g1[idx(pr[i])])));
add(s[j],f(pr[i],e+1));
}
}
}
vec.resize(wtot);
for(int i=0;i<wtot;++i) vec[i]=adc(s[i+1],1);
std::sort(vec.begin(),vec.end());
vec.erase(std::unique(vec.begin(),vec.end()),vec.end());
int ans=0;
for(auto v:vec) ans^=v;
printf("%d",ans);
return 0;
}
小试牛刀
例 1
\[\sum \limits_{T=1}^n S_k\left(\dfrac nT\right)^2 T^k \sum \limits_{t \mid T} \mu(t) \cdot t^k \](杜教筛 例 3)
求
\[\sum \limits_{i = 1}^{n} \sum \limits_{j =1}^{n} \operatorname{lcm}(i,j)^k \]\(1 \leq n \leq 10^9,1 \leq k \leq 1000\)
快进到这里。
那么我们只要瞪眼法看出后面是积性函数,直接上 Min_25 筛就完事了,不需要观察到可以卷 \(id_{2k}\)。
Typora 崩掉了,不写了。
标签:right,函数,limits,leq,数论,sum,mu,求和,left From: https://www.cnblogs.com/Meatherm/p/18314262