From 3daee1d6c3783e32a226d145f85a887007f5754b Mon Sep 17 00:00:00 2001 From: David Schultz Date: Mon, 17 Oct 2011 05:41:03 +0000 Subject: [PATCH] Add c{cos,sin,tan}{,h}{,f} functions. This is joint work with bde and kargl. --- include/complex.h | 14 +++- lib/msun/Makefile | 12 +++- lib/msun/Symbol.map | 12 ++++ lib/msun/man/ccos.3 | 80 +++++++++++++++++++++++ lib/msun/man/ccosh.3 | 80 +++++++++++++++++++++++ lib/msun/man/complex.3 | 22 ++++--- lib/msun/man/cos.3 | 1 + lib/msun/man/cosh.3 | 1 + lib/msun/man/sin.3 | 1 + lib/msun/man/sinh.3 | 1 + lib/msun/man/tan.3 | 1 + lib/msun/man/tanh.3 | 1 + lib/msun/src/s_ccosh.c | 138 +++++++++++++++++++++++++++++++++++++++ lib/msun/src/s_ccoshf.c | 87 +++++++++++++++++++++++++ lib/msun/src/s_csinh.c | 140 ++++++++++++++++++++++++++++++++++++++++ lib/msun/src/s_csinhf.c | 88 +++++++++++++++++++++++++ lib/msun/src/s_ctanh.c | 137 +++++++++++++++++++++++++++++++++++++++ lib/msun/src/s_ctanhf.c | 81 +++++++++++++++++++++++ 18 files changed, 884 insertions(+), 13 deletions(-) create mode 100644 lib/msun/man/ccos.3 create mode 100644 lib/msun/man/ccosh.3 create mode 100644 lib/msun/src/s_ccosh.c create mode 100644 lib/msun/src/s_ccoshf.c create mode 100644 lib/msun/src/s_csinh.c create mode 100644 lib/msun/src/s_csinhf.c create mode 100644 lib/msun/src/s_ctanh.c create mode 100644 lib/msun/src/s_ctanhf.c diff --git a/include/complex.h b/include/complex.h index f37b003e0d50..414230472d2b 100644 --- a/include/complex.h +++ b/include/complex.h @@ -1,5 +1,5 @@ /*- - * Copyright (c) 2001-2008 The FreeBSD Project. + * Copyright (c) 2001-2011 The FreeBSD Project. * All rights reserved. * * Redistribution and use in source and binary forms, with or without @@ -49,6 +49,10 @@ long double cabsl(long double complex); double carg(double complex); float cargf(float complex); long double cargl(long double complex); +double complex ccos(double complex); +float complex ccosf(float complex); +double complex ccosh(double complex); +float complex ccoshf(float complex); double complex cexp(double complex); float complex cexpf(float complex); double cimag(double complex) __pure2; @@ -65,10 +69,18 @@ long double complex double creal(double complex) __pure2; float crealf(float complex) __pure2; long double creall(long double complex) __pure2; +double complex csin(double complex); +float complex csinf(float complex); +double complex csinh(double complex); +float complex csinhf(float complex); double complex csqrt(double complex); float complex csqrtf(float complex); long double complex csqrtl(long double complex); +double complex ctan(double complex); +float complex ctanf(float complex); +double complex ctanh(double complex); +float complex ctanhf(float complex); __END_DECLS diff --git a/lib/msun/Makefile b/lib/msun/Makefile index bb99e45a25b2..7542f639457e 100644 --- a/lib/msun/Makefile +++ b/lib/msun/Makefile @@ -101,9 +101,11 @@ COMMON_SRCS+= e_acosl.c e_asinl.c e_atan2l.c e_fmodl.c \ .endif # C99 complex functions -COMMON_SRCS+= s_cexp.c s_cexpf.c s_cimag.c s_cimagf.c s_cimagl.c \ +COMMON_SRCS+= s_ccosh.c s_ccoshf.c s_cexp.c s_cexpf.c \ + s_cimag.c s_cimagf.c s_cimagl.c \ s_conj.c s_conjf.c s_conjl.c \ - s_cproj.c s_cprojf.c s_creal.c s_crealf.c s_creall.c + s_cproj.c s_cprojf.c s_creal.c s_crealf.c s_creall.c \ + s_csinh.c s_csinhf.c s_ctanh.c s_ctanhf.c # FreeBSD's C library supplies these functions: #COMMON_SRCS+= s_fabs.c s_frexp.c s_isnan.c s_ldexp.c s_modf.c @@ -125,7 +127,8 @@ SRCS= ${COMMON_SRCS} ${ARCH_SRCS} INCS= fenv.h math.h -MAN= acos.3 acosh.3 asin.3 asinh.3 atan.3 atan2.3 atanh.3 ceil.3 cexp.3 \ +MAN= acos.3 acosh.3 asin.3 asinh.3 atan.3 atan2.3 atanh.3 \ + ceil.3 ccos.3 ccosh.3 cexp.3 \ cimag.3 copysign.3 cos.3 cosh.3 csqrt.3 erf.3 exp.3 fabs.3 fdim.3 \ feclearexcept.3 feenableexcept.3 fegetenv.3 \ fegetround.3 fenv.3 floor.3 \ @@ -143,6 +146,9 @@ MLINKS+=atan.3 atanf.3 atan.3 atanl.3 MLINKS+=atanh.3 atanhf.3 MLINKS+=atan2.3 atan2f.3 atan2.3 atan2l.3 \ atan2.3 carg.3 atan2.3 cargf.3 atan2.3 cargl.3 +MLINKS+=ccos.3 ccosf.3 ccos.3 csin.3 ccos.3 csinf.3 ccos.3 ctan.3 ccos.3 ctanf.3 +MLINKS+=ccosh.3 ccoshf.3 ccosh.3 csinh.3 ccosh.3 csinhf.3 \ + ccosh.3 ctanh.3 ccosh.3 ctanhf.3 MLINKS+=ceil.3 ceilf.3 ceil.3 ceill.3 MLINKS+=cexp.3 cexpf.3 MLINKS+=cimag.3 cimagf.3 cimag.3 cimagl.3 \ diff --git a/lib/msun/Symbol.map b/lib/msun/Symbol.map index cd823ae73b4b..f3d30003b33a 100644 --- a/lib/msun/Symbol.map +++ b/lib/msun/Symbol.map @@ -237,4 +237,16 @@ FBSD_1.3 { fegetround; fesetround; fesetenv; + csin; + csinf; + csinh; + csinhf; + ccos; + ccosf; + ccosh; + ccoshf; + ctan; + ctanf; + ctanh; + ctanhf; }; diff --git a/lib/msun/man/ccos.3 b/lib/msun/man/ccos.3 new file mode 100644 index 000000000000..cf708c15f20f --- /dev/null +++ b/lib/msun/man/ccos.3 @@ -0,0 +1,80 @@ +.\" Copyright (c) 2011 David Schultz +.\" 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, 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 AND CONTRIBUTORS ``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 OR CONTRIBUTORS 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. +.\" +.\" $FreeBSD$ +.\" +.Dd October 17, 2011 +.Dt CCOS 3 +.Os +.Sh NAME +.Nm ccos , +.Nm ccosf , +.Nm csin , +.Nm csinf +.Nm ctan , +.Nm ctanf +.Nd complex trigonometric functions +.Sh LIBRARY +.Lb libm +.Sh SYNOPSIS +.In complex.h +.Ft double complex +.Fn ccos "double complex z" +.Ft float complex +.Fn ccosf "float complex z" +.Ft double complex +.Fn csin "double complex z" +.Ft float complex +.Fn csinf "float complex z" +.Ft double complex +.Fn ctan "double complex z" +.Ft float complex +.Fn ctanf "float complex z" +.Sh DESCRIPTION +The +.Fn ccos , +.Fn csin , +and +.Fn ctan +functions compute the cosine, sine, and tangent of the complex number +.Fa z , +respectively. +The +.Fn ccosf , +.Fn csinf , +and +.Fn ctanf +functions perform the same operations in +.Fa float +precision. +.Sh SEE ALSO +.Xr ccosh 3 , +.Xr complex 3 , +.Xr cos 3 , +.Xr math 3 , +.Xr sin 3 , +.Xr tan 3 +.Sh STANDARDS +These functions conform to +.St -isoC-99 . diff --git a/lib/msun/man/ccosh.3 b/lib/msun/man/ccosh.3 new file mode 100644 index 000000000000..01688b5bdb8c --- /dev/null +++ b/lib/msun/man/ccosh.3 @@ -0,0 +1,80 @@ +.\" Copyright (c) 2011 David Schultz +.\" 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, 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 AND CONTRIBUTORS ``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 OR CONTRIBUTORS 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. +.\" +.\" $FreeBSD$ +.\" +.Dd October 17, 2011 +.Dt CCOSH 3 +.Os +.Sh NAME +.Nm ccosh , +.Nm ccoshf , +.Nm csinh , +.Nm csinhf +.Nm ctanh , +.Nm ctanhf +.Nd complex hyperbolic functions +.Sh LIBRARY +.Lb libm +.Sh SYNOPSIS +.In complex.h +.Ft double complex +.Fn ccosh "double complex z" +.Ft float complex +.Fn ccoshf "float complex z" +.Ft double complex +.Fn csinh "double complex z" +.Ft float complex +.Fn csinhf "float complex z" +.Ft double complex +.Fn ctanh "double complex z" +.Ft float complex +.Fn ctanhf "float complex z" +.Sh DESCRIPTION +The +.Fn ccosh , +.Fn csinh , +and +.Fn ctanh +functions compute the hyperbolic cosine, sine, and tangent of the complex number +.Fa z , +respectively. +The +.Fn ccoshf , +.Fn csinhf , +and +.Fn ctanhf +functions perform the same operations in +.Fa float +precision. +.Sh SEE ALSO +.Xr ccos 3 , +.Xr complex 3 , +.Xr cosh 3 , +.Xr math 3 , +.Xr sinh 3 , +.Xr tanh 3 +.Sh STANDARDS +These functions conform to +.St -isoC-99 . diff --git a/lib/msun/man/complex.3 b/lib/msun/man/complex.3 index 3f7c18c6f77e..a8d9f950e4e5 100644 --- a/lib/msun/man/complex.3 +++ b/lib/msun/man/complex.3 @@ -24,7 +24,7 @@ .\" .\" $FreeBSD$ .\" -.Dd March 7, 2011 +.Dd October 17, 2011 .Dt COMPLEX 3 .Os .Sh NAME @@ -86,6 +86,16 @@ conj compute the complex conjugate cproj compute projection onto Riemann sphere creal compute the real part .El +.\" Section 7.3.5-6 of ISO C99 standard +.Ss Trigonometric and Hyperbolic Functions +.Cl +ccos cosine +ccosh hyperbolic cosine +csin sine +csinh hyperbolic sine +ctan tangent +ctanh hyperbolic tangent +.El .Sh SEE ALSO .Xr math 3 , .Xr fenv 3 , @@ -101,20 +111,14 @@ The functions described here conform to .St -isoC-99 . .Sh BUGS -The trigonmetric and hyperbolic functions +The inverse trigonmetric and hyperbolic functions .Fn cacos , .Fn cacosh , .Fn casin , .Fn casinh , .Fn catan , -.Fn catanh , -.Fn ccos , -.Fn ccosh , -.Fn csin , -.Fn csinh , -.Fn ctan , and -.Fn ctanh +.Fn catanh are not implemented. .Pp The logarithmic functions diff --git a/lib/msun/man/cos.3 b/lib/msun/man/cos.3 index d0cdf68bf4b4..a430f81e165c 100644 --- a/lib/msun/man/cos.3 +++ b/lib/msun/man/cos.3 @@ -71,6 +71,7 @@ functions return the cosine value. .Xr asin 3 , .Xr atan 3 , .Xr atan2 3 , +.Xr ccos 3 , .Xr cosh 3 , .Xr math 3 , .Xr sin 3 , diff --git a/lib/msun/man/cosh.3 b/lib/msun/man/cosh.3 index 517ec61395dc..96abc5aa7c9f 100644 --- a/lib/msun/man/cosh.3 +++ b/lib/msun/man/cosh.3 @@ -55,6 +55,7 @@ functions compute the hyperbolic cosine of .Xr asin 3 , .Xr atan 3 , .Xr atan2 3 , +.Xr ccosh 3 , .Xr cos 3 , .Xr math 3 , .Xr sin 3 , diff --git a/lib/msun/man/sin.3 b/lib/msun/man/sin.3 index 2a34281d150e..c7daf094dd08 100644 --- a/lib/msun/man/sin.3 +++ b/lib/msun/man/sin.3 @@ -70,6 +70,7 @@ functions return the sine value. .Xr asin 3 , .Xr atan 3 , .Xr atan2 3 , +.Xr csin 3 , .Xr cos 3 , .Xr cosh 3 , .Xr math 3 , diff --git a/lib/msun/man/sinh.3 b/lib/msun/man/sinh.3 index a2a2d8ce4439..02944ccc3a6a 100644 --- a/lib/msun/man/sinh.3 +++ b/lib/msun/man/sinh.3 @@ -56,6 +56,7 @@ functions compute the hyperbolic sine of .Xr atan2 3 , .Xr cos 3 , .Xr cosh 3 , +.Xr csinh 3 , .Xr math 3 , .Xr sin 3 , .Xr tan 3 , diff --git a/lib/msun/man/tan.3 b/lib/msun/man/tan.3 index 1577815ad961..f73c61208fa4 100644 --- a/lib/msun/man/tan.3 +++ b/lib/msun/man/tan.3 @@ -73,6 +73,7 @@ functions return the tangent value. .Xr atan2 3 , .Xr cos 3 , .Xr cosh 3 , +.Xr ctan 3 , .Xr math 3 , .Xr sin 3 , .Xr sinh 3 , diff --git a/lib/msun/man/tanh.3 b/lib/msun/man/tanh.3 index 699abe5d94dc..6fb185c418b6 100644 --- a/lib/msun/man/tanh.3 +++ b/lib/msun/man/tanh.3 @@ -65,6 +65,7 @@ functions return the hyperbolic tangent value. .Xr atan2 3 , .Xr cos 3 , .Xr cosh 3 , +.Xr ctanh 3 , .Xr math 3 , .Xr sin 3 , .Xr sinh 3 , diff --git a/lib/msun/src/s_ccosh.c b/lib/msun/src/s_ccosh.c new file mode 100644 index 000000000000..40863947fb3e --- /dev/null +++ b/lib/msun/src/s_ccosh.c @@ -0,0 +1,138 @@ +/*- + * Copyright (c) 2005 Bruce D. Evans and 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. + * + * cosh(z) = cosh(x+iy) + * = cosh(x) cos(y) + i sinh(x) sin(y). + * + * Exceptional values are noted in the comments within the source code. + * These values and the return value were taken from n1124.pdf. + */ + +#include +__FBSDID("$FreeBSD$"); + +#include +#include + +#include "math_private.h" + +double complex +ccosh(double complex z) +{ + double x, y; + int32_t hx, hy, ix, iy, lx, ly; + + x = creal(z); + y = cimag(z); + + EXTRACT_WORDS(hx, lx, x); + EXTRACT_WORDS(hy, ly, y); + + ix = 0x7fffffff & hx; + iy = 0x7fffffff & hy; + + /* Handle the nearly-non-exceptional cases where x and y are finite. */ + if (ix < 0x7ff00000 && iy < 0x7ff00000) { + if ((iy | ly) == 0) + return (cpack(cosh(x), x * y)); + /* XXX We don't handle |x| > DBL_MAX ln(2) yet. */ + return (cpack(cosh(x) * cos(y), sinh(x) * sin(y))); + } + + /* + * cosh(+-0 +- I Inf) = dNaN + I sign(d(+-0, dNaN))0. + * The sign of 0 in the result is unspecified. Choice = normally + * the same as dNaN. Raise the invalid floating-point exception. + * + * cosh(+-0 +- I NaN) = d(NaN) + I sign(d(+-0, NaN))0. + * The sign of 0 in the result is unspecified. Choice = normally + * the same as d(NaN). + */ + if ((ix | lx) == 0 && iy >= 0x7ff00000) + return (cpack(y - y, copysign(0, x * (y - y)))); + + /* + * cosh(+-Inf +- I 0) = +Inf + I (+-)(+-)0. + * + * cosh(NaN +- I 0) = d(NaN) + I sign(d(NaN, +-0))0. + * The sign of 0 in the result is unspecified. + */ + if ((iy | ly) == 0 && ix >= 0x7ff00000) { + if (((hx & 0xfffff) | lx) == 0) + return (cpack(x * x, copysign(0, x) * y)); + return (cpack(x * x, copysign(0, (x + x) * y))); + } + + /* + * cosh(x +- I Inf) = dNaN + I dNaN. + * Raise the invalid floating-point exception for finite nonzero x. + * + * cosh(x + I NaN) = d(NaN) + I d(NaN). + * Optionally raises the invalid floating-point exception for finite + * nonzero x. Choice = don't raise (except for signaling NaNs). + */ + if (ix < 0x7ff00000 && iy >= 0x7ff00000) + return (cpack(y - y, x * (y - y))); + + /* + * cosh(+-Inf + I NaN) = +Inf + I d(NaN). + * + * cosh(+-Inf +- I Inf) = +Inf + I dNaN. + * The sign of Inf in the result is unspecified. Choice = always +. + * Raise the invalid floating-point exception. + * + * cosh(+-Inf + I y) = +Inf cos(y) +- I Inf sin(y) + */ + if (ix >= 0x7ff00000 && ((hx & 0xfffff) | lx) == 0) { + if (iy >= 0x7ff00000) + return (cpack(x * x, x * (y - y))); + return (cpack((x * x) * cos(y), x * sin(y))); + } + + /* + * cosh(NaN + I NaN) = d(NaN) + I d(NaN). + * + * cosh(NaN +- I Inf) = d(NaN) + I d(NaN). + * Optionally raises the invalid floating-point exception. + * Choice = raise. + * + * cosh(NaN + I y) = d(NaN) + I d(NaN). + * Optionally raises the invalid floating-point exception for finite + * nonzero y. Choice = don't raise (except for signaling NaNs). + */ + return (cpack((x * x) * (y - y), (x + x) * (y - y))); +} + +double complex +ccos(double complex z) +{ + + /* ccos(z) = ccosh(I * z) */ + return (ccosh(cpack(-cimag(z), creal(z)))); +} diff --git a/lib/msun/src/s_ccoshf.c b/lib/msun/src/s_ccoshf.c new file mode 100644 index 000000000000..bf9b78ffe022 --- /dev/null +++ b/lib/msun/src/s_ccoshf.c @@ -0,0 +1,87 @@ +/*- + * Copyright (c) 2005 Bruce D. Evans and 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. See s_ccosh.c for details. + */ + +#include +__FBSDID("$FreeBSD$"); + +#include +#include + +#include "math_private.h" + +float complex +ccoshf(float complex z) +{ + float x, y; + int32_t hx, hy, ix, iy; + + x = crealf(z); + y = cimagf(z); + + GET_FLOAT_WORD(hx, x); + GET_FLOAT_WORD(hy, y); + + ix = 0x7fffffff & hx; + iy = 0x7fffffff & hy; + + if (ix < 0x7f800000 && iy < 0x7f800000) { + if (iy == 0) + return (cpackf(coshf(x), x * y)); + /* XXX We don't handle |x| > FLT_MAX ln(2) yet. */ + return (cpackf(coshf(x) * cosf(y), sinhf(x) * sinf(y))); + } + + if (ix == 0 && iy >= 0x7f800000) + return (cpackf(y - y, copysignf(0, x * (y - y)))); + + if (iy == 0 && ix >= 0x7f800000) { + if ((hx & 0x7fffff) == 0) + return (cpackf(x * x, copysignf(0, x) * y)); + return (cpackf(x * x, copysignf(0, (x + x) * y))); + } + + if (ix < 0x7f800000 && iy >= 0x7f800000) + return (cpackf(y - y, x * (y - y))); + + if (ix >= 0x7f800000 && (hx & 0x7fffff) == 0) { + if (iy >= 0x7f800000) + return (cpackf(x * x, x * (y - y))); + return (cpackf((x * x) * cosf(y), x * sinf(y))); + } + + return (cpackf((x * x) * (y - y), (x + x) * (y - y))); +} + +float complex +ccosf(float complex z) +{ + + return (ccoshf(cpackf(-cimagf(z), crealf(z)))); +} diff --git a/lib/msun/src/s_csinh.c b/lib/msun/src/s_csinh.c new file mode 100644 index 000000000000..74a0aff0abbd --- /dev/null +++ b/lib/msun/src/s_csinh.c @@ -0,0 +1,140 @@ +/*- + * Copyright (c) 2005 Bruce D. Evans and 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 sine of a complex argument z = x + i y. + * + * sinh(z) = sinh(x+iy) + * = sinh(x) cos(y) + i cosh(x) sin(y). + * + * Exceptional values are noted in the comments within the source code. + * These values and the return value were taken from n1124.pdf. + */ + +#include +__FBSDID("$FreeBSD$"); + +#include +#include + +#include "math_private.h" + +double complex +csinh(double complex z) +{ + double x, y; + int32_t hx, hy, ix, iy, lx, ly; + + x = creal(z); + y = cimag(z); + + EXTRACT_WORDS(hx, lx, x); + EXTRACT_WORDS(hy, ly, y); + + ix = 0x7fffffff & hx; + iy = 0x7fffffff & hy; + + /* Handle the nearly-non-exceptional cases where x and y are finite. */ + if (ix < 0x7ff00000 && iy < 0x7ff00000) { + if ((iy | ly) == 0) + return (cpack(sinh(x), y)); + /* XXX We don't handle |x| > DBL_MAX ln(2) yet. */ + return (cpack(sinh(x) * cos(y), cosh(x) * sin(y))); + } + + /* + * sinh(+-0 +- I Inf) = sign(d(+-0, dNaN))0 + I dNaN. + * The sign of 0 in the result is unspecified. Choice = normally + * the same as dNaN. Raise the invalid floating-point exception. + * + * sinh(+-0 +- I NaN) = sign(d(+-0, NaN))0 + I d(NaN). + * The sign of 0 in the result is unspecified. Choice = normally + * the same as d(NaN). + */ + if ((ix | lx) == 0 && iy >= 0x7ff00000) + return (cpack(copysign(0, x * (y - y)), y - y)); + + /* + * sinh(+-Inf +- I 0) = +-Inf + I (+-)(+-)0. + * + * sinh(NaN +- I 0) = d(NaN) + I +-0. + */ + if ((iy | ly) == 0 && ix >= 0x7ff00000) { + if (((hx & 0xfffff) | lx) == 0) + return (cpack(x, copysign(0, x) * y)); + return (cpack(x, copysign(0, y))); + } + + /* + * sinh(x +- I Inf) = dNaN + I dNaN. + * Raise the invalid floating-point exception for finite nonzero x. + * + * sinh(x + I NaN) = d(NaN) + I d(NaN). + * Optionally raises the invalid floating-point exception for finite + * nonzero x. Choice = don't raise (except for signaling NaNs). + */ + if (ix < 0x7ff00000 && iy >= 0x7ff00000) + return (cpack(y - y, x * (y - y))); + + /* + * sinh(+-Inf + I NaN) = +-Inf + I d(NaN). + * The sign of Inf in the result is unspecified. Choice = normally + * the same as d(NaN). + * + * sinh(+-Inf +- I Inf) = +Inf + I dNaN. + * The sign of Inf in the result is unspecified. Choice = always +. + * Raise the invalid floating-point exception. + * + * sinh(+-Inf + I y) = +-Inf cos(y) + I Inf sin(y) + */ + if (ix >= 0x7ff00000 && ((hx & 0xfffff) | lx) == 0) { + if (iy >= 0x7ff00000) + return (cpack(x * x, x * (y - y))); + return (cpack(x * cos(y), INFINITY * sin(y))); + } + + /* + * sinh(NaN + I NaN) = d(NaN) + I d(NaN). + * + * sinh(NaN +- I Inf) = d(NaN) + I d(NaN). + * Optionally raises the invalid floating-point exception. + * Choice = raise. + * + * sinh(NaN + I y) = d(NaN) + I d(NaN). + * Optionally raises the invalid floating-point exception for finite + * nonzero y. Choice = don't raise (except for signaling NaNs). + */ + return (cpack((x * x) * (y - y), (x + x) * (y - y))); +} + +double complex +csin(double complex z) +{ + + /* csin(z) = -I * csinh(I * z) */ + z = csinh(cpack(-cimag(z), creal(z))); + return (cpack(cimag(z), -creal(z))); +} diff --git a/lib/msun/src/s_csinhf.c b/lib/msun/src/s_csinhf.c new file mode 100644 index 000000000000..e1478acac809 --- /dev/null +++ b/lib/msun/src/s_csinhf.c @@ -0,0 +1,88 @@ +/*- + * Copyright (c) 2005 Bruce D. Evans and 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 sine of a complex argument z. See s_csinh.c for details. + */ + +#include +__FBSDID("$FreeBSD$"); + +#include +#include + +#include "math_private.h" + +float complex +csinhf(float complex z) +{ + float x, y; + int32_t hx, hy, ix, iy; + + x = crealf(z); + y = cimagf(z); + + GET_FLOAT_WORD(hx, x); + GET_FLOAT_WORD(hy, y); + + ix = 0x7fffffff & hx; + iy = 0x7fffffff & hy; + + if (ix < 0x7f800000 && iy < 0x7f800000) { + if (iy == 0) + return (cpackf(sinhf(x), y)); + /* XXX We don't handle |x| > FLT_MAX ln(2) yet. */ + return (cpackf(sinhf(x) * cosf(y), coshf(x) * sinf(y))); + } + + if (ix == 0 && iy >= 0x7f800000) + return (cpackf(copysignf(0, x * (y - y)), y - y)); + + if (iy == 0 && ix >= 0x7f800000) { + if ((hx & 0x7fffff) == 0) + return (cpackf(x, copysignf(0, x) * y)); + return (cpackf(x, copysignf(0, y))); + } + + if (ix < 0x7f800000 && iy >= 0x7f800000) + return (cpackf(y - y, x * (y - y))); + + if (ix >= 0x7f800000 && (hx & 0x7fffff) == 0) { + if (iy >= 0x7f800000) + return (cpackf(x * x, x * (y - y))); + return (cpackf(x * cosf(y), INFINITY * sinf(y))); + } + + return (cpackf((x * x) * (y - y), (x + x) * (y - y))); +} + +float complex +csinf(float complex z) +{ + + z = csinhf(cpackf(-cimagf(z), crealf(z))); + return (cpackf(cimagf(z), -crealf(z))); +} diff --git a/lib/msun/src/s_ctanh.c b/lib/msun/src/s_ctanh.c new file mode 100644 index 000000000000..1705e909e4dd --- /dev/null +++ b/lib/msun/src/s_ctanh.c @@ -0,0 +1,137 @@ +/*- + * Copyright (c) 2011 David Schultz + * 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 tangent of a complex argument z = x + i y. + * + * The algorithm is from: + * + * W. Kahan. Branch Cuts for Complex Elementary Functions or Much + * Ado About Nothing's Sign Bit. In The State of the Art in + * Numerical Analysis, pp. 165 ff. Iserles and Powell, eds., 1987. + * + * Method: + * + * Let t = tan(x) + * beta = 1/cos^2(y) + * s = sinh(x) + * rho = cosh(x) + * + * We have: + * + * tanh(z) = sinh(z) / cosh(z) + * + * sinh(x) cos(y) + i cosh(x) sin(y) + * = --------------------------------- + * cosh(x) cos(y) + i sinh(x) sin(y) + * + * cosh(x) sinh(x) / cos^2(y) + i tan(y) + * = ------------------------------------- + * 1 + sinh^2(x) / cos^2(y) + * + * beta rho s + i t + * = ---------------- + * 1 + beta s^2 + * + * Modifications: + * + * I omitted the original algorithm's handling of overflow in tan(x) after + * verifying with nearpi.c that this can't happen in IEEE single or double + * precision. I also handle large x differently. + */ + +#include +__FBSDID("$FreeBSD$"); + +#include +#include + +#include "math_private.h" + +double complex +ctanh(double complex z) +{ + double x, y; + double t, beta, s, rho, denom; + uint32_t hx, ix, lx; + + x = creal(z); + y = cimag(z); + + EXTRACT_WORDS(hx, lx, x); + ix = hx & 0x7fffffff; + + /* + * ctanh(NaN + i 0) = NaN + i 0 + * + * ctanh(NaN + i y) = NaN + i NaN for y != 0 + * + * The imaginary part has the sign of x*sin(2*y), but there's no + * special effort to get this right. + * + * ctanh(+-Inf +- i Inf) = +-1 +- 0 + * + * ctanh(+-Inf + i y) = +-1 + 0 sin(2y) for y finite + * + * The imaginary part of the sign is unspecified. This special + * case is only needed to avoid a spurious invalid exception when + * y is infinite. + */ + if (ix >= 0x7ff00000) { + if ((ix & 0xfffff) | lx) /* x is NaN */ + return (cpack(x, (y == 0 ? y : x * y))); + SET_HIGH_WORD(x, hx - 0x40000000); /* x = copysign(1, x) */ + return (cpack(x, copysign(0, isinf(y) ? y : sin(y) * cos(y)))); + } + + /* + * ctanh(+-huge + i +-y) ~= +-1 +- i 2sin(2y)/exp(2x), using the + * approximation sinh^2(huge) ~= exp(2*huge) / 4. + * We use a modified formula to avoid spurious overflow. + */ + if (ix >= 0x40360000) { /* x >= 22 */ + double exp_mx = exp(-fabs(x)); + return (cpack(copysign(1, x), + 4 * sin(y) * cos(y) * exp_mx * exp_mx)); + } + + /* Kahan's algorithm */ + t = tan(y); + beta = 1.0 + t * t; /* = 1 / cos^2(y) */ + s = sinh(x); + rho = sqrt(1 + s * s); /* = cosh(x) */ + denom = 1 + beta * s * s; + return (cpack((beta * rho * s) / denom, t / denom)); +} + +double complex +ctan(double complex z) +{ + + /* ctan(z) = -I * ctanh(I * z) */ + z = ctanh(cpack(-cimag(z), creal(z))); + return (cpack(cimag(z), -creal(z))); +} diff --git a/lib/msun/src/s_ctanhf.c b/lib/msun/src/s_ctanhf.c new file mode 100644 index 000000000000..ae833dfa952d --- /dev/null +++ b/lib/msun/src/s_ctanhf.c @@ -0,0 +1,81 @@ +/*- + * Copyright (c) 2011 David Schultz + * 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 tangent of a complex argument z. See s_ctanh.c for details. + */ + +#include +__FBSDID("$FreeBSD$"); + +#include +#include + +#include "math_private.h" + +float complex +ctanhf(float complex z) +{ + float x, y; + float t, beta, s, rho, denom; + uint32_t hx, ix; + + x = crealf(z); + y = cimagf(z); + + GET_FLOAT_WORD(hx, x); + ix = hx & 0x7fffffff; + + if (ix >= 0x7f800000) { + if (ix & 0x7fffff) + return (cpackf(x, (y == 0 ? y : x * y))); + SET_FLOAT_WORD(x, hx - 0x40000000); + return (cpackf(x, + copysignf(0, isinf(y) ? y : sinf(y) * cosf(y)))); + } + + if (ix >= 0x41300000) { /* x >= 11 */ + float exp_mx = expf(-fabsf(x)); + return (cpackf(copysignf(1, x), + 4 * sinf(y) * cosf(y) * exp_mx * exp_mx)); + } + + t = tanf(y); + beta = 1.0 + t * t; + s = sinhf(x); + rho = sqrtf(1 + s * s); + denom = 1 + beta * s * s; + return (cpackf((beta * rho * s) / denom, t / denom)); +} + +float complex +ctanf(float complex z) +{ + + z = ctanhf(cpackf(-cimagf(z), crealf(z))); + return (cpackf(cimagf(z), -crealf(z))); +} +