SuFac.cin

From TORI
Jump to: navigation, search

// SuFac.cin defines z_type SuFac(z_type) that evaluates // Superfunction of Factorial built up at the fixed pont 2. // Of order of 14 correct decimal digits are expected in the result.

z_type superfac0(z_type z){ int n; z_type s;
//      DB K=1.8455686701969342788;
DB k=0.61278745233070836381366079016859252;  //k=log(K);
DB u[21]={2.,1., //0,1
.798731835172434541585621072345730147,  // 2
.577880975476483235803807592348110833,  // 3
.393978809662971757177848639852917378,  // 4
.257533958032332679820773329133486586,  // 5
.162901958103705249541496101752195514,  // 6
.100282419171352371943554511785342142,  // 7
.0603184725913977494512136774562415014, // 8
.0355544582258061836048059212969418417, // 9
.0205859954874424134686332481358935023, //10
.0117302279624549548734823541033644211, //11
.00658835541777254650743317221091667507,//12
.00365218351418374834372649788987162842,//13
.00200039479760669665711545138631474960,//14
.00108362752868222808502286098449166985,//15
.000581036636299227699924018045799185045,//16
.000308601963223618214714523083268563975,//17
.000162 ,.000084, 0.000043      //18,19,20
};
z_type e=exp(k*z);   
       s=u[20]; for(n=19;n>=0;n--){s*=e; s+=u[n];}
//      s=u[15]; for(n=14;n>=0;n--){s*=e; s+=u[n];}
return s;}
z_type SuFac(z_type z){
       if(Re(z)>-2.)   return fac(SuFac(z-1.));
                       return superfac0(z);
       }

//