1 /* 2 * CDDL HEADER START 3 * 4 * The contents of this file are subject to the terms of the 5 * Common Development and Distribution License (the "License"). 6 * You may not use this file except in compliance with the License. 7 * 8 * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE 9 * or http://www.opensolaris.org/os/licensing. 10 * See the License for the specific language governing permissions 11 * and limitations under the License. 12 * 13 * When distributing Covered Code, include this CDDL HEADER in each 14 * file and include the License file at usr/src/OPENSOLARIS.LICENSE. 15 * If applicable, add the following below this CDDL HEADER, with the 16 * fields enclosed by brackets "[]" replaced with your own identifying 17 * information: Portions Copyright [yyyy] [name of copyright owner] 18 * 19 * CDDL HEADER END 20 */ 21 /* 22 * Copyright 2011 Nexenta Systems, Inc. All rights reserved. 23 */ 24 /* 25 * Copyright 2005 Sun Microsystems, Inc. All rights reserved. 26 * Use is subject to license terms. 27 */ 28 29 #pragma weak __csinhf = csinhf 30 31 #include "libm.h" 32 #include "complex_wrapper.h" 33 34 #if defined(__i386) && !defined(__amd64) 35 extern int __swapRP(int); 36 #endif 37 38 static const float zero = 0.0F, half = 0.5F; 39 40 fcomplex 41 csinhf(fcomplex z) { 42 float x, y, S, C; 43 double t; 44 int hx, ix, hy, iy, n; 45 fcomplex ans; 46 47 x = F_RE(z); 48 y = F_IM(z); 49 hx = THE_WORD(x); 50 ix = hx & 0x7fffffff; 51 hy = THE_WORD(y); 52 iy = hy & 0x7fffffff; 53 x = fabsf(x); 54 y = fabsf(y); 55 56 sincosf(y, &S, &C); 57 if (ix >= 0x41600000) { /* |x| > 14 = prec/2 (14,28,34,60) */ 58 if (ix >= 0x42B171AA) { /* |x| > 88.722... ~ log(2**128) */ 59 if (ix >= 0x7f800000) { /* |x| is inf or NaN */ 60 if (iy == 0) { 61 F_RE(ans) = x; 62 F_IM(ans) = y; 63 } else if (iy >= 0x7f800000) { 64 F_RE(ans) = x; 65 F_IM(ans) = x - y; 66 } else { 67 F_RE(ans) = C * x; 68 F_IM(ans) = S * x; 69 } 70 } else { 71 #if defined(__i386) && !defined(__amd64) 72 int rp = __swapRP(fp_extended); 73 #endif 74 /* return (C, S) * exp(x) / 2 */ 75 t = __k_cexp((double)x, &n); 76 F_RE(ans) = (float)scalbn(C * t, n - 1); 77 F_IM(ans) = (float)scalbn(S * t, n - 1); 78 #if defined(__i386) && !defined(__amd64) 79 if (rp != fp_extended) 80 (void) __swapRP(rp); 81 #endif 82 } 83 } else { 84 t = expf(x) * half; 85 F_RE(ans) = C * t; 86 F_IM(ans) = S * t; 87 } 88 } else { 89 if (ix == 0) { /* x = 0, return (0,S) */ 90 F_RE(ans) = zero; 91 F_IM(ans) = S; 92 } else { 93 F_RE(ans) = C * sinhf(x); 94 F_IM(ans) = S * coshf(x); 95 } 96 } 97 if (hx < 0) 98 F_RE(ans) = -F_RE(ans); 99 if (hy < 0) 100 F_IM(ans) = -F_IM(ans); 101 return (ans); 102 } 103