Arctran.cin

// arctran.cin is routine in C++ that evaluates function ArcTra of complex argument.

// The same function can be expressed also through the Tania function, but the direct implementation below seems to be faster.

//

z_type tra(z_type z){ return z+exp(z);}

z_type arctra1(z_type z){ static DB c[22]={0.,0.5,-0.0625,.005208333333333333, // 0 - 3
0.0003255208333333333, -0.00021158854166666667, //5
0.00003187391493055556, 1.7680819072420636e-6, //7
-1.8520960732111855e-6, 3.5344398000673434e-7, //9
 8.173825669330684e-9, -2.0624513960198102e-8, //11
 4.727254469326131e-9, -5.254838295493782e-11, //13
-2.5497322696797093e-10,6.897078914225712e-11, //15
-3.004636375311634e-12,-3.313983848679892e-12, //17
 1.054773386200646e-12,-7.811289918947213e-14, //19
-4.391790189701172e-14, 1.6573692169992586e-14}; //21
z_type t=z-1.; z_type s=c[21]*z; int n;
for(n=20;n>0;n--){ s+=c[n]; s*=t;}
DO(n,3) {s+=(z-s-exp(s))/(1.+exp(s));} return s;}

z_type arctra2(z_type z) { z_type s,d,t; int n;
//t=z+z_type(1.,-M_PI); t/=2.; t=sqrt(t);
d=z+z_type(1.,-M_PI); t=.5*d; t=sqrt(t);
s= I*M_PI
+t*( z_type(0.,-2.)
+t*( 2./3
+t*( z_type(0.,2./9.)
+t*( -8./135.
+t*( z_type(0.,-1./135.)
+t*( -32./8505.
+t*( z_type(0.,-139./42525.)
+t*( 32./25515.
+t*( z_type(0.,571./4592700.)
+t*( 35968./189448875.
))))))))));
if(abs(d)<.02) return s;
DO(n,3) {s+=(z-s-exp(s))/(1.+exp(s));}
return s;
}


z_type arctra3(z_type z){z_type c[5]; z_type L=log(z); z_type s; int n,M=4;
c[0]= -1.;
c[1]= 1. + L*(-1./2. );
c[2]= -1. + L*( 1.5 + L*( -1./3.)) ;
c[3]= 1. + L*(-3. + L*( 11./6. + L*(-1./4.))) ;
c[4]= -1. + L*( 5. + L*(-35./6. + L*(25./12. + L*(-1./5.))));
s=c[M]; for(n=M-1;n>=0;n--){ s/=z; s+=c[n];}
s=L*(1.+s/z);
DO(n,3) {s+=(z-s-exp(s))/(1.+exp(s));} return s;
}

z_type arctra4(z_type z) { z_type e=exp(z),s;
s= z +e*(-1. +e*( 1. +e*(-0.5)));
int n; DO(n,3) {s+=(z-s-exp(s))/(1.+exp(s));} return s;
}

z_type arctran(z_type z) { DB x=Re(z), y=Im(z);
if( x>2.) return arctra3(z);
DB Y=fabs(y);
if(Y<M_PI){ if(x<-1.5) return arctra4(z);
            if(Y<2.) return arctra1(z); }
if( Y>5. || x<-4. ) return arctra3(z);
if( y>0. ) return arctra2(z);
return conj( arctra2(conj(z)) );
}
//

//