Skip site navigation (1)Skip section navigation (2)
Date:      Sat, 22 Sep 2012 05:05:10 +1000 (EST)
From:      Bruce Evans <brde@optusnet.com.au>
To:        Stephen Montgomery-Smith <stephen@missouri.edu>
Cc:        freebsd-numerics@FreeBSD.org
Subject:   Re: Complex arg-trig functions
Message-ID:  <20120922042112.E3044@besplex.bde.org>
In-Reply-To: <505C7490.90600@missouri.edu>
References:  <5017111E.6060003@missouri.edu> <50526050.2070303@missouri.edu> <20120914212403.H1983@besplex.bde.org> <50538E28.6050400@missouri.edu> <20120915231032.C2669@besplex.bde.org> <50548E15.3010405@missouri.edu> <5054C027.2040008@missouri.edu> <5054C200.7090307@missouri.edu> <20120916041132.D6344@besplex.bde.org> <50553424.2080902@missouri.edu> <20120916134730.Y957@besplex.bde.org> <5055ECA8.2080008@missouri.edu> <20120917022614.R2943@besplex.bde.org> <50562213.9020400@missouri.edu> <20120917060116.G3825@besplex.bde.org> <50563C57.60806@missouri.edu> <20120918012459.V5094@besplex.bde.org> <5057A932.3000603@missouri.edu> <5057F24B.7020605@missouri.edu> <20120918162105.U991@besplex.bde.org> <20120918232850.N2144@besplex.bde.org> <20120919010613.T2493@besplex.bde.org> <505BD9B4.8020801@missouri.edu> <20120921172402.W945@besplex.bde.org> <20120921212525.W1732@besplex.bde.org> <505C7490.90600@missouri.edu>

next in thread | previous in thread | raw e-mail | index | archive | help
On Fri, 21 Sep 2012, Stephen Montgomery-Smith wrote:

> I will keep the inexact raising up-front (mostly because I forgot how I did 
> it earlier).

It is working very well.  I have minor cleanups to it using a 1-line
raise_inexact() macro.

> I will still use atanh(fl), and rely on someone else to fix it.  (If it is 
> inexact near 0, it is only a few ULP, and that is good enough for me.)

I see you simplified the SQRT_[34]_EPSILON thresholds and fixed their
initialization by removing them (the function calls in the static
initializers didn't compile).  I spent more time on them and found the
best values in practice (take the sqrt's accurately and divide them by
2 or 4 instead of your 8), but they are painful to initialize, especiallu
for long doubles.  A few points of more general interest turned up while
debugging this,
- atan()'s series is alternating, while the others are not.  Alternation
   causes more cancelation errors.
- FOO_EPSILON only applies to 1 side of an addition.  E.g., it applies
   to 1.0 + x where x > 0, for 1.0 + x where x < 0 the size of the
   corresponding epsilon is half as much.  Non-alternation means that
   the FOO_EPSILON side applies.
- the general approximation in cacos(z) and casin(z) is quite good for
   small z, so larger thresholds for using the special approximations
   don't affect accuracy much.  However, the general approximation in
   catanh(z) is not so good for small z, so using larger threshold for
   the special approximation affects accuracy significantly.  I found
   the approximate best point to switch the approximations.
- some combination of the previous 3 points means that the switching
   point is about twice as large relative to the SQRT_N_EPSILON threshold
   for catan() as for the others (divide by 2 instead of 4).

> I'll go ahead and see about the pio2h and pio2l.

Please wait for my patch for this.  It has all the details for all
precisions including pretty-printing the declarations.  Or you
can do a nearly-global substition of m_pi_2 by pio2_hi + pio2_lo.

My other mostly-complete changes:
- avoid all the scalb() and related calls.  This makes do_hard_work()
   a bit faster and simpler and real_value_reciprocal() much faster and
   a bit more complex.
- make real_value_reciprocal() handle signs (everything is automatic
   except for x = inf), and avoid a copysign() after it
- a few improvements in comments

My other unfinished changes:
- figure out if the up-front things in catanh() are best placed there.
- decide whether to handle pure real and pure imaginary args specially
   like I do for both in catanhf().  This interacts with the previous
   point.
- decide whether my old change to remove unnecessary accuracy for the
   case where ax == 1, ay < FLT_EPSILON in catanh() is correct (you
   didn't accept it, and maybe other accuracy changes make it extra
   accuracy more interesting).

Bruce



Want to link to this message? Use this URL: <https://mail-archive.FreeBSD.org/cgi/mid.cgi?20120922042112.E3044>