求和





内容

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

等比数列的和

给定整数 ，正整数 ，求 ，模 。

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

我们要求的就是 。

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

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

注意到

而

于是有

是偶数

当 是奇数时，

。对于 ，有

是偶数是奇数

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

例题 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

给定长为 的整数序列 。求

表示 和 的最小公倍数。

由于答案可能很大，输出它除以 的余数。

我们知道 。

对于每个 ，我们计算

那么答案就是

一如前两题的套路，我们有

且

现在考虑如何计算 且。有不止一种计算方法。
不难看出它等于

且

即

对每个 ，定义

  • ： 中有多少项等于 。

那么

而

int main() {
    int n; cin >> n;
    const int maxv = 1e6;
    vector<int> cnt(maxv + 1);
    for (int i = 0; i < n; i++) {
        int x; cin >> x;
        cnt[x]++;
    }
    const int mod = 998244353;
    vector<int> inv(maxv + 1); //递推法计算逆元
    inv[1] = 1;
    for (int i = 2; i <= maxv; i++)
        inv[i] = (long long) (mod - mod / i) * inv[mod % i] % mod;
    long long ans = 0;
    vector<long long> f(maxv + 1);
    for (int d = maxv; d >= 1; d--) {
        long long sum = 0;
        long long sum2 = 0;
        for (int x = d; x <= maxv; x += d) {
            sum += (long long) cnt[x] * x;
            sum2 += (long long) cnt[x] * x * x;
        }
        sum %= mod;
        f[d] = (sum * sum - sum2) % mod * inv[2] % mod; 
        for (int i = 2 * d; i <= maxv; i += d)
            f[d] -= f[i];
        f[d] %= mod;
        ans += f[d] * inv[d] % mod;
    }
    ans %= mod;
    if (ans < 0) ans += mod;
    cout << ans << '\n';
}

对子集求和

设 是一个 元集合。 有 个子集，它们构成的集合称为 的幂集， 记作 。

给定函数 ，定义函数 如下

对每个 ，求 。

不失一般性，取 。

对于每个 ，我们把 的子集进行分类。对每个 ，定义

和的差别只在前个元素

所谓前 个元素指的是 。

根据定义，我们有 而 。
对于 ，有递推式

若若

  • 表示无交并，也就是说 暗含着集合 不相交。
  • 表示把从 中去掉 所得的集合。
  • 等于 。

对于每个 和每个 ，定义

边界条件：。
递推式

若若

就是 。

为了编程的方便，我们通常取 ，用整数 表示 的子集

的二进制第位是

对应

for (int s = 0; s < 1 << n; s++)
    F[s] = f[s];

for (int i = 0; i < n; i++)
  for (int s = 0; s < 1 << n; s++)
    if (s >> i & 1)
      F[s] += F[s ^ (1 << i)];

时间：。

例题 Compatible Numbers

Codeforces 165E

称整数 和 相容，如果 . 例如 （）和 （）相容，因为 ，而 （）和 （）不相容，因为 。

给你 个正整数 ，对于每个 ，在 中找出一个与 相容的数，若不存在这样的数，输出 。

分析

注意到 。

与 相容也就是说 是 的子集。

对于每个 ，定义

若存在使得若不存在使得

和

某个若存在满足若不存在满足

用上述对子集求和的方法，可以从 算出 。

int main() {
    const int L = 22;
    vector<int> f(1 << L, -1);
    int n;
    cin >> n;
    vector<int> a(n);
    for (int i = 0; i < n; i++) {
        cin >> a[i];
        f[a[i]] = a[i];
    }
    for (int i = 0; i < L; i++) {
        for (int j = 0; j < 1 << L; j++) {
            if (j >> i & 1 && f[j] == -1)
                f[j] = f[j ^ (1 << i)];
        }
    }
    for (int i = 0; i < n; i++) {
        cout << f[((1 << L) - 1) ^ a[i]] << ' ';
    }
    cout << '\n';
}

例题 Vowels

Codeforces 383E

集合 由前 24 个小写英文字母，a 到 x，构成。
给你 个长度为 3 的单词，每个单词也由 a 到 x 构成。

对于 的每个子集 （一共 个子集），定义

  • ： 个单词中，有多少个单词至少包含一个 中的字母。
    约定 。

求 的异或和。

。

解法

对于 的每个子集 ，定义

  • ： 个单词中，有多少个单词只含有 中的字母。

那么

又定义

  • ： 个单词中，有多少个单词恰含有 中的字母。

那么

int main() {
    const int L = 24;
    vector<int> f(1 << L);
    int n; cin >> n;
    for (int i = 0; i < n; i++) {
        string s; cin >> s;
        int mask = 0;
        for (char c : s)
            mask |= 1 << (c - 'a');
        f[mask]++;
    }
    for (int i = 0; i < L; i++)
        for (int j = 0; j < 1 << L; j++)
            if (j >> i & 1)
                f[j] += f[j ^ 1 << i];
    int ans = 0;
    for (int x : f)
        ans ^= (n - x) * (n - x);
    cout << ans << '\n';
}

例题 Or Plus Max

arc100c

给你一个长为 的整数序列 。

对每个整数 ，解决下述问题：

  • 设 是整数。求 且，
    其中 表示按位或运算。


约定

对每个 ，我们把

  • 设 是整数。求 且，

称为问题 。

指标 和 的取值范围总是 到 ，以下不再写出。

解法

对每个 ，定义

  • 且

于是问题 的答案就是 .

注意到

  • 条件 也就是 且 .

更进一步，可见

  • 就是集合 中第一大和第二大的两个数之和。
int main() {
    int n; cin >> n;
    vector<pair<int,int>> f(1 << n);
    for (int i = 0; i < 1 << n; i++)
        cin >> f[i].first;

    for (int i = 0; i < n; i++)
        for (int s = 0; s < 1 << n; s++)
            if (s >> i & 1) {
                auto [v1, v2] = f[s];
                auto [v3, v4] = f[s ^ (1 << i)];
                vector<int> t = {v1, v2, v3, v4};
                sort(t.rbegin(), t.rend());
                f[s] = {t[0], t[1]};
            }

    int ans = 0;
    for (int k = 1; k < 1 << n; k++) {
        ans = max(ans, f[k].first + f[k].second);
        cout << ans << '\n';
    }
}

可见 $0 \le s < 2^n$。这样，函数 $f$ 和 $F$ 就是长为 $2^n$ 的数组。 对于每个 $s = 0, \dots, 2^{n} - 1$ 和每个 $i = 0, \dots, n$，定义 $$ G(s, i) := \set{t \bitand s = t : t 和 s只在前i个数位上不同} $$ 其中 $\bitand$ 表示按位与。 根据定义，有 $G(s, 0) = \set{s}$。对于 $i = 0, \dots, n-1$，有递推式 $$ G(s, i+1) = \begin{cases} G(s,i), & 若 s 的第 i 位是 0 \\ G(s,i) \sqcup G(s \oplus 2^i, i+1), & 若 s 的第 i 位是 1 \end{cases} $$ 其中 $\oplus$ 表示按位异或。