Ausin.cin

From TORI
Jump to: navigation, search

// ausin.cin is complex double implementation of function AuSin, which is inverse of SuSin and the Abel function for sin.

//


z_type ausin0(z_type z){ z_type b,d,L; int n,N; DB c[32]=
{ 0.075238095238095238,0.011047619047619048, 0.0025161272589844018,0.00064214089261708309,
 0.00015466661256774635,0.000027370601208859475, -2.7882266241028593e-7,-2.6398530640882673e-6,
 -8.4436657963874166e-7,2.1258409455641387e-7, 2.9625812021061085e-7,5.5084490596214015e-8,
 -8.505934497279864e-8,-5.7880498792647227e-8, 1.931195823640328e-8,4.1375071763898608e-8,
 3.8895251916390044e-9,-2.9644790466835329e-8, -1.5302998576174588e-8,2.2013832459456463e-8,
 2.5112957261097041e-8,-1.5610882396047076e-8, -3.8703209933396839e-8,6.3863715618218346e-9,
 6.2295996364195277e-8,1.4512488433891138e-8, -1.0741668102638471e-7,-7.2006309159592681e-8,
 1.9825395033112708e-7,2.4402848447224951e-7, -3.8457536962314424e-7,-7.917263057277091e-7};
N=21; L=log(z); b=z*z; //d=b*(c[N]*.5);
d=b*c[N]*.5;
for(n=N;n>=0;n--) { d+=c[n]; d*=b;} //d=b*(79./1050.+ b*(.29/2635.));
d+=1.2*L+3./b;
//return d- 2.08962271972955094;
//return d- 2.089622719729546;
return d- 2.089622719729524;
}

z_type ausin(z_type z){ int m,M;
if(Im(z)<0) return conj(ausin(conj(z)));
M=50; DO(m,M) z=sin(z);
return ausin0(z)-(0.+M);
}

//