PE 521【DP】

题目大意：设ms(i)代表i的最小质因数，求∑ni=2ms(i),n=1012mod109

对于pi>n−√,pi只会出现一次，所以可以先处理<=n−√的素数
设 F[i,j] 表示不含p1..i的因数，范围在[1,j]的数的个数
设 G[i,j] 表示不含p1..i的因数，范围在[1,j]的数的和
对于p2i>j，有

F[i,j]=F[i−1,j]−1G[i,j]=G[i−1,j]−pi

对于p2i<j，令k=j/pi有

F[i,j]=F[i−1,j]−F[i−1,k]G[i,j]=G[i−1,j]−G[i−1,k]∗pi

因为要用到的只有f[i,n],g[tot,n]，所以第二维可以简化成只有2n−√个

素数分布比较均匀，积分算下来，复杂度O(n34lnn)

开始想多了是我sb。。。
线筛写错也是醉了→_→

【答案】44389811

#include<iostream>#include<algorithm>#include<cstdlib>#include<cstdio>#include<cstring>#include<string>#include<ctime>#include<cmath>#define N 2000005#define M 1000000000#define LL long longusing namespace std;int tot,m;LL n = (LL)1e12,S = (LL)1e6,ans;LL f[N],g[N],f_[N],g_[N],a[N],sum[N],prime[N];bool v[N];void get_prime(){    for (int i = 2;i <= S;i ++)    {        if (!v[i]) prime[++ tot] = i;        for (int j = 1;j <= tot && i * prime[j] <= S;j ++)        {            v[i * prime[j]] = true;            if (i % prime[j] == 0) break;        }    }}LL mul2(LL x){    LL y = x + 1;    (x & 1) ? y >>= 1 : x >>= 1;    x %= M,y %= M;    return x * y % M;}int main(){    get_prime();    for (int i = 1;i <= S;i ++) a[++ m] = i;    for (int i = S;i;i --) if (n / i != a[m]) a[++ m] = n / i;    for (int i = 1;i <= m;i ++)    {        g[i] = mul2(a[i]);        f[i] = a[i];    }    for (int i = tot;i;i --)        sum[i] = (sum[i + 1] + prime[i]) % M;    for (int i = 1,last = 1,j;i <= tot;i ++)    {        for (j = last;a[j] < prime[i] * prime[i];j ++)        {            (f[j] -= tot - i + 1) %= M;            (g[j] -= sum[i]) %= M;        }        last = j;        for (j = last;j <= m;j ++)        {            f_[j] = f[j];            g_[j] = g[j];        }        for (j = last;j <= m;j ++)        {            LL k = a[j] / prime[i],t = f[j];            if (k > S) k = m + 1 - n / k;            if (k < last)            {                f[j] -= f[k] + tot - i + 1;                g[j] -= (g[k] + sum[i]) * prime[i];            }            else            {                f[j] -= f_[k];                g[j] -= g_[k] * prime[i];            }            f[j] %= M,g[j] %= M;            if (j == m) (ans += (t - f[j]) % M * prime[i]) %= M;        }    }    (ans += g[m] - 1) %= M;    if (ans < 0) ans += M;    cout << ans << endl;    return 0;}