专栏文章

2024.11.26

LLeNotFound
个人记录参与者 1已保存评论 0

文章操作

快速查看文章及其快照的属性，并进行相关操作。

当前评论
0 条
当前快照
1 份
快照标识符
@miqzng1f
此快照首次捕获于
2025/12/04 13:21
3 个月前
此快照最后确认于
2025/12/04 13:21
3 个月前
查看原文
今天上课突然发现有人作业写了 std::tgamma 这个函数，去查了一下，好像可以 O(1)O(1) 计算 n−1n-1 的阶乘，C++11 引入的，ICPC 可用
在 GCC 中，std::tgamma 的实现通常依赖于底层的数学库（libm）。在 GNU libm 库中，tgamma 函数的实现使用了一些经典的方法，如 Lanczos 逼近法 和 反射公式。以下是具体的实现细节概述：

关于 tgamma

1. 基本结构

tgamma 函数的实现可以分为几个部分：
  1. 特殊值处理：
    • 检测输入是否为负整数或零（伽马函数在这些点未定义）。
    • 对极端值（如 NaN 或 Inf）返回适当的结果。
  2. 正数的伽马值计算：
    • 对于大于一定范围的正数，直接使用 斯特林近似。
    • 对于中等范围的正数，使用 Lanczos 逼近法。
  3. 负数和小数的处理：
    • 使用 反射公式： [ \Gamma(1-x)\Gamma(x) = \frac{\pi}{\sin(\pi x)} ] 将负数或小数转化为正数计算。

2. 特殊值处理

在伽马函数中，以下特殊情况需要单独处理：
  • 负整数或零：
    伽马函数在这些点具有奇点（无穷大），因此返回 HUGE_VAL 并设置 errno 为 EDOM。
  • NaN 或无穷：
    对输入 NaN 返回 NaN，对正无穷返回正无穷。
这部分通常用简单的分支语句处理。

3. 中等范围正数：Lanczos 逼近法

在 GNU 的数学库中，tgamma 的核心计算部分大多基于 Lanczos 逼近法。这是一个快速且高精度的伽马函数近似算法。
Lanczos 逼近公式 Lanczos 方法通过一个展开式逼近伽马函数： [ \Gamma(x) \approx \sqrt{2\pi} \cdot \left(x + g - \frac{1}{2}\right)^{x-\frac{1}{2}} \cdot e^{-(x+g-\frac{1}{2})} \cdot A(x) ] 其中：
  • ( g ) 是一个常数（通常选择 ( g \approx 7 )）。
  • ( A(x) ) 是一个由系数 ( a_k ) 定义的级数： [ A(x) = \frac{a_0}{x+g} + \frac{a_1}{x+g+1} + \cdots ]
  • 系数 ( a_k ) 是预先计算并存储的常量。
这种方法在保证精度的同时计算复杂度较低。

4. 大范围正数：斯特林近似

对于较大的 ( x )，伽马函数值增长迅速，因此采用 斯特林近似： [ \Gamma(x) \sim \sqrt{2\pi} x^{x-\frac{1}{2}} e^{-x} ] 该公式适用于 ( x > 10 ) 或更大时，计算简单且精度足够。

5. 负数和小数：反射公式

对于负数输入 ( x ) 或小数输入 ( x < 1 )，使用反射公式： [ \Gamma(x) = \frac{\pi}{\sin(\pi x) \cdot \Gamma(1-x)} ] 这将计算转化为对 ( 1-x ) 的正数处理，避免直接处理负数。

6. GCC 的具体实现（简化伪代码）

以下是 tgamma 在 libm 中的一个典型实现结构：
C
double tgamma(double x) {
    if (x <= 0 && x == (int)x) {
        // 负整数或零，伽马函数未定义
        errno = EDOM;
        return HUGE_VAL;
    }

    if (x < 0) {
        // 负数：使用反射公式
        double sin_pi_x = sin(M_PI * x);
        if (sin_pi_x == 0) {
            errno = ERANGE;
            return HUGE_VAL;
        }
        return M_PI / (sin_pi_x * tgamma(1 - x));
    }

    if (x > THRESHOLD) {
        // 大数：斯特林近似
        return stirling_approximation(x);
    }

    // 中等范围：Lanczos 逼近
    return lanczos_approximation(x);
}

7. 性能优化

  • 分段处理：根据输入的范围选择最优算法（反射公式、Lanczos 逼近或斯特林公式）。
  • 查表优化：预存 Lanczos 系数或某些特殊值以加速计算。
  • 特殊输入快速返回：提前处理简单输入（如正整数）。

总结

GCC 中的 std::tgamma 主要依赖 Lanczos 逼近法 和 反射公式，并结合 斯特林近似 来高效处理伽马函数的不同输入范围。它的实现经过高度优化，既保证了高精度又能适应各种数值计算场景。如果需要具体实现代码，可以查阅 GNU libm 的源码，通常位于 glibc 仓库的 math 模块中。

关于 midpoint

GPT 结束，然后刚想去逛逛菜园，发现他们在讨论 std::midpoint，更好的计算区间中点，向第一个参数靠拢，可以避免溢出以及精度问题。
更逆天的是还能 midpoint(r,l)

关于 bits/////////stdc++.h

帮舍友调代码的时候发现了一个骚操作，原来编译器早就想到了可能的错误。
GPT FS 的结果
后续测试了 Linux 下的编译，无论多少个 / 都能过编，反向的话只有两条的情况 \\\\ 可以。

讨论区 可能有用的东西

关于重载[]运算符
芝士++ 大量空格可以用来卡 getchar() 快读
快读快写为什么会负优化？

评论

共 0 条评论，欢迎与作者交流。

正在加载评论...