[Mpc-discuss] asin, atan, etc

Stephen Montgomery-Smith stephen at missouri.edu
Sun Aug 12 19:58:06 CEST 2012


I was using mpc to check a casinh that I am writing for FreeBSD.  I 
noticed that your asinh functions can take a very long time to compute 
when for certain data.  I was wondering if you guys would be interested 
in using the same algorithm that I used.

You guys effectively use asinh(z) = log(z + sqrt(z^2+1)), and if z=x+I*y 
where x is very small and |y| <= 1, then this computation may have to be 
performed in double the precision you actually want. And my sample data 
was something like x=1e-3000, y=1-1e-3000.

I used the algorithm that can be found in "Implementing the complex 
arcsine and arccosine functions using exception handling" by T. E. Hull, 
Thomas F. Fairgrieve, and Ping Tak Peter Tang, published in ACM 
Transactions on Mathematical Software, Volume 23 Issue 3, 1997, Pages 
299-335, http://dl.acm.org/citation.cfm?id=275324.

It really is a very nice algorithm, and because you guys don't have 
problems with underflow or overflow, you don't need to mess with the 
exception handling.

I have described the algorithm at 
http://people.freebsd.org/~stephen/catrig.c, but the part you need is 
completely described in the opening comments, because underflow and 
overflow is not a problem for you guys.

Also, regarding atan, instead of using the formula:
        Im(atan(x+I*y)) = 1/4 * [log(x^2+(1+y)^2) - log (x^2 +(1-y)^2)]
you might like to use the formula
        Im(atan(x+I*y)) = 1/4 * [log1p (4y / (x^2 +(1-y)^2))]
when y > 0, and then use the fact that atan is odd for the case y<0.

Stephen




More information about the Mpc-discuss mailing list