【题解】洛谷 P1829 [国家集训队]Crash的数字表格 / JZPTAB
题解(公式推导)
\[\begin{aligned}
\sum_{i=1}^{n}\sum_{j=1}^{m}\operatorname{lcm}(i,j)
&=\sum_{i=1}^{n}\sum_{j=1}^{m}\frac{ij}{\gcd(i,j)}\\
&=\sum_{i=1}^{n}\sum_{j=1}^{m}\sum_{g=1}^{\min(i,j)}[\gcd(i,j)=g]\frac{ij}{g}\\
&=\sum_{i=1}^{n}\sum_{j=1}^{m}\sum_{g=1}^{\min(i,j)}[\gcd(i,j)=g]\frac{i}{g}\frac{j}{g}g\\
&=\sum_{g=1}^{\min(n,m)}g\sum_{i=1}^{n}\sum_{j=1}^{m}[\gcd(i,j)=g]\frac{i}{g}\frac{j}{g}\\
&=\sum_{g=1}^{\min(n,m)}g\sum_{i=1}^{\lfloor\frac{n}{g}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{g}\rfloor}[\gcd(i,j)=1]ij\\
&=\sum_{g=1}^{\min(n,m)}g\sum_{i=1}^{\lfloor\frac{n}{g}\rfloor}i\sum_{j=1}^{\lfloor\frac{m}{g}\rfloor}j[\gcd(i,j)=1]\\
&=\sum_{g=1}^{\min(n,m)}g\sum_{i=1}^{\lfloor\frac{n}{g}\rfloor}i\sum_{j=1}^{\lfloor\frac{m}{g}\rfloor}j\sum_{g|i}\sum_{g|j}\mu(g)\\
\end{aligned}
\]
设 \(\mathbf{F}(n,m)=\sum\limits_{i=1}^{n}i\sum\limits_{j=1}^{m}j\sum\limits_{g|i}\sum\limits_{g|j}\mu(g)\),则原式 \(=\sum\limits_{g=1}^{\min(n,m)}g\mathbf{F}(\frac{n}{g},\frac{m}{g})\)
考虑化简 \(\mathbf{F}(n,m)\):
\[\begin{aligned}
\mathbf{F}(n,m)
&=\sum_{i=1}^{n}i\sum_{j=1}^{m}j\sum_{g=1}^{\min(n, m)}\sum_{g|i}\sum_{g|j}\mu(g)\\
&=\sum_{g=1}^{\min(n, m)}\mu(g)\sum_{i=1}^{n}i\sum_{j=1}^{m}j[g|i][g|j]\\
&=\sum_{g=1}^{\min(n, m)}\mu(g)\sum_{i=1}^{\lfloor\frac{n}{g}\rfloor}id\sum_{j=1}^{\lfloor\frac{m}{g}\rfloor}jd\\
&=\sum_{g=1}^{\min(n, m)}\mu(g)g^2\sum_{i=1}^{\lfloor\frac{n}{g}\rfloor}i\sum_{j=1}^{\lfloor\frac{m}{g}\rfloor}j\\
\end{aligned}
\]
记 \(\mathbf{S}(l, r)=\frac{(l+r)(r-l+1)}{2}\) (等差数列求和公式)
那么
\[\begin{aligned}
\mathbf{F}(n,m)
&=\sum_{g=1}^{\min(n, m)}\mu(g)g^2\mathbf{S}(1,\lfloor\frac{n}{g}\rfloor)\mathbf{S}(1,\lfloor\frac{m}{g}\rfloor)\\
\end{aligned}
\]
至此,\(\mathbf{F}(n,m)\) 已经可以用数论分块在 \(\Theta(\sqrt{\min(n, m)})\) 的时间内求出。
而对原式 \(\sum\limits_{g=1}^{\min(n,m)}g\mathbf{F}(\frac{n}{g},\frac{m}{g})\) 也可用用一次数论分块(\(\Theta(\sqrt{\min(n, m)})\))
这样的总时间复杂度是 \(\Theta(\min(n,m))=\Theta(n)\) 的(认为 \(n,m\) 同阶)
代码
const int N = 1e7 + 7, P = 20101009, inv6 = 16750841;
bool IsntPrime[N];
int primes[N], pcnt;
int mu[N];
int64 sum[N];
void init(int maxx) {
sum[1] = mu[1] = 1;
for (int i = 2; i <= maxx; ++i) {
if (!IsntPrime[i]) primes[++pcnt] = i, mu[i] = -1;
for (int p = 1; p <= pcnt && i * primes[p] <= maxx; ++p) {
IsntPrime[i * primes[p]] = true;
if (i % primes[p] == 0) break;
mu[i * primes[p]] = -mu[i];
}
sum[i] = (sum[i - 1] + (int64)i * i * mu[i]) % P;
}
}
inline int64 sum_i1(int64 n) { return (n * (n + 1) >> 1) % P; }
#define query_sum(l, r) (sum[r] - sum[(l) - 1])
inline int64 solve(int n, int m) {
int64 ans = 0;
for (int l = 1, r; l <= min(n, m); l = r + 1) {
r = min(n / (n / l), m / (m / l));
ans += query_sum(l, r) * sum_i1(n / l) % P * sum_i1(m / l) % P;
}
return ans % P;
}
signed main() {
int n, m;
read >> n >> m;
init(min(n, m));
int64 ans = 0;
for (int l = 1, r; l <= min(n, m); l = r + 1) {
r = min(n / (n / l), m / (m / l));
ans += ((int64)(l + r) * (r - l + 1) >> 1) * solve(n / l, m / l) % P;
}
write << (ans % P + P) % P << '\n';
return 0;
}