complex acos is occasionally wrong

[Stephen Montgomery-Smith <stephen@…>]

I am finding that the acos function is getting some of them wrong. For example, your program evaluates acos(1.00000002785941 + I*5.72464869028403e-200) as 0 - I*0.000236048349018331 whereas it should be 2.42520172707401e-196 - I*0.000236048349018331.

Looking at http://www.boost.org/doc/libs/1_51_0/boost/math/complex/acos.hpp, I am somewhat sure that the mistake is in the last line of this code fragment:

This is the Hull et al exception handling code from Fig 6 of their paper: if(y <= (std::numeric_limits<T>::epsilon() * std::fabs(xm1))) {

if(x < one) {

real = std::acos(x); imag = y / std::sqrt(xp1*(one-x));

} else {

real = 0;

For asin, setting real = half_pi does just fine. But for acos, real should be something extremely small, but definitely not 0.

[Stephen Montgomery-Smith <stephen@…>]

It looks like the mistake is in Figure 6 of the paper by Hull, Fairgrieve and Tang. It looks like this patch fixes it.

diff -u boost-trunk/boost/math/complex/acos.hpp boost-trunk-new/boost/math/complex/acos.hpp
--- boost-trunk/boost/math/complex/acos.hpp     2012-08-26 18:24:02.000000000 +0000
+++ boost-trunk-new/boost/math/complex/acos.hpp 2012-08-26 18:25:48.000000000 +0000
@@ -172,14 +172,15 @@
             }
             else
             {
-               real = 0;
                if(((std::numeric_limits<T>::max)() / xp1) > xm1)
                {
                   // xp1 * xm1 won't overflow:
+                  real = y / std::sqrt(xm1*xp1);
                   imag = boost::math::log1p(xm1 + std::sqrt(xp1*xm1));
                }
                else
                {
+                  real = y / x;
                   imag = log_two + std::log(x);
                }
             }

[Stephen Montgomery-Smith <stephen@…>]

Looking at page 327 of the paper by Hull et al, they say that the real part can be set to zero, because the absolute error will be negligible compared to the absolute error in the imaginary part. But everywhere else in the paper, they manage to control the relative error of the real and imaginary parts separately.

Looking at page 304, it seems that their goal is only to control the relative error of the real and imaginary parts together. So this wasn't really a mistake in their paper.

But I think the "higher goal" of controlling the relative error of the real and imaginary parts separately is worthy, so I still commend my patch to you.

[Stephen Montgomery-Smith <stephen@…>]

Another thing - through experiment I have found that setting a_crossover to 10 instead of 1.5 works a little better. This is for asin as well as acos.

[John Maddock]

(In [81624]) Lots of small patches. Update and regenerate docs. Fixes #7183. Fixes #7290. Fixes #7291. Fixes #7649. Refs #7492.