经过一系列有理有据的底层优化以后,此代码目前在bzoj用时
2640ms
排名第一。
这题是【bzoj2154 Crash的数字表格】的强化版。接着上一题推导。
记
S(n)=∑i=1n
那么
Ans==∑d=1min(n,m)d∑x=1min(⌊nd⌋,⌊md⌋)x2μ(x)S(⌊ndx⌋)S(⌊mdx⌋)∑p=1min(n,m)S(⌊np⌋)S(⌊mp⌋)∑x|ppxx2μ(x)
记函数
f(p)=p∑x|pxμ(x)=id∗(1×(id∗μ))
是积性函数,可以线性筛。当 i=pmq 时
f(pmq∗p)=f(pm+1)f(q)=pm+1(1−p)f(q)=pf(pmq)
#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
#define LL long long
int rd()
{
int x=0;
char c=getchar();
while (c<'0'||c>'9') c=getchar();
while (c>='0'&&c<='9')
{
x=x*10+c-'0';
c=getchar();
}
return x;
}
const int maxn=10000000,p=100000009;
int prm[10000010],sid[10000010],f[10000010],have[10000010],sf[10000010],tot;
void init()
{
f[1]=1;
for (int i=2;i<=maxn;i++)
{
if (!have[i])
{
prm[++tot]=i;
f[i]=(LL)i*(1-i+p)%p;
}
for (int j=1;j<=tot&&(LL)i*prm[j]<=maxn;j++)
{
have[i*prm[j]]=1;
if (i%prm[j]) f[i*prm[j]]=(LL)f[i]*f[prm[j]]%p;
else
{
f[i*prm[j]]=(LL)f[i]*prm[j]%p;
break;
}
}
}
for (int i=1;i<=maxn;i++)
{
sid[i]=(sid[i-1]+i)%p;
sf[i]=(sf[i-1]+f[i])%p;
}
}
int main()
{
int T,m,n,x;
LL ans;
init();
T=rd();
while (T--)
{
n=rd();
m=rd();
if (n>m) swap(m,n);
ans=0;
for (int i=1;i<=n;i=x+1)
{
x=min(n/(n/i),m/(m/i));
ans+=(LL)sid[n/i]*sid[m/i]%p*(sf[x]-sf[i-1]+p)%p;
}
printf("%lld\n",ans%p);
}
}