bzoj千题计划253：bzoj2154: Crash的数字表格

http://www.lydsy.com/JudgeOnline/problem.php?id=2154

#include<cstdio>
#include<algorithm>

using namespace std;

#define N 10000001

const int mod=20101009;

int cnt,p[N];
int miu[N];
bool vis[N];

int sum[N];

void pre()
{
    miu[1]=1;
    for(int i=2;i<N;++i)
    {
        if(!vis[i])
        {
            p[++cnt]=i;
            miu[i]=-1;
        }
        for(int j=1;j<=cnt;++j)
        {
            if(p[j]*i>=N) break;
            vis[p[j]*i]=true;
            if(!(i%p[j])) break;
            miu[i*p[j]]=-miu[i];
        }
    }
    for(int i=1;i<N;++i) sum[i]=(sum[i-1]+1LL*i*i*miu[i]%mod)%mod;
}

int get_S(int n,int m)
{
    return 1LL*n*(n+1)/2%mod*(1LL*m*(m+1)/2%mod)%mod;
}

int get_F(int n,int m)
{
    if(n>m) swap(n,m);
    int ans=0,j;
    for(int i=1;i<=n;i=j+1)
    {
        j=min(n/(n/i),m/(m/i));
        ans=(ans+1LL*(sum[j]-sum[i-1]+mod)%mod*get_S(n/i,m/i)%mod)%mod;
    }
    return ans;
}

int main()
{
    pre();
    int n,m;
    scanf("%d%d",&n,&m);
    if(n>m) swap(n,m);
    int ans=0,j;
    for(int i=1;i<=n;i=j+1)
    {
        j=min(n/(n/i),m/(m/i));
        ans=(ans+(1LL*(j+i)*(j-i+1)/2%mod*get_F(n/i,m/i)%mod))%mod;
    }
    printf("%d",ans);
}