From: pws AT s-a DOT amtp DOT liv DOT ac DOT uk (Peter Stephenson)
Subject: Re: asin (& atan) ok?
To: djgpp AT sun DOT soe DOT clarkson DOT edu (The DJGPP mailing list)
Date: Sat, 21 Nov 1992 19:28:01 +0000 (GMT)

Christoph Kukulies writes:
> As we were having some trouble with the exp() function (or was it pow())
> in previous versions of djgpp (--> 1.08),
> I recently had doubts whether asin() is OK over the defined range because
> a big fortran (f2c) program crashed in asin(). The fault may be elsewhere
> but just in case ... 

If you're using emu387, the atan() emulation is certainly a bit off for
large arguments: it gives 4*atan(1.0) as 3.05 instead of pi (does the
state of Tennessee know?  :-) (historical joke)).  The asin() routine
calls fpatan so suffers from this. wmemu387 seems to have the same problem.

It looks like it's computing atan as a power series for all arguments
(execept negatives are properly handled), so it doesn't converge for
>0.1 or so.  The first guess at what to do would be (this is a bit out
of my range of expertise)
1) for arguments >1, take the inverse (and subtract from pi/2 at the end),
2) for arguments 0.1 (say) to 1, use the half-angle formulae to reduce
the argument until it's less than 0.1 (this needs a square root each
time and so is not fast), perhaps with 1.0 special-cased.
3) use the power series as before (fewer iterations should be needed)
4) double the answer as many times as necessary.

I have this in straight C  (which I'll gladly pass on), but my low-level
DOS knowledge isn't currently good enough to translate it.  Perhaps
someone out there knows the algorithm the 387 actutally uses.  If no-one
more knowledgeable comes forward I might give it a go.

Hope this is a bit useful.
-- 
Peter Stephenson     <pws AT s-a DOT amtp DOT liv DOT ac DOT uk> or <PWS AT UKACRL>  (BITNET)
Dept. of Applied Mathematics and Theoretical Physics
University of Liverpool, P.O. Box 147, Liverpool, L69 3BX, U.K.
Tel. +44 51 794 3784 or 3791