[bzoj2154]Crash的数字表格

http://www.lydsy.com/JudgeOnline/problem.php?id=2154
Orz w_yqts 一遍过的大佬
首先,lcm(i,j)=i*j/gcd(i,j)
显然,ans=ΣΣlcm(i,j)
为了方便，我们令n< m
于是……

PS:%为mod
PSS:sum(x)=1+2+…+x
PSSS:有多个*的地方请多加mod和long long
嗯..应该没什么问题吧……

提供个样例
10000000 10000000
4328779
代码如下

#include <bits/stdc++.h>
using namespace std;
#define N 10000005
#define MOD 20101009
int sum[N],miu[N],pr[N],pn;
bool isprime[N];
int n,m;
inline int calc(int n,int m)
{
    int res=0;
    for (int i=1,pos;i<=n;i=pos+1)
    {
        pos=min(n/(n/i),m/(m/i));
        res=(res+(long long)(miu[pos]-miu[i-1])*sum[n/i]%MOD*sum[m/i]%MOD)%MOD;
    }
    return res;
}
int main()
{
    cin>>n>>m;
    if (n>m) swap(n,m);
    isprime[1]=miu[1]=1;
    sum[1]=1;
    for (int i=2;i<=m;++i) sum[i]=sum[i-1]+i,sum[i]-=MOD*(sum[i]>=MOD);
    for (int i=2;i<=m;++i)
    {
        if (!isprime[i])
        {
            pr[++pn]=i;
            miu[i]=-1;
        }
        for (int j=1;j<=pn && (long long)(pr[j])*i<=m;++j)
        {
            isprime[i*pr[j]]=1;
            if (i%pr[j]==0)
            {
                miu[i*pr[j]]=0;
                break;
            }
            miu[i*pr[j]]=-miu[i];
        }
        miu[i]=(long long)miu[i]*i*i%MOD;
        miu[i]=(miu[i]+miu[i-1])%MOD;
        if (miu[i]<0) miu[i]+=MOD;
    }
    int ans=0;
    for (int d=1;d<=n;++d) ans=(ans+(long long)(d)*calc(n/d,m/d)%MOD)%MOD;
    cout<<(ans+MOD)%MOD<<endl;
}