数学求和公式c++

1. 等差数列求和

// 首项a1，末项an，项数n
ll arithmetic_sum(ll a1, ll an, ll n) {
    return n * (a1 + an) / 2;
}

// 首项a1，公差d，项数n
ll arithmetic_sum(ll a1, ll d, ll n) {
    return n * a1 + n * (n - 1) / 2 * d;
}

2. 等比数列求和

// 首项a，公比r，项数n
ll geometric_sum(ll a, ll r, ll n) {
    if (r == 1) return a * n;
    return a * (pow(r, n) - 1) / (r - 1);
}

// 模意义下的等比数列求和
ll geometric_sum_mod(ll a, ll r, ll n, ll mod) {
    if (r == 1) return (a % mod) * (n % mod) % mod;
    
    ll numerator = (a % mod) * (pow_mod(r, n, mod) - 1 + mod) % mod;
    ll denominator = (r - 1 + mod) % mod;
    return numerator * mod_inverse(denominator, mod) % mod;
}

3. 平方和公式

// 1² + 2² + ... + n² = n(n+1)(2n+1)/6
ll sum_of_squares(ll n) {
    return n * (n + 1) * (2 * n + 1) / 6;
}

4. 立方和公式

// 1³ + 2³ + ... + n³ = [n(n+1)/2]²
ll sum_of_cubes(ll n) {
    ll s = n * (n + 1) / 2;
    return s * s;
}

5. 调和级数近似

// H(n) = 1 + 1/2 + 1/3 + ... + 1/n ≈ ln(n) + γ
double harmonic_approx(ll n) {
    const double gamma = 0.5772156649; // 欧拉常数
    return log(n) + gamma + 1.0 / (2 * n) - 1.0 / (12 * n * n);
}

6. 组合数求和公式

// C(n,0) + C(n,1) + ... + C(n,n) = 2^n
ll sum_of_combinations(ll n) {
    return 1LL << n; // 2^n
}

// C(n,0) + C(n,2) + C(n,4) + ... = 2^(n-1)
ll sum_even_combinations(ll n) {
    return 1LL << (n - 1);
}

7. 前缀和技巧

// 一维前缀和
vector<ll> build_prefix(const vector<ll>& arr) {
    int n = arr.size();
    vector<ll> prefix(n + 1, 0);
    for (int i = 0; i < n; i++) {
        prefix[i + 1] = prefix[i] + arr[i];
    }
    return prefix;
}

// 二维前缀和
vector<vector<ll>> build_2d_prefix(const vector<vector<ll>>& mat) {
    int n = mat.size(), m = mat[0].size();
    vector<vector<ll>> prefix(n + 1, vector<ll>(m + 1, 0));
    
    for (int i = 0; i < n; i++) {
        for (int j = 0; j < m; j++) {
            prefix[i + 1][j + 1] = mat[i][j] + prefix[i][j + 1] 
                                 + prefix[i + 1][j] - prefix[i][j];
        }
    }
    return prefix;
}

8. 数论分块（整除求和）

// ∑_{i=1}^n ⌊n/i⌋ 的高效计算
ll floor_sum(ll n) {
    ll res = 0;
    for (ll l = 1, r; l <= n; l = r + 1) {
        r = n / (n / l);
        res += (n / l) * (r - l + 1);
    }
    return res;
}

// ∑_{i=1}^n ⌊n/i⌋ * i
ll floor_sum_weighted(ll n) {
    ll res = 0;
    for (ll l = 1, r; l <= n; l = r + 1) {
        r = n / (n / l);
        ll sum_i = (l + r) * (r - l + 1) / 2;
        res += (n / l) * sum_i;
    }
    return res;
}

9. 斐波那契数列求和

// F(1) + F(2) + ... + F(n) = F(n+2) - 1
ll fibonacci_sum(ll n) {
    if (n <= 0) return 0;
    
    // 矩阵快速幂计算F(n+2)
    // 实现矩阵快速幂...
    return fib(n + 2) - 1;
}

10. 欧拉函数求和

// ∑_{i=1}^n φ(i) 使用杜教筛
unordered_map<ll, ll> phi_sum_map;

ll phi_sum(ll n) {
    if (n < maxn) return pre[n]; // pre为预处理的小范围值
    if (phi_sum_map.count(n)) return phi_sum_map[n];
    
    ll res = n * (n + 1) / 2;
    for (ll l = 2, r; l <= n; l = r + 1) {
        r = n / (n / l);
        res -= (r - l + 1) * phi_sum(n / l);
    }
    return phi_sum_map[n] = res;
}

11. 莫比乌斯函数求和

// ∑_{i=1}^n μ(i) 使用杜教筛
unordered_map<ll, ll> mu_sum_map;

ll mu_sum(ll n) {
    if (n < maxn) return pre_mu[n]; // pre_mu为预处理的小范围值
    if (mu_sum_map.count(n)) return mu_sum_map[n];
    
    ll res = 1;
    for (ll l = 2, r; l <= n; l = r + 1) {
        r = n / (n / l);
        res -= (r - l + 1) * mu_sum(n / l);
    }
    return mu_sum_map[n] = res;
}

辅助函数

// 快速幂
ll pow_mod(ll a, ll b, ll mod) {
    ll res = 1;
    while (b) {
        if (b & 1) res = res * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return res;
}

// 模逆元
ll mod_inverse(ll a, ll mod) {
    return pow_mod(a, mod - 2, mod);
}

// 斐波那契数列矩阵快速幂
pair<ll, ll> fib(ll n) {
    if (n == 0) return {0, 1};
    auto p = fib(n >> 1);
    ll c = p.first * (2 * p.second - p.first);
    ll d = p.first * p.first + p.second * p.second;
    if (n & 1) return {d, c + d};
    else return {c, d};
}

12. 伯努利数求幂和

// 1^k + 2^k + ... + n^k 使用伯努利数
const int MAXK = 20;
ll bernoulli[MAXK + 1], comb[MAXK + 1][MAXK + 1];

void init_bernoulli() {
    // 预处理组合数
    for (int i = 0; i <= MAXK; i++) {
        comb[i][0] = comb[i][i] = 1;
        for (int j = 1; j < i; j++) {
            comb[i][j] = comb[i-1][j-1] + comb[i-1][j];
        }
    }
    
    // 计算伯努利数
    bernoulli[0] = 1;
    for (int m = 1; m <= MAXK; m++) {
        bernoulli[m] = 0;
        for (int k = 0; k < m; k++) {
            bernoulli[m] -= comb[m+1][k] * bernoulli[k];
        }
        bernoulli[m] /= (m + 1);
    }
}

ll power_sum(ll n, int k) {
    ll res = 0;
    for (int j = 0; j <= k; j++) {
        res += comb[k+1][j] * bernoulli[j] * pow(n+1, k+1-j);
    }
    return res / (k+1);
}

13. 拉格朗日插值求幂和

// 使用拉格朗日插值求幂和
ll lagrange_power_sum(ll n, int k, vector<ll>& pre) {
    // pre[i] = 1^k + 2^k + ... + i^k, i=0..k+2
    if (n <= k + 2) return pre[n];
    
    ll res = 0;
    for (int i = 0; i <= k + 2; i++) {
        ll numerator = pre[i];
        ll denominator = 1;
        
        for (int j = 0; j <= k + 2; j++) {
            if (i == j) continue;
            numerator = numerator * (n - j) % MOD;
            denominator = denominator * (i - j + MOD) % MOD;
        }
        res = (res + numerator * mod_inverse(denominator, MOD)) % MOD;
    }
    return res;
}

14. 斯特林数求和

// 第一类斯特林数 ∑S1(n,k)
void stirling_first(int n, vector<vector<ll>>& s1) {
    s1.resize(n+1, vector<ll>(n+1, 0));
    s1[0][0] = 1;
    for (int i = 1; i <= n; i++) {
        for (int j = 1; j <= i; j++) {
            s1[i][j] = (s1[i-1][j-1] + (i-1) * s1[i-1][j]) % MOD;
        }
    }
}

// 第二类斯特林数 ∑S2(n,k)
void stirling_second(int n, vector<vector<ll>>& s2) {
    s2.resize(n+1, vector<ll>(n+1, 0));
    s2[0][0] = 1;
    for (int i = 1; i <= n; i++) {
        for (int j = 1; j <= i; j++) {
            s2[i][j] = (s2[i-1][j-1] + j * s2[i-1][j]) % MOD;
        }
    }
}

15. 二项式系数求和恒等式

// ∑C(i,k) = C(n+1, k+1)
ll sum_comb_up_to_n(ll n, ll k) {
    if (k < 0 || k > n) return 0;
    // 需要实现组合数计算
    return comb(n + 1, k + 1);
}

// ∑C(n,2k) 和 ∑C(n,2k+1)
ll sum_even_binomial(ll n) {
    return (1LL << (n - 1));
}

ll sum_odd_binomial(ll n) {
    return (1LL << (n - 1));
}

// 范德蒙恒等式 ∑C(a,i)*C(b,n-i) = C(a+b,n)
ll vandermonde(ll a, ll b, ll n) {
    return comb(a + b, n);
}

16. 容斥原理求和

// 计算1..n中与m互质的数的和
ll sum_coprime(ll n, ll m) {
    vector<ll> factors;
    ll temp = m;
    
    // 质因数分解
    for (ll i = 2; i * i <= temp; i++) {
        if (temp % i == 0) {
            factors.push_back(i);
            while (temp % i == 0) temp /= i;
        }
    }
    if (temp > 1) factors.push_back(temp);
    
    // 容斥原理
    ll res = 0;
    int k = factors.size();
    for (int mask = 0; mask < (1 << k); mask++) {
        ll product = 1;
        int bits = 0;
        for (int i = 0; i < k; i++) {
            if (mask >> i & 1) {
                product *= factors[i];
                bits++;
            }
        }
        ll cnt = n / product;
        ll sum = product * cnt * (cnt + 1) / 2;
        if (bits % 2 == 0) res += sum;
        else res -= sum;
    }
    return res;
}

17. 生成函数相关求和

// 分拆数 p(n) - 将n写成正整数之和的方案数
void partition_numbers(int n, vector<ll>& p) {
    p.resize(n + 1, 0);
    p[0] = 1;
    for (int i = 1; i <= n; i++) {
        for (int j = 1, k = 1; i >= (3*j*j - j)/2; j++, k *= -1) {
            int idx1 = i - (3*j*j - j)/2;
            int idx2 = i - (3*j*j + j)/2;
            p[i] += k * (p[idx1] + (idx2 >= 0 ? p[idx2] : 0));
        }
    }
}

18. 多项式求和

// 多项式在等差点上的求和 ∑P(a + i*d)
ll polynomial_sum(const vector<ll>& poly, ll a, ll d, ll n) {
    // 使用有限差分法
    int deg = poly.size() - 1;
    vector<ll> diff(deg + 1);
    
    // 构建差分表
    for (int i = 0; i <= deg; i++) {
        ll x = a;
        ll sum = 0;
        for (int j = 0; j <= i; j++) {
            ll term = poly[j];
            for (int k = 0; k < j; k++) term *= (x - (a + k * d));
            sum += term;
        }
        diff[i] = sum;
    }
    
    // 使用组合数计算总和
    ll res = 0;
    for (int i = 0; i <= deg; i++) {
        res += diff[i] * comb(n, i + 1);
    }
    return res;
}

19. 三角函数求和

// ∑sin(kx) 和 ∑cos(kx)
pair<double, double> trig_sum(double x, int n) {
    if (sin(x/2) == 0) return {0, n}; // 特殊情况
    
    double sin_sum = sin(n*x/2) * sin((n+1)*x/2) / sin(x/2);
    double cos_sum = sin(n*x/2) * cos((n+1)*x/2) / sin(x/2);
    
    return {sin_sum, cos_sum};
}

20. 矩阵幂求和

// 计算 S = A + A² + A³ + ... + A^k
vector<vector<ll>> matrix_power_sum(vector<vector<ll>> A, ll k, ll mod) {
    int n = A.size();
    
    // 构建分块矩阵 [A, A; 0, I]
    vector<vector<ll>> M(2*n, vector<ll>(2*n, 0));
    for (int i = 0; i < n; i++) {
        for (int j = 0; j < n; j++) {
            M[i][j] = A[i][j];
            M[i][j+n] = A[i][j];
        }
        M[i+n][i+n] = 1; // 单位矩阵
    }
    
    // 计算 M^(k)
    auto Mk = matrix_power(M, k, mod);
    
    // 提取结果
    vector<vector<ll>> res(n, vector<ll>(n, 0));
    for (int i = 0; i < n; i++) {
        for (int j = 0; j < n; j++) {
            res[i][j] = Mk[i][j+n];
        }
    }
    return res;
}

21. 子集和问题相关

// 子集和计数 - 使用生成函数
vector<ll> subset_sum_count(const vector<int>& nums, int max_sum) {
    vector<ll> dp(max_sum + 1, 0);
    dp[0] = 1;
    
    for (int num : nums) {
        for (int s = max_sum; s >= num; s--) {
            dp[s] += dp[s - num];
        }
    }
    return dp;
}

// 子集异或和计数 - 使用线性基
ll xor_subset_count(const vector<int>& nums, int target) {
    // 构建线性基
    vector<int> basis(32, 0);
    int free_vars = 0;
    
    for (int x : nums) {
        for (int i = 31; i >= 0; i--) {
            if ((x >> i) & 1) {
                if (basis[i] == 0) {
                    basis[i] = x;
                    break;
                } else {
                    x ^= basis[i];
                }
            }
        }
        if (x == 0) free_vars++;
    }
    
    // 检查target是否能被表示
    for (int i = 31; i >= 0; i--) {
        if ((target >> i) & 1) {
            if (basis[i] == 0) return 0;
            target ^= basis[i];
        }
    }
    
    return target == 0 ? (1LL << free_vars) : 0;
}

22. 快速傅里叶变换求多项式乘积

// FFT实现（简化版）
void fft(vector<complex<double>>& a, bool invert) {
    int n = a.size();
    if (n == 1) return;
    
    vector<complex<double>> a0(n/2), a1(n/2);
    for (int i = 0; 2*i < n; i++) {
        a0[i] = a[2*i];
        a1[i] = a[2*i+1];
    }
    
    fft(a0, invert);
    fft(a1, invert);
    
    double angle = 2 * M_PI / n * (invert ? -1 : 1);
    complex<double> w(1), wn(cos(angle), sin(angle));
    
    for (int i = 0; 2*i < n; i++) {
        a[i] = a0[i] + w * a1[i];
        a[i + n/2] = a0[i] - w * a1[i];
        if (invert) {
            a[i] /= 2;
            a[i + n/2] /= 2;
        }
        w *= wn;
    }
}

vector<ll> multiply_poly(const vector<ll>& a, const vector<ll>& b) {
    vector<complex<double>> fa(a.begin(), a.end()), fb(b.begin(), b.end());
    int n = 1;
    while (n < a.size() + b.size()) n <<= 1;
    fa.resize(n); fb.resize(n);
    
    fft(fa, false);
    fft(fb, false);
    
    for (int i = 0; i < n; i++) fa[i] *= fb[i];
    fft(fa, true);
    
    vector<ll> res(n);
    for (int i = 0; i < n; i++) {
        res[i] = round(fa[i].real());
    }
    return res;
}