标签：前缀 limits 狄利克 int 题解 sum varphi ij
先撇开这题不谈，通过这题我倒是发现了推了十几道莫反后连狄利克雷前缀和都没学扎实。
狄利克雷前缀和定义：\(f*1=g\)，我们称 \(g\) 是 \(f\) 的狄利克雷前缀和。
给定 \(f_1,f_2...f_n\)，考虑如何在 \(O(n\log\log n)\) 内求出 \(g\) 函数。
沿用高维前缀和的思路，我们可以把每一个质因子看作一维，然后套用高维前缀和的求法。
对于二位前缀和，我们是先累加其中一维，然后再累加另一维，形象地来说，如果把二位前缀和理解为一个二维矩阵，就是先对每行求出前缀和。
那么我们可以沿用一样的思路求出狄利克雷前缀和：
for (int i = 1; i <= cnt; ++ i)
for (int j = 1; j * Prime[i] <= n; ++ j) f[j * Prime[i]] += f[j];
狄利克雷后缀和：
\(g_i=\sum\limits^{d\le n}_{i\vert d}f_d\)，我们称 \(g\) 是 \(f\) 的狄利克雷后缀和。
求法：与狄利克雷前缀和不能说是非常相像，只能说是一模一样
for (int i = 1; i <= cnt; ++ i)
for (int j = n / Prime[i]; j; -- j) f[j] += f[j * Prime[i]];
倒推狄利克雷前缀和：
如果我们已知 \(f*1=g\)，如何快速求出 \(f\) 函数呢？
解法：前缀和的差分数组就是原数组，我们还是拿二维前缀和举例。首先对行做差分，得到的数组是对应列的前缀和。再做一遍差分即可。
但注意差分要从大到小做。
for (int i = 1; i <= cnt; ++ i)
for (int j = n / Prime[i]; j; -- j) f[j * Prime[i]] -= f[j * Prime[i]];
然后我们回到这道题：求
\[\sum\limits^n_{i=1}\sum\limits^m_{j=1}\text{gcd}^n(i,j)\sum\limits_{k=1}^{ij}[i\perp k][j\perp k]k \]\(1\le n,m\le 10^7\)。
将原式变为
\[\sum\limits^n_{i=1}\sum\limits^m_{j=1}\text{gcd}^n(i,j)\sum\limits_{k=1}^{ij}[ij\perp k]k \]然后我们发现这个 \(\sum\limits_{k=1}^{ij}[ij\perp k]k\) 很不好动，但如果 \(ij\perp k\)，则 \(ij\perp ij-k\)，所以这样的 \(k\) 总是成对出现，因此 \(\sum\limits_{k=1}^{ij}[ij\perp k]k=\frac{\varphi(ij)ij}{2}\)，\(ij=1\) 例外，需要特判。
我怎么会告诉你我是打表发现的呢
于是原式化为（我们忽略那个 \(\div 2\) 和 \(ij=1\) 的特例，因为这无关紧要）
\[\sum\limits^n_{i=1}\sum\limits^m_{j=1}\text{gcd}^n(i,j)\varphi(ij)ij \]仍然不好动，因为 \(ij\) 很大。这个时候就只有利用 \(\varphi(ij)=\frac{\varphi(i)\varphi(j)\gcd(i,j)}{\varphi(\gcd(i,j))}\)。
然后变成
\[\sum\limits^n_{i=1}\sum\limits^m_{j=1}\frac{\text{gcd}^{n+1}(i,j)}{\varphi(\gcd(i,j))}\varphi(i)\varphi(j)ij\\ =\sum\limits^n_{d=1}\frac{d^{n+1}}{\varphi(d)}\sum_{k=1}\mu(k)\sum^{\lfloor\frac{n}{dk}\rfloor}_{i=1}\varphi(idk)idk\sum^{\lfloor\frac{m}{dk}\rfloor}_{j=1}\varphi(jdk)jdk \]令 \(f_1(i)=\sum\limits_{j=1}^{\lfloor\frac{n}{i}\rfloor}\varphi(ij)ij,f_2(i)=\sum\limits_{j=1}^{\lfloor\frac{m}{i}\rfloor}\varphi(ij)ij\)，这俩玩意儿就是个狄利克雷后缀和，可以快速预处理。
枚举 \(dk\)，化为
\[\sum\limits^n_{i=1}f_1(i)f_2(i)\sum_{d\vert i}\frac{d^{n+1}}{\varphi(d)}\mu(\frac{i}{d}) \]令 \(f(i)=\frac{i^{n+1}}{\varphi(i)},g(i)=\sum\limits_{d\vert i}\frac{d^{n+1}}{\varphi(d)}\mu(\frac{i}{d})\)，显然 \(g=f*\mu\)。
狄利克雷卷积要 \(O(n\log n)\) 求，无法通过本题。但 \(g=f*\mu\) 意味着 \(f=g*1\)，也就是倒推狄利克雷前缀和。综上，总时间复杂度 \(O(n\log\log n)\)。
然后别忘了最后有个 \(+1\) 再 \(\div 2\)。
#include <cstdio>
const int mod = 998244353, N = 1e7 + 5;
int qpow(int a, int b) {
int ret = 1ll;
while (b) {
if (b & 1) ret = 1ll * ret * a % mod;
a = 1ll * a * a % mod, b >>= 1;
}
return ret;
}
int phi[N], inv[N], Prime[N / 10], cnt, f1[N], f2[N], g[N], n, m;
bool mark[10000005];
inline void add(int &x, const int y) {
if ((x += y) >= mod) x -= mod;
}
inline void mns(int &x, const int y) {
if ((x -= y) < 0) x += mod;
}
void init(int N) {
g[1] = inv[1] = phi[1] = 1;
for (int i = 2; i <= N; ++ i) {
if (!mark[i]) Prime[++ cnt] = i, phi[i] = i - 1, g[i] = qpow(i, n + 1), inv[i] = qpow(i, mod - 2);
for (int j = 1; i * Prime[j] <= N; ++ j) {
mark[i * Prime[j]] = true, g[i * Prime[j]] = 1ll * g[i] * g[Prime[j]] % mod;
inv[i * Prime[j]] = 1ll * inv[i] * inv[Prime[j]] % mod;
if (i % Prime[j] == 0) {phi[i * Prime[j]] = phi[i] * Prime[j]; break;}
phi[i * Prime[j]] = phi[i] * (Prime[j] - 1);
}
}
for (int i = 1; i <= n; ++ i) g[i] = 1ll * g[i] * inv[phi[i]] % mod;
for (int i = 1; i <= n; ++ i) f1[i] = 1ll * phi[i] * i % mod;
for (int i = 1; i <= m; ++ i) f2[i] = 1ll * phi[i] * i % mod;
for (register int i = cnt; i; -- i)
for (register int j = n / Prime[i]; j; -- j) mns(g[j * Prime[i]], g[j]);
for (register int i = 1; i <= cnt && Prime[i] <= n; ++ i)
for (register int j = n / Prime[i]; j; -- j) add(f1[j], f1[j * Prime[i]]);
for (register int i = 1; i <= cnt && Prime[i] <= m; ++ i)
for (register int j = m / Prime[i]; j; -- j) add(f2[j], f2[j * Prime[i]]);
}
int main() {
scanf("%d%d", &n, &m);
init(n > m ? n : m);
int ans = 0;
for (int i = 1; i <= n; ++ i) add(ans, 1ll * f1[i] * f2[i] % mod * g[i] % mod);
printf("%d", ((int)((1ll * ans + 1) * qpow(2, mod - 2) % mod) + mod) % mod);
return 0;
}
标签：前缀,limits,狄利克,int,题解,sum,varphi,ij 来源： https://www.cnblogs.com/stinger/p/15778463.html
本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享； 2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关； 3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关； 4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除； 5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。