求和

内容

• 等比数列求和及其扩展
• 数论分块
• 对倍数求和

等比数列的和

我们要求的就是 。

• ，如果 和 互素，那么 在模 下有乘法逆元。
• 以下我们考虑一般情况，即 和 未必互素。

当 是偶数时，我们把偶次项和奇次项分开，写成

注意到

而

于是有

当 是奇数时，

• 这跟计算 的快速幂算法有点像。

例题 Geometric Progression

abc293_e

给定整数 ，求 ，模 。

递归写法

int mod;

int power_sum(long long x, long long n) {
    if (n == 0)
        return 0;

    if (n & 1)
        return (1 + x * power_sum(x, n - 1)) % mod;

    return (1 + x) * power_sum(x * x % mod, n >> 1) % mod;
}

非递归写法

int power_sum(long long x, long long n, int mod) {
    long long ans = 0;
    long long prod = 1;
    // 答案 == ans + prod * power_sum(x, n)
    while (n) {
        if (n & 1) {
            ans = (ans + prod) % mod;
            prod = (prod * x) % mod;
        }
        prod = prod * (1 + x) % mod;
        x = x * x % mod;
        n >>= 1;
    }
    return ans;
}

习题 LCM 111

arc050_c

整数 的十进制写法是 个 ，整数 的十进制写法是 个 。求 和 的最小公倍数除以 的余数。

分析

我们知道 ，考虑 。

我们把十进制写法是 个 的整数记作 ，那么 。

注意到

设 ，又设 。

注意到

所以问题就归结为计算两个等比数列的和。

int power(int x, long long n, int mod) {
    int ans = 1;
    while (n > 0) {
        if (n & 1)
            ans = (long long) ans * x % mod;
        x = (long long) x * x % mod;
        n >>= 1;
    }
    return ans;
}

int power_sum(int x, long long n, int mod) {
   int ans = 0;
   int prod = 1;
   while (n > 0) {
        if (n & 1) {
            ans = (ans + prod) % mod;
            prod = (long long) prod * x % mod;
        }
        prod = (long long) prod * (1 + x) % mod;
        x = (long long) x * x % mod;
        n >>= 1;
   }
   return ans;
}

int main() {
    long long a, b;
    int m;
    cin >> a >> b >> m;
    long long d = __gcd(a, b);
    int ans = (long long) power_sum(power(10, d, m), a / d, m) * power_sum(10, b, m) % m;
    cout << ans << '\n';
}

扩展

对于非负整数 ，考虑多项式 。

我们要求的就是 。

当 是偶数时，我们把偶次项和奇次项分开，写成

当 是奇数时

我们可以把上面的式子看作一个递推式：

• 若 是偶数，我们由 推出 。
• 若 是奇数，我们由 推出 。

更进一步，我们也可以说

• 若 是偶数，我们由 推出 。
• 若 是奇数，我们由 推出 。

总之，利用上面的式子，我们可以在 的时间内，由一个长为 的序列推出另一个长为 的序列。

例题 数列求和

洛谷P4948

给定整数 ，求 ，模 。

递归写法

const int mod = 1e9 + 7;
const int maxk = 2005;
int C[maxk][maxk]; // 组合数
int p2[maxk]; // p2[i]：2的i次方
int k;

vector<long long> solve(int x, long long n) { //返回一个长度是 k + 1 的数组
    vector<long long> ans(k + 1);
    if (n == 0)
        return ans;
    if (n & 1) {
        auto f = solve(x, n - 1);
        for (int i = 0; i <= k; i++) {
            for (int j = 0; j <= i; j++)
                ans[i] += C[i][j] * f[j] % mod;
            ans[i] = ans[i] % mod * x;
        }
        ans[0]++;
        for (int i = 0; i <= k; i++)
            ans[i] %= mod;
        return ans;
    }
    auto f = solve((long long) x * x % mod, n / 2);
    for (int i = 0; i <= k; i++) {
        for (int j = 0; j < i; j++)
            ans[i] += f[j] * C[i][j] % mod * p2[j] % mod;
        ans[i] = (ans[i] % mod * x + (f[i] * p2[i] % mod * (1 + x) % mod)) % mod;
    }
    return ans;
}

int main() {
    long long n;
    int a;
    cin >> n >> a >> k;
    for (int i = 0; i <= k; i++) {
        C[i][0] = 1;
        for (int j = 1; j <= i; j++)
            C[i][j] = (C[i - 1][j] + C[i - 1][j - 1]) % mod;
    }

    p2[0] = 1;
    for (int i = 1; i <= k; i++)
        p2[i] = p2[i - 1] * 2 % mod;

    long long ans = solve(a, n + 1)[k];

    if (k == 0) {
        ans--;
        if (ans < 0)
            ans += mod;
    }
    cout << ans << '\n';
}

非递归写法

仿照前面的等比数列求和的非递归写法，我们维护一列系数 以及三个变量 、、，使得等式

总是成立。

一开始 ，，，， 都是 。

最后， 变成 而 就是所求的 。

当 是奇数时，我们有

当 是偶数时，我们有

const int mod = 1e9 + 7;
const int maxk = 2005;
int C[maxk][maxk]; // 组合数
int p2[maxk]; // 二的幂

int solve(int x, long long n, int k) {
    vector<long long> w(k + 1);
    w[k] = 1;
    long long s = 0;
    while (n > 0) {
        if (n & 1) {
            s += w[0];
            for (int i = 0; i <= k; i++) {
                for (int j = i + 1; j <= k; j++)
                    w[i] += w[j] * C[j][i] % mod;
                w[i] = w[i] % mod * x % mod;
            }
        }
        for (int i = 0; i <= k; i++) {
            long long t = 0;
            for (int j = i; j <= k; j++)
                t += C[j][i] * w[j] % mod;
            w[i] = (w[i] + t % mod * x) % mod * p2[i] % mod;
        }
        x = (long long) x * x % mod;
        n >>= 1;
    }
    return s % mod;
}

int main() {
    long long n;
    int a, k;
    cin >> n >> a >> k;

    // 计算组合数和2的幂

    int ans = solve(a, n + 1, k);
    if (k == 0) {
        ans--;
        if (ans < 0) ans += mod;
    }
    cout << ans << '\n';
}

例题 数学作业

对于正整数 ，定义 Concatenate() 为，把 依次连接起来得到的数。

例如 Concatenate() = 。

给定正整数 ，求 Concatenate() 。

• ，。

对于正整数 （），我们把依次连接 所得的数记作 。那么 Concatenate() 就是 。

如果 在十进制下都是 位数，那么 可表为

做指标变换，令 ，把上式写成

模板

int k_power_sum(int x, long long n, int k, int mod) {
    vector<vector<int>> C(k + 1, vector<int>(k + 1));
    for (int i = 0; i <= k; i++) {
        C[i][0] = 1;
        for (int j = 1; j <= i; j++)
            C[i][j] = (C[i - 1][j] + C[i - 1][j - 1]) % mod;
    }
    vector<int> p2(k + 1);
    p2[0] = 1;
    for (int i = 1; i <= k; i++)
        p2[i] = p2[i - 1] * 2 % mod;

    vector<long long> w(k + 1);
    w[k] = 1;
    long long s = 0;
    // 循环不变量：答案 == ans + \sum_{i=0..k} w[i] * f_{n,k}(x)
    while (n > 0) {
        if (n & 1) {
            s += w[0];
            for (int i = 0; i <= k; i++) {
                for (int j = i + 1; j <= k; j++)
                    w[i] += w[j] * C[j][i] % mod;
                w[i] = w[i] % mod * x % mod;
            }
        }
        for (int i = 0; i <= k; i++) {
            long long t = 0;
            for (int j = i; j <= k; j++)
                t += C[j][i] * w[j] % mod;
            w[i] = (w[i] + t % mod * x) % mod * p2[i] % mod;
        }
        x = (long long) x * x % mod;
        n >>= 1;
    }
    return s % mod;
}

int power(long long x, unsigned long long n, int mod) {
    long long ans = 1;
    while (n > 0) {
        if (n & 1)
            ans = ans * x % mod;
        x = x * x % mod;
        n >>= 1;
    }
    return (int) ans;
}

int main() {
    unsigned long long n;
    int m;
    cin >> n >> m;
    n++;
    long long ans = 0;
    int len = 1;
    for (unsigned long long l = 1; l < n; l *= 10) {
        unsigned long long r = min(l * 10, n);
        long long t = (r - 1) % m * k_power_sum(l * 10 % m, r - l, 0, m) % m - k_power_sum(l * 10 % m, r - l, 1, m);
        ans = ans * power(10, len * (r - l), m) % m + t;
        len++;
    }
    ans %= m;
    if (ans < 0)
        ans += m;
    cout << ans << '\n';
}

数论分块

long long ans = 0;
for (int i = 1; i <= N; ) {
    int j = N / i;
    int ni = N / j + 1;
    ans += (ni - i) * f(j);
    i = ni;
}

long long ans = 0;
for (int i = 1; i <= N; ) {
    int j = N / i;
    int ni = N / j + 1;
    ans += (i + ni - 1) * (ni - i) / 2 * f(j);
    i = ni;
}

例题 余数求和

洛谷P2261

求 .

• .

分析

注意到

于是

int main() {
    int n, k;
    cin >> n >> k;
    long long ans = (long long) n * k;
    for (int i = 1; i <= n;) {
        if (i > k)
            break;
        int j = k / i;
        int ni = min(n, k / j) + 1;
        ans -= (long long) (i + ni - 1) * (ni - i) / 2 * j;
        i = ni;
    }
    cout << ans - sum << '\n';
}

扩展

对倍数求和

例题 Couting Rhyme

codeforces 1884D

给你整数序列 。
求满足下列条件的整数对 的数量。

• 不存在整数 （）使得 整除 且 整除 。

，。
多组数据， 的总和不超过 。

对每个 定义

• 满足 的 有多少对。

对每个 定义

• 中有多少项是 的因数。

答案就是

现在考虑如何计算 。
先放松限制，对每个 定义

• 满足 且 的 有多少对。

也就是说，把条件从「 是 和 的最大公因数」放松为「 是 和 的公因数」。

我们有

对每个 定义

• 中有多少项是 的倍数。

于是有

对每个 定义

• 中有多少项等于 。

那么

而

二者都可以通过枚举倍数来计算。

int main() {
    int T; cin >> T;
    while (T--) {
        int n; cin >> n;
        vector<int> cnt(n + 1);
        vector<int> mul(n + 1);
        vector<int> div(n + 1);
        for (int i = 0; i < n; i++) {
            int x; cin >> x;
            cnt[x]++;
        }
        for (int i = 1; i <= n; i++)
            for (int j = i; j <= n; j += i)
                mul[i] += cnt[j];

        for (int i = 1; i <= n; i++)
            for (int j = i; j <= n; j += i)
                div[j] += cnt[i];

        vector<long long> f(n + 1);

        long long ans = 0;
        for (int  i = n; i >= 1; i--) {
            f[i] = (long long) mul[i] * (mul[i] - 1) / 2;
            for (int j = 2 * i; j <= n; j += i)
                f[i] -= f[j];
            if (div[i] == 0)
                ans += f[i];
        }
        cout << ans << '\n';
    }
}

我们知道

所以上述解法的时间是 。

例题 GCD 卷积

https://judge.yosupo.jp/problem/gcd_convolution

给定整数序列 和 ，计算序列 ，定义如下

一如上一题的思路，我们有

而

int main() {
    ios::sync_with_stdio(0);
    cin.tie(0);
    int n; cin >> n;
    vector<int> a(n + 1), b(n + 1);
    for (int i = 1; i <= n; i++)
        cin >> a[i];
    for (int i = 1; i <= n; i++)
        cin >> b[i];

    const int mod = 998244353;
    vector<long long> c(n + 1);
    for (int i = n; i >= 1; i--) {
        long long sa = 0, sb = 0;
        for (int j = i; j <= n; j += i) {
            sa += a[j];
            sb += b[j];
        }
        c[i] = (sa % mod) * (sb % mod);
        for (int j = 2 * i; j <= n; j += i) {
            c[i] -= c[j];
        }
        c[i] %= mod;
        if (c[i] < 0)
            c[i] += mod;
    }
    for (int i = 1; i <= n; i++)
        cout << c[i] << ' ';
    cout << '\n';
}

例题 LCMs

agc038c

给定长为 的整数序列 。求

表示 和 的最小公倍数。