Use of C99 extra long double math functions after r236148

[Bruce Evans]

On Thu, 19 Jul 2012, Stephen Montgomery-Smith wrote:

> On 07/19/2012 01:37 AM, Bruce Evans wrote:
>> ...
>> problem.  Complex functions should have only poles and zeros, with
>> projective infinity and "projective zero" (= inverse of projective
>> infinity).  Real functions can and do have affine infinities and zeros
>> (+-Inf and +-0), with more detailed special cases.  It's just impossible
>> to have useful, detailed special cases for all the ways of approaching
>> complex (projective) infinity and 0.
>> ...
>> sign functions).  Hopefully, the specification of imag(clog()) is
>> that it has the same sign behaviour as atan2(), so you can just use
>> atan2().  The sign conventions for both are arbitrary, but they
>> shouldn't be gratuitously different.  You still have to check that
>> they aren't non-gratuitously different, because different conventions
>> became established.
>
> I checked.  Actually the sign conventions are not that arbitrary.  But as a 
> mathematician I would say they are a bit useless, e.g.
> atan(infinity,infinity) = pi/4 = 45 degrees
> How do you know that the two infinities are the same?  One could be double 
> the other.

It should be NaN with projective infinity.

> If it had been up to me, there would have been finite numbers, and nan.  And 
> none of this -0.

I think Kahan is a mathematician, and is primarily responsible for +-0.
+-0 give poor man's branch cuts for real functions.

>> I don't know what happens with zeros of complex inverse trig functions.
>> I think they don't have many (like log()), but their real and imaginary
>> parts do, and they are too general for accurate behaviour of the real
>> and imaginary parts relative to themselves to fall out.
>
> casinh(z) is zero only when z=0, and near that point I could use Taylor's 
> series (but a lot of terms would be needed because the Taylot series 
> converges quite slowly).

To get accuracy near zeros, you only need to use a series method in a
very small radius.  Even a linear approximation may be enough, and the
main difficulty is the linear term:

 	f() ~= f'(z0) * (z-z0)

where z0 typically needs to be known to hundreds or thousands of bits of
precision and the subtraction must be done in this precision.  f'(z0) and
the multiplkication only need a couple of extra bits.  This is only easy
when z0 is a dyadic rational.

> I can now see that the separate cases of the real part and imaginary parts of 
> casinh being zero is going to be hard.

I won't ask for that and will measure errors relative to the absolute value
of the result.

Bruce