PE 439 【莫比乌斯反演】【杜教筛】

来源：互联网发布：手机淘宝买家秀在哪看编辑：程序博客网时间：2024/04/28 13:26

~~垃圾题，坑了我一天QAQ~~

题目大意：求∑Ni=1∑Nj=1σ(ij),N=1011

这貌似跟SDOI2015的那道题差不多，但是这里有个jy,就半天没有想出来QAQ
对于∑Ni=1σ(i)，前N23项可以线筛预处理，后面可以暴力分块求

然后是处理∑Ni=1iμ(i),然后又卡了TAT。。。
这就是一个杜教筛啊QAQ

g(1)S(n)=∑ni=1(f∗g)(i)−∑i=2g(i)S(⌊ni⌋)
这里 f(n)=iμ(i),g(n)=n
就有 S(n)=1−∑ni=2iS(⌊ni⌋)
前N23项预处理，后面记忆化搜索

复杂度O(N23)

【答案】968697378

#include<iostream>#include<algorithm>#include<ctime>#include<cmath>#include<cstring>#include<string>#include<cstdlib>#include<cstdio>#define N 21543983#define M 1000000000#define LL long longusing namespace std;LL n = (LL)1e11,tot,siz;LL s0[N],pr[N],num[N],f[N],miu[N];int prime[N];bool v[N];LL first[N],to[N],next[N],S[N];void get_prime(){    num[1] = 1;    miu[1] = 1;    f[1] = 1;    for (int i = 2;i < N;i ++)    {        if (!v[i])        {            prime[tot ++] = i;            miu[i] = -1;            num[i] = 1 + i;            s0[i] = 1 + i;            pr[i] = 1;        }        for (int j = 0;j < tot && i * prime[j] < N;j ++)        {            v[i * prime[j]] = true;            if (i % prime[j] == 0)            {                pr[i * prime[j]] = pr[i];                s0[i * prime[j]] = (s0[i] * prime[j] + 1) % M;                num[i * prime[j]] = pr[i] * s0[i * prime[j]] % M;                break;            }            else             {                miu[i * prime[j]] = -miu[i];                pr[i * prime[j]] = num[i];                s0[i * prime[j]] = prime[j] + 1;                num[i * prime[j]] = num[i] * s0[i * prime[j]] % M;            }        }        f[i] = (f[i - 1] + i * miu[i]) % M;    }    for (int i = 2;i < N;i ++)        (num[i] += num[i - 1]) %= M;}LL find(LL n){    LL x = n % N;    for (LL i = first[x];i;i = next[i])        if (to[i] == n) return S[i];    return -1;}void inser(LL n,LL w){    if (w < 0) w += M;    LL x = n % N;    next[++ siz] = first[x];    first[x] = siz;    to[siz] = n;    S[siz] = w;}LL mul2(LL a,LL b){    (~a & 1) ? a >>= 1 : b >>= 1;    return ((a % M) * (b % M)) % M;}LL cal(LL n){    if (n < N) return f[n];    LL S = find(n);    if (~ S) return S;    LL ret = 1;    for (LL i = 2,j;i <= n;i = j + 1)    {        j = min(n,n / (n / i));        (ret -= mul2(j + i,j - i + 1) * cal(n / i)) %= M;    }    inser(n,ret);    return ret;}LL sigma(LL n){    if (n < N) return num[n];    LL ret = 0;    for (LL i = 1,j;i <= n;i = j + 1)    {        j = min(n,n / (n / i));        (ret += mul2(j + i,j - i + 1) * (n / i)) %= M;    }    return ret;}LL pf(LL x){    return x * x % M;}int main(){    get_prime();    LL ans = 0;    for (LL i = 1,j;i <= n;i = j + 1)    {        j = min(n,n / (n / i));        (ans += (cal(j) - cal(i - 1)) * pf(sigma(n / i))) %= M;    }    if (ans < 0) ans += M;    cout << ans << endl;    return 0;}

0 0