125c28e83SPiotr Jasiukajtis /*
225c28e83SPiotr Jasiukajtis * CDDL HEADER START
325c28e83SPiotr Jasiukajtis *
425c28e83SPiotr Jasiukajtis * The contents of this file are subject to the terms of the
525c28e83SPiotr Jasiukajtis * Common Development and Distribution License (the "License").
625c28e83SPiotr Jasiukajtis * You may not use this file except in compliance with the License.
725c28e83SPiotr Jasiukajtis *
825c28e83SPiotr Jasiukajtis * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
925c28e83SPiotr Jasiukajtis * or http://www.opensolaris.org/os/licensing.
1025c28e83SPiotr Jasiukajtis * See the License for the specific language governing permissions
1125c28e83SPiotr Jasiukajtis * and limitations under the License.
1225c28e83SPiotr Jasiukajtis *
1325c28e83SPiotr Jasiukajtis * When distributing Covered Code, include this CDDL HEADER in each
1425c28e83SPiotr Jasiukajtis * file and include the License file at usr/src/OPENSOLARIS.LICENSE.
1525c28e83SPiotr Jasiukajtis * If applicable, add the following below this CDDL HEADER, with the
1625c28e83SPiotr Jasiukajtis * fields enclosed by brackets "[]" replaced with your own identifying
1725c28e83SPiotr Jasiukajtis * information: Portions Copyright [yyyy] [name of copyright owner]
1825c28e83SPiotr Jasiukajtis *
1925c28e83SPiotr Jasiukajtis * CDDL HEADER END
2025c28e83SPiotr Jasiukajtis */
2125c28e83SPiotr Jasiukajtis
2225c28e83SPiotr Jasiukajtis /*
2325c28e83SPiotr Jasiukajtis * Copyright 2011 Nexenta Systems, Inc. All rights reserved.
2425c28e83SPiotr Jasiukajtis */
2525c28e83SPiotr Jasiukajtis /*
2625c28e83SPiotr Jasiukajtis * Copyright 2006 Sun Microsystems, Inc. All rights reserved.
2725c28e83SPiotr Jasiukajtis * Use is subject to license terms.
2825c28e83SPiotr Jasiukajtis */
2925c28e83SPiotr Jasiukajtis
30*ddc0e0b5SRichard Lowe #pragma weak __ctanhl = ctanhl
3125c28e83SPiotr Jasiukajtis
3225c28e83SPiotr Jasiukajtis #include "libm.h" /* expl/expm1l/fabsl/isinfl/isnanl/sincosl/sinl/tanhl */
3325c28e83SPiotr Jasiukajtis #include "complex_wrapper.h"
3425c28e83SPiotr Jasiukajtis #include "longdouble.h"
3525c28e83SPiotr Jasiukajtis
3625c28e83SPiotr Jasiukajtis /* INDENT OFF */
3725c28e83SPiotr Jasiukajtis static const long double four = 4.0L, two = 2.0L, one = 1.0L, zero = 0.0L;
3825c28e83SPiotr Jasiukajtis /* INDENT ON */
3925c28e83SPiotr Jasiukajtis
4025c28e83SPiotr Jasiukajtis ldcomplex
ctanhl(ldcomplex z)4125c28e83SPiotr Jasiukajtis ctanhl(ldcomplex z) {
4225c28e83SPiotr Jasiukajtis long double r, u, v, t, x, y, S, C;
4325c28e83SPiotr Jasiukajtis int hx, ix, hy, iy;
4425c28e83SPiotr Jasiukajtis ldcomplex ans;
4525c28e83SPiotr Jasiukajtis
4625c28e83SPiotr Jasiukajtis x = LD_RE(z);
4725c28e83SPiotr Jasiukajtis y = LD_IM(z);
4825c28e83SPiotr Jasiukajtis hx = HI_XWORD(x);
4925c28e83SPiotr Jasiukajtis ix = hx & 0x7fffffff;
5025c28e83SPiotr Jasiukajtis hy = HI_XWORD(y);
5125c28e83SPiotr Jasiukajtis iy = hy & 0x7fffffff;
5225c28e83SPiotr Jasiukajtis x = fabsl(x);
5325c28e83SPiotr Jasiukajtis y = fabsl(y);
5425c28e83SPiotr Jasiukajtis
5525c28e83SPiotr Jasiukajtis if (y == zero) { /* ctanh(x,0) = (x,0) for x = 0 or NaN */
5625c28e83SPiotr Jasiukajtis LD_RE(ans) = tanhl(x);
5725c28e83SPiotr Jasiukajtis LD_IM(ans) = zero;
5825c28e83SPiotr Jasiukajtis } else if (iy >= 0x7fff0000) { /* y is inf or NaN */
5925c28e83SPiotr Jasiukajtis if (ix < 0x7fff0000) /* catanh(finite x,inf/nan) is nan */
6025c28e83SPiotr Jasiukajtis LD_RE(ans) = LD_IM(ans) = y - y;
6125c28e83SPiotr Jasiukajtis else if (isinfl(x)) { /* x is inf */
6225c28e83SPiotr Jasiukajtis LD_RE(ans) = one;
6325c28e83SPiotr Jasiukajtis LD_IM(ans) = zero;
6425c28e83SPiotr Jasiukajtis } else {
6525c28e83SPiotr Jasiukajtis LD_RE(ans) = x + y;
6625c28e83SPiotr Jasiukajtis LD_IM(ans) = y - y;
6725c28e83SPiotr Jasiukajtis }
6825c28e83SPiotr Jasiukajtis } else if (ix >= 0x4004e000) {
6925c28e83SPiotr Jasiukajtis /* INDENT OFF */
7025c28e83SPiotr Jasiukajtis /*
7125c28e83SPiotr Jasiukajtis * |x| > 60 = prec/2 (14,28,34,60)
7225c28e83SPiotr Jasiukajtis * ctanh z ~ 1 + i (sin2y)/(exp(2x))
7325c28e83SPiotr Jasiukajtis */
7425c28e83SPiotr Jasiukajtis /* INDENT ON */
7525c28e83SPiotr Jasiukajtis LD_RE(ans) = one;
7625c28e83SPiotr Jasiukajtis if (iy < 0x7ffe0000) /* t = sin(2y) */
7725c28e83SPiotr Jasiukajtis S = sinl(y + y);
7825c28e83SPiotr Jasiukajtis else {
7925c28e83SPiotr Jasiukajtis (void) sincosl(y, &S, &C);
8025c28e83SPiotr Jasiukajtis S = (S + S) * C;
8125c28e83SPiotr Jasiukajtis }
8225c28e83SPiotr Jasiukajtis if (ix >= 0x7ffe0000) { /* |x| > max/2 */
8325c28e83SPiotr Jasiukajtis if (ix >= 0x7fff0000) { /* |x| is inf or NaN */
8425c28e83SPiotr Jasiukajtis if (isnanl(x)) /* x is NaN */
8525c28e83SPiotr Jasiukajtis LD_RE(ans) = LD_IM(ans) = x + y;
8625c28e83SPiotr Jasiukajtis else
8725c28e83SPiotr Jasiukajtis LD_IM(ans) = zero * S; /* x is inf */
8825c28e83SPiotr Jasiukajtis } else
8925c28e83SPiotr Jasiukajtis LD_IM(ans) = S * expl(-x); /* underflow */
9025c28e83SPiotr Jasiukajtis } else
9125c28e83SPiotr Jasiukajtis LD_IM(ans) = (S + S) * expl(-(x + x));
9225c28e83SPiotr Jasiukajtis /* 2 sin 2y / exp(2x) */
9325c28e83SPiotr Jasiukajtis } else {
9425c28e83SPiotr Jasiukajtis /* INDENT OFF */
9525c28e83SPiotr Jasiukajtis /*
9625c28e83SPiotr Jasiukajtis * t*t+2t
9725c28e83SPiotr Jasiukajtis * ctanh z = ---------------------------
9825c28e83SPiotr Jasiukajtis * t*t+[4(t+1)(cos y)](cos y)
9925c28e83SPiotr Jasiukajtis *
10025c28e83SPiotr Jasiukajtis * [4(t+1)(cos y)]*(sin y)
10125c28e83SPiotr Jasiukajtis * i --------------------------
10225c28e83SPiotr Jasiukajtis * t*t+[4(t+1)(cos y)](cos y)
10325c28e83SPiotr Jasiukajtis */
10425c28e83SPiotr Jasiukajtis /* INDENT ON */
10525c28e83SPiotr Jasiukajtis sincosl(y, &S, &C);
10625c28e83SPiotr Jasiukajtis t = expm1l(x + x);
10725c28e83SPiotr Jasiukajtis r = (four * C) * (t + one);
10825c28e83SPiotr Jasiukajtis u = t * t;
10925c28e83SPiotr Jasiukajtis v = one / (u + r * C);
11025c28e83SPiotr Jasiukajtis LD_RE(ans) = (u + two * t) * v;
11125c28e83SPiotr Jasiukajtis LD_IM(ans) = (r * S) * v;
11225c28e83SPiotr Jasiukajtis }
11325c28e83SPiotr Jasiukajtis if (hx < 0)
11425c28e83SPiotr Jasiukajtis LD_RE(ans) = -LD_RE(ans);
11525c28e83SPiotr Jasiukajtis if (hy < 0)
11625c28e83SPiotr Jasiukajtis LD_IM(ans) = -LD_IM(ans);
11725c28e83SPiotr Jasiukajtis return (ans);
11825c28e83SPiotr Jasiukajtis }
119