From owner-freebsd-standards@FreeBSD.ORG Sun Oct 2 19:19:02 2005 Return-Path: X-Original-To: freebsd-standards@freebsd.org Delivered-To: freebsd-standards@freebsd.org Received: from mx1.FreeBSD.org (mx1.freebsd.org [216.136.204.125]) by hub.freebsd.org (Postfix) with ESMTP id E81FA16A420 for ; Sun, 2 Oct 2005 19:19:02 +0000 (GMT) (envelope-from sgk@troutmask.apl.washington.edu) Received: from troutmask.apl.washington.edu (troutmask.apl.washington.edu [128.208.78.105]) by mx1.FreeBSD.org (Postfix) with ESMTP id 39D8E43D5E for ; Sun, 2 Oct 2005 19:18:49 +0000 (GMT) (envelope-from sgk@troutmask.apl.washington.edu) Received: from troutmask.apl.washington.edu (localhost [127.0.0.1]) by troutmask.apl.washington.edu (8.13.4/8.13.4) with ESMTP id j92JIl3H040651; Sun, 2 Oct 2005 12:18:47 -0700 (PDT) (envelope-from sgk@troutmask.apl.washington.edu) Received: (from sgk@localhost) by troutmask.apl.washington.edu (8.13.4/8.13.1/Submit) id j92JIfen040650; Sun, 2 Oct 2005 12:18:41 -0700 (PDT) (envelope-from sgk) Date: Sun, 2 Oct 2005 12:18:41 -0700 From: Steve Kargl To: Bruce Evans Message-ID: <20051002191841.GA40568@troutmask.apl.washington.edu> References: <20050929195552.GA14982@troutmask.apl.washington.edu> <20050930221846.T34595@delplex.bde.org> <20050930152734.GA20655@troutmask.apl.washington.edu> <20051001204005.N37704@delplex.bde.org> Mime-Version: 1.0 Content-Type: text/plain; charset=us-ascii Content-Disposition: inline In-Reply-To: <20051001204005.N37704@delplex.bde.org> User-Agent: Mutt/1.4.2.1i Cc: freebsd-standards@freebsd.org Subject: Re: complex.h math functions X-BeenThere: freebsd-standards@freebsd.org X-Mailman-Version: 2.1.5 Precedence: list List-Id: Standards compliance List-Unsubscribe: , List-Archive: List-Post: List-Help: List-Subscribe: , X-List-Received-Date: Sun, 02 Oct 2005 19:19:03 -0000 On Sat, Oct 01, 2005 at 09:38:44PM +1000, Bruce Evans wrote: > >A slightly better implementation than > >the above would be > > > > if (y == 0.) return (cos(x)); > > if (x == 0.) return (cosh(y)); > > return (cosh(y) * (cos(x) - I * sin(x) * tanh(y))); > > I'm not sure if this is better. The special cases may be so rare that > checking for them slows down the average case. I think special handling > is unnecessary since cos(0) is exactly 1.0, etc, so the checks are just > an optimization of the special args. OTOH, the type-generic version > should have special checks using something like > "if (__builtin_const_p(y) && y == 0.)" to classify the pure real case, > etc. Inline function args are never constant, so the type-generic > version must be implemented as a macro. Type-generic stuff is way above my C competency level. I'll leave that to bde, das, of stefanf. :-) > > >... > >OK, I'll use the implementation details of ccosh to infer the > >handling of exceptional values. > > At the risk of embarassing myself once again on freebsd-standards, here's an implementation of ccosh. If this is acceptable, I'll proceed with implementations for csinh, ctanh, ccos, csin, and ctan. One question: Can we use a #define for complex float or do we need an trivial C program? That is, in math.h can we do #define ccoshf(z) ccosh((float complex)(z)) -- Steve /*- * Copyright (c) 2005, Steven G. Kargl * All rights reserved. * * Redistribution and use in source and binary forms, with or without * modification, are permitted provided that the following conditions * are met: * 1. Redistributions of source code must retain the above copyright * notice unmodified, this list of conditions, and the following * disclaimer. * 2. Redistributions in binary form must reproduce the above copyright * notice, this list of conditions and the following disclaimer in the * documentation and/or other materials provided with the distribution. * * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. * IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT, * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. */ /* * Hyperbolic cosine of a complex argument z = x + i y such that i = sqrt(-1). * * cosh(z) = cosh(x+iy) * = cosh(x) if y == 0. * = cos(y) if x == 0. * = cosh(x) [cos(y) + i tanh(x) sin(y)] otherwise. * * Exception values are noted in the comments within the source code. * These values and the return result were taken from n869.pdf. */ #include #include #include "math_private.h" double complex ccosh(double complex z) { int con = 1; double x, y; int32_t hx, hy, ix, iy, lx, ly; x = creal(z); y = cimag(z); if (x < 0.) { /* Even function: ccosh(-z) = ccosh(z). */ x = -x; y = -y; /* Required identity: ccosh(conj(z)) = conj(ccosh(z)). */ if (y < 0.) { con = -1; y = -y; } } EXTRACT_WORDS(hx, lx, x); EXTRACT_WORDS(hy, ly, y); ix = 0x7fffffff&hx; /* if ix >= 0x7ff00000, then inf or NaN */ iy = 0x7fffffff&hy; /* if iy >= 0x7ff00000, then inf or NaN */ /* Filter out the non-exceptional cases. x and y are finite. */ if (ix < 0x7ff00000 && iy < 0x7ff00000) { /* cosh(+0 + I 0) = 1. */ if ((ix|lx) == 0 && (iy|ly) == 0) return 1.; return (cosh(x) * (cos(y) + I * con * tanh(x) * sin(y))); } /* * cosh(+0 + I Inf) = NaN +- I 0, sign of 0 is unspecified, * invalid exception. * cosh(+0 + I NaN) = NaN +- I 0, sign of 0 is unspecified. */ if ((ix|lx) == 0 && iy >= 0x7ff00000) return (y - y) / (y - y); /* * cosh(x + I Inf) = NaN + I NaN, finite x != 0, invalid exception. * cosh(x + I NaN) = NaN + I NaN, finite x != 0, opt. inval. except. */ if (ix < 0x7ff00000 && iy >= 0x7ff00000) return (x - x) / (x - x) + I * (x - x) / (x - x); /* * cosh(+Inf + I 0) = +Inf + I 0 * cosh(+Inf + I Inf) = +Inf + I NaN, invalid exception. * cosh(+Inf + I y) = +Inf [cos(y) + I sin(y)], y != 0 and finite. * cosh(+Inf + I NaN) = +Inf + I NaN. */ if (ix >= 0x7ff00000 && ((hx&0xfffff)|lx) == 0) { if ((iy|ly) == 0) return (x * x); else if (iy >= 0x7ff00000) return (x * x) + I * (y - y) / (y - y); return (x * x) * (cos(y) + I * con * sin(y)); } /* * cosh(NaN + I 0) = NaN +- I 0, sign of 0 unspecified. * cosh(NaN + I y) = NaN + I NaN, nonzero y, opt. inval. except. * cosh(NaN + I NaN) = NaN + I NaN. */ if (ix >= 0x7ff00000 && ((hx&0xfffff)|lx) != 0) { if ((iy|ly) == 0) return x * x; return x * x + I * x * x; } }