Java实现伽马函数gamma()

实习的项目里遇到了第一类贝塞尔函数，C++里有直接的贝塞尔函数可以使用，但是是实数域的而且仅C++17支持。项目里是要用到复数域的半整阶的第一类贝塞尔函数和第二类汉克尔函数，所以就直接根据原始定义实现了半整阶的第一类贝塞尔函数。贝塞尔函数定义公式可以用伽马函数求解，就顺便把伽马函数也实现了。写完了发现C++的cmath里有现成的伽马函数tgamma()有float、double、long double三种形式。

思路分析

几个特殊的伽马函数值

Γ(0) = 1

Γ(1) = 1

Γ(1/2) = √π

将这些特殊值作为递归结束的条件。

在x > 1时，使用递推公式Γ(x + 1) = x * Γ(x)推至0 < x < 1

在x < 0时，使用递推公式Γ(x) = Γ(x + 1) / x推至0 < x < 1

在(0,1)区间内的伽马值若无法递推到特殊值，则使用伽马函数原始定义公式实现。

注*：因为程序在计算时，会损失精度，导致计算的结果不精确。故在运行时设置setAbsRelaErr=0.00001相对误差绝对值来保护精度。

递归实现伽马函数

public static double gamma(double x, double setAbsRelaErr) {
	//setAbsRelaErr 相对误差绝对值
	//递归结束条件
	if(x < 0) {
		return gamma(x + 1, setAbsRelaErr) / x;
	}
	if(Math.abs(1.0 - x) < 0.00001) {
		return 1;
	}
	if(Math.abs(0.5 - x) < 0.00001) {
		return Math.sqrt(3.1415926);
	}
	
	if(x > 1.0) {
		return (x - 1) * gamma(x - 1, setAbsRelaErr);
	}
	
	double res = 0.0;
	double temp = 1.0;
	double check = 0.0;
	int i = 1;
	while(Math.abs((check - temp) / temp) > setAbsRelaErr){
		check = temp;
		temp *= i / (x - 1 + i);
		i++;
	}
	res = temp * Math.pow(i, x - 1);
	return res;
}

public static void main(String[] args) {
	System.out.println(gamma(-0.5, 0.00001));
}
           

循环实现伽马函数

使用伽马函数的递推公式Γ(x + 1) = x * Γ(x)实现，传进去的形参是x，根据公式在函数里先减1，再根据公式计算。while循环的判断条件是前后两个数的差值不超过0.00001，以此来限制精度。跟递归实现相比，时间复杂度相当高。

public static double gamma(double x, double setAbsRelaErr) {
	x -= 1;
	double res = 0.0;
	int i = 1;
	double temp = 1.0;
	double temp0 = 1.0;
	double temp1 = 1.0;
	temp = 1 / (x + 1);
	do {
		temp0 = temp * Math.pow(i, x);
		i++;
		temp *= i / (x + i);
		temp1 = temp * Math.pow(i, x);
	}while(Math.abs(temp1 - temp0) > 0.00001);
	res = temp1;
	return res;
}