C++中大型二项式系数的稳定计算方案问询

针对你的需求——支持n=100（最大系数约1e29），同时兼顾精度或精确性，以下是几种可行的实现方案：

1. 高精度整数实现（完全精确） #

n=100时的最大组合数C(100,50)有30位，远超64位整数的存储上限，要得到完全精确的结果，需要自己实现高精度整数运算。利用组合数计算时每一步乘除均为整除的特性，我们可以用数组存储数字的每一位，逐步完成乘除操作：

#include <vector>
#include <algorithm>
#include <string>

using namespace std;

// 高精度整数（低位存于数组前端）乘小整数
vector<int> multiply(const vector<int>& num, int factor) {
    vector<int> res;
    int carry = 0;
    for (int digit : num) {
        int product = digit * factor + carry;
        res.push_back(product % 10);
        carry = product / 10;
    }
    while (carry > 0) {
        res.push_back(carry % 10);
        carry /= 10;
    }
    return res;
}

// 高精度整数（低位存于数组前端）除以小整数（保证整除）
vector<int> divide(const vector<int>& num, int divisor) {
    vector<int> res;
    int remainder = 0;
    // 从高位到低位处理
    for (auto it = num.rbegin(); it != num.rend(); ++it) {
        int current = remainder * 10 + *it;
        res.push_back(current / divisor);
        remainder = current % divisor;
    }
    // 移除前导零
    while (!res.empty() && res.front() == 0) {
        res.erase(res.begin());
    }
    // 转回低位在前的格式
    reverse(res.begin(), res.end());
    return res.empty() ? vector<int>{0} : res;
}

// 计算C(n,r)，返回精确结果的字符串表示
string comb_precise(long long n, long long r) {
    if (r < 0 || r > n) return "0";
    r = min(r, n - r); // 利用对称性减少计算量
    vector<int> result = {1};
    for (int i = 0; i < r; ++i) {
        result = multiply(result, n - i);
        result = divide(result, i + 1);
    }
    // 转换为字符串输出
    string s;
    for (auto it = result.rbegin(); it != result.rend(); ++it) {
        s += to_string(*it);
    }
    return s;
}

优点：结果完全精确，无精度损失；缺点：需要手动实现高精度运算逻辑，代码量稍大。

2. long double近似实现（保留15+有效数字） #

多数平台的long double是80位扩展精度类型，拥有64位尾数，可提供约18-19位有效数字，完全满足你保留15位有效数字的要求。可以直接沿用你原有的逐步乘除逻辑，将类型替换为long double：

#include <cmath>
#include <algorithm>

long double comb_approx(long long n, long long r) {
    if (r < 0 || r > n) return 0.0L;
    r = min(r, n - r);
    long double result = 1.0L;
    for (int i = 0; i < r; ++i) {
        result = result * (n - i) / (i + 1);
    }
    return result;
}

如果担心大数直接乘除的累积误差，也可以用对数求和的方式计算，稳定性更强：

#include <cmath>
#include <algorithm>

long double comb_approx_log(long long n, long long r) {
    if (r < 0 || r > n) return 0.0L;
    r = min(r, n - r);
    long double log_sum = 0.0L;
    for (int i = 1; i <= r; ++i) {
        log_sum += logl(n - r + i) - logl(i);
    }
    return expl(log_sum);
}

优点：代码极简，无需额外逻辑；缺点：结果是近似值，但精度完全符合需求。

3. __int128实现（精确且代码简洁） #

GCC、Clang等编译器支持__int128类型，这是128位整数，最大可存储约1e38的数值，完全能容纳n=100时的所有组合数。计算逻辑和你的原代码几乎一致，只需将类型替换为__int128，再手动实现字符串转换（标准库不支持__int128的流输出）：

#include <algorithm>
#include <string>

using namespace std;

// 计算C(n,r)，返回精确结果的字符串表示
string comb_int128(long long n, long long r) {
    if (r < 0 || r > n) return "0";
    r = min(r, n - r);
    __int128 result = 1;
    for (int i = 0; i < r; ++i) {
        result = result * (n - i) / (i + 1);
    }
    // 转换__int128为字符串
    if (result == 0) return "0";
    string s;
    while (result > 0) {
        s += char('0' + (result % 10));
        result /= 10;
    }
    reverse(s.begin(), s.end());
    return s;
}

优点：结果完全精确，代码简洁；缺点：依赖非标准的__int128类型，仅部分编译器支持。

内容的提问来源于stack exchange，提问作者Alice Schwarze