BZOJ3944: Sum 杜教筛

杜教筛，用于在n的2/3次方的时间复杂度内求某些积性函数的前缀和。

记目标函数为f(i)，前缀和为F(i)

令g(i)=sigma(d|i) f(d)

G(i)为g(i)前缀和

G(n)=sigma(i=1~n) g(i)

=sigma(i=1~n) sigma(d|i) f(d)

=sigma(d=1~n) f(d)*(n/d)

=sigma(i=1~n) F(n/i)

由此得 F(n)=G(n)-sigma(i=2~n) F(n/i)

若G(i)可以O(1)求，那么预处理前n^(2/3)个F(i)，用记忆化搜索的方式求F，总时间复杂度为n^(2/3)。

对于欧拉函数，g(i)=i，所以G(i)=(i+1)*i/2

对于梦比优斯函数，g(1)=1,g(n)=0(n!=1)，所以G(i)=1

于是就可以求了。

#include<cstdio>
#include<algorithm>
#include<cmath>
#include<ext/pb_ds/assoc_container.hpp>
#define gm 2000000
using namespace std;
typedef long long ll;
int T,size;
int a[10];
ll b[10];
typedef __gnu_pbds::gp_hash_table<int,ll> hat;
hat h;
int stk[gm],top;
bool np[gm];
ll phi[gm],mu[gm];
void liner()
{
phi[1]=mu[1]=1;
for(int i=2;i<=size;++i)
{
if(!np[i])
{
stk[++top]=i;
phi[i]=i-1;
mu[i]=-1;
}
for(int j=1;j<=top;++j)
{
int x=stk[j];
if(i*x>size) break;
np[i*x]=1;
if(i%x)
{
phi[i*x]=phi[i]*phi[x];
mu[i*x]=mu[i]*mu[x];
}
else
{
phi[i*x]=phi[i]*x;
mu[i*x]=0;
break;
}
}
}
for(int i=1;i<=size;++i)
{
phi[i]+=phi[i-1];
mu[i]+=mu[i-1];
}
}
struct Phi
{
ll* A;
Phi():A(phi){}
ll G(int n)
{
return (n+1ll)*n>>1;
}
};
struct Mu
{
ll *A;
Mu():A(mu){}
ll G(int n)
{
return 1ll;
}
};
template<typename T>
struct F:public T
{
ll operator() (int n)
{
if(n<=size) return T::A
;
hat::point_iterator it=h.find(n);
if(it!=h.end()) return it->second;
ll res=T::G(n);
for(int i=2,j;j!=n;i=j+1)
{
j=n/(n/i);
ll x=operator()(n/i);
res-=(j-i+1)*x;
}
return h
=res;
}
};
F<Phi> F_phi;
F<Mu> F_mu;
int main()
{
int maxn=0;
scanf("%d",&T);
for(int i=0;i<T;++i)
{
scanf("%d",a+i);
maxn=max(maxn,a[i]);
}
size=pow(maxn,2.0/3.0)+0.5;
liner();
for(int i=0;i<T;++i)
{
b[i]=F_phi(a[i]);
}
h.clear();
for(int i=0;i<T;++i)
{
printf("%lld %lld\n",b[i],F_mu(a[i]));
}
return 0;
}