Gauss's blog

Gauss's blog

I’d go back to December turn around and make it alright.

题解「湖北省队互测2014」一个人的数论

posted on 2021-04-06 15:59:19 | under 题解 |

今天心情不好,来写一道简单数论题调节一下~


题目描述

给出$n$的质因数分解和$m$,求$$\sum_{i=1}^n[i\perp n]i^m$$

题解

$$\begin{aligned}\sum_{i=1}^n[i\perp n]i^m &=\sum_{i=1}^mi^m\sum_{k|\gcd(i,n)}\mu(k)\\&=\sum_{k|n}\mu(k)\sum_{i=1}^{n/k}(ik)^m\\&=\sum_{k|n}\mu(k)k^m\sum_{i=1}^{n/k}i^m\end{aligned}$$ 看到后面那个自然数幂和,直接拿伯努利数淦它就好了.

现在令$$S_m(n)=\sum_{i=1}^ni^m=\sum_{i=1}^{m+1}a_in^i$$ 带入得$$\begin{aligned}\sum_{i=1}^{n}[i\perp n]i^m&=\sum_{k|n}\mu(k)k^m\sum_{i=1}^{m+1}a_i\left(\frac{n}{k}\right)^i\\&=\sum_{i=1}^{m+1}a_i\sum_{k|n}\mu(k)k^m\left(\frac{n}{k}\right)^i\\&=\sum_{i=1}^{m+1}a_ig_i(n)\end{aligned}$$ 显然,$\mu(n)n^m$与$n^i$均为积性函数,而积性函数的狄利克雷卷积也为积性函数,所以$g_i(n)$为积性函数.

所以有$$g_i(n)=\prod_{j=1}^{\omega(n)}g(p_j^{c_j})$$ 考虑到$\mu$只在free square number处有值,所以$$g_i(n)=\prod_{j=1}^{\omega(n)}(p_j^{ic_j}-p^{m+i(c_j-1)})$$ 时间复杂度$O(m\omega(n))$.

代码

考虑了下要不要写多项式求逆

发现模数不对(懒得写MTT了),于是没写

$m$太小了怎么搞都行.

#include <cstdio>
#include <iostream>

using namespace std;
const int M = 100 + 1;
const int mod = 1e9 + 7;
int B[M + 1], inv[M + 2], fac[M + 1], facinv[M + 1];

int pow(int a, int b)
{
    int ans = 1;
    for(; b; b >>= 1, a = 1ll * a * 1ll * a % mod)
        if(b & 1) ans = 1ll * ans * 1ll * a % mod;
    return ans;
}
int Comb(int n, int m)
{
    return 1ll * fac[n] * 1ll * facinv[m] % mod * 1ll * facinv[n - m] % mod;
}
void init()
{
    fac[0] = 1;
    for(int i = 1; i <= M; i++)
        fac[i] = 1ll * fac[i - 1] * 1ll * i % mod;
    facinv[M] = pow(fac[M], mod - 2);
    for(int i = M; i >= 1; i--)
        facinv[i - 1] = 1ll * facinv[i] * 1ll * i % mod;
    inv[1] = 1;
    for(int i = 2; i <= M + 1; i++)
        inv[i] = 1ll * inv[mod % i] * 1ll * (mod - mod / i) % mod;
    B[0] = 1;
    for(int i = 1; i <= M; i++)
    {
        for(int j = 0; j < i; j++)
            B[i] = (1ll * B[i] - 1ll * Comb(i + 1, j) * 1ll * B[j] % mod + 1ll * mod) % mod;
        B[i] = 1ll * B[i] * 1ll * inv[i + 1] % mod;
    }
}
const int W = 1e3;
int w, m, ans, a[M + 1], p[W + 1], c[W + 1];
int read()
{
    int ans = 0; char ch;
    for(; !isdigit(ch = getchar()); );
    for(; isdigit(ch); ans = (ans << 1) + (ans << 3) + (ch ^ 48), ch = getchar());
    return ans;
}
int g(int i)
{
    int ans = 1;
    for(int j = 1; j <= w; j++)
        ans = 1ll * ans * (1ll * pow(p[j], 1ll * c[j] * 1ll * i % (mod - 1)) - 1ll * pow(p[j], (1ll * m + 1ll * i * 1ll * (c[j] - 1)) % (mod - 1)) + 1ll * mod) % mod;
    return ans;
}
int main()
{
    init();
    m = read(), w = read();
    for(int i = 1; i <= w; i++)
        p[i] = read(), c[i] = read();
    for(int i = 1; i <= m + 1; i++)
    {
        for(int j = 0; j <= i; j++)
            a[j] = (1ll * a[j] + 1ll * Comb(m + 1, i) * 1ll * B[m + 1 - i] % mod * 1ll * Comb(i, j)) % mod;
    }
    for(int i = 1; i <= m + 1; i++)
        a[i] = 1ll * a[i] * 1ll * pow(m + 1, mod - 2) % mod;
    for(int i = 1; i <= m + 1; i++)
        ans = (1ll * ans + 1ll * a[i] * 1ll * g(i)) % mod;
    printf("%d\n", ans);
    return 0;
}