杜教筛

积性函数

1、$f(n)$ 定义域为正整数域，值域为复数，称为数论函数

2、$f(n)$ 为数论函数，$f(1)=1$，对于互质的正整数有$f(p * q)=f(p) * f(q)$，称为积性函数

3、$f(n)$ 为积性函数，对于任意 $f(p * q) = f(p) * f(q)$，称为完全积性函数

常见的积性函数

1、除数函数 $σ_k(n)=\sum_{d|n}d^k$，表示 $n$ 的约数$k$次幂和

2、欧拉函数 $φ(n)=\sum_{i=1}^n[(n, i)=1]*1$，即不大于$n$ 且与$n$ 互质的整数个数。有性质：对于正整数$n>2$ 有$φ(n)$ 是偶数。

3、莫比乌斯函数$μ(n)$，在迪利克雷卷积中与恒等函数互为逆元。$μ(1)=1$, 对无平方因子$n = \prod_{i=1}^tp_i,有μ(n)=(-1)^t$，对有平方因子数$μ(n)=0$

4、元函数$e(n)=[n==1]$,迪利克雷卷积乘法单元，完全积性

5、恒等函数$I(n)=1$，完全积性

6、单位函数$id(n)=n$，完全积性

7、幂函数${id}^k(n) = n^k$，完全积性

经典公式

\[\sum_{i=1}^n[(n, i)=1]*i=\frac {n*φ(n)+[n==1]} 2 \]

\[[n=1] = \sum_{d|n}μ(d) \]

将$μ(d)$ 看作容斥系数可证明

\[n=\sum_{d|n}φ(d) \]

将$\frac i n$ 化为最简分数统计个数即可

4、若$f(n)$ 为积性函数，那么对于整数$n=\prod_{i=1}^t p_i^{k_i}$, 有$f(n)=\prod_{i=1}^tf(p_i^{k_i})$

若$f(n)$ 为完全积性函数，那么对于整数$n$ 有$f(n)=\prod_{i=1}^t f(p_i)^{k_i}$

迪利克雷卷积

规定两个数论函数$(f(x),g(x))$的迪利克雷卷积为$(f*g)(n)=\sum_{d|n}f(d)*g(\frac n d)$

满足交换律，结合律，对加法满足分配律,满足单位元函数$e(n)$,$f*e=f=e * f$，若$f(x), g(x)$为积性函数，则$f*g$也是积性函数。

迪利克雷常用技巧是对积性函数$f$，和恒等函数$I$的卷积处理，

\[g(n)=(f*I)(n)=\sum_{d|n}f(d) \]

那么有

\[g(n)=\prod_{i=1}^t\sum_{j=0}^{k_i}f(p_i^j) \]

就是每个子集的乘积就组成了所有因子

$g(n)$为完全积性函数满足,当只是积性函数时，只有$[(i, j)=1]$情况下满足

\[F(i * j) = \sum_{d|i*j}g(d) = \sum_{a|i}\sum_{b|j}[(a, b)=1]g(a)g(\frac j b) \]

莫比乌斯反演

对于$f$ 是数论函数，$g(n)=\sum_{d|n}f(d)$，已知$g(n)$ 求$f(n)$

由于$I * μ = e$，则$g * μ = f * I * μ = f * e = f$即

\[f(n)=\sum_{d|n}g(d)*μ(\frac n d) \]

\[g(n)=\sum_{n|d}f(d)⇒f(n)=\sum_{n|d}g(d)⋅μ(\frac d n) \]

$\sum_{d|n}φ(d) = id(n)$，所以$φ(d)=\sum_{d|n}μ(d)*\frac n d$

基本运用实例

整除分块

求除数函数和如求$\sum_{i=1}^nσ_1(i)$

\[\sum_{i=1}^nσ_1(i)=\sum_{i=1}^n\sum_{d|i}d=\sum_{d=1}^n⌊\frac n d⌋ = \sum_{d=1}^n\frac {⌊\frac n d⌋ + (⌊\frac n d⌋ + 1)} 2 \]

所以所有的除数函数都可以用$\lfloor \frac n d \rfloor$ 代替$i$ 来整除分块直接带入公式求和，求区间就是求两次减一减。

    for(int i = 1; i <= n; ++ i) {
        for(int j = 1; j <= i; ++ j) {
            if(i % j == 0) {
                ans += j * j * j;
            }
        }
    }
    dbg(ans);   ans = 0;
    for(int l = 1, r; l <= n; l = r + 1) {
        r = n / (n / l);
        ans += (r - l + 1) * (n / l) * (n / l + 1) * (n / l) * ((n / l) + 1) / 4;
    }

杜教筛

求欧拉函数和 $\sum_{i=1}^nφ(i) $

前面已经知道了$\sum_{d|n}φ(d)=n$

所以$\sum_{i=1}^n\sum_{d|i}φ(d) = \sum_{i=1}^ni = \frac {n*(n+1)} 2$

考虑转化这个式子，

\[\sum_{i=1}^n\sum_{d|i}φ(d)=\sum_{d=1}^n\lfloor\frac n d \rfloorφ(d)=\sum_{\frac i d =1}^n\sum_{d=1}^{\lfloor\frac n {\frac i d} \rfloor}φ(d) \]

这三个是等价的    
for(int i = x; i <= n; ++ i) 
    for(int j = x; j <= i; ++ j) 
        if(i % j == 0) 
            ans += phi[j];

for(int i = 1; i <= n; ++ i)
    for(int j = x; j <= n / i; ++ j)
        ans += phi[j];

for(int l = x, r; l <= n; l = r + 1) {
    r = n / (n / l);
    ans += (n / l) * (prephi[r] - prephi[l - 1]);
 }

对于$\sum_{i=x}^n\sum_{d|n}f(d)$，就是

\[\sum_{i=x}^n\sum_{d|i,d>=x}f(d)=\sum_{d=x}^n\lfloor\frac n d \rfloor f(d)=\sum_{i=1}^n\sum_{d=x}^{\lfloor\frac n i \rfloor}f(j) \]

所以可得

\[\sum_{i=1}^nφ(i)=\frac {n*(n+1)} 2 - \sum_{i=2}^n\sum_{d=1}^{\lfloor \frac n i \rfloor}φ(d) \]

我们令$\phi(x) = \sum_{i=1}^xφ(i)$,可得递推式

\[\phi(n)=\frac {n*(n+1)} 2-\sum_{d=2}^n\phi(\lfloor \frac n d \rfloor) \]

$dfs + $记忆化，总共求大约 $\sqrt n$ 个就可以解决了，对于前面的前缀和，可以预处理出来，后面的就用递归加记忆化解决。

预处理后复杂度差不多是$log(n) * \sqrt n -> O(n^{\frac 2 3})$

预处理的数量也在$O(n^{\frac 2 3})$ 较优,记忆化用$unorderd\_map$或者手写$hash$。

求梅滕斯函数$M(n)=\sum_{i=1}^nμ(i)$

已知$[1 = n] = \sum_{d|n}μ(d)$,那么

\[\sum_{i=1}^n\sum_{d|i}μ(d)=\sum_{d=1}^n\lfloor\frac n d \rfloorμ(d)=\sum_{i=1}^n\sum_{d=1}^{\lfloor\frac n i \rfloor}μ(d)=1 \]

所以得到了

\[M(n)=1-\sum_{d=2}^nM(\lfloor\frac n d \rfloor) \]

void init() {
    phi[1] = 1; mu[1] = 1; prephi[1] = 1;  premu[1] = 1;
    for(int i = 2; i < maxn; ++ i) {
        if(!ispr[i])    pri[cnt ++] = i, phi[i] = i - 1, mu[i] = -1;
        for(int j = 0; j < cnt && i * pri[j] < maxn; ++ j) {
            ispr[i * pri[j]] = 1;
            if(i % pri[j] == 0) {
                phi[i * pri[j]] = phi[i] * pri[j];
                break;
            }
            phi[i * pri[j]] = phi[i] * phi[pri[j]];
            mu[i * pri[j]] = -mu[i];
        }
        prephi[i] = prephi[i - 1] +  phi[i];
        premu[i] = premu[i - 1] + mu[i];
    }
}
ll eular(int n) {
    if(n < maxn)    return prephi[n];
    ll x = eu.find(n);
    if(~x) return x;
    ll ans = n * (n + 1ll) / 2;
    for(ll l = 2, r; l <= n; l = r + 1) {
        int tmp = n / l;
        r = n / tmp;
        ans -= (r - l + 1) * eular(tmp);
    };
    eu.insert(n, ans);
    return ans;
}

ll mub(int n) {
    if(n < maxn)    return premu[n];
    ll x = mb.find(n);
    if(~x) return x;
    ll ans = 1;
    for(ll l = 2, r; l <= n; l = r + 1) {
        int tmp = n / l;
        r = n / tmp;
        ans -= (r - l + 1) * mub(tmp);
    }
    mb.insert(n, ans);
    return ans;
}

再来一个$A(n)=\sum_{i=1}^n\frac i {(n, i)}, F(n)=\sum_{i=1}^nA(i)$

\[A(n)=\sum_{i=1}^n\frac i {(n, i)} = \\\sum_{i=1}^n\sum_{d|n}[(n, i)=d]\frac i d=\\ \sum_{d|n}\sum_{i=1}^n[(n, i)=d]\frac i d=\\ \sum_{d|n}\sum_{\cfrac i d = 1}^{\cfrac n d}[(\frac n d, \frac i d) = 1]\frac i d = \\ \sum_{d|n}\cfrac {d * φ(d) + [1 = d]} 2 = \\ \frac 1 2 * (1 + \sum_{d|n}d*φ(d)) \]

$f(n)=\sum_{i=1}^nA(i)$，令 $g(n)=2*F(n)-n$,则

\[g(n)=\sum_{i=1}^n\sum_{d|n}d*φ(d) =\\ \sum_{i=1}^n\sum_{d=1}^{\lfloor\frac n i\rfloor}d *φ(d) \]

令$\phi(x) = \sum_{i=1}^xi *φ(i) $，现在就需要求出这个

考虑杜教筛的形式$\sum_{i=1}^n\sum_{d|i}d*φ(d)$,尝试展开, $n = \prod_{i=1}^tp_i^{k_i}$

$\sum_{d|n}d*φ(d) = \prod_{i=1}^t\sum_{j=0}^{k_i}p_i^jφ(p_i^j)$

$φ(p_i^j)=p_i^{j-1}*(p_i-1)$ 得到

\[\prod_{i=1}^t\cfrac {p_i^{2*k_i+1}+1} {p_i+1} \]

这个不好求,发现对原式乘上$\frac i d$，得到

\[\sum_{\frac i d=1}^n\sum_{d=1}^{\lfloor \frac n {\frac i d} \rfloor}d * φ(d) * \frac i d = \sum_{i=1}^ni^2=\frac {n*(n+1)*(2n+1)} 6 \]

所以$\phi(n)=\frac {n*(n+1)*(2n+1)} 6 - \sum_{i=2}^ni * \phi(\lfloor \frac n i\rfloor)$

可以套两次求出来即可

 prephi[i] = qmod(prephi[i - 1] + (ll)i * phi[i] % mod);
int getf(int n) {
    if(n < maxn)    return prephi[n];
    int x = hs.find(n);
    if(~x)  return x;
    int res = n * (n + 1) % mod * (2ll * n + 1) % mod * inv6 % mod;
    for(int l = 2, r; l <= n; l = r + 1) {
        r = n / (n / l);
        res -= (r - l + 1) * (l + r) % mod * inv2 % mod * getf(n / l) % mod;
        if(res < 0) res += mod;
    }
    hs.insert(n, res);
    return res;
}

void run() {
    int n; scanf("%d", &n);
    int ans = 0;
    for(int l = 1, r; l <= n; l = r + 1) {
        r = n / (n / l);
        ans = qmod(ans + (r - l + 1ll) * getf(n / l) % mod);
    }
    //f(n) = (g(n) + n) / 2
    printf("%d\n", ans + n >> 1);
    test(n);
    return ;
}

常用公式：

$i^3 = n^2 * (n+1)^2 / 4$

$i^4 = n*(1+n)*(1+2*n)*(-1+3 * n+3 * n^2) / 30$

$i^5=(n) = n^2*(n+1)^2*(2*n^2+2*n-1)/12$

$i^6=a(n) = n*(n+1)*(2*n+1)*(3*n^4+6*n^3-3*n+1)/42$

1、$M(n)=\sum_{i=1}^nμ(i)=1-\sum_{i=2}^nM(\frac n i)$

2、$M(n)=\sum_{i=1}^nμ(i) * i=1-\sum_{i=2}^ni * M(\frac n i)$

3、$M(n)=\sum_{i=1}^nμ(i)*i^2=1-\sum_{i=2}^ni^2*M(\frac n i)$

4、$F(n)=\sum_{i=1}^nφ(i)=\frac {n*(n+1)} 2-\sum_{i=2}^nF(\frac n i)$

5、$F(n)=\sum_{i=1}^nφ(i) * i=\frac {n*(n+1)*(2n + 1)} 6-\sum_{i=2}^ni * F(\frac n i)$

6、$F(n)=\sum_{i=1}^nφ(i)*i^2=\frac {n^2*(n+1)^2} 4-\sum_{i=2}^ni^2*F(\frac n i)$

差不多的规律.

下面来说说公式该怎么推，变成前缀和的形式

\[F(n) = \sum_{i=1}^n\sum_{d|i}i^2d\\ =\sum_{d=1}^n\sum_{d|i}i^2d\\ =\sum_{d=1}^n\sum_{j=1}^{\lfloor \frac n d\rfloor}(d * j)^2d\\ =\sum_{d=1}^nd^3\sum_{j=1}^{\lfloor \frac n d\rfloor}j^2 \]

难点的,约数是完全积性函数，前面有公式的

\[F(n)=\sum_{i=1}^n\sum_{j=1}^nσ_1(i * j)\\ =\sum_{i=1}^n\sum_{j=1}^n\sum_{a|i}\sum_{b|j}[(a, b)=1]\frac {a * j} b\\ =\sum_{a=1}^n\sum_{i=1}^{\lfloor \frac n a\rfloor}\sum_{b=1}^n\sum_{j=1}^{\lfloor \frac n b \rfloor}[(a, b)=1]\frac {a * (b * j)} b= aj\\ =\sum_{a=1}^na\sum_{b=1}^n[(a, b)=1]\sum_{i=1}^{\lfloor \frac n a\rfloor}\sum_{j=1}^{\lfloor \frac n b \rfloor}j\\ =\sum_{a=1}^na\lfloor\frac n a\rfloor\sum_{b=1}^n[(a, b)=1]\cfrac {\lfloor \frac n b \rfloor + (\lfloor \frac n b \rfloor + 1)} 2 \]

剩下的化简$gcd$，有公式$[(a, b)=1] = \sum_{d=1}^{min(a, b)}μ(d)[d|gcd(a, b)]$，莫比乌斯反演的容斥，带进去继续推 $d|gcd(a, b)=[d|a,d|b]$

\[\sum_{a=1}^na\lfloor\frac n a\rfloor\sum_{b=1}^n[(a, b)=1]\cfrac {\lfloor \frac n b \rfloor + (\lfloor \frac n b \rfloor + 1)} 2\\ =\sum_{a=1}^na\lfloor\frac n a\rfloor\sum_{b=1}^n\cfrac {\lfloor \frac n b \rfloor + (\lfloor \frac n b \rfloor + 1)} 2\sum_{d=1}^nμ(d)[d|gcd(a, b)]\\ =\sum_{d=1}^n\sum_{d|a}\sum_{d|b}a*\lfloor\frac n a\rfloor*\cfrac {\lfloor \frac n b \rfloor + (\lfloor \frac n b \rfloor + 1)} 2*μ(d)\\ =\sum_{d=1}^n\sum_{p=1}^{\lfloor \frac n d\rfloor}\sum_{q=1}^{\lfloor \frac n d\rfloor}(pd)*\lfloor\frac n {pd}\rfloor*\cfrac {\lfloor \frac n {dq} \rfloor + (\lfloor \frac n {dq} \rfloor + 1)} 2*μ(d)\\ =\sum_{d=1}^ndμ(d)\sum_{p=1}^{\lfloor \frac n d\rfloor}p*\lfloor\frac n {pd}\rfloor\sum_{q=1}^{\lfloor \frac n d\rfloor}\cfrac {\lfloor \frac n {dq} \rfloor + (\lfloor \frac n {dq} \rfloor + 1)} 2 \]

现在几个公式都已经独立了

\[\sum_{q=1}^{\lfloor \frac n d\rfloor}\cfrac {\lfloor \frac n {dq} \rfloor + (\lfloor \frac n {dq} \rfloor + 1)} 2 = \sum_{i=1}^n\cfrac {\lfloor \frac n {i} \rfloor + (\lfloor \frac n {i} \rfloor + 1)} 2\\ \sum_{p=1}^{\lfloor \frac n d\rfloor}p*\lfloor\frac n {pd}\rfloor=\sum_{i=1}^ni*\lfloor\frac n i\rfloor=\sum_{i=1}^n\cfrac {\lfloor \frac n {i} \rfloor + (\lfloor \frac n {i} \rfloor + 1)} 2 \]

显然这两个都是约数和的公式

\[F(n)=\sum_{d=1}^ndμ(d)\sum_{p=1}^{\lfloor \frac n d\rfloor}p*\lfloor\frac n {pd}\rfloor\sum_{q=1}^{\lfloor \frac n d\rfloor}\cfrac {\lfloor \frac n {dq} \rfloor + (\lfloor \frac n {dq} \rfloor + 1)} 2\\ =\sum_{i=1}^niμ(i)(σ_1(\lfloor\frac n i\rfloor))^2 \]

这题是题目-约数之和 (51nod.com)

代码实现：

const int maxn = 1e6 + 10;
const int hsh = 774547;
const int M = 6e5 + 10;
struct Hsh{
    int head[hsh], nxt[M], to[M], ecnt;
    int val[M];
    inline int find(int x) {
        int id = x % hsh;
        for(int i = head[id]; i; i = nxt[i]) if(to[i] == x) return val[i];
        return -1;
    }
    inline void insert(int x, int v) {
        int id = x % hsh;
        to[++ecnt] = x; nxt[ecnt] = head[id]; head[id] = ecnt; val[ecnt] = v;
    }
    inline void clear() {
        memset(head, 0, sizeof head);
        ecnt = 0;
    }
}hs1, hs2;
int pri[maxn], cnt, phi[maxn], premu[maxn], ssum[maxn], dsum[maxn], tsum[maxn];
short mu[maxn];
bool ispr[maxn];
const int mod = 1e9 + 7;
inline int qmod(int x) {
    return x >= mod ? x - mod : x;
}
void init() {//premu[i] = premu[i - 1] + i * mu[i], ssum[i] = ssum[i - 1] + dsum[i]
    phi[1] = 1;  mu[1] = 1; premu[1] = 1; dsum[1] = ssum[1] = 1;
    for(int i = 2; i < maxn; ++ i) {
        if(!ispr[i])    pri[cnt ++] = i, phi[i] = i - 1, mu[i] = -1, dsum[i] = tsum[i] = i + 1;
        for(int j = 0; j < cnt && i * pri[j] < maxn; ++ j) {
            ispr[i * pri[j]] = 1;
            if(i % pri[j] == 0) {
                phi[i * pri[j]] = phi[i] * pri[j];
                dsum[i * pri[j]] = dsum[i] / tsum[i] * 1ll * (tsum[i] * pri[j] + 1) % mod;
                tsum[i * pri[j]] = qmod((ll)tsum[i] * pri[j] % mod + 1);
                break;
            }
            phi[i * pri[j]] = phi[i] * phi[pri[j]];
            mu[i * pri[j]] = -mu[i];
            tsum[i * pri[j]] = tsum[pri[j]];
            dsum[i * pri[j]] = (ll)dsum[i] * tsum[pri[j]] % mod;
        }
        premu[i] = qmod(premu[i - 1] + (ll)i * mu[i] % mod);
        if(premu[i] < 0)    premu[i] += mod;
        ssum[i] = qmod(ssum[i - 1] + dsum[i]);
    }
}

inline int gcd(int a, int b) {
    return !b ? a : gcd(b, a % b);
}
inline int ksm(int a, int b) {
    int res = 1;
    while(b) {
        if(b & 1)   res = 1ll * res * a % mod;
        a = 1ll * a * a % mod;
        b >>= 1;
    }
    return res;
}
const int inv6 = ksm(6, mod - 2);
const int inv2 = ksm(2, mod - 2);

int getd(int n) {
    if(n < maxn)    return ssum[n];
    int x = hs1.find(n);
    if(~x)  return x;
    int res = 0;
    for(int l = 1, r; l <= n; l = r + 1) {
        int sub = n / l;
        r = n / sub;
        res = qmod(res + (r - l + 1ll) * sub % mod * (sub + 1) % mod * inv2 % mod);
    }
    hs1.insert(n, res);
    return res;
}
//3331232

int getmb(int n) {
    if(n < maxn)    return premu[n];
    int x = hs2.find(n);
    if(~x)  return x;
    int res = 1;
    for(int l = 2, r; l <= n; l = r + 1) {
        int sub = n / l;
        r = n / sub;
        res = qmod(res - (r - l + 1ll) * (r + l) % mod * inv2 % mod * getmb(sub) % mod + mod);
        while(res >= mod)  res -= mod;
    }
    hs2.insert(n, res);
    return res;
}


void test(int n) {
    //f(n) = (1 to n)i*mu(i)*d(n/i)^2
    int ans = 0;
    for(int i = 1; i <= n; ++ i)    {
        ans = qmod(ans + (ll)i * mu[i] * getd(n / i) % mod * getd(n / i) % mod);
        if(ans < 0) ans += mod;
    }
    dbg(ans);
}

void run() {
    int n; scanf("%d", &n);
    int ans = 0;
    //f(n) = (1 to n)i*mu(i)*d(n/i)^2
    for(int l = 1, r; l <= n; l = r + 1) {
        int sub = n / l;
        r = n / sub;
        int dd = getd(sub);
        ans = qmod(ans + qmod(getmb(r) - getmb(l - 1) + mod) * 1ll * dd % mod * dd % mod);
    }
    printf("%d\n", ans);
    //test(n);
    return ;
}

signed mian() {
    init(); int t = 1;
    while(t--) run();
    return 0;
}

posted @ 2021-08-08 17:18 wlhp 阅读(70) 评论(0) 收藏举报

刷新页面返回顶部

wlhp

杜教筛