Description
Solution
里面随便变换一下
原式化为
后面的式子对于每个数单独考虑,原式可以化为
后面的式子可以杜教筛弄出来
考虑前面
根据容斥原理有
设
Fd(n)=∑i=1nfd(i),λ(n)=f∞(n)
这是本题的关键,解释一下这个式子
一开始i=1,那就是所有的总和,
μ(1)=1
对于i=2,应该将所有2的指数至少是d+1的数的贡献去掉,枚举它的倍数,因此
μ(2)=−1
对于i=3与i=2同理
但是会减重复,因为
i=6
的情况被减了两次,又要加回来
因此
μ(6)=1
此处 μ 刚好是容斥系数
继续
可以发现
λ
是完全积性函数
设
Λ(n)=∑i=1nλ(i)
原式
λ(id+1)
提出
i明显到不了n,修改上界
那前面很小,可以暴力枚举
考虑如何算
Λ
根据定义可以发现 ∑i|nλ(i)=[n 为完全平方数 ] <script type="math/tex" id="MathJax-Element-4669">]</script>,n的每个质因子指数均为偶数
单独考虑某一个质因子,如果是奇数它的和刚好为0,偶数才有答案,全为偶数则为1
于是也可以杜教筛
Code
本人被卡常,TLE80分
#include <cstdio>
#include <algorithm>
#include <iostream>
#include <cmath>
#include <cstdlib>
#include <cstring>
#define fo(i,a,b) for(int i=a;i<=b;i++)
#define fod(i,a,b) for(int i=a;i>=b;i--)
#define N 10000005
#define mo 1073741824
#define mo1 1073741823
#define H 30000007
#define LL long long
using namespace std;
LL n,m;
int pr[N],l,mu[N],s[N],phi[N],sp[N],sum[N],f[N];
LL h[H][2],hs[H][2];
int hash(LL k)
{
int k1=k%H;
while(h[k1][0]!=0&&h[k1][0]!=k) (k1+=1)%=H;
return k1;
}
int has(LL k)
{
int k1=k%H;
while(hs[k1][0]!=0&&hs[k1][0]!=k) (k1+=1)%=H;
return k1;
}
bool bz[N];
void prp()
{
s[1]=mu[1]=1;
f[1]=sum[1]=sp[1]=1;
phi[1]=1;
fo(i,2,10000000)
{
if(!bz[i])
{
mu[i]=-1;
pr[++l]=i;
phi[i]=i-1;
f[i]=-1;
}
for(int j=1;j<=l&&pr[j]*i<=10000000;j++)
{
bz[i*pr[j]]=1;
f[i*pr[j]]=-f[i];
if(i%pr[j]==0)
{
mu[i*pr[j]]=0;
phi[i*pr[j]]=phi[i]*pr[j];
break;
}
mu[i*pr[j]]=-mu[i];
phi[i*pr[j]]=phi[i]*phi[pr[j]];
}
s[i]=s[i-1]+mu[i];
sp[i]=(sp[i-1]+phi[i])%mo;
sum[i]=sum[i-1]+f[i];
}
}
LL get(LL k)
{
if(k<=10000000) return sp[k];
int p=hash(k);
if(h[p][0]!=0) return h[p][1];
h[p][0]=k;
h[p][1]=((k*(k+1))/2)&mo1;
LL i=2;
while(i<=k)
{
LL i1=k/(k/i);
h[p][1]=(h[p][1]-(get(k/i)*(i1-i+1)&mo1)+mo)&mo1;
i=i1+1;
}
return h[p][1];
}
LL ksm(LL k,LL n)
{
if(!n) return 1;
LL s=ksm(k,n/2);
return (n%2)?s*s*k:s*s;
}
LL fd(LL k)
{
if(k<=10000000) return sum[k];
int p=has(k);
if(hs[p][0]==k) return hs[p][1];
hs[p][0]=k,hs[p][1]=sqrt(k);
LL d=2;
while(d<=k)
{
LL d1=k/(k/d);
hs[p][1]=(hs[p][1]-(fd(k/d)*(d1-d+1)&mo1)+mo)&mo1;
d=d1+1;
}
return hs[p][1];
}
LL gt(LL n,LL d)
{
LL s=0;
int lm=(int)(pow(n,1.0/(d+1))+0.001);
fo(i,1,lm)
{
LL v=1;
if((d+1)%2) v*=f[i];
s=(s+v*mu[i]*fd(n/ksm(i,d+1)))&mo1;
}
return s;
}
int main()
{
freopen("sum.in","r",stdin);
//freopen("sum.out","w",stdout);
cin>>n>>m;
prp();
LL ans=0;
fo(ld,1,m)
{
LL p=1;
while(p<=n)
{
LL m1=n/p,p1=n/m1,s1=0,i=1;
(ans+=(LL)(gt(p1,ld)-gt(p-1,ld)+mo)%mo*(((LL)2*get(m1)-1+mo)%mo)+mo)%=mo;
p=p1+1;
}
}
printf("%lld\n",ans);
}