Skip site navigation (1)Skip section navigation (2)
Date:      Sun, 12 Aug 2012 23:06:23 -0000
From:      Bruce Evans <brde@optusnet.com.au>
To:        Bruce Evans <brde@optusnet.com.au>
Cc:        Diane Bruce <db@db.net>, John Baldwin <jhb@freebsd.org>, David Chisnall <theraven@freebsd.org>, Stephen Montgomery-Smith <stephen@missouri.edu>, 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:  <20120717211625.E6848@besplex.bde.org>
Resent-Message-ID: <20120812230616.GX20453@server.rulingia.com>
In-Reply-To: <20120717200931.U6624@besplex.bde.org>
References:  <20120529045612.GB4445@server.rulingia.com> <20120711223247.GA9964@troutmask.apl.washington.edu> <20120713114100.GB83006@server.rulingia.com> <201207130818.38535.jhb@freebsd.org> <9EB2DA4F-19D7-4BA5-8811-D9451CB1D907@theravensnest.org> <C527B388-3537-406F-BA6D-2FA45B9EAA3B@FreeBSD.org> <20120713155805.GC81965@zim.MIT.EDU> <20120714120432.GA70706@server.rulingia.com> <20120717084457.U3890@besplex.bde.org> <5004A5C7.1040405@missouri.edu> <5004DEA9.1050001@missouri.edu> <20120717200931.U6624@besplex.bde.org>

next in thread | previous in thread | raw e-mail | index | archive | help
On Tue, 17 Jul 2012, Bruce Evans wrote:

> On Mon, 16 Jul 2012, Stephen Montgomery-Smith wrote:
>
>> On 07/16/2012 06:37 PM, Stephen Montgomery-Smith wrote:
>>> ...
>>> We might get lucky, and find that the definitions of csqrt and clog in
>>> the C99 standard are already set up so that the naive formulas for
>>> cacosh, etc, just work.  But whether they do or whether they don't, I
>>> think I can do it.  (As a first guess, I think that catanh and casinh
>>> will work "out of the box" but cacosh is going to take a bit more work.)
>
> See below what happened for naive formulars for ccosh.

I forgot the below before.

The following is ccosh() using mainly the naive formula.

% /*
%  * Hyperbolic cosine of a double complex argument.
%  *
%  * Most exceptional values are automatically correctly handled by the
%  * standard formula.  See n1124.pdf for details.
%  */
% 
% #include <complex.h>
% #include <math.h>
% 
% #include "math_private.h"
% 
% double complex
% ccosh1(double complex z)
% {
% 	double x, y;
% 
% 	x = creal(z);
% 	y = cimag(z);
% 
% 	/*
% 	 * This is subtler than it looks.  We handle the cases x == 0 and
% 	 * y == 0 directly not for efficiency, but to avoid multiplications
% 	 * that don't work like we need.  In these cases, the result must
% 	 * be almost pure real or a pure imaginary, except it has type
% 	 * float complex and its impure part may be -0.  Multiplication of
% 	 * +-0 by an infinity or a NaN would give a NaN for the impure part,
% 	 * and would raise an unwanted exception.
% 	 *
% 	 * We depend on cos(y) and sin(y) never being precisely +-0 except
% 	 * when y is +-0 to prevent getting NaNs from other cases of
% 	 * +-Inf*+-0.  This is true in infinite precision (since pi is
% 	 * irrational), and checking shows that it is also true after
% 	 * rounding to float precision.
% 	 */
% 	if (x == 0 && !isfinite(y))
% 		return (cpack(y - y, copysign(0, x * (y - y))));
% 	if (y == 0)
% 		return (cpack(cosh(x), isnan(x) ? copysign(0, (x + x) * y) :
% 		    copysign(0, x) * y));
% 	if (isinf(x) && !isfinite(y))
% 		return (cpack(x * x, x * (y - y)));
% 	if (fabs(x) > 710 && fabs(x) < 1455) {
% 		z = __ldexp_cexp(cpack(fabs(x), y), -1);
% 		return (cpack(creal(z), cimag(z) * copysign(1, x)));
% 	}
% 	return (cpack(cosh(x) * cos(y), sinh(x) * sin(y)));
% }

This was compared with the committed version and changed minimally to have
the same behaviour for exceptional args.  (Actually, more the reverse --
try the naive version and see what behaviour falls out of it for
exceptional args.  Then change both to match the behaviour specified in
C99.  This version was changed minimally to keep it as short as possible
and the committed version was changed maximally to have at least a comment
for all the exceptional cases.)

One problem with this is that it is hard to see even what the exceptional
cases are.  Another is that the classification macros are very slow and
only slightly easier to use than the bit-based classifications in the
committed version (the latter could probably be macro-ized to make them
almost as asy to use as the classification macros).

I added the __ldexp_cexp() code later to maintain the binary compatibility
of this with the comiitted version.  So it's not using the naive formula
for large |x|.

Bruce



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