求
\[\sum_{i=1}^n\sum_{j=1}^mLCM(i, j)\; mod\; 20101009\],
\(n, m\le {10}^7\) 以下
\(n \le m\)\[ans=\sum_{i=1}^n\sum_{j=1}^mLCM(i, j)\]\[=\sum_{i=1}^n\sum_{j=1}^m\frac{ij}{GCD(i, j)}\]\[=\sum_{d=1}^n\sum_{i=1}^n\sum_{j=1}^m\frac{ij[GCD(i, j) == d]}{d}\]\[=\sum_{d=1}^nd\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}ij[GCD(i, j) == 1]\] 令
\[f(d)=\sum_{i=1}^x\sum_{j=1}^yij[GCD(i,j) == d]\]\[g(d)=\sum_{d|n}f(n)\] 则
\[g(d)=\sum_{i=1}^x\sum_{j=1}^yij[d|GCD(i,j)]\]\[=d^2\sum_{i=1}^{\lfloor\frac{x}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{y}{d}\rfloor}ij[1|GCD(i,j)]\] 令
\(sum(x)=\sum_{i=1}^xi\) 则
\[g(d)=d^2\sum_{i=1}^{\lfloor\frac{x}{d}\rfloor}i*sum(\lfloor\frac{y}{d}\rfloor)\]\[=d^2*sum(\lfloor\frac{x}{d}\rfloor)*sum(\lfloor\frac{y}{d}\rfloor)\] 则
\[f(n)=\sum_{n|d}\mu(\frac{d}{n})g(d)\]\[f(1)=\sum_{d=1}^x\mu(d)g(d)\]\[=\sum_{d=1}^x\mu(d)d^2*sum(\lfloor\frac{x}{d}\rfloor)*sum(\lfloor\frac{y}{d}\rfloor)\]\[ans=\sum_{d=1}^nd\sum_{x=1}^{\lfloor\frac{n}{d}\rfloor}\mu(x)x^2*sum(\lfloor\frac{n}{xd}\rfloor)*sum(\lfloor\frac{m}{xd}\rfloor)\] 好像到这样就可以了？
void init(){ miu[1]=1; for(int i=2; i < Maxn; i++){ if(!pr[i]) pr[++ptot]=i, miu[i]=-1; for(int j=1, x; j <= ptot && (x=i*pr[j]) < Maxn; j++){ pr[x]=1; if(i%pr[j] == 0) break; miu[x]=-miu[i]; } } for(ll i=1; i < Maxn; i++) miu[i]=(miu[i-1]+miu[i]*i*i%p)%p, sum[i]=(sum[i-1]+i)%p;}ll cal(ll n, ll m){ if(n > m) swap(n, m); ll ans=0; for(ll l=1, r=0; r < n; l=r+1){ r=min(n/(n/l), m/(m/l)); ans=(ans+(miu[r]-miu[l-1])*sum[n/l]%p*sum[m/l]%p)%p; } return ans;}void solve(){ init(); n=read(), m=read(); if(n > m) swap(n, m); ll ans=0; for(ll l=1, r=0; r < n; l=r+1){ r=min(n/(n/l), m/(m/l)); ans=(ans+1ll*(r+l)*(r-l+1)/2%p*cal(n/l, m/l)%p)%p; } printf("%lld\n", (ans+p)%p);}