伽玛函数（Gamma函数）C++实现

最新推荐文章于 2023-06-27 17:44:05 发布

东风吹柳

最新推荐文章于 2023-06-27 17:44:05 发布

阅读量6.1k

点赞数 6

分类专栏：分析基础

本文链接：https://blog.csdn.net/A_Pointer/article/details/108021352

版权

分析基础专栏收录该内容

7 篇文章 4 订阅

订阅专栏

前几天做项目中的一个算法，需要用到贝塞尔函数（半整阶，复数域），网上的代码一般都要收费，索性自己撸了一个，其中涉及伽马函数值求取，下面是我的思路与实现，供大家免费使用。

第一步：首先介绍一下Gamma函数。

（1）在实数域上伽玛函数定义为：

$\Gamma (x)=\int_{0}^{\infty }t^{x-1}e^{-t}dt(x>0)$

（2）在复数域上伽玛函数定义为：

$\Gamma (z)=\int_{0}^{\infty }t^{z-1}e^{-t}dt (Re(z)>0)$

上述两个定义在百度百科上可以查到。

第二步：设计思路简介

本实现的设计思路如下：

当 0<x<1 时，Gamma函数为采用下列极限求值：

$\Gamma (x+1)=\lim_{n\rightarrow \infty } \frac{n!n^{x}}{\prod_{m=1}^{n}(x+m)}(x\in R^{*})$

当 x> 1 时，Gamma函数为采用 $\Gamma (s)=(s-1)\Gamma(s-1)$ 递推至 0<s<1 。

当 x=1 时， $\Gamma (1)=0!=1$ 。

当 x< 0 时，Gamma函数为采用 $\Gamma (s)=\Gamma(s+1)/s$ 递推至 0<s<1 。

第三步：实现Gamma函数。

#include <cmath>

//伽马函数（辅助函数，用于求解贝塞尔函数,任意阶数）
//x 变量值
//setAbsRelaErr 相对误差绝对值
double gamma(const double x, const double setAbsRelaErr){
	//初始条件判断
	if (abs(1.0 - x) < 0.00001) {
		return 1.0;
	}else if (abs(x - 0.5) < 0.00001) {
		return sqrt(3.1415926);
	}

    //递归求取伽马函数值
	if (x > 1.0){
		return (x-1)*gamma(x - 1, setAbsRelaErr);
	}else if (x < 0.0) {
		return gamma(x + 1, setAbsRelaErr)/ x ;
	}

	double res = 0.0;
	double temp = 1.0;
	double check = 0.0;
	int i = 1;
	while (abs((check - temp)/ temp)> setAbsRelaErr)
	{
		check = temp;
		temp *= i / (x - 1 + i);
		i++;
	}
	res = temp * pow(i, x - 1);

	return res;
}

第四步：验证