Rでガウスの超幾何関数を計算する。
hypergeo 関数は一部パラメータでバグがあるのでgslの超幾何関数でおきかえる。
ただし、gslの超幾何関数は複素数に解析接続されていないので、そこはhypergeoを使う。
library("hypergeo")
library(gsl)
library("BAS")
###2017.03.09
#beta5<-0.4
Gauss2F1 <- function(a,b,c,x){
if(abs(Re(x))>=0 & abs(Re(x))<1){
if(x<0){
##For a==b and a<0, are there the bugs in hyperg_2F1, in which sign is reversed?
#-hyperg_2F1(a,c-b,c,1-1/(1-x))/(1-x)^a
if(a==b & a<0){
hyperg_2F1(a,b,c,x)
}else{
hyperg_2F1(a,b,c,x)
}
}else{
#hypergeo(c-a,b,c,1-1/(1-x))/(1-x)^b
hypergeo(a,b,c,x)
}
}else{
if(Re(1-1/(1-x))>0){
# hyperg_2F1(c-a,b,c,1-1/(1-x))/(1-x)^b
if(a!=b){
#hyperg_2F1(a,c-b,c,1-1/(1-x))/(1-x)^a
hyperg_2F1(c-a,b,c,1-1/(1-x))/(1-x)^b
}else{
##For a==b and a<0, are there the bugs in hyperg_2F1, in which sign is reversed?
#-hyperg_2F1(a,c-b,c,1-1/(1-x))/(1-x)^a
if(a<0){
-hyperg_2F1(c-a,b,c,1-1/(1-x))/(1-x)^b
}else{
hyperg_2F1(c-a,b,c,1-1/(1-x))/(1-x)^b
}
}
#coef1=gamma(c)*gamma(b-a)*(1-x)^-a/(gamma(b)*gamma(c-a));
#coef2=gamma(c)*gamma(a-b)*(1-x)^-b/(gamma(a)*gamma(c-b));
#coef1*hyperg_2F1(a,c-b,a-b+1,1/(1-x))+coef2*hyperg_2F1(b,c-a,b-a+1,1/(1-x));
}else{
#hypergeo(c-a,b,c,1-1/(1-x))/(1-x)^b
hypergeo(a,b,c,x)
}
}
}