1. 程式人生 > >[SPOJ20174]DIVCNT3 - Counting Divisors (cube):Min_25篩

[SPOJ20174]DIVCNT3 - Counting Divisors (cube):Min_25篩

分析

首先,STO ywy OTZ,ywy TQL%%%!

說一下這道題用min_25篩怎麼做。

容易發現,對於所有質數\(p\),都滿足\(f(p)=4\),於是我們就可以直接通過\([1,x]\)內的質數的個數\(h(x)\)來求出\(g(x)=\sum_{i=1}^{x}f(i) \times [i \in prime]\)了,即\(g(x)\)可以等價地表示為\(g(x)=4 \times h(x)\)。如何求\(h(x)\)是min_25篩的基本操作就不過多贅述了。而且進一步分析我們可以得出,\(f(p^q)=3 \times q+1,p \in prime\),這個結論在實現min_25篩時也很重要。

程式碼裡的g[x]其實是我上面寫的\(h(x)\)

說一句題外話,min_25篩的話,使用雜湊表實現更直觀,但是如果用id1,id2兩個陣列實現的話常數會顯著減小(雜湊表的地方我註釋掉了)。

程式碼

#include <bits/stdc++.h>
#define rin(i,a,b) for(register int i=(a);i<=(b);++i)
#define irin(i,a,b) for(register int i=(a);i>=(b);--i)
#define trav(i,a) for(register int i=head[a];i;i=e[i].nxt)
typedef long long LL;
using std::cin;
using std::cout;
using std::endl;

inline LL read(){
    LL x=0,f=1;char ch=getchar();
    while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
    while(ch>='0'&&ch<='9'){x=(x<<3)+(x<<1)+ch-'0';ch=getchar();}
    return x*f;
}

const LL MAXN=1e11+5;
const int MAXR=1e6+5;
const int HMOD=3e6-1;
LL n,num[MAXR],g[MAXR];
int cnt,tot,prm[MAXR],id1[MAXR],id2[MAXR];
bool vis[MAXR];

/*
struct Hash_Map{
    int siz,head[HMOD+5],nxt[MAXR],to[MAXR],sta[MAXR],top;
    LL val[MAXR];
    inline void insert(LL x,int y){
        ++siz;
        nxt[siz]=head[x%HMOD];
        val[siz]=x;
        to[siz]=y;
        head[x%HMOD]=siz;
        sta[++top]=x%HMOD;
    }
    inline int operator [] (LL x){
        for(register int i=head[x%HMOD];i;i=nxt[i])
            if(val[i]==x) return to[i];
    }
    void clear(){
        siz=0;
        while(top) head[sta[top--]]=0;
    }
}mp;
*/

LL solve(LL x,int y){
//  int k=mp[x];
    int k=x<=MAXR-5?id1[x]:id2[n/x];
    LL ret=(g[k]<<2)-((y-1)<<2);
    rin(i,y,cnt){
        if(1ll*prm[i]*prm[i]>x) break;
        LL temp1=prm[i],temp2=1ll*prm[i]*prm[i];
        for(register int j=1;temp2<=x;++j,temp1=temp2,temp2*=prm[i])
            ret+=(3*j+1)*solve(x/temp1,i+1)+3*j+4;
    }
    return ret;
}

void pre_process(int n){
    rin(i,2,n){
        if(!vis[i]) prm[++cnt]=i;
        rin(j,1,cnt){
            if(i*prm[j]>n) break;
            vis[i*prm[j]]=true;
            if(i%prm[j]==0) break;
        }
    }
}

int main(){
    int T=read();
    pre_process(500000);
    while(T--){
        tot=0;
//      mp.clear();
        n=read();
        LL nxti;
        for(register LL i=1;i<=n;i=nxti){
            nxti=n/(n/i)+1;
            num[++tot]=n/i;
//          mp.insert(num[tot],tot);
            if(num[tot]<=MAXR-5) id1[num[tot]]=tot;
            else id2[nxti-1]=tot;
            g[tot]=num[tot]-1;
        }
        rin(i,1,cnt){
            if(1ll*prm[i]*prm[i]>n) break;// 這個一定要加,否則會T第一個點!
            rin(j,1,tot){
                if(1ll*prm[i]*prm[i]>num[j]) break;
//              int k=mp[num[j]/prm[i]];
                int k=num[j]/prm[i]<=MAXR-5?id1[num[j]/prm[i]]:id2[n/(num[j]/prm[i])];
                g[j]-=g[k]-(i-1);
            }
        }
        printf("%lld\n",solve(n,1)+1);
    }
    return 0;
}