mxqz PN 筛
查看原帖
mxqz PN 筛
114082
Sya_Resory楼主2022/1/29 09:51

rt.

7 个 MLE,大概是 dfs 的时候爆栈了/kk

code:

#include <cstdio>
#include <map>

typedef long long ll;
const int N = 2e6,maxn = N + 5,maxm = 35;
const ll P = 1e9 + 7,inv3 = (P + 1) / 3;

ll n,cnt,ans;
ll pri[maxn]; bool isp[maxn];
ll phi[maxn],g[maxn],G[maxn];
ll h[maxn][maxm]; bool vis[maxn][maxm];
std::map < ll,ll > memG;

inline void Init() {
    phi[1] = g[1] = G[1] = 1;
    for(ll i = 2;i <= N;i ++) {
        if(!isp[i]) pri[++ cnt] = i,phi[i] = i - 1;
        for(ll t,j = 1;j <= cnt && (t = i * pri[j]) <= N;j ++) {
            isp[t] = true;
            if(!(i % pri[j])) { phi[t] = phi[i] * pri[j] % P; break; }
            phi[t] = phi[i] * phi[pri[j]];
        }
    }
    for(ll i = 2;i <= N;i ++) G[i] = (G[i - 1] + (g[i] = i * phi[i] % P)) % P;
    for(ll i = 1;i <= cnt;i ++) h[i][0] = 1,h[i][1] = 0,vis[i][0] = vis[i][1] = true;
}

inline ll fpow(ll a,ll b) {
    ll res = 1;
    for(;b;a = a * a % P,b >>= 1) if(b & 1) res = res * a % P;
    return res;
}

inline ll S1(ll n) { return n * (n + 1) / 2 % P; }
inline ll S2(ll n) { return S1(n) * (2 * n + 1) % P * inv3 % P; }

ll SumG(ll n) {
    if(n <= N) return G[n];
    if(memG.count(n)) return memG[n];
    ll res = S2(n);
    for(ll r,l = 1;l <= n;l = r + 1) {
        r = n / (n / l);
        res = (res - (S1(r) - S1(l - 1)) * SumG(n / l) % P + P) % P;
    }
    return memG[n] = res;
}

void dfs(ll nd,ll hd,ll priId) {
    if(!nd || !hd) return ;
    ans = (ans + hd * SumG(nd) % P) % P;
    for(ll i = priId,p = pri[i];i <= cnt;p = pri[++ i]) {
        if(nd < p * p) break;
        for(ll e = 2,pe = p * p;pe <= nd;pe *= p,e ++)
            dfs(nd / pe,hd * (e - 1) % P * (p - 1) % P * pe % P,i + 1);
    }
}

inline ll solve() {
    ans = 0;
    dfs(n,1,1);
    return ans;
}

int main() {
    scanf("%lld",&n); Init();
    printf("%lld\n",solve());
    return 0;
}
2022/1/29 09:51
加载中...