X-Git-Url: http://nsz.repo.hu/git/?a=blobdiff_plain;f=src%2Fmath%2Fcosl.c;h=79c41c77fa6a0addc57e33c62f31d879392b36e6;hb=a7dbcf5c8ca7edb7a46eb276031ad1df4131135b;hp=25d9005ac5b98545ab44f5e2bcb3c65c51dc4117;hpb=2e8c8fbe7d65ba0026cb084dc8570d94cbc908ff;p=musl diff --git a/src/math/cosl.c b/src/math/cosl.c index 25d9005a..79c41c77 100644 --- a/src/math/cosl.c +++ b/src/math/cosl.c @@ -1,34 +1,3 @@ -/* origin: FreeBSD /usr/src/lib/msun/src/s_cosl.c */ -/*- - * Copyright (c) 2007 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. - */ -/* - * Limited testing on pseudorandom numbers drawn within [-2e8:4e8] shows - * an accuracy of <= 0.7412 ULP. - */ - #include "libm.h" #if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024 @@ -38,44 +7,33 @@ long double cosl(long double x) { #elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384 long double cosl(long double x) { - union IEEEl2bits z; - int e0; - long double y[2]; - long double hi, lo; - - z.e = x; - z.bits.sign = 0; - - /* If x = +-0 or x is a subnormal number, then cos(x) = 1 */ - if (z.bits.exp == 0) - return 1.0; - - /* If x = NaN or Inf, then cos(x) = NaN. */ - if (z.bits.exp == 32767) - return (x - x) / (x - x); - - /* Optimize the case where x is already within range. */ - if (z.e < M_PI_4) - return __cosl(z.e, 0); - - e0 = __rem_pio2l(x, y); + union ldshape u = {x}; + unsigned n; + long double y[2], hi, lo; + + u.i.se &= 0x7fff; + if (u.i.se == 0x7fff) + return x - x; + x = u.f; + if (x < M_PI_4) { + if (u.i.se < 0x3fff - LDBL_MANT_DIG) + /* raise inexact if x!=0 */ + return 1.0 + x; + return __cosl(x, 0); + } + n = __rem_pio2l(x, y); hi = y[0]; lo = y[1]; - - switch (e0 & 3) { + switch (n & 3) { case 0: - hi = __cosl(hi, lo); - break; + return __cosl(hi, lo); case 1: - hi = -__sinl(hi, lo, 1); - break; + return -__sinl(hi, lo, 1); case 2: - hi = -__cosl(hi, lo); - break; + return -__cosl(hi, lo); case 3: - hi = __sinl(hi, lo, 1); - break; + default: + return __sinl(hi, lo, 1); } - return hi; } #endif