BZOJ3944: Sum 杜教筛

来源:互联网 发布:吉他知乎 编辑:程序博客网 时间:2024/06/18 14:17

这里写图片描述
杜教筛,用于在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 2000000using 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[n];        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[n]=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;}
0 0
原创粉丝点击