题意
求
∑ni=1∑mj=1lcm(i,j)
n,m<=10000000
分析
题解见Po神的论文
或者看这篇比较详细的博客
恩全OJ倒数第二……
代码
#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<algorithm>
#define MOD 20101009
#define ll long long
using namespace std;
const int MAXN=10000000;
const int N=10000005;
int prime[1000005],tot,mu[N],n,m;
ll s[N];
bool not_prime[N];
void getmu(int n)
{
mu[1]=1;
for (int i=2;i<=n;i++)
{
if (!not_prime[i])
{
prime[++tot]=i;
mu[i]=-1;
}
for (int j=1;j<=tot&&i*prime[j]<=n;j++)
{
not_prime[i*prime[j]]=1;
if (i%prime[j]==0)
{
mu[i*prime[j]]=0;
break;
}
mu[i*prime[j]]=-mu[i];
}
}
for (int i=1;i<=n;i++)
s[i]=(s[i-1]+(ll)mu[i]*i*i%MOD)%MOD;
}
ll getsum(int n)
{
return (ll)n*(n+1)%MOD*10050505%MOD;
}
ll F(int n,int m)
{
if (n>m) swap(n,m);
ll ans=0;
for (int i=1,last;i<=n;i=last+1)
{
last=min(n/(n/i),m/(m/i));
ans=(ans+getsum(n/i)*getsum(m/i)%MOD*(s[last]-s[i-1])%MOD)%MOD;
}
return ans;
}
int main()
{
getmu(MAXN);
scanf("%d%d",&n,&m);
if (n>m) swap(n,m);
ll ans=0;
for (int i=1,last;i<=n;i=last+1)
{
last=min(n/(n/i),m/(m/i));
ans=(ans+F(n/i,m/i)*(getsum(last)-getsum(i-1))%MOD)%MOD;
}
printf("%lld",(ans+MOD)%MOD);
return 0;
}