Skip site navigation (1)Skip section navigation (2)
Date:      Sun, 12 Aug 2012 23:02:54 -0000
From:      Bruce Evans <brde@optusnet.com.au>
To:        Stephen Montgomery-Smith <stephen@missouri.edu>
Cc:        Diane Bruce <db@db.net>, Bruce Evans <brde@optusnet.com.au>, John Baldwin <jhb@freebsd.org>, David Chisnall <theraven@freebsd.org>, Bruce Evans <bde@freebsd.org>, Steve Kargl <sgk@troutmask.apl.washington.edu>, David Schultz <das@freebsd.org>, Peter Jeremy <peter@rulingia.com>, Warner Losh <imp@bsdimp.com>
Subject:   Re: Use of C99 extra long double math functions after r236148
Message-ID:  <20120724113214.G934@besplex.bde.org>
Resent-Message-ID: <20120812230247.GD20453@server.rulingia.com>
In-Reply-To: <500DAD41.5030104@missouri.edu>
References:  <20120717084457.U3890@besplex.bde.org> <5004A5C7.1040405@missouri.edu> <5004DEA9.1050001@missouri.edu> <20120717040118.GA86840@troutmask.apl.washington.edu> <20120717042125.GF66913@server.rulingia.com> <20120717043848.GB87001@troutmask.apl.washington.edu> <20120717225328.GA86902@server.rulingia.com> <20120717232740.GA95026@troutmask.apl.washington.edu> <20120718001337.GA87817@server.rulingia.com> <20120718123627.D1575@besplex.bde.org> <20120722121219.GC73662@server.rulingia.com> <500DAD41.5030104@missouri.edu>

next in thread | previous in thread | raw e-mail | index | archive | help
On Mon, 23 Jul 2012, Stephen Montgomery-Smith wrote:

> I just realized that catan(z) = reverse( catanh(reverse (z))), just like 
> casin relates to casinh (remember reverse(x+I*y) = y+I*x).  This is a 
> consequence of catan and catanh being odd functions, as well as the standard 
> relation catan(z) = -I*catanh(I*z).
>
> So I would modify Peter's code by taking out the minus signs.
>
> Maybe it would make a difference if the answers involved -0.

Quite likely.  Note that keeping the real and imaginary parts separate
is giving perfect reflection across some axes.  Reflecting 0 is for
other axes and it doesn't fall out so accidentally.

> I feel that I am done with these functions for now.  I tried to change my 
> comments to conform to the style given to me by Bruce.  However spacing 
> inside mathematical expressions is something where I am inconsistent.

Run it through indent to find all the inconsistencies.

> The functions still need a lot of work to handle -0, infs and NaNs correctly. 
> I will leave that to you guys, because you seem so much better at it than me. 
> I still don't understand why the proper test is "if (x!=x) return(x+x)" 
> rather than "if (isnan(x)) return(NAN)".

First, signaling NaNs are required to signal, in a strong way:
   (IEEE doesn't specify complex functions of course, but it requires this
   for all the functions and most of the operations that it specifies.  The
   main exceptions (in IEEE-854-1987) are:
   - copying a signaling NaN without changing its precision.  Whether this
     triggers is implementation-defined.  So i387 conforms by triggering
     for changing the precision on loads iff the load does a non-null
     change to long double precision
   - the sign of a NaN is not specified.  Triggering for signaling NaNs is
     still required "for every operation listed in Section 5".  These
     operations don't seem to include absolute values or negation.  Thus
     the i387 fabs and fchs conform because they are extensions, but they
     should never be used!?, since any use doesn't conform to the spirit
     of the standard.)
Every "operation" involving a signaling NaN or invalid operation shall,
if no "trap" occurs and a floating point result is to be delivered,
deliver "a" quiet NaN as its result.

"if (isnan(x)) return(NAN)" does extra work to fail to conform to this.
First it uses the slow isnan() classification to be sure that no "trap"
occurs.
   (Whether (x != x) triggers signaling NaNs, and whether the triggering
   involves a trap or just setting FE_INVALID is delicate.  This was
   broken in the original i8087.  Actually, on checking the details, I
   found that it was for quiet NaNs that was broken.  (x != x) should
   and does set FE_INVALID for signaling NaNs on all x87.  This traps
   iff it is unmasked.  Comparison doesn't deliver an FP result, so
   quietening of x is irrelevant.  So (x != x) works right for signaling
   NaNs on all x87.  The i8087 bug was that its comparison operators
   all set FE_INVALID, etc for quiet NaNs too.  Later x87 have comparison
   operators that work right on all types of NaNs, and compilers use them.
   Thus you can write isnan() as (x != x) on at least x87 without getting
   a spurious FE_INVALID for quiet NaNs.  I think this is unportable, so
   isnan() should be used if it is necessary to avoid such FE_INVALID.
   But where we want FE_INVALID, and we don't care if we get it from
   both the comparison and (x+x).

Second, returning x wouldn't conform, since the returned x might not be
quiet if the original x isn't.  Returning the `NAN' does conform,
provided NaN is quiet, since the standard allows "a" (that is, any
quiet NaN to be returned).  BTW, C has very bad names for NaNs and Infs.
`NaN' should be spelled with its 'a' in lower case, as in the standard.
Here the NaN returned should be spelled qNaN or as my dNaN (default NaN)
to emphasize that it is a quiet default NaN.

The standard then recommends which NaN "should" be returned, and it isn't
the C default NaN:
   Every operation [...] shall [...] if a floating point result is to be
   delivered, [when this result is specified to be "a" quiet NaN,] then
   this NaN _should_ be _one_ of its input NaNs.  Note that precision
   conversions might be unable to deliver the same NaN.

This is only a "should", but returning the default C NaN fails to
conform to at least its spirit, by making no attempt to return the
original NaN.  Returning the an original one is especially simple when
there is only one arg.  This arg just needs to be quietened and undergo
any necessary conversions.  (x+x) does this (modulo excessive
conversions due to C's promotion rules).  When the arg is a signaling
NaN, it is impossible to return exactly it, so returning an unrelated
default would be more reasonable, but doing this would take more code.
It is best to depend on the hardware doing the right thing for (x+x)
on signaling NaN x.  The hardware is specified to give a quiet result
and would need much the same complications as us return a default
instead of x slightly modified for signaling NaNs only.

This also specifies NaN mixing: when there are _two_ args, _one_ of them
_should_ be returned.  Unfortunately, which one is not specified.  But
hardware operations should pick one, and we can do about as well as the
hardware by mixing (= choosing) using (x+y).

Bruce



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