对于非质数p的分解来求组合(完全不明白)

#include <iostream>
#include <vector>
using namespace std;

// 预处理组合数(适用于较小的数)
const int mx = 5;
int c[mx][mx];
auto init = []() {
    for (int i = 0; i < mx; i++) {
        c[i][0] = c[i][i] = 1;
        for (int j = 1; j < i; j++) {
            c[i][j] = c[i - 1][j] + c[i - 1][j - 1];
        }
    }
    return 0;
}();

// 通用的 Lucas 定理,用于计算 C(m, n) % p
long long lucas(long long m, long long n, long long p) {
    if (n == 0) return 1;
    return c[m % p][n % p] * lucas(m / p, n / p, p) % p;
}

// 快速幂算法:计算 (base^exp) % MOD
long long power(long long base, long long exp, long long mod) {
    long long result = 1;
    while (exp > 0) {
        if (exp % 2 == 1) result = (result * base) % mod;
        base = (base * base) % mod;
        exp /= 2;
    }
    return result;
}

// 扩展欧几里得算法,用于求 ax + by = gcd(a, b) 的解
long long extended_gcd(long long a, long long b, long long &x, long long &y) {
    if (b == 0) {
        x = 1;
        y = 0;
        return a;
    }
    long long g = extended_gcd(b, a % b, y, x);
    y -= a / b * x;
    return g;
}

// 求 a 在模 m 下的逆元
long long mod_inverse(long long a, long long m) {
    long long x, y;
    long long g = extended_gcd(a, m, x, y);
    if (g != 1) return -1; // 无逆元
    return (x % m + m) % m;
}

// 分解质因数及其幂次
vector<pair<long long, long long>> factorize(long long p) {
    vector<pair<long long, long long>> factors;
    for (long long i = 2; i * i <= p; i++) {
        if (p % i == 0) {
            long long count = 0;
            while (p % i == 0) {
                p /= i;
                count++;
            }
            factors.emplace_back(i, count);
        }
    }
    if (p > 1) {
        factors.emplace_back(p, 1);
    }
    return factors;
}

// 中国剩余定理合并结果
long long chinese_remainder_theorem(const vector<long long>& remainders, const vector<long long>& moduli) {
    long long M = 1;
    for (long long mod : moduli) {
        M *= mod;
    }
    long long result = 0;
    for (int i = 0; i < remainders.size(); i++) {
        long long Mi = M / moduli[i];
        long long inv_Mi = mod_inverse(Mi, moduli[i]);
        result = (result + remainders[i] * Mi * inv_Mi) % M;
    }
    return result;
}

// 计算组合数 C(n, m) % p,p 为非质数
long long cmb(long long m, long long n, long long p) {
    auto factors = factorize(p);
    vector<long long> remainders;
    vector<long long> moduli;
    for (const auto& [prime, exponent] : factors) {
        long long mod = 1;
        for (int i = 0; i < exponent; i++) {
            mod *= prime;
        }
        moduli.push_back(mod);
        remainders.push_back(lucas(m, n, mod));
    }
    return chinese_remainder_theorem(remainders, moduli);
}

int main() {
    // 示例:计算 C(100000, 50000) % 12(12 为非质数)
    long long n = 100000;
    long long m = 50000;
    long long p = 12;
    long long result = cmb(m, n, p);
    cout << "C(" << m << ", " << n << ") % " << p << " = " << result << endl;
    return 0;
}
posted @ 2025-02-25 22:01  Qacter  阅读(12)  评论(0)    收藏  举报