【2011集训队出题】Crash的数字表格

Description

求

$$\sum_{i=1}^n\sum_{j=1}^mlcm(i,j)\mod 20101009$$

Solution

对于一个式子的做法思路，有思路以后应该注意的东西

对于一个式子必须化简

$$\sum_{i=1}^n\sum_{j=1}^mi*j{1\over gcd(i,j)}\mod 20101009$$

对于有gcd的题目，常见思路

我们设n< m
我们可以枚举最小公倍数d
我们设一个

$$f(d)=\sum_{i=1}^{n}\sum_{j=1}^{m}i*j(gcd(i,j)=d)$$
所以答案是
$$ans=\sum_{d=1}^{n}{1\over d}*f(d)$$
那么f(d)的求解自然就可以想到莫比乌斯反演
$$g(d)=\sum_{d=1}^{\lfloor {n\over d}\rfloor}f(d)$$
莫比乌斯反演后
$$f(d)=\sum_{i=1}^{\lfloor {n\over d}\rfloor}g(id)*\mu(i)$$
g(d)其实可以化简
那么g(d)就是等于
$$(d+2d+3d+......+\lfloor {n\over d}\rfloor)*(d+2d+3d+......+\lfloor {m\over d}\rfloor)$$
因为这些都是d的倍数，然后在矩阵中再两两相乘。
提取供应式d，然后在用等差数列求和，得到
$$d^2\lfloor {n\over d}\rfloor(\lfloor {n\over d}\rfloor+1)*\lfloor {m\over d}\rfloor(\lfloor {m\over d}\rfloor+1)$$
把这个带入答案求解里面
$$ans=\sum_{d=1}^{n}{1\over d}\sum_{i=1}^{\lfloor {n\over d}\rfloor}d^2i^2\lfloor {n\over di}\rfloor(\lfloor {n\over di}\rfloor+1)*\lfloor {m\over di}\rfloor(\lfloor {m\over di}\rfloor+1)*\mu(i)$$
然后发现可以消掉一个d
$$ans=\sum_{d=1}^{n}\sum_{i=1}^{\lfloor {n\over d}\rfloor}i^2d\lfloor {n\over di}\rfloor(\lfloor {n\over di}\rfloor+1)*\lfloor {m\over di}\rfloor(\lfloor {m\over di}\rfloor+1)*\mu(i)$$

分块大法

有整除的式子题，就用分块大法

预处理

用一个前缀和预处理$mu(i)*i^2$

注意会爆空间

坑死人

Code

#include<iostream>
#include<algorithm>
#include<cmath>
#include<cstring>
#include<cstdio>
#define fo(i,a,b) for(i=a;i<=b;i++)
#define ll long long
const int mo=20101009;
const int maxn=10000005;
using namespace std;
ll i,j,k,l,t,n,m,r,l1,r1,o,q;
ll sum,sum1,ni,ans2,ans;
int qian[maxn];
int miu[maxn],zhi[maxn];
bool bz[maxn];
ll chu2(ll x){
    if(x%2==0)return x/2*(x+1)%mo;
    else return (x+1)/2*x%mo;
}
int main(){
    scanf("%d%d",&n,&m);
    if(n>m)swap(n,m);
    miu[1]=1;qian[1]=1;
    fo(i,2,n){
        if(!bz[i]){
            zhi[++zhi[0]]=i;
            miu[i]=-1;
        }
        fo(j,1,zhi[0]){
            t=zhi[j]*i;
            if(t>n)break;
            bz[t]=1;
            if(i%zhi[j]==0)break;
            miu[t]=-miu[i];
        }
        qian[i]=(qian[i-1]+miu[i]*i%mo*i%mo)%mo;
    }
    l1=1;
    while(l1<=n){
        r1=min(n / (n / l1),m / (m / l1));
        t=n / l1;
        k=m / l1;
        l=1;
        ans2=0;
        while(l<=t){
            r=min(t/(t/l),k/(k/l));o=t/l,q=k/l;
            ans2=(ans2+(qian[r]-qian[l-1]+mo)*chu2(o)%mo*chu2(q))%mo;    
            l=r+1;
        }
        ans=(ans+(chu2(r1)-chu2(l1-1)+mo)%mo*ans2%mo)%mo;
        l1=r1+1;
    }
    printf("%lld\n",ans);
}