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[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;
}