diff options
| author | joseph <joseph@7b3dc134-2b1b-0410-93df-9e9f96275f8d> | 2014-01-03 17:51:28 +0000 |
|---|---|---|
| committer | joseph <joseph@7b3dc134-2b1b-0410-93df-9e9f96275f8d> | 2014-01-03 17:51:28 +0000 |
| commit | f71e9a28759c1f5d6e8dc0baced440e72507b433 (patch) | |
| tree | d1ef7f1f84de67d0446b60c04f10a47c84ac1ddc /libc/sysdeps/ieee754 | |
| parent | 036dfd44bcc8ed8e97dcd9c763f61d3b59796941 (diff) | |
| download | eglibc2-f71e9a28759c1f5d6e8dc0baced440e72507b433.tar.gz | |
Merge changes between r24468 and r24941 from /fsf/trunk.
git-svn-id: svn://svn.eglibc.org/trunk@24942 7b3dc134-2b1b-0410-93df-9e9f96275f8d
Diffstat (limited to 'libc/sysdeps/ieee754')
242 files changed, 696 insertions, 1756 deletions
diff --git a/libc/sysdeps/ieee754/bits/huge_val.h b/libc/sysdeps/ieee754/bits/huge_val.h index 7dcdefe61..111922a28 100644 --- a/libc/sysdeps/ieee754/bits/huge_val.h +++ b/libc/sysdeps/ieee754/bits/huge_val.h @@ -1,6 +1,6 @@ /* `HUGE_VAL' constant for IEEE 754 machines (where it is infinity). Used by <stdlib.h> and <math.h> functions for overflow. - Copyright (C) 1992-2013 Free Software Foundation, Inc. + Copyright (C) 1992-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/bits/huge_valf.h b/libc/sysdeps/ieee754/bits/huge_valf.h index 7974820ac..a8025a903 100644 --- a/libc/sysdeps/ieee754/bits/huge_valf.h +++ b/libc/sysdeps/ieee754/bits/huge_valf.h @@ -1,6 +1,6 @@ /* `HUGE_VALF' constant for IEEE 754 machines (where it is infinity). Used by <stdlib.h> and <math.h> functions for overflow. - Copyright (C) 1992-2013 Free Software Foundation, Inc. + Copyright (C) 1992-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/bits/inf.h b/libc/sysdeps/ieee754/bits/inf.h index cf66bf0fe..07d080eb3 100644 --- a/libc/sysdeps/ieee754/bits/inf.h +++ b/libc/sysdeps/ieee754/bits/inf.h @@ -1,5 +1,5 @@ /* `INFINITY' constant for IEEE 754 machines. - Copyright (C) 2004-2013 Free Software Foundation, Inc. + Copyright (C) 2004-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/bits/nan.h b/libc/sysdeps/ieee754/bits/nan.h index 41f47ba09..f0704da5c 100644 --- a/libc/sysdeps/ieee754/bits/nan.h +++ b/libc/sysdeps/ieee754/bits/nan.h @@ -1,5 +1,5 @@ /* `NAN' constant for IEEE 754 machines. - Copyright (C) 1992-2013 Free Software Foundation, Inc. + Copyright (C) 1992-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/dbl-64/Makefile b/libc/sysdeps/ieee754/dbl-64/Makefile index 1a7b31158..35f545ff8 100644 --- a/libc/sysdeps/ieee754/dbl-64/Makefile +++ b/libc/sysdeps/ieee754/dbl-64/Makefile @@ -1,4 +1,5 @@ ifeq ($(subdir),math) # branred depends on precise IEEE double rounding CFLAGS-branred.c = $(config-cflags-nofma) +CFLAGS-e_sqrt.c = $(config-cflags-nofma) endif diff --git a/libc/sysdeps/ieee754/dbl-64/MathLib.h b/libc/sysdeps/ieee754/dbl-64/MathLib.h index cf56606e3..395a57669 100644 --- a/libc/sysdeps/ieee754/dbl-64/MathLib.h +++ b/libc/sysdeps/ieee754/dbl-64/MathLib.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/asincos.tbl b/libc/sysdeps/ieee754/dbl-64/asincos.tbl index b50a366d8..840b61852 100644 --- a/libc/sysdeps/ieee754/dbl-64/asincos.tbl +++ b/libc/sysdeps/ieee754/dbl-64/asincos.tbl @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/atnat.h b/libc/sysdeps/ieee754/dbl-64/atnat.h index a1a3572f7..de612b09c 100644 --- a/libc/sysdeps/ieee754/dbl-64/atnat.h +++ b/libc/sysdeps/ieee754/dbl-64/atnat.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/atnat2.h b/libc/sysdeps/ieee754/dbl-64/atnat2.h index f12498bf3..d5a141deb 100644 --- a/libc/sysdeps/ieee754/dbl-64/atnat2.h +++ b/libc/sysdeps/ieee754/dbl-64/atnat2.h @@ -2,7 +2,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/branred.c b/libc/sysdeps/ieee754/dbl-64/branred.c index 524d091dc..9d15fdf4d 100644 --- a/libc/sysdeps/ieee754/dbl-64/branred.c +++ b/libc/sysdeps/ieee754/dbl-64/branred.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/branred.h b/libc/sysdeps/ieee754/dbl-64/branred.h index 4e27c37b8..2306e2e02 100644 --- a/libc/sysdeps/ieee754/dbl-64/branred.h +++ b/libc/sysdeps/ieee754/dbl-64/branred.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/dbl2mpn.c b/libc/sysdeps/ieee754/dbl-64/dbl2mpn.c index 087d64377..f7d818609 100644 --- a/libc/sysdeps/ieee754/dbl-64/dbl2mpn.c +++ b/libc/sysdeps/ieee754/dbl-64/dbl2mpn.c @@ -1,4 +1,4 @@ -/* Copyright (C) 1993-2013 Free Software Foundation, Inc. +/* Copyright (C) 1993-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/dbl-64/dla.h b/libc/sysdeps/ieee754/dbl-64/dla.h index 6666982bf..f90a6f8d2 100644 --- a/libc/sysdeps/ieee754/dbl-64/dla.h +++ b/libc/sysdeps/ieee754/dbl-64/dla.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/doasin.c b/libc/sysdeps/ieee754/dbl-64/doasin.c index ff444d5ca..132a3f011 100644 --- a/libc/sysdeps/ieee754/dbl-64/doasin.c +++ b/libc/sysdeps/ieee754/dbl-64/doasin.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/doasin.h b/libc/sysdeps/ieee754/dbl-64/doasin.h index 554821078..39cfcf4ae 100644 --- a/libc/sysdeps/ieee754/dbl-64/doasin.h +++ b/libc/sysdeps/ieee754/dbl-64/doasin.h @@ -2,7 +2,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/dosincos.c b/libc/sysdeps/ieee754/dbl-64/dosincos.c index e1c8836b7..de6f57bff 100644 --- a/libc/sysdeps/ieee754/dbl-64/dosincos.c +++ b/libc/sysdeps/ieee754/dbl-64/dosincos.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/dosincos.h b/libc/sysdeps/ieee754/dbl-64/dosincos.h index 49a22c8d7..516afd397 100644 --- a/libc/sysdeps/ieee754/dbl-64/dosincos.h +++ b/libc/sysdeps/ieee754/dbl-64/dosincos.h @@ -2,7 +2,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/e_asin.c b/libc/sysdeps/ieee754/dbl-64/e_asin.c index 27cd2a10d..5bb5aeb07 100644 --- a/libc/sysdeps/ieee754/dbl-64/e_asin.c +++ b/libc/sysdeps/ieee754/dbl-64/e_asin.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/e_atan2.c b/libc/sysdeps/ieee754/dbl-64/e_atan2.c index e36305cda..a287ca665 100644 --- a/libc/sysdeps/ieee754/dbl-64/e_atan2.c +++ b/libc/sysdeps/ieee754/dbl-64/e_atan2.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/e_atanh.c b/libc/sysdeps/ieee754/dbl-64/e_atanh.c index 958e9375c..21bb9908d 100644 --- a/libc/sysdeps/ieee754/dbl-64/e_atanh.c +++ b/libc/sysdeps/ieee754/dbl-64/e_atanh.c @@ -1,4 +1,4 @@ -/* Copyright (C) 2011-2013 Free Software Foundation, Inc. +/* Copyright (C) 2011-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@gmail.com>, 2011. diff --git a/libc/sysdeps/ieee754/dbl-64/e_exp.c b/libc/sysdeps/ieee754/dbl-64/e_exp.c index df3aa5efa..100fc397a 100644 --- a/libc/sysdeps/ieee754/dbl-64/e_exp.c +++ b/libc/sysdeps/ieee754/dbl-64/e_exp.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by @@ -39,6 +39,7 @@ #include "uexp.tbl" #include <math_private.h> #include <fenv.h> +#include <float.h> #ifndef SECTION # define SECTION @@ -169,7 +170,7 @@ __ieee754_exp (double x) else { retval = __slowexp (x); - goto ret; + goto check_uflow_ret; } /*if error is over bound */ } ex = -(1022 + ex); @@ -185,13 +186,23 @@ __ieee754_exp (double x) { binexp.i[HIGH_HALF] = 0x00100000; retval = (res - 1.0) * binexp.x; - goto ret; + goto check_uflow_ret; } else { retval = __slowexp (x); - goto ret; + goto check_uflow_ret; } /* if error is over bound */ + check_uflow_ret: + if (retval < DBL_MIN) + { +#if FLT_EVAL_METHOD != 0 + volatile +#endif + double force_underflow = tiny * tiny; + math_force_eval (force_underflow); + } + goto ret; } else { diff --git a/libc/sysdeps/ieee754/dbl-64/e_exp10.c b/libc/sysdeps/ieee754/dbl-64/e_exp10.c index 50b97f17a..eac609b71 100644 --- a/libc/sysdeps/ieee754/dbl-64/e_exp10.c +++ b/libc/sysdeps/ieee754/dbl-64/e_exp10.c @@ -1,4 +1,4 @@ -/* Copyright (C) 2012-2013 Free Software Foundation, Inc. +/* Copyright (C) 2012-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/dbl-64/e_exp2.c b/libc/sysdeps/ieee754/dbl-64/e_exp2.c index e1ba940e6..10e23e221 100644 --- a/libc/sysdeps/ieee754/dbl-64/e_exp2.c +++ b/libc/sysdeps/ieee754/dbl-64/e_exp2.c @@ -1,5 +1,5 @@ /* Double-precision floating point 2^x. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Geoffrey Keating <geoffk@ozemail.com.au> diff --git a/libc/sysdeps/ieee754/dbl-64/e_gamma_r.c b/libc/sysdeps/ieee754/dbl-64/e_gamma_r.c index 13e389d7c..1c427556b 100644 --- a/libc/sysdeps/ieee754/dbl-64/e_gamma_r.c +++ b/libc/sysdeps/ieee754/dbl-64/e_gamma_r.c @@ -1,5 +1,5 @@ /* Implementation of gamma function according to ISO C. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/e_hypot.c b/libc/sysdeps/ieee754/dbl-64/e_hypot.c index 9f4321efe..88242bc4f 100644 --- a/libc/sysdeps/ieee754/dbl-64/e_hypot.c +++ b/libc/sysdeps/ieee754/dbl-64/e_hypot.c @@ -89,7 +89,7 @@ __ieee754_hypot (double x, double y) SET_HIGH_WORD (a, ha); SET_HIGH_WORD (b, hb); } - if (__builtin_expect (hb < 0x20b00000, 0)) /* b < 2**-500 */ + if (__builtin_expect (hb < 0x23d00000, 0)) /* b < 2**-450 */ { if (hb <= 0x000fffff) /* subnormal b or 0 */ { @@ -102,6 +102,17 @@ __ieee754_hypot (double x, double y) b *= t1; a *= t1; k -= 1022; + GET_HIGH_WORD (ha, a); + GET_HIGH_WORD (hb, b); + if (hb > ha) + { + t1 = a; + a = b; + b = t1; + j = ha; + ha = hb; + hb = j; + } } else /* scale a and b by 2^600 */ { diff --git a/libc/sysdeps/ieee754/dbl-64/e_log.c b/libc/sysdeps/ieee754/dbl-64/e_log.c index a7ab54435..0b2889cb3 100644 --- a/libc/sysdeps/ieee754/dbl-64/e_log.c +++ b/libc/sysdeps/ieee754/dbl-64/e_log.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/e_pow.c b/libc/sysdeps/ieee754/dbl-64/e_pow.c index 9cf230917..1c5f4b311 100644 --- a/libc/sysdeps/ieee754/dbl-64/e_pow.c +++ b/libc/sysdeps/ieee754/dbl-64/e_pow.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/e_remainder.c b/libc/sysdeps/ieee754/dbl-64/e_remainder.c index c6a1901b1..28689614f 100644 --- a/libc/sysdeps/ieee754/dbl-64/e_remainder.c +++ b/libc/sysdeps/ieee754/dbl-64/e_remainder.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/e_sqrt.c b/libc/sysdeps/ieee754/dbl-64/e_sqrt.c index 854ae38c4..f095d6cb3 100644 --- a/libc/sysdeps/ieee754/dbl-64/e_sqrt.c +++ b/libc/sysdeps/ieee754/dbl-64/e_sqrt.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by @@ -66,8 +66,9 @@ __ieee754_sqrt (double x) /*----------------- 2^-1022 <= | x |< 2^1024 -----------------*/ if (k > 0x000fffff && k < 0x7ff00000) { + int rm = fegetround (); fenv_t env; - libc_feholdexcept (&env); + libc_feholdexcept_setround (&env, FE_TONEAREST); double ret; y = 1.0 - t * (t * s); t = t * (rt0 + y * (rt1 + y * (rt2 + y * rt3))); @@ -82,15 +83,44 @@ __ieee754_sqrt (double x) { res1 = res + 1.5 * ((y - res) + del); EMULV (res, res1, z, zz, p, hx, tx, hy, ty); /* (z+zz)=res*res1 */ - ret = ((((z - s) + zz) < 0) ? max (res, res1) : - min (res, res1)) * c.x; + res = ((((z - s) + zz) < 0) ? max (res, res1) : + min (res, res1)); + ret = res * c.x; } math_force_eval (ret); libc_fesetenv (&env); - if (x / ret != ret) + double dret = x / ret; + if (dret != ret) { double force_inexact = 1.0 / 3.0; math_force_eval (force_inexact); + /* The square root is inexact, ret is the round-to-nearest + value which may need adjusting for other rounding + modes. */ + switch (rm) + { +#ifdef FE_UPWARD + case FE_UPWARD: + if (dret > ret) + ret = (res + 0x1p-1022) * c.x; + break; +#endif + +#ifdef FE_DOWNWARD + case FE_DOWNWARD: +#endif +#ifdef FE_TOWARDZERO + case FE_TOWARDZERO: +#endif +#if defined FE_DOWNWARD || defined FE_TOWARDZERO + if (dret < ret) + ret = (res - 0x1p-1022) * c.x; + break; +#endif + + default: + break; + } } /* Otherwise (x / ret == ret), either the square root was exact or the division was inexact. */ diff --git a/libc/sysdeps/ieee754/dbl-64/gamma_product.c b/libc/sysdeps/ieee754/dbl-64/gamma_product.c index 2a3fc1aff..511cb3356 100644 --- a/libc/sysdeps/ieee754/dbl-64/gamma_product.c +++ b/libc/sysdeps/ieee754/dbl-64/gamma_product.c @@ -1,5 +1,5 @@ /* Compute a product of X, X+1, ..., with an error estimate. - Copyright (C) 2013 Free Software Foundation, Inc. + Copyright (C) 2013-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/dbl-64/gamma_productf.c b/libc/sysdeps/ieee754/dbl-64/gamma_productf.c index 46072f16e..60c6f7bce 100644 --- a/libc/sysdeps/ieee754/dbl-64/gamma_productf.c +++ b/libc/sysdeps/ieee754/dbl-64/gamma_productf.c @@ -1,5 +1,5 @@ /* Compute a product of X, X+1, ..., with an error estimate. - Copyright (C) 2013 Free Software Foundation, Inc. + Copyright (C) 2013-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/dbl-64/halfulp.c b/libc/sysdeps/ieee754/dbl-64/halfulp.c index 382ad7aad..68d613412 100644 --- a/libc/sysdeps/ieee754/dbl-64/halfulp.c +++ b/libc/sysdeps/ieee754/dbl-64/halfulp.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/mpa-arch.h b/libc/sysdeps/ieee754/dbl-64/mpa-arch.h index 779a43d9b..4a9f45ecc 100644 --- a/libc/sysdeps/ieee754/dbl-64/mpa-arch.h +++ b/libc/sysdeps/ieee754/dbl-64/mpa-arch.h @@ -1,5 +1,5 @@ /* Overridable constants and operations. - Copyright (C) 2013 Free Software Foundation, Inc. + Copyright (C) 2013-2014 Free Software Foundation, Inc. This program is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/mpa.c b/libc/sysdeps/ieee754/dbl-64/mpa.c index 190e8a637..9f7f44fa4 100644 --- a/libc/sysdeps/ieee754/dbl-64/mpa.c +++ b/libc/sysdeps/ieee754/dbl-64/mpa.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/mpa.h b/libc/sysdeps/ieee754/dbl-64/mpa.h index 5fad58439..bf1ad873d 100644 --- a/libc/sysdeps/ieee754/dbl-64/mpa.h +++ b/libc/sysdeps/ieee754/dbl-64/mpa.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/mpatan.c b/libc/sysdeps/ieee754/dbl-64/mpatan.c index a6ae61195..49ecdd2d8 100644 --- a/libc/sysdeps/ieee754/dbl-64/mpatan.c +++ b/libc/sysdeps/ieee754/dbl-64/mpatan.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/mpatan.h b/libc/sysdeps/ieee754/dbl-64/mpatan.h index 743a1b98c..4ac4541ff 100644 --- a/libc/sysdeps/ieee754/dbl-64/mpatan.h +++ b/libc/sysdeps/ieee754/dbl-64/mpatan.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/mpatan2.c b/libc/sysdeps/ieee754/dbl-64/mpatan2.c index f66f9eb88..524949266 100644 --- a/libc/sysdeps/ieee754/dbl-64/mpatan2.c +++ b/libc/sysdeps/ieee754/dbl-64/mpatan2.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/mpexp.c b/libc/sysdeps/ieee754/dbl-64/mpexp.c index 06347b74b..0096afb83 100644 --- a/libc/sysdeps/ieee754/dbl-64/mpexp.c +++ b/libc/sysdeps/ieee754/dbl-64/mpexp.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/mplog.c b/libc/sysdeps/ieee754/dbl-64/mplog.c index f8d5c1095..75adac9df 100644 --- a/libc/sysdeps/ieee754/dbl-64/mplog.c +++ b/libc/sysdeps/ieee754/dbl-64/mplog.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/mpn2dbl.c b/libc/sysdeps/ieee754/dbl-64/mpn2dbl.c index 0fada79a0..4e8a76432 100644 --- a/libc/sysdeps/ieee754/dbl-64/mpn2dbl.c +++ b/libc/sysdeps/ieee754/dbl-64/mpn2dbl.c @@ -1,4 +1,4 @@ -/* Copyright (C) 1995-2013 Free Software Foundation, Inc. +/* Copyright (C) 1995-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/dbl-64/mpsqrt.c b/libc/sysdeps/ieee754/dbl-64/mpsqrt.c index 230d1f36e..86418237b 100644 --- a/libc/sysdeps/ieee754/dbl-64/mpsqrt.c +++ b/libc/sysdeps/ieee754/dbl-64/mpsqrt.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/mpsqrt.h b/libc/sysdeps/ieee754/dbl-64/mpsqrt.h index 2b83c4cbf..9cf725b02 100644 --- a/libc/sysdeps/ieee754/dbl-64/mpsqrt.h +++ b/libc/sysdeps/ieee754/dbl-64/mpsqrt.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/mptan.c b/libc/sysdeps/ieee754/dbl-64/mptan.c index a5c3fb36c..aa16782c2 100644 --- a/libc/sysdeps/ieee754/dbl-64/mptan.c +++ b/libc/sysdeps/ieee754/dbl-64/mptan.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/mydefs.h b/libc/sysdeps/ieee754/dbl-64/mydefs.h index a430397c6..b53dd9476 100644 --- a/libc/sysdeps/ieee754/dbl-64/mydefs.h +++ b/libc/sysdeps/ieee754/dbl-64/mydefs.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/powtwo.tbl b/libc/sysdeps/ieee754/dbl-64/powtwo.tbl index 2963d14e3..d1d207155 100644 --- a/libc/sysdeps/ieee754/dbl-64/powtwo.tbl +++ b/libc/sysdeps/ieee754/dbl-64/powtwo.tbl @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/root.tbl b/libc/sysdeps/ieee754/dbl-64/root.tbl index 6ceddee39..f52ab96c2 100644 --- a/libc/sysdeps/ieee754/dbl-64/root.tbl +++ b/libc/sysdeps/ieee754/dbl-64/root.tbl @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/s_atan.c b/libc/sysdeps/ieee754/dbl-64/s_atan.c index 744e96146..0aa145cad 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_atan.c +++ b/libc/sysdeps/ieee754/dbl-64/s_atan.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/s_cbrt.c b/libc/sysdeps/ieee754/dbl-64/s_cbrt.c index 208a369a6..ada34b663 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_cbrt.c +++ b/libc/sysdeps/ieee754/dbl-64/s_cbrt.c @@ -1,5 +1,5 @@ /* Compute cubic root of double value. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Dirk Alboth <dirka@uni-paderborn.de> and Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/s_erf.c b/libc/sysdeps/ieee754/dbl-64/s_erf.c index aab4ed593..3f3739722 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_erf.c +++ b/libc/sysdeps/ieee754/dbl-64/s_erf.c @@ -112,6 +112,8 @@ static char rcsid[] = "$NetBSD: s_erf.c,v 1.8 1995/05/10 20:47:05 jtc Exp $"; */ +#include <errno.h> +#include <float.h> #include <math.h> #include <math_private.h> @@ -391,14 +393,25 @@ __erfc (double x) r = __ieee754_exp (-z * z - 0.5625) * __ieee754_exp ((z - x) * (z + x) + R / S); if (hx > 0) - return r / x; + { +#if FLT_EVAL_METHOD != 0 + volatile +#endif + double ret = r / x; + if (ret == 0) + __set_errno (ERANGE); + return ret; + } else return two - r / x; } else { if (hx > 0) - return tiny * tiny; + { + __set_errno (ERANGE); + return tiny * tiny; + } else return two - tiny; } diff --git a/libc/sysdeps/ieee754/dbl-64/s_finite.c b/libc/sysdeps/ieee754/dbl-64/s_finite.c index 5dc398efc..49986bbde 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_finite.c +++ b/libc/sysdeps/ieee754/dbl-64/s_finite.c @@ -23,8 +23,12 @@ static char rcsid[] = "$NetBSD: s_finite.c,v 1.8 1995/05/10 20:47:17 jtc Exp $"; #include <math_private.h> #undef __finite -int -__finite (double x) + +#ifndef FINITE +# define FINITE __finite +#endif + +int FINITE(double x) { int32_t hx; GET_HIGH_WORD (hx, x); diff --git a/libc/sysdeps/ieee754/dbl-64/s_fma.c b/libc/sysdeps/ieee754/dbl-64/s_fma.c index 06c0d2ae6..cfaa22d18 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_fma.c +++ b/libc/sysdeps/ieee754/dbl-64/s_fma.c @@ -1,5 +1,5 @@ /* Compute x * y + z as ternary operation. - Copyright (C) 2010-2013 Free Software Foundation, Inc. + Copyright (C) 2010-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jakub@redhat.com>, 2010. diff --git a/libc/sysdeps/ieee754/dbl-64/s_fmaf.c b/libc/sysdeps/ieee754/dbl-64/s_fmaf.c index e92145054..9a062f0bb 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_fmaf.c +++ b/libc/sysdeps/ieee754/dbl-64/s_fmaf.c @@ -1,5 +1,5 @@ /* Compute x * y + z as ternary operation. - Copyright (C) 2010-2013 Free Software Foundation, Inc. + Copyright (C) 2010-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jakub@redhat.com>, 2010. diff --git a/libc/sysdeps/ieee754/dbl-64/s_fpclassify.c b/libc/sysdeps/ieee754/dbl-64/s_fpclassify.c index 0f72454f1..d875ca4db 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_fpclassify.c +++ b/libc/sysdeps/ieee754/dbl-64/s_fpclassify.c @@ -1,5 +1,5 @@ /* Return classification value corresponding to argument. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/s_issignaling.c b/libc/sysdeps/ieee754/dbl-64/s_issignaling.c index f475bc7bd..4631ff9c5 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_issignaling.c +++ b/libc/sysdeps/ieee754/dbl-64/s_issignaling.c @@ -1,5 +1,5 @@ /* Test for signaling NaN. - Copyright (C) 2013 Free Software Foundation, Inc. + Copyright (C) 2013-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/dbl-64/s_llrint.c b/libc/sysdeps/ieee754/dbl-64/s_llrint.c index 0cbc39a55..3bcca2898 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_llrint.c +++ b/libc/sysdeps/ieee754/dbl-64/s_llrint.c @@ -1,6 +1,6 @@ /* Round argument to nearest integral value according to current rounding direction. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/s_llround.c b/libc/sysdeps/ieee754/dbl-64/s_llround.c index 7d50a57dd..e1cc1770c 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_llround.c +++ b/libc/sysdeps/ieee754/dbl-64/s_llround.c @@ -1,5 +1,5 @@ /* Round double value to long long int. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/s_lrint.c b/libc/sysdeps/ieee754/dbl-64/s_lrint.c index 1c3037364..3ab988473 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_lrint.c +++ b/libc/sysdeps/ieee754/dbl-64/s_lrint.c @@ -1,6 +1,6 @@ /* Round argument to nearest integral value according to current rounding direction. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/s_lround.c b/libc/sysdeps/ieee754/dbl-64/s_lround.c index e9dfcbaab..ae7f178b9 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_lround.c +++ b/libc/sysdeps/ieee754/dbl-64/s_lround.c @@ -1,5 +1,5 @@ /* Round double value to long int. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/s_remquo.c b/libc/sysdeps/ieee754/dbl-64/s_remquo.c index c40615004..04b1d81bb 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_remquo.c +++ b/libc/sysdeps/ieee754/dbl-64/s_remquo.c @@ -1,5 +1,5 @@ /* Compute remainder and a congruent to the quotient. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/s_round.c b/libc/sysdeps/ieee754/dbl-64/s_round.c index 7913ac061..7eacdf3dc 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_round.c +++ b/libc/sysdeps/ieee754/dbl-64/s_round.c @@ -1,5 +1,5 @@ /* Round double to integer away from zero. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/s_signbit.c b/libc/sysdeps/ieee754/dbl-64/s_signbit.c index bc1edac48..89bc8b6d6 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_signbit.c +++ b/libc/sysdeps/ieee754/dbl-64/s_signbit.c @@ -1,5 +1,5 @@ /* Return nonzero value if number is negative. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/s_sin.c b/libc/sysdeps/ieee754/dbl-64/s_sin.c index 53eef6002..6105e9fbd 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_sin.c +++ b/libc/sysdeps/ieee754/dbl-64/s_sin.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by @@ -56,10 +56,9 @@ #include <fenv.h> /* Helper macros to compute sin of the input values. */ -#define POLYNOMIAL2(xx) ((((s5.x * (xx) + s4.x) * (xx) + s3.x) * (xx) + s2.x) \ - * (xx)) +#define POLYNOMIAL2(xx) ((((s5 * (xx) + s4) * (xx) + s3) * (xx) + s2) * (xx)) -#define POLYNOMIAL(xx) (POLYNOMIAL2 (xx) + s1.x) +#define POLYNOMIAL(xx) (POLYNOMIAL2 (xx) + s1) /* The computed polynomial is a variation of the Taylor series expansion for sin(a): @@ -68,7 +67,7 @@ The constants s1, s2, s3, etc. are pre-computed values of 1/3!, 1/5! and so on. The result is returned to LHS and correction in COR. */ -#define TAYLOR_SINCOS(xx, a, da, cor) \ +#define TAYLOR_SIN(xx, a, da, cor) \ ({ \ double t = ((POLYNOMIAL (xx) * (a) - 0.5 * (da)) * (xx) + (da)); \ double res = (a) + t; \ @@ -88,11 +87,11 @@ static const double th2_36 = 206158430208.0; /* 1.5*2**37 */ \ double xx = (x0) * (x0); \ double x1 = ((x0) + th2_36) - th2_36; \ - double y = aa.x * x1 * x1 * x1; \ + double y = aa * x1 * x1 * x1; \ double r = (x0) + y; \ double x2 = ((x0) - x1) + (dx); \ - double t = (((POLYNOMIAL2 (xx) + bb.x) * xx + 3.0 * aa.x * x1 * x2) \ - * (x0) + aa.x * x2 * x2 * x2 + (dx)); \ + double t = (((POLYNOMIAL2 (xx) + bb) * xx + 3.0 * aa * x1 * x2) \ + * (x0) + aa * x2 * x2 * x2 + (dx)); \ t = (((x0) - r) + y) + t; \ double res = r + t; \ (cor) = (r - res) + t; \ @@ -125,6 +124,8 @@ static const double cs4 = -4.16666666666664434524222570944589E-02, cs6 = 1.38888874007937613028114285595617E-03; +static const double t22 = 0x1.8p22; + void __dubsin (double x, double dx, double w[]); void __docos (double x, double dx, double w[]); double __mpsin (double x, double dx, bool reduce_range); @@ -133,7 +134,7 @@ static double slow (double x); static double slow1 (double x); static double slow2 (double x); static double sloww (double x, double dx, double orig); -static double sloww1 (double x, double dx, double orig); +static double sloww1 (double x, double dx, double orig, int m); static double sloww2 (double x, double dx, double orig, int n); static double bsloww (double x, double dx, double orig, int n); static double bsloww1 (double x, double dx, double orig, int n); @@ -141,17 +142,113 @@ static double bsloww2 (double x, double dx, double orig, int n); int __branred (double x, double *a, double *aa); static double cslow2 (double x); static double csloww (double x, double dx, double orig); -static double csloww1 (double x, double dx, double orig); +static double csloww1 (double x, double dx, double orig, int m); static double csloww2 (double x, double dx, double orig, int n); +/* Given a number partitioned into U and X such that U is an index into the + sin/cos table, this macro computes the cosine of the number by combining + the sin and cos of X (as computed by a variation of the Taylor series) with + the values looked up from the sin/cos table to get the result in RES and a + correction value in COR. */ +static double +do_cos (mynumber u, double x, double *corp) +{ + double xx, s, sn, ssn, c, cs, ccs, res, cor; + xx = x * x; + s = x + x * xx * (sn3 + xx * sn5); + c = xx * (cs2 + xx * (cs4 + xx * cs6)); + SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); + cor = (ccs - s * ssn - cs * c) - sn * s; + res = cs + cor; + cor = (cs - res) + cor; + *corp = cor; + return res; +} + +/* A more precise variant of DO_COS where the number is partitioned into U, X + and DX. EPS is the adjustment to the correction COR. */ +static double +do_cos_slow (mynumber u, double x, double dx, double eps, double *corp) +{ + double xx, y, x1, x2, e1, e2, res, cor; + double s, sn, ssn, c, cs, ccs; + xx = x * x; + s = x * xx * (sn3 + xx * sn5); + c = x * dx + xx * (cs2 + xx * (cs4 + xx * cs6)); + SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); + x1 = (x + t22) - t22; + x2 = (x - x1) + dx; + e1 = (sn + t22) - t22; + e2 = (sn - e1) + ssn; + cor = (ccs - cs * c - e1 * x2 - e2 * x) - sn * s; + y = cs - e1 * x1; + cor = cor + ((cs - y) - e1 * x1); + res = y + cor; + cor = (y - res) + cor; + if (cor > 0) + cor = 1.0005 * cor + eps; + else + cor = 1.0005 * cor - eps; + *corp = cor; + return res; +} + +/* Given a number partitioned into U and X and DX such that U is an index into + the sin/cos table, this macro computes the sine of the number by combining + the sin and cos of X (as computed by a variation of the Taylor series) with + the values looked up from the sin/cos table to get the result in RES and a + correction value in COR. */ +static double +do_sin (mynumber u, double x, double dx, double *corp) +{ + double xx, s, sn, ssn, c, cs, ccs, cor, res; + xx = x * x; + s = x + (dx + x * xx * (sn3 + xx * sn5)); + c = x * dx + xx * (cs2 + xx * (cs4 + xx * cs6)); + SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); + cor = (ssn + s * ccs - sn * c) + cs * s; + res = sn + cor; + cor = (sn - res) + cor; + *corp = cor; + return res; +} + +/* A more precise variant of res = do_sin where the number is partitioned into U, X + and DX. EPS is the adjustment to the correction COR. */ +static double +do_sin_slow (mynumber u, double x, double dx, double eps, double *corp) +{ + double xx, y, x1, x2, c1, c2, res, cor; + double s, sn, ssn, c, cs, ccs; + xx = x * x; + s = x * xx * (sn3 + xx * sn5); + c = xx * (cs2 + xx * (cs4 + xx * cs6)); + SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); + x1 = (x + t22) - t22; + x2 = (x - x1) + dx; + c1 = (cs + t22) - t22; + c2 = (cs - c1) + ccs; + cor = (ssn + s * ccs + cs * s + c2 * x + c1 * x2 - sn * x * dx) - sn * c; + y = sn + c1 * x1; + cor = cor + ((sn - y) + c1 * x1); + res = y + cor; + cor = (y - res) + cor; + if (cor > 0) + cor = 1.0005 * cor + eps; + else + cor = 1.0005 * cor - eps; + *corp = cor; + return res; +} + /* Reduce range of X and compute sin of a + da. K is the amount by which to rotate the quadrants. This allows us to use the same routine to compute cos by simply rotating the quadrants by 1. */ static inline double __always_inline -reduce_and_compute (double x, double a, double da, unsigned int k) +reduce_and_compute (double x, unsigned int k) { - double retval = 0; + double retval = 0, a, da; unsigned int n = __branred (x, &a, &da); k = (n + k) % 4; switch (k) @@ -211,8 +308,8 @@ __sin (double x) /*---------------------------- 0.25<|x|< 0.855469---------------------- */ else if (k < 0x3feb6000) { - u.x = (m > 0) ? big.x + x : big.x - x; - y = (m > 0) ? x - (u.x - big.x) : x + (u.x - big.x); + u.x = (m > 0) ? big + x : big - x; + y = (m > 0) ? x - (u.x - big) : x + (u.x - big); xx = y * y; s = y + y * xx * (sn3 + xx * sn5); c = xx * (cs2 + xx * (cs4 + xx * cs6)); @@ -232,36 +329,30 @@ __sin (double x) else if (k < 0x400368fd) { - y = (m > 0) ? hp0.x - x : hp0.x + x; + y = (m > 0) ? hp0 - x : hp0 + x; if (y >= 0) { - u.x = big.x + y; - y = (y - (u.x - big.x)) + hp1.x; + u.x = big + y; + y = (y - (u.x - big)) + hp1; } else { - u.x = big.x - y; - y = (-hp1.x) - (y + (u.x - big.x)); + u.x = big - y; + y = (-hp1) - (y + (u.x - big)); } - xx = y * y; - s = y + y * xx * (sn3 + xx * sn5); - c = xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - cor = (ccs - s * ssn - cs * c) - sn * s; - res = cs + cor; - cor = (cs - res) + cor; + res = do_cos (u, y, &cor); retval = (res == res + 1.020 * cor) ? ((m > 0) ? res : -res) : slow2 (x); } /* else if (k < 0x400368fd) */ /*-------------------------- 2.426265<|x|< 105414350 ----------------------*/ else if (k < 0x419921FB) { - t = (x * hpinv.x + toint.x); - xn = t - toint.x; + t = (x * hpinv + toint); + xn = t - toint; v.x = t; - y = (x - xn * mp1.x) - xn * mp2.x; + y = (x - xn * mp1) - xn * mp2; n = v.i[LOW_HALF] & 3; - da = xn * mp3.x; + da = xn * mp3; a = y - da; da = (y - a) - da; eps = ABS (x) * 1.2e-30; @@ -279,36 +370,26 @@ __sin (double x) if (xx < 0.01588) { /* Taylor series. */ - res = TAYLOR_SINCOS (xx, a, da, cor); + res = TAYLOR_SIN (xx, a, da, cor); cor = (cor > 0) ? 1.02 * cor + eps : 1.02 * cor - eps; retval = (res == res + cor) ? res : sloww (a, da, x); } else { if (a > 0) - { - m = 1; - t = a; - db = da; - } + m = 1; else { m = 0; - t = -a; - db = -da; + a = -a; + da = -da; } - u.x = big.x + t; - y = t - (u.x - big.x); - xx = y * y; - s = y + (db + y * xx * (sn3 + xx * sn5)); - c = y * db + xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - cor = (ssn + s * ccs - sn * c) + cs * s; - res = sn + cor; - cor = (sn - res) + cor; + u.x = big + a; + y = a - (u.x - big); + res = do_sin (u, y, da, &cor); cor = (cor > 0) ? 1.035 * cor + eps : 1.035 * cor - eps; retval = ((res == res + cor) ? ((m) ? res : -res) - : sloww1 (a, da, x)); + : sloww1 (a, da, x, m)); } break; @@ -319,15 +400,9 @@ __sin (double x) a = -a; da = -da; } - u.x = big.x + a; - y = a - (u.x - big.x) + da; - xx = y * y; - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - s = y + y * xx * (sn3 + xx * sn5); - c = xx * (cs2 + xx * (cs4 + xx * cs6)); - cor = (ccs - s * ssn - cs * c) - sn * s; - res = cs + cor; - cor = (cs - res) + cor; + u.x = big + a; + y = a - (u.x - big) + da; + res = do_cos (u, y, &cor); cor = (cor > 0) ? 1.025 * cor + eps : 1.025 * cor - eps; retval = ((res == res + cor) ? ((n & 2) ? -res : res) : sloww2 (a, da, x, n)); @@ -338,17 +413,17 @@ __sin (double x) /*---------------------105414350 <|x|< 281474976710656 --------------------*/ else if (k < 0x42F00000) { - t = (x * hpinv.x + toint.x); - xn = t - toint.x; + t = (x * hpinv + toint); + xn = t - toint; v.x = t; xn1 = (xn + 8.0e22) - 8.0e22; xn2 = xn - xn1; - y = ((((x - xn1 * mp1.x) - xn1 * mp2.x) - xn2 * mp1.x) - xn2 * mp2.x); + y = ((((x - xn1 * mp1) - xn1 * mp2) - xn2 * mp1) - xn2 * mp2); n = v.i[LOW_HALF] & 3; - da = xn1 * pp3.x; + da = xn1 * pp3; t = y - da; da = (y - t) - da; - da = (da - xn2 * pp3.x) - xn * pp4.x; + da = (da - xn2 * pp3) - xn * pp4; a = t + da; da = (t - a) + da; eps = 1.0e-24; @@ -366,7 +441,7 @@ __sin (double x) if (xx < 0.01588) { /* Taylor series. */ - res = TAYLOR_SINCOS (xx, a, da, cor); + res = TAYLOR_SIN (xx, a, da, cor); cor = (cor > 0) ? 1.02 * cor + eps : 1.02 * cor - eps; retval = (res == res + cor) ? res : bsloww (a, da, x, n); } @@ -375,24 +450,17 @@ __sin (double x) if (a > 0) { m = 1; - t = a; db = da; } else { m = 0; - t = -a; + a = -a; db = -da; } - u.x = big.x + t; - y = t - (u.x - big.x); - xx = y * y; - s = y + (db + y * xx * (sn3 + xx * sn5)); - c = y * db + xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - cor = (ssn + s * ccs - sn * c) + cs * s; - res = sn + cor; - cor = (sn - res) + cor; + u.x = big + a; + y = a - (u.x - big); + res = do_sin (u, y, db, &cor); cor = (cor > 0) ? 1.035 * cor + eps : 1.035 * cor - eps; retval = ((res == res + cor) ? ((m) ? res : -res) : bsloww1 (a, da, x, n)); @@ -406,15 +474,9 @@ __sin (double x) a = -a; da = -da; } - u.x = big.x + a; - y = a - (u.x - big.x) + da; - xx = y * y; - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - s = y + y * xx * (sn3 + xx * sn5); - c = xx * (cs2 + xx * (cs4 + xx * cs6)); - cor = (ccs - s * ssn - cs * c) - sn * s; - res = cs + cor; - cor = (cs - res) + cor; + u.x = big + a; + y = a - (u.x - big) + da; + res = do_cos (u, y, &cor); cor = (cor > 0) ? 1.025 * cor + eps : 1.025 * cor - eps; retval = ((res == res + cor) ? ((n & 2) ? -res : res) : bsloww2 (a, da, x, n)); @@ -424,7 +486,7 @@ __sin (double x) /* -----------------281474976710656 <|x| <2^1024----------------------------*/ else if (k < 0x7ff00000) - retval = reduce_and_compute (x, a, da, 0); + retval = reduce_and_compute (x, 0); /*--------------------- |x| > 2^1024 ----------------------------------*/ else @@ -447,7 +509,7 @@ double SECTION __cos (double x) { - double y, xx, res, t, cor, s, c, sn, ssn, cs, ccs, xn, a, da, db, eps, xn1, + double y, xx, res, t, cor, xn, a, da, db, eps, xn1, xn2; mynumber u, v; int4 k, m, n; @@ -467,27 +529,21 @@ __cos (double x) else if (k < 0x3feb6000) { /* 2^-27 < |x| < 0.855469 */ y = ABS (x); - u.x = big.x + y; - y = y - (u.x - big.x); - xx = y * y; - s = y + y * xx * (sn3 + xx * sn5); - c = xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - cor = (ccs - s * ssn - cs * c) - sn * s; - res = cs + cor; - cor = (cs - res) + cor; + u.x = big + y; + y = y - (u.x - big); + res = do_cos (u, y, &cor); retval = (res == res + 1.020 * cor) ? res : cslow2 (x); } /* else if (k < 0x3feb6000) */ else if (k < 0x400368fd) { /* 0.855469 <|x|<2.426265 */ ; - y = hp0.x - ABS (x); - a = y + hp1.x; - da = (y - a) + hp1.x; + y = hp0 - ABS (x); + a = y + hp1; + da = (y - a) + hp1; xx = a * a; if (xx < 0.01588) { - res = TAYLOR_SINCOS (xx, a, da, cor); + res = TAYLOR_SIN (xx, a, da, cor); cor = (cor > 0) ? 1.02 * cor + 1.0e-31 : 1.02 * cor - 1.0e-31; retval = (res == res + cor) ? res : csloww (a, da, x); } @@ -496,27 +552,19 @@ __cos (double x) if (a > 0) { m = 1; - t = a; - db = da; } else { m = 0; - t = -a; - db = -da; + a = -a; + da = -da; } - u.x = big.x + t; - y = t - (u.x - big.x); - xx = y * y; - s = y + (db + y * xx * (sn3 + xx * sn5)); - c = y * db + xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - cor = (ssn + s * ccs - sn * c) + cs * s; - res = sn + cor; - cor = (sn - res) + cor; + u.x = big + a; + y = a - (u.x - big); + res = do_sin (u, y, da, &cor); cor = (cor > 0) ? 1.035 * cor + 1.0e-31 : 1.035 * cor - 1.0e-31; retval = ((res == res + cor) ? ((m) ? res : -res) - : csloww1 (a, da, x)); + : csloww1 (a, da, x, m)); } } /* else if (k < 0x400368fd) */ @@ -524,12 +572,12 @@ __cos (double x) else if (k < 0x419921FB) { /* 2.426265<|x|< 105414350 */ - t = (x * hpinv.x + toint.x); - xn = t - toint.x; + t = (x * hpinv + toint); + xn = t - toint; v.x = t; - y = (x - xn * mp1.x) - xn * mp2.x; + y = (x - xn * mp1) - xn * mp2; n = v.i[LOW_HALF] & 3; - da = xn * mp3.x; + da = xn * mp3; a = y - da; da = (y - a) - da; eps = ABS (x) * 1.2e-30; @@ -546,7 +594,7 @@ __cos (double x) } if (xx < 0.01588) { - res = TAYLOR_SINCOS (xx, a, da, cor); + res = TAYLOR_SIN (xx, a, da, cor); cor = (cor > 0) ? 1.02 * cor + eps : 1.02 * cor - eps; retval = (res == res + cor) ? res : csloww (a, da, x); } @@ -555,27 +603,19 @@ __cos (double x) if (a > 0) { m = 1; - t = a; - db = da; } else { m = 0; - t = -a; - db = -da; + a = -a; + da = -da; } - u.x = big.x + t; - y = t - (u.x - big.x); - xx = y * y; - s = y + (db + y * xx * (sn3 + xx * sn5)); - c = y * db + xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - cor = (ssn + s * ccs - sn * c) + cs * s; - res = sn + cor; - cor = (sn - res) + cor; + u.x = big + a; + y = a - (u.x - big); + res = do_sin (u, y, da, &cor); cor = (cor > 0) ? 1.035 * cor + eps : 1.035 * cor - eps; retval = ((res == res + cor) ? ((m) ? res : -res) - : csloww1 (a, da, x)); + : csloww1 (a, da, x, m)); } break; @@ -586,15 +626,9 @@ __cos (double x) a = -a; da = -da; } - u.x = big.x + a; - y = a - (u.x - big.x) + da; - xx = y * y; - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - s = y + y * xx * (sn3 + xx * sn5); - c = xx * (cs2 + xx * (cs4 + xx * cs6)); - cor = (ccs - s * ssn - cs * c) - sn * s; - res = cs + cor; - cor = (cs - res) + cor; + u.x = big + a; + y = a - (u.x - big) + da; + res = do_cos (u, y, &cor); cor = (cor > 0) ? 1.025 * cor + eps : 1.025 * cor - eps; retval = ((res == res + cor) ? ((n) ? -res : res) : csloww2 (a, da, x, n)); @@ -604,17 +638,17 @@ __cos (double x) else if (k < 0x42F00000) { - t = (x * hpinv.x + toint.x); - xn = t - toint.x; + t = (x * hpinv + toint); + xn = t - toint; v.x = t; xn1 = (xn + 8.0e22) - 8.0e22; xn2 = xn - xn1; - y = ((((x - xn1 * mp1.x) - xn1 * mp2.x) - xn2 * mp1.x) - xn2 * mp2.x); + y = ((((x - xn1 * mp1) - xn1 * mp2) - xn2 * mp1) - xn2 * mp2); n = v.i[LOW_HALF] & 3; - da = xn1 * pp3.x; + da = xn1 * pp3; t = y - da; da = (y - t) - da; - da = (da - xn2 * pp3.x) - xn * pp4.x; + da = (da - xn2 * pp3) - xn * pp4; a = t + da; da = (t - a) + da; eps = 1.0e-24; @@ -631,7 +665,7 @@ __cos (double x) } if (xx < 0.01588) { - res = TAYLOR_SINCOS (xx, a, da, cor); + res = TAYLOR_SIN (xx, a, da, cor); cor = (cor > 0) ? 1.02 * cor + eps : 1.02 * cor - eps; retval = (res == res + cor) ? res : bsloww (a, da, x, n); } @@ -640,24 +674,17 @@ __cos (double x) if (a > 0) { m = 1; - t = a; db = da; } else { m = 0; - t = -a; + a = -a; db = -da; } - u.x = big.x + t; - y = t - (u.x - big.x); - xx = y * y; - s = y + (db + y * xx * (sn3 + xx * sn5)); - c = y * db + xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - cor = (ssn + s * ccs - sn * c) + cs * s; - res = sn + cor; - cor = (sn - res) + cor; + u.x = big + a; + y = a - (u.x - big); + res = do_sin (u, y, db, &cor); cor = (cor > 0) ? 1.035 * cor + eps : 1.035 * cor - eps; retval = ((res == res + cor) ? ((m) ? res : -res) : bsloww1 (a, da, x, n)); @@ -671,15 +698,9 @@ __cos (double x) a = -a; da = -da; } - u.x = big.x + a; - y = a - (u.x - big.x) + da; - xx = y * y; - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - s = y + y * xx * (sn3 + xx * sn5); - c = xx * (cs2 + xx * (cs4 + xx * cs6)); - cor = (ccs - s * ssn - cs * c) - sn * s; - res = cs + cor; - cor = (cs - res) + cor; + u.x = big + a; + y = a - (u.x - big) + da; + res = do_cos (u, y, &cor); cor = (cor > 0) ? 1.025 * cor + eps : 1.025 * cor - eps; retval = ((res == res + cor) ? ((n) ? -res : res) : bsloww2 (a, da, x, n)); @@ -689,7 +710,7 @@ __cos (double x) /* 281474976710656 <|x| <2^1024 */ else if (k < 0x7ff00000) - retval = reduce_and_compute (x, a, da, 1); + retval = reduce_and_compute (x, 1); else { @@ -734,25 +755,12 @@ SECTION slow1 (double x) { mynumber u; - double sn, ssn, cs, ccs, s, c, w[2], y, y1, y2, c1, c2, xx, cor, res; - static const double t22 = 6291456.0; + double w[2], y, cor, res; y = ABS (x); - u.x = big.x + y; - y = y - (u.x - big.x); - xx = y * y; - s = y * xx * (sn3 + xx * sn5); - c = xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - y1 = (y + t22) - t22; - y2 = y - y1; - c1 = (cs + t22) - t22; - c2 = (cs - c1) + ccs; - cor = (ssn + s * ccs + cs * s + c2 * y + c1 * y2) - sn * c; - y = sn + c1 * y1; - cor = cor + ((sn - y) + c1 * y1); - res = y + cor; - cor = (y - res) + cor; - if (res == res + 1.0005 * cor) + u.x = big + y; + y = y - (u.x - big); + res = do_sin_slow (u, y, 0, 0, &cor); + if (res == res + cor) return (x > 0) ? res : -res; else { @@ -773,42 +781,30 @@ SECTION slow2 (double x) { mynumber u; - double sn, ssn, cs, ccs, s, c, w[2], y, y1, y2, e1, e2, xx, cor, res, del; - static const double t22 = 6291456.0; + double w[2], y, y1, y2, cor, res, del; + y = ABS (x); - y = hp0.x - y; + y = hp0 - y; if (y >= 0) { - u.x = big.x + y; - y = y - (u.x - big.x); - del = hp1.x; + u.x = big + y; + y = y - (u.x - big); + del = hp1; } else { - u.x = big.x - y; - y = -(y + (u.x - big.x)); - del = -hp1.x; + u.x = big - y; + y = -(y + (u.x - big)); + del = -hp1; } - xx = y * y; - s = y * xx * (sn3 + xx * sn5); - c = y * del + xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - y1 = (y + t22) - t22; - y2 = (y - y1) + del; - e1 = (sn + t22) - t22; - e2 = (sn - e1) + ssn; - cor = (ccs - cs * c - e1 * y2 - e2 * y) - sn * s; - y = cs - e1 * y1; - cor = cor + ((cs - y) - e1 * y1); - res = y + cor; - cor = (y - res) + cor; - if (res == res + 1.0005 * cor) + res = do_cos_slow (u, y, del, 0, &cor); + if (res == res + cor) return (x > 0) ? res : -res; else { - y = ABS (x) - hp0.x; - y1 = y - hp1.x; - y2 = (y - y1) - hp1.x; + y = ABS (x) - hp0; + y1 = y - hp1; + y2 = (y - y1) - hp1; __docos (y1, y2, w); if (w[0] == w[0] + 1.000000005 * w[1]) return (x > 0) ? w[0] : -w[0]; @@ -830,17 +826,14 @@ SECTION sloww (double x, double dx, double orig) { double y, t, res, cor, w[2], a, da, xn; - union - { - int4 i[2]; - double x; - } v; + mynumber v; int4 n; res = TAYLOR_SLOW (x, dx, cor); - cor = - (cor > - 0) ? 1.0005 * cor + ABS (orig) * 3.1e-30 : 1.0005 * cor - - ABS (orig) * 3.1e-30; + if (cor > 0) + cor = 1.0005 * cor + ABS (orig) * 3.1e-30; + else + cor = 1.0005 * cor - ABS (orig) * 3.1e-30; + if (res == res + cor) return res; else @@ -855,15 +848,15 @@ sloww (double x, double dx, double orig) return (x > 0) ? w[0] : -w[0]; else { - t = (orig * hpinv.x + toint.x); - xn = t - toint.x; + t = (orig * hpinv + toint); + xn = t - toint; v.x = t; - y = (orig - xn * mp1.x) - xn * mp2.x; + y = (orig - xn * mp1) - xn * mp2; n = v.i[LOW_HALF] & 3; - da = xn * pp3.x; + da = xn * pp3; t = y - da; da = (y - t) - da; - y = xn * pp4.x; + y = xn * pp4; a = t - y; da = ((t - a) - y) + da; if (n & 2) @@ -894,40 +887,20 @@ sloww (double x, double dx, double orig) static double SECTION -sloww1 (double x, double dx, double orig) +sloww1 (double x, double dx, double orig, int m) { mynumber u; - double sn, ssn, cs, ccs, s, c, w[2], y, y1, y2, c1, c2, xx, cor, res; - static const double t22 = 6291456.0; - - y = ABS (x); - u.x = big.x + y; - y = y - (u.x - big.x); - dx = (x > 0) ? dx : -dx; - xx = y * y; - s = y * xx * (sn3 + xx * sn5); - c = xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - y1 = (y + t22) - t22; - y2 = (y - y1) + dx; - c1 = (cs + t22) - t22; - c2 = (cs - c1) + ccs; - cor = (ssn + s * ccs + cs * s + c2 * y + c1 * y2 - sn * y * dx) - sn * c; - y = sn + c1 * y1; - cor = cor + ((sn - y) + c1 * y1); - res = y + cor; - cor = (y - res) + cor; + double w[2], y, cor, res; - if (cor > 0) - cor = 1.0005 * cor + 3.1e-30 * ABS (orig); - else - cor = 1.0005 * cor - 3.1e-30 * ABS (orig); + u.x = big + x; + y = x - (u.x - big); + res = do_sin_slow (u, y, dx, 3.1e-30 * ABS (orig), &cor); if (res == res + cor) - return (x > 0) ? res : -res; + return (m > 0) ? res : -res; else { - __dubsin (ABS (x), dx, w); + __dubsin (x, dx, w); if (w[1] > 0) cor = 1.000000005 * w[1] + 1.1e-30 * ABS (orig); @@ -935,7 +908,7 @@ sloww1 (double x, double dx, double orig) cor = 1.000000005 * w[1] - 1.1e-30 * ABS (orig); if (w[0] == w[0] + cor) - return (x > 0) ? w[0] : -w[0]; + return (m > 0) ? w[0] : -w[0]; else return __mpsin (orig, 0, true); } @@ -953,38 +926,17 @@ SECTION sloww2 (double x, double dx, double orig, int n) { mynumber u; - double sn, ssn, cs, ccs, s, c, w[2], y, y1, y2, e1, e2, xx, cor, res; - static const double t22 = 6291456.0; + double w[2], y, cor, res; - y = ABS (x); - u.x = big.x + y; - y = y - (u.x - big.x); - dx = (x > 0) ? dx : -dx; - xx = y * y; - s = y * xx * (sn3 + xx * sn5); - c = y * dx + xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - - y1 = (y + t22) - t22; - y2 = (y - y1) + dx; - e1 = (sn + t22) - t22; - e2 = (sn - e1) + ssn; - cor = (ccs - cs * c - e1 * y2 - e2 * y) - sn * s; - y = cs - e1 * y1; - cor = cor + ((cs - y) - e1 * y1); - res = y + cor; - cor = (y - res) + cor; - - if (cor > 0) - cor = 1.0005 * cor + 3.1e-30 * ABS (orig); - else - cor = 1.0005 * cor - 3.1e-30 * ABS (orig); + u.x = big + x; + y = x - (u.x - big); + res = do_cos_slow (u, y, dx, 3.1e-30 * ABS (orig), &cor); if (res == res + cor) return (n & 2) ? -res : res; else { - __docos (ABS (x), dx, w); + __docos (x, dx, w); if (w[1] > 0) cor = 1.000000005 * w[1] + 1.1e-30 * ABS (orig); @@ -1042,27 +994,13 @@ SECTION bsloww1 (double x, double dx, double orig, int n) { mynumber u; - double sn, ssn, cs, ccs, s, c, w[2], y, y1, y2, c1, c2, xx, cor, res; - static const double t22 = 6291456.0; + double w[2], y, cor, res; y = ABS (x); - u.x = big.x + y; - y = y - (u.x - big.x); + u.x = big + y; + y = y - (u.x - big); dx = (x > 0) ? dx : -dx; - xx = y * y; - s = y * xx * (sn3 + xx * sn5); - c = xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - y1 = (y + t22) - t22; - y2 = (y - y1) + dx; - c1 = (cs + t22) - t22; - c2 = (cs - c1) + ccs; - cor = (ssn + s * ccs + cs * s + c2 * y + c1 * y2 - sn * y * dx) - sn * c; - y = sn + c1 * y1; - cor = cor + ((sn - y) + c1 * y1); - res = y + cor; - cor = (y - res) + cor; - cor = (cor > 0) ? 1.0005 * cor + 1.1e-24 : 1.0005 * cor - 1.1e-24; + res = do_sin_slow (u, y, dx, 1.1e-24, &cor); if (res == res + cor) return (x > 0) ? res : -res; else @@ -1093,28 +1031,13 @@ SECTION bsloww2 (double x, double dx, double orig, int n) { mynumber u; - double sn, ssn, cs, ccs, s, c, w[2], y, y1, y2, e1, e2, xx, cor, res; - static const double t22 = 6291456.0; + double w[2], y, cor, res; y = ABS (x); - u.x = big.x + y; - y = y - (u.x - big.x); + u.x = big + y; + y = y - (u.x - big); dx = (x > 0) ? dx : -dx; - xx = y * y; - s = y * xx * (sn3 + xx * sn5); - c = y * dx + xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - - y1 = (y + t22) - t22; - y2 = (y - y1) + dx; - e1 = (sn + t22) - t22; - e2 = (sn - e1) + ssn; - cor = (ccs - cs * c - e1 * y2 - e2 * y) - sn * s; - y = cs - e1 * y1; - cor = cor + ((cs - y) - e1 * y1); - res = y + cor; - cor = (y - res) + cor; - cor = (cor > 0) ? 1.0005 * cor + 1.1e-24 : 1.0005 * cor - 1.1e-24; + res = do_cos_slow (u, y, dx, 1.1e-24, &cor); if (res == res + cor) return (n & 2) ? -res : res; else @@ -1143,26 +1066,13 @@ SECTION cslow2 (double x) { mynumber u; - double sn, ssn, cs, ccs, s, c, w[2], y, y1, y2, e1, e2, xx, cor, res; - static const double t22 = 6291456.0; + double w[2], y, cor, res; y = ABS (x); - u.x = big.x + y; - y = y - (u.x - big.x); - xx = y * y; - s = y * xx * (sn3 + xx * sn5); - c = xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - y1 = (y + t22) - t22; - y2 = y - y1; - e1 = (sn + t22) - t22; - e2 = (sn - e1) + ssn; - cor = (ccs - cs * c - e1 * y2 - e2 * y) - sn * s; - y = cs - e1 * y1; - cor = cor + ((cs - y) - e1 * y1); - res = y + cor; - cor = (y - res) + cor; - if (res == res + 1.0005 * cor) + u.x = big + y; + y = y - (u.x - big); + res = do_cos_slow (u, y, 0, 0, &cor); + if (res == res + cor) return res; else { @@ -1188,15 +1098,11 @@ SECTION csloww (double x, double dx, double orig) { double y, t, res, cor, w[2], a, da, xn; - union - { - int4 i[2]; - double x; - } v; + mynumber v; int4 n; /* Taylor series */ - t = TAYLOR_SLOW (x, dx, cor); + res = TAYLOR_SLOW (x, dx, cor); if (cor > 0) cor = 1.0005 * cor + ABS (orig) * 3.1e-30; @@ -1218,15 +1124,15 @@ csloww (double x, double dx, double orig) return (x > 0) ? w[0] : -w[0]; else { - t = (orig * hpinv.x + toint.x); - xn = t - toint.x; + t = (orig * hpinv + toint); + xn = t - toint; v.x = t; - y = (orig - xn * mp1.x) - xn * mp2.x; + y = (orig - xn * mp1) - xn * mp2; n = v.i[LOW_HALF] & 3; - da = xn * pp3.x; + da = xn * pp3; t = y - da; da = (y - t) - da; - y = xn * pp4.x; + y = xn * pp4; a = t - y; da = ((t - a) - y) + da; if (n == 1) @@ -1258,46 +1164,26 @@ csloww (double x, double dx, double orig) static double SECTION -csloww1 (double x, double dx, double orig) +csloww1 (double x, double dx, double orig, int m) { mynumber u; - double sn, ssn, cs, ccs, s, c, w[2], y, y1, y2, c1, c2, xx, cor, res; - static const double t22 = 6291456.0; - - y = ABS (x); - u.x = big.x + y; - y = y - (u.x - big.x); - dx = (x > 0) ? dx : -dx; - xx = y * y; - s = y * xx * (sn3 + xx * sn5); - c = xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); - y1 = (y + t22) - t22; - y2 = (y - y1) + dx; - c1 = (cs + t22) - t22; - c2 = (cs - c1) + ccs; - cor = (ssn + s * ccs + cs * s + c2 * y + c1 * y2 - sn * y * dx) - sn * c; - y = sn + c1 * y1; - cor = cor + ((sn - y) + c1 * y1); - res = y + cor; - cor = (y - res) + cor; + double w[2], y, cor, res; - if (cor > 0) - cor = 1.0005 * cor + 3.1e-30 * ABS (orig); - else - cor = 1.0005 * cor - 3.1e-30 * ABS (orig); + u.x = big + x; + y = x - (u.x - big); + res = do_sin_slow (u, y, dx, 3.1e-30 * ABS (orig), &cor); if (res == res + cor) - return (x > 0) ? res : -res; + return (m > 0) ? res : -res; else { - __dubsin (ABS (x), dx, w); + __dubsin (x, dx, w); if (w[1] > 0) cor = 1.000000005 * w[1] + 1.1e-30 * ABS (orig); else cor = 1.000000005 * w[1] - 1.1e-30 * ABS (orig); if (w[0] == w[0] + cor) - return (x > 0) ? w[0] : -w[0]; + return (m > 0) ? w[0] : -w[0]; else return __mpcos (orig, 0, true); } @@ -1316,38 +1202,17 @@ SECTION csloww2 (double x, double dx, double orig, int n) { mynumber u; - double sn, ssn, cs, ccs, s, c, w[2], y, y1, y2, e1, e2, xx, cor, res; - static const double t22 = 6291456.0; - - y = ABS (x); - u.x = big.x + y; - y = y - (u.x - big.x); - dx = (x > 0) ? dx : -dx; - xx = y * y; - s = y * xx * (sn3 + xx * sn5); - c = y * dx + xx * (cs2 + xx * (cs4 + xx * cs6)); - SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); + double w[2], y, cor, res; - y1 = (y + t22) - t22; - y2 = (y - y1) + dx; - e1 = (sn + t22) - t22; - e2 = (sn - e1) + ssn; - cor = (ccs - cs * c - e1 * y2 - e2 * y) - sn * s; - y = cs - e1 * y1; - cor = cor + ((cs - y) - e1 * y1); - res = y + cor; - cor = (y - res) + cor; - - if (cor > 0) - cor = 1.0005 * cor + 3.1e-30 * ABS (orig); - else - cor = 1.0005 * cor - 3.1e-30 * ABS (orig); + u.x = big + x; + y = x - (u.x - big); + res = do_cos_slow (u, y, dx, 3.1e-30 * ABS (orig), &cor); if (res == res + cor) return (n) ? -res : res; else { - __docos (ABS (x), dx, w); + __docos (x, dx, w); if (w[1] > 0) cor = 1.000000005 * w[1] + 1.1e-30 * ABS (orig); else diff --git a/libc/sysdeps/ieee754/dbl-64/s_sincos.c b/libc/sysdeps/ieee754/dbl-64/s_sincos.c index 01d1bdf97..e4371069d 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_sincos.c +++ b/libc/sysdeps/ieee754/dbl-64/s_sincos.c @@ -1,5 +1,5 @@ /* Compute sine and cosine of argument. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/s_tan.c b/libc/sysdeps/ieee754/dbl-64/s_tan.c index 09db096d0..a72530dad 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_tan.c +++ b/libc/sysdeps/ieee754/dbl-64/s_tan.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/s_trunc.c b/libc/sysdeps/ieee754/dbl-64/s_trunc.c index e23447721..bd777f686 100644 --- a/libc/sysdeps/ieee754/dbl-64/s_trunc.c +++ b/libc/sysdeps/ieee754/dbl-64/s_trunc.c @@ -1,5 +1,5 @@ /* Truncate argument to nearest integral value not larger than the argument. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/sincos32.c b/libc/sysdeps/ieee754/dbl-64/sincos32.c index 49aa14895..6b2fa878a 100644 --- a/libc/sysdeps/ieee754/dbl-64/sincos32.c +++ b/libc/sysdeps/ieee754/dbl-64/sincos32.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by @@ -43,6 +43,7 @@ #include "mpa.h" #include "sincos32.h" #include <math_private.h> +#include <stap-probe.h> #ifndef SECTION # define SECTION @@ -149,6 +150,7 @@ __sin32 (double x, double res, double res1) /* if a > 0 return min (res, res1), otherwise return max (res, res1). */ if ((a.d[0] > 0 && res >= res1) || (a.d[0] <= 0 && res <= res1)) res = res1; + LIBC_PROBE (slowasin, 2, &res, &x); return res; } @@ -182,6 +184,7 @@ __cos32 (double x, double res, double res1) /* if a > 0 return max (res, res1), otherwise return min (res, res1). */ if ((a.d[0] > 0 && res <= res1) || (a.d[0] <= 0 && res >= res1)) res = res1; + LIBC_PROBE (slowacos, 2, &res, &x); return res; } @@ -240,6 +243,7 @@ __mpsin (double x, double dx, bool reduce_range) default: __mp_dbl (&s, &y, p); } + LIBC_PROBE (slowsin, 3, &x, &dx, &y); return y; } @@ -298,6 +302,7 @@ __mpcos (double x, double dx, bool reduce_range) default: __mp_dbl (&c, &y, p); } + LIBC_PROBE (slowcos, 3, &x, &dx, &y); return y; } diff --git a/libc/sysdeps/ieee754/dbl-64/sincos32.h b/libc/sysdeps/ieee754/dbl-64/sincos32.h index 1e5ac8bcf..e289b165c 100644 --- a/libc/sysdeps/ieee754/dbl-64/sincos32.h +++ b/libc/sysdeps/ieee754/dbl-64/sincos32.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/sincostab.c b/libc/sysdeps/ieee754/dbl-64/sincostab.c index 43833f0b5..e3e24e0f8 100644 --- a/libc/sysdeps/ieee754/dbl-64/sincostab.c +++ b/libc/sysdeps/ieee754/dbl-64/sincostab.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/slowexp.c b/libc/sysdeps/ieee754/dbl-64/slowexp.c index 525224f44..e32742d5d 100644 --- a/libc/sysdeps/ieee754/dbl-64/slowexp.c +++ b/libc/sysdeps/ieee754/dbl-64/slowexp.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/slowpow.c b/libc/sysdeps/ieee754/dbl-64/slowpow.c index d200c39e5..b2e403308 100644 --- a/libc/sysdeps/ieee754/dbl-64/slowpow.c +++ b/libc/sysdeps/ieee754/dbl-64/slowpow.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/t_exp.c b/libc/sysdeps/ieee754/dbl-64/t_exp.c index 7466e2997..f32fa57a4 100644 --- a/libc/sysdeps/ieee754/dbl-64/t_exp.c +++ b/libc/sysdeps/ieee754/dbl-64/t_exp.c @@ -1,5 +1,5 @@ /* Accurate tables for exp(). - Copyright (C) 1998-2013 Free Software Foundation, Inc. + Copyright (C) 1998-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Geoffrey Keating <geoffk@ozemail.com.au> diff --git a/libc/sysdeps/ieee754/dbl-64/uasncs.h b/libc/sysdeps/ieee754/dbl-64/uasncs.h index 66030d1c6..db56a5354 100644 --- a/libc/sysdeps/ieee754/dbl-64/uasncs.h +++ b/libc/sysdeps/ieee754/dbl-64/uasncs.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/uatan.tbl b/libc/sysdeps/ieee754/dbl-64/uatan.tbl index d259587b7..24e72e710 100644 --- a/libc/sysdeps/ieee754/dbl-64/uatan.tbl +++ b/libc/sysdeps/ieee754/dbl-64/uatan.tbl @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/uexp.h b/libc/sysdeps/ieee754/dbl-64/uexp.h index 073c0ee42..dd59451e3 100644 --- a/libc/sysdeps/ieee754/dbl-64/uexp.h +++ b/libc/sysdeps/ieee754/dbl-64/uexp.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/uexp.tbl b/libc/sysdeps/ieee754/dbl-64/uexp.tbl index 68396b79d..90ff55de4 100644 --- a/libc/sysdeps/ieee754/dbl-64/uexp.tbl +++ b/libc/sysdeps/ieee754/dbl-64/uexp.tbl @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/ulog.h b/libc/sysdeps/ieee754/dbl-64/ulog.h index eec1eef67..70389b4cd 100644 --- a/libc/sysdeps/ieee754/dbl-64/ulog.h +++ b/libc/sysdeps/ieee754/dbl-64/ulog.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/ulog.tbl b/libc/sysdeps/ieee754/dbl-64/ulog.tbl index 32f384635..454fb61b2 100644 --- a/libc/sysdeps/ieee754/dbl-64/ulog.tbl +++ b/libc/sysdeps/ieee754/dbl-64/ulog.tbl @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/upow.h b/libc/sysdeps/ieee754/dbl-64/upow.h index 130fb3296..f26240652 100644 --- a/libc/sysdeps/ieee754/dbl-64/upow.h +++ b/libc/sysdeps/ieee754/dbl-64/upow.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/upow.tbl b/libc/sysdeps/ieee754/dbl-64/upow.tbl index f836b9b7b..f5e7660d9 100644 --- a/libc/sysdeps/ieee754/dbl-64/upow.tbl +++ b/libc/sysdeps/ieee754/dbl-64/upow.tbl @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/urem.h b/libc/sysdeps/ieee754/dbl-64/urem.h index b576691ce..3bd2883be 100644 --- a/libc/sysdeps/ieee754/dbl-64/urem.h +++ b/libc/sysdeps/ieee754/dbl-64/urem.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/uroot.h b/libc/sysdeps/ieee754/dbl-64/uroot.h index f2484994c..2db841ffa 100644 --- a/libc/sysdeps/ieee754/dbl-64/uroot.h +++ b/libc/sysdeps/ieee754/dbl-64/uroot.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/usncs.h b/libc/sysdeps/ieee754/dbl-64/usncs.h index 11c189ce3..b3296ca9c 100644 --- a/libc/sysdeps/ieee754/dbl-64/usncs.h +++ b/libc/sysdeps/ieee754/dbl-64/usncs.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by @@ -27,51 +27,22 @@ #ifndef USNCS_H #define USNCS_H -#ifdef BIG_ENDI -static const mynumber - -/**/ s1 = {{0xBFC55555, 0x55555555 }}, /* -0.16666666666666666 */ -/**/ s2 = {{0x3F811111, 0x11110ECE }}, /* 0.0083333333333323288 */ -/**/ s3 = {{0xBF2A01A0, 0x19DB08B8 }}, /* -0.00019841269834414642 */ -/**/ s4 = {{0x3EC71DE2, 0x7B9A7ED9 }}, /* 2.755729806860771e-06 */ -/**/ s5 = {{0xBE5ADDFF, 0xC2FCDF59 }}, /* -2.5022014848318398e-08 */ -/**/ aa = {{0xBFC55580, 0x00000000 }}, /* -0.1666717529296875 */ -/**/ bb = {{0x3ED55555, 0x55556E24 }}, /* 5.0862630208387126e-06 */ -/**/ big = {{0x42c80000, 0x00000000 }}, /* 52776558133248 */ -/**/ hp0 = {{0x3FF921FB, 0x54442D18 }}, /* 1.5707963267948966 */ -/**/ hp1 = {{0x3C91A626, 0x33145C07 }}, /* 6.123233995736766e-17 */ -/**/ mp1 = {{0x3FF921FB, 0x58000000 }}, /* 1.5707963407039642 */ -/**/ mp2 = {{0xBE4DDE97, 0x3C000000 }}, /* -1.3909067564377153e-08 */ -/**/ mp3 = {{0xBC8CB3B3, 0x99D747F2 }}, /* -4.9789962505147994e-17 */ -/**/ pp3 = {{0xBC8CB3B3, 0x98000000 }}, /* -4.9789962314799099e-17 */ -/**/ pp4 = {{0xbacd747f, 0x23e32ed7 }}, /* -1.9034889620193266e-25 */ -/**/ hpinv = {{0x3FE45F30, 0x6DC9C883 }}, /* 0.63661977236758138 */ -/**/ toint = {{0x43380000, 0x00000000 }}; /* 6755399441055744 */ - -#else -#ifdef LITTLE_ENDI -static const mynumber - -/**/ s1 = {{0x55555555, 0xBFC55555 }},/* -0.16666666666666666 */ -/**/ s2 = {{0x11110ECE, 0x3F811111 }},/* 0.0083333333333323288 */ -/**/ s3 = {{0x19DB08B8, 0xBF2A01A0 }},/* -0.00019841269834414642 */ -/**/ s4 = {{0x7B9A7ED9, 0x3EC71DE2 }},/* 2.755729806860771e-06 */ -/**/ s5 = {{0xC2FCDF59, 0xBE5ADDFF }},/* -2.5022014848318398e-08 */ -/**/ aa = {{0x00000000, 0xBFC55580 }},/* -0.1666717529296875 */ -/**/ bb = {{0x55556E24, 0x3ED55555 }},/* 5.0862630208387126e-06 */ -/**/ big = {{0x00000000, 0x42c80000 }},/* 52776558133248 */ -/**/ hp0 = {{0x54442D18, 0x3FF921FB }},/* 1.5707963267948966 */ -/**/ hp1 = {{0x33145C07, 0x3C91A626 }},/* 6.123233995736766e-17 */ -/**/ mp1 = {{0x58000000, 0x3FF921FB }},/* 1.5707963407039642 */ -/**/ mp2 = {{0x3C000000, 0xBE4DDE97 }},/* -1.3909067564377153e-08 */ -/**/ mp3 = {{0x99D747F2, 0xBC8CB3B3 }},/* -4.9789962505147994e-17 */ -/**/ pp3 = {{0x98000000, 0xBC8CB3B3 }},/* -4.9789962314799099e-17 */ -/**/ pp4 = {{0x23e32ed7, 0xbacd747f }},/* -1.9034889620193266e-25 */ -/**/ hpinv = {{0x6DC9C883, 0x3FE45F30 }},/* 0.63661977236758138 */ -/**/ toint = {{0x00000000, 0x43380000 }};/* 6755399441055744 */ - - -#endif -#endif +static const double s1 = -0x1.5555555555555p-3; /* -0.16666666666666666 */ +static const double s2 = 0x1.1111111110ECEp-7; /* 0.0083333333333323288 */ +static const double s3 = -0x1.A01A019DB08B8p-13; /* -0.00019841269834414642 */ +static const double s4 = 0x1.71DE27B9A7ED9p-19; /* 2.755729806860771e-06 */ +static const double s5 = -0x1.ADDFFC2FCDF59p-26; /* -2.5022014848318398e-08 */ +static const double aa = -0x1.5558000000000p-3; /* -0.1666717529296875 */ +static const double bb = 0x1.5555555556E24p-18; /* 5.0862630208387126e-06 */ +static const double big = 0x1.8000000000000p45; /* 52776558133248 */ +static const double hp0 = 0x1.921FB54442D18p0; /* 1.5707963267948966 */ +static const double hp1 = 0x1.1A62633145C07p-54; /* 6.123233995736766e-17 */ +static const double mp1 = 0x1.921FB58000000p0; /* 1.5707963407039642 */ +static const double mp2 = -0x1.DDE973C000000p-27; /* -1.3909067564377153e-08 */ +static const double mp3 = -0x1.CB3B399D747F2p-55; /* -4.9789962505147994e-17 */ +static const double pp3 = -0x1.CB3B398000000p-55; /* -4.9789962314799099e-17 */ +static const double pp4 = -0x1.d747f23e32ed7p-83; /* -1.9034889620193266e-25 */ +static const double hpinv = 0x1.45F306DC9C883p-1; /* 0.63661977236758138 */ +static const double toint = 0x1.8000000000000p52; /* 6755399441055744 */ #endif diff --git a/libc/sysdeps/ieee754/dbl-64/utan.h b/libc/sysdeps/ieee754/dbl-64/utan.h index 3bdeee1c4..374132be5 100644 --- a/libc/sysdeps/ieee754/dbl-64/utan.h +++ b/libc/sysdeps/ieee754/dbl-64/utan.h @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/utan.tbl b/libc/sysdeps/ieee754/dbl-64/utan.tbl index 271f6fbeb..25c59c275 100644 --- a/libc/sysdeps/ieee754/dbl-64/utan.tbl +++ b/libc/sysdeps/ieee754/dbl-64/utan.tbl @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * Written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/dbl-64/w_exp.c b/libc/sysdeps/ieee754/dbl-64/w_exp.c index 4f6da6c15..3c2c3e5f4 100644 --- a/libc/sysdeps/ieee754/dbl-64/w_exp.c +++ b/libc/sysdeps/ieee754/dbl-64/w_exp.c @@ -1,4 +1,4 @@ -/* Copyright (C) 2011-2013 Free Software Foundation, Inc. +/* Copyright (C) 2011-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@gmail.com>, 2011. diff --git a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_floor.c b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_floor.c index 914a3c823..3823a53ad 100644 --- a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_floor.c +++ b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_floor.c @@ -1,5 +1,5 @@ /* Round double to integer away from zero. - Copyright (C) 2011-2013 Free Software Foundation, Inc. + Copyright (C) 2011-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 2011. diff --git a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_frexp.c b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_frexp.c index 59dafe8e1..d03e33ee1 100644 --- a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_frexp.c +++ b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_frexp.c @@ -1,4 +1,4 @@ -/* Copyright (C) 2011-2013 Free Software Foundation, Inc. +/* Copyright (C) 2011-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@gmail.com>, 2011. diff --git a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_issignaling.c b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_issignaling.c index 3fa025437..3b8be8430 100644 --- a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_issignaling.c +++ b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_issignaling.c @@ -1,5 +1,5 @@ /* Test for signaling NaN. - Copyright (C) 2013 Free Software Foundation, Inc. + Copyright (C) 2013-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_llround.c b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_llround.c index d9214298a..da180fef5 100644 --- a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_llround.c +++ b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_llround.c @@ -1,5 +1,5 @@ /* Round double value to long long int. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_logb.c b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_logb.c index 42349bd90..e51c84959 100644 --- a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_logb.c +++ b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_logb.c @@ -1,5 +1,5 @@ /* Compute radix independent exponent. - Copyright (C) 2011-2013 Free Software Foundation, Inc. + Copyright (C) 2011-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@gmail.com>, 2011. diff --git a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_remquo.c b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_remquo.c index 29e62874b..b22503f5c 100644 --- a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_remquo.c +++ b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_remquo.c @@ -1,5 +1,5 @@ /* Compute remainder and a congruent to the quotient. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_round.c b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_round.c index bea796083..684858c06 100644 --- a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_round.c +++ b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_round.c @@ -1,5 +1,5 @@ /* Round double to integer away from zero. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_trunc.c b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_trunc.c index dcdf569d6..768e36fdd 100644 --- a/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_trunc.c +++ b/libc/sysdeps/ieee754/dbl-64/wordsize-64/s_trunc.c @@ -1,5 +1,5 @@ /* Truncate argument to nearest integral value not larger than the argument. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/dbl-64/x2y2m1.c b/libc/sysdeps/ieee754/dbl-64/x2y2m1.c index d36a950e3..d39738fc2 100644 --- a/libc/sysdeps/ieee754/dbl-64/x2y2m1.c +++ b/libc/sysdeps/ieee754/dbl-64/x2y2m1.c @@ -1,5 +1,5 @@ /* Compute x^2 + y^2 - 1, without large cancellation error. - Copyright (C) 2012-2013 Free Software Foundation, Inc. + Copyright (C) 2012-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/dbl-64/x2y2m1f.c b/libc/sysdeps/ieee754/dbl-64/x2y2m1f.c index bd1444dad..fd80b9f51 100644 --- a/libc/sysdeps/ieee754/dbl-64/x2y2m1f.c +++ b/libc/sysdeps/ieee754/dbl-64/x2y2m1f.c @@ -1,5 +1,5 @@ /* Compute x^2 + y^2 - 1, without large cancellation error. - Copyright (C) 2012-2013 Free Software Foundation, Inc. + Copyright (C) 2012-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/flt-32/e_atanhf.c b/libc/sysdeps/ieee754/flt-32/e_atanhf.c index e87e34ec7..dbd3fe4be 100644 --- a/libc/sysdeps/ieee754/flt-32/e_atanhf.c +++ b/libc/sysdeps/ieee754/flt-32/e_atanhf.c @@ -1,4 +1,4 @@ -/* Copyright (C) 2011-2013 Free Software Foundation, Inc. +/* Copyright (C) 2011-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@gmail.com>, 2011. diff --git a/libc/sysdeps/ieee754/flt-32/e_exp2f.c b/libc/sysdeps/ieee754/flt-32/e_exp2f.c index d586979f6..544226a63 100644 --- a/libc/sysdeps/ieee754/flt-32/e_exp2f.c +++ b/libc/sysdeps/ieee754/flt-32/e_exp2f.c @@ -1,5 +1,5 @@ /* Single-precision floating point 2^x. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Geoffrey Keating <geoffk@ozemail.com.au> diff --git a/libc/sysdeps/ieee754/flt-32/e_expf.c b/libc/sysdeps/ieee754/flt-32/e_expf.c index 2dc42606f..678b97cf4 100644 --- a/libc/sysdeps/ieee754/flt-32/e_expf.c +++ b/libc/sysdeps/ieee754/flt-32/e_expf.c @@ -1,5 +1,5 @@ /* Single-precision floating point e^x. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Geoffrey Keating <geoffk@ozemail.com.au> diff --git a/libc/sysdeps/ieee754/flt-32/e_gammaf_r.c b/libc/sysdeps/ieee754/flt-32/e_gammaf_r.c index f58f4c805..e8da51a42 100644 --- a/libc/sysdeps/ieee754/flt-32/e_gammaf_r.c +++ b/libc/sysdeps/ieee754/flt-32/e_gammaf_r.c @@ -1,5 +1,5 @@ /* Implementation of gamma function according to ISO C. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/flt-32/mpn2flt.c b/libc/sysdeps/ieee754/flt-32/mpn2flt.c index 9e7365dba..90fbfa1ac 100644 --- a/libc/sysdeps/ieee754/flt-32/mpn2flt.c +++ b/libc/sysdeps/ieee754/flt-32/mpn2flt.c @@ -1,4 +1,4 @@ -/* Copyright (C) 1995-2013 Free Software Foundation, Inc. +/* Copyright (C) 1995-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/flt-32/s_cbrtf.c b/libc/sysdeps/ieee754/flt-32/s_cbrtf.c index d0dcba08b..7d94dfe1f 100644 --- a/libc/sysdeps/ieee754/flt-32/s_cbrtf.c +++ b/libc/sysdeps/ieee754/flt-32/s_cbrtf.c @@ -1,5 +1,5 @@ /* Compute cubic root of float value. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Dirk Alboth <dirka@uni-paderborn.de> and Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/flt-32/s_erff.c b/libc/sysdeps/ieee754/flt-32/s_erff.c index 7d17f426e..7c0958964 100644 --- a/libc/sysdeps/ieee754/flt-32/s_erff.c +++ b/libc/sysdeps/ieee754/flt-32/s_erff.c @@ -17,6 +17,8 @@ static char rcsid[] = "$NetBSD: s_erff.c,v 1.4 1995/05/10 20:47:07 jtc Exp $"; #endif +#include <errno.h> +#include <float.h> #include <math.h> #include <math_private.h> @@ -203,9 +205,22 @@ float __erfcf(float x) SET_FLOAT_WORD(z,ix&0xffffe000); r = __ieee754_expf(-z*z-(float)0.5625)* __ieee754_expf((z-x)*(z+x)+R/S); - if(hx>0) return r/x; else return two-r/x; + if(hx>0) { +#if FLT_EVAL_METHOD != 0 + volatile +#endif + float ret = r/x; + if (ret == 0) + __set_errno (ERANGE); + return ret; + } else + return two-r/x; } else { - if(hx>0) return tiny*tiny; else return two-tiny; + if(hx>0) { + __set_errno (ERANGE); + return tiny*tiny; + } else + return two-tiny; } } weak_alias (__erfcf, erfcf) diff --git a/libc/sysdeps/ieee754/flt-32/s_finitef.c b/libc/sysdeps/ieee754/flt-32/s_finitef.c index dfdf4ad0b..4ea270ae0 100644 --- a/libc/sysdeps/ieee754/flt-32/s_finitef.c +++ b/libc/sysdeps/ieee754/flt-32/s_finitef.c @@ -26,7 +26,12 @@ static char rcsid[] = "$NetBSD: s_finitef.c,v 1.4 1995/05/10 20:47:18 jtc Exp $" #include <math_private.h> #undef __finitef -int __finitef(float x) + +#ifndef FINITEF +# define FINITEF __finitef +#endif + +int FINITEF(float x) { int32_t ix; GET_FLOAT_WORD(ix,x); diff --git a/libc/sysdeps/ieee754/flt-32/s_fpclassifyf.c b/libc/sysdeps/ieee754/flt-32/s_fpclassifyf.c index dd5c74c74..5d374c250 100644 --- a/libc/sysdeps/ieee754/flt-32/s_fpclassifyf.c +++ b/libc/sysdeps/ieee754/flt-32/s_fpclassifyf.c @@ -1,5 +1,5 @@ /* Return classification value corresponding to argument. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/flt-32/s_issignalingf.c b/libc/sysdeps/ieee754/flt-32/s_issignalingf.c index 59a189252..ac5e7292b 100644 --- a/libc/sysdeps/ieee754/flt-32/s_issignalingf.c +++ b/libc/sysdeps/ieee754/flt-32/s_issignalingf.c @@ -1,5 +1,5 @@ /* Test for signaling NaN. - Copyright (C) 2013 Free Software Foundation, Inc. + Copyright (C) 2013-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/flt-32/s_llrintf.c b/libc/sysdeps/ieee754/flt-32/s_llrintf.c index fee0b1a6c..d0625da6d 100644 --- a/libc/sysdeps/ieee754/flt-32/s_llrintf.c +++ b/libc/sysdeps/ieee754/flt-32/s_llrintf.c @@ -1,6 +1,6 @@ /* Round argument to nearest integral value according to current rounding direction. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/flt-32/s_llroundf.c b/libc/sysdeps/ieee754/flt-32/s_llroundf.c index 4ea7b38ba..36c66dc76 100644 --- a/libc/sysdeps/ieee754/flt-32/s_llroundf.c +++ b/libc/sysdeps/ieee754/flt-32/s_llroundf.c @@ -1,5 +1,5 @@ /* Round float value to long long int. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/flt-32/s_lrintf.c b/libc/sysdeps/ieee754/flt-32/s_lrintf.c index 9f3d5f084..77f150556 100644 --- a/libc/sysdeps/ieee754/flt-32/s_lrintf.c +++ b/libc/sysdeps/ieee754/flt-32/s_lrintf.c @@ -1,6 +1,6 @@ /* Round argument to nearest integral value according to current rounding direction. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/flt-32/s_lroundf.c b/libc/sysdeps/ieee754/flt-32/s_lroundf.c index 00b874fcb..b5e063898 100644 --- a/libc/sysdeps/ieee754/flt-32/s_lroundf.c +++ b/libc/sysdeps/ieee754/flt-32/s_lroundf.c @@ -1,5 +1,5 @@ /* Round float value to long int. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/flt-32/s_remquof.c b/libc/sysdeps/ieee754/flt-32/s_remquof.c index d55704c44..c98f9aa22 100644 --- a/libc/sysdeps/ieee754/flt-32/s_remquof.c +++ b/libc/sysdeps/ieee754/flt-32/s_remquof.c @@ -1,5 +1,5 @@ /* Compute remainder and a congruent to the quotient. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/flt-32/s_roundf.c b/libc/sysdeps/ieee754/flt-32/s_roundf.c index 5ee3d29d2..58ec5bcbf 100644 --- a/libc/sysdeps/ieee754/flt-32/s_roundf.c +++ b/libc/sysdeps/ieee754/flt-32/s_roundf.c @@ -1,5 +1,5 @@ /* Round float to integer away from zero. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/flt-32/s_signbitf.c b/libc/sysdeps/ieee754/flt-32/s_signbitf.c index a42bfb976..1cbb2d029 100644 --- a/libc/sysdeps/ieee754/flt-32/s_signbitf.c +++ b/libc/sysdeps/ieee754/flt-32/s_signbitf.c @@ -1,5 +1,5 @@ /* Return nonzero value if number is negative. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/flt-32/s_sincosf.c b/libc/sysdeps/ieee754/flt-32/s_sincosf.c index 1b3fd9fd5..c86fd5721 100644 --- a/libc/sysdeps/ieee754/flt-32/s_sincosf.c +++ b/libc/sysdeps/ieee754/flt-32/s_sincosf.c @@ -1,5 +1,5 @@ /* Compute sine and cosine of argument. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/flt-32/s_truncf.c b/libc/sysdeps/ieee754/flt-32/s_truncf.c index cc8e5f3af..f1acf179c 100644 --- a/libc/sysdeps/ieee754/flt-32/s_truncf.c +++ b/libc/sysdeps/ieee754/flt-32/s_truncf.c @@ -1,5 +1,5 @@ /* Truncate argument to nearest integral value not larger than the argument. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/flt-32/t_exp2f.h b/libc/sysdeps/ieee754/flt-32/t_exp2f.h index 51ac6f5ec..6f84cfce0 100644 --- a/libc/sysdeps/ieee754/flt-32/t_exp2f.h +++ b/libc/sysdeps/ieee754/flt-32/t_exp2f.h @@ -1,5 +1,5 @@ /* Accurate tables for exp2f(). - Copyright (C) 1998-2013 Free Software Foundation, Inc. + Copyright (C) 1998-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Geoffrey Keating <geoffk@ozemail.com.au> diff --git a/libc/sysdeps/ieee754/flt-32/w_expf.c b/libc/sysdeps/ieee754/flt-32/w_expf.c index 1fbfe6a96..cf5cae5f9 100644 --- a/libc/sysdeps/ieee754/flt-32/w_expf.c +++ b/libc/sysdeps/ieee754/flt-32/w_expf.c @@ -1,4 +1,4 @@ -/* Copyright (C) 2011-2013 Free Software Foundation, Inc. +/* Copyright (C) 2011-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@gmail.com>, 2011. diff --git a/libc/sysdeps/ieee754/ieee754.h b/libc/sysdeps/ieee754/ieee754.h index 4c694bb95..48d9360ef 100644 --- a/libc/sysdeps/ieee754/ieee754.h +++ b/libc/sysdeps/ieee754/ieee754.h @@ -1,4 +1,4 @@ -/* Copyright (C) 1992-2013 Free Software Foundation, Inc. +/* Copyright (C) 1992-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/k_standard.c b/libc/sysdeps/ieee754/k_standard.c index 4a0d82d1a..5399c6682 100644 --- a/libc/sysdeps/ieee754/k_standard.c +++ b/libc/sysdeps/ieee754/k_standard.c @@ -248,7 +248,7 @@ __kernel_standard(double x, double y, int type) else exc.retval = -HUGE_VAL; if (_LIB_VERSION == _POSIX_) - __set_errno (EDOM); + __set_errno (ERANGE); else if (!matherr(&exc)) { if (_LIB_VERSION == _SVID_) { (void) WRITE2("y0: DOMAIN error\n", 17); @@ -265,7 +265,7 @@ __kernel_standard(double x, double y, int type) if (_LIB_VERSION == _SVID_) exc.retval = -HUGE; else - exc.retval = -HUGE_VAL; + exc.retval = NAN; if (_LIB_VERSION == _POSIX_) __set_errno (EDOM); else if (!matherr(&exc)) { @@ -286,7 +286,7 @@ __kernel_standard(double x, double y, int type) else exc.retval = -HUGE_VAL; if (_LIB_VERSION == _POSIX_) - __set_errno (EDOM); + __set_errno (ERANGE); else if (!matherr(&exc)) { if (_LIB_VERSION == _SVID_) { (void) WRITE2("y1: DOMAIN error\n", 17); @@ -303,7 +303,7 @@ __kernel_standard(double x, double y, int type) if (_LIB_VERSION == _SVID_) exc.retval = -HUGE; else - exc.retval = -HUGE_VAL; + exc.retval = NAN; if (_LIB_VERSION == _POSIX_) __set_errno (EDOM); else if (!matherr(&exc)) { @@ -322,9 +322,11 @@ __kernel_standard(double x, double y, int type) if (_LIB_VERSION == _SVID_) exc.retval = -HUGE; else - exc.retval = -HUGE_VAL; + exc.retval = ((x < 0 && ((int) x & 1) != 0) + ? HUGE_VAL + : -HUGE_VAL); if (_LIB_VERSION == _POSIX_) - __set_errno (EDOM); + __set_errno (ERANGE); else if (!matherr(&exc)) { if (_LIB_VERSION == _SVID_) { (void) WRITE2("yn: DOMAIN error\n", 17); @@ -341,7 +343,7 @@ __kernel_standard(double x, double y, int type) if (_LIB_VERSION == _SVID_) exc.retval = -HUGE; else - exc.retval = -HUGE_VAL; + exc.retval = NAN; if (_LIB_VERSION == _POSIX_) __set_errno (EDOM); else if (!matherr(&exc)) { diff --git a/libc/sysdeps/ieee754/ldbl-128/e_exp10l.c b/libc/sysdeps/ieee754/ldbl-128/e_exp10l.c index e6f4a8679..3c5096379 100644 --- a/libc/sysdeps/ieee754/ldbl-128/e_exp10l.c +++ b/libc/sysdeps/ieee754/ldbl-128/e_exp10l.c @@ -1,4 +1,4 @@ -/* Copyright (C) 2012-2013 Free Software Foundation, Inc. +/* Copyright (C) 2012-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128/e_expl.c b/libc/sysdeps/ieee754/ldbl-128/e_expl.c index 23c204ce0..8259758cf 100644 --- a/libc/sysdeps/ieee754/ldbl-128/e_expl.c +++ b/libc/sysdeps/ieee754/ldbl-128/e_expl.c @@ -1,5 +1,5 @@ /* Quad-precision floating point e^x. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jj@ultra.linux.cz> Partly based on double-precision code diff --git a/libc/sysdeps/ieee754/ldbl-128/e_gammal_r.c b/libc/sysdeps/ieee754/ldbl-128/e_gammal_r.c index e8d49e987..ffa27bc2b 100644 --- a/libc/sysdeps/ieee754/ldbl-128/e_gammal_r.c +++ b/libc/sysdeps/ieee754/ldbl-128/e_gammal_r.c @@ -1,5 +1,5 @@ /* Implementation of gamma function according to ISO C. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and Jakub Jelinek <jj@ultra.linux.cz, 1999. diff --git a/libc/sysdeps/ieee754/ldbl-128/e_hypotl.c b/libc/sysdeps/ieee754/ldbl-128/e_hypotl.c index f5ab901e6..01444cfb4 100644 --- a/libc/sysdeps/ieee754/ldbl-128/e_hypotl.c +++ b/libc/sysdeps/ieee754/ldbl-128/e_hypotl.c @@ -89,6 +89,17 @@ __ieee754_hypotl(long double x, long double y) b *= t1; a *= t1; k -= 16382; + GET_LDOUBLE_MSW64 (ha, a); + GET_LDOUBLE_MSW64 (hb, b); + if (hb > ha) + { + t1 = a; + a = b; + b = t1; + j = ha; + ha = hb; + hb = j; + } } else { /* scale a and b by 2^9600 */ ha += 0x2580000000000000LL; /* a *= 2^9600 */ hb += 0x2580000000000000LL; /* b *= 2^9600 */ diff --git a/libc/sysdeps/ieee754/ldbl-128/e_jnl.c b/libc/sysdeps/ieee754/ldbl-128/e_jnl.c index 70d5672fd..c2a49235c 100644 --- a/libc/sysdeps/ieee754/ldbl-128/e_jnl.c +++ b/libc/sysdeps/ieee754/ldbl-128/e_jnl.c @@ -316,7 +316,7 @@ __ieee754_ynl (int n, long double x) if (x <= 0.0L) { if (x == 0.0L) - return -HUGE_VALL + x; + return ((n < 0 && (n & 1) != 0) ? 1.0L : -1.0L) / 0.0L; if (se & 0x80000000) return zero / (zero * x); } diff --git a/libc/sysdeps/ieee754/ldbl-128/e_lgammal_r.c b/libc/sysdeps/ieee754/ldbl-128/e_lgammal_r.c index 2b44afb75..23ab9b9a4 100644 --- a/libc/sysdeps/ieee754/ldbl-128/e_lgammal_r.c +++ b/libc/sysdeps/ieee754/ldbl-128/e_lgammal_r.c @@ -782,6 +782,8 @@ __ieee754_lgammal_r (long double x, int *signgamp) *signgamp = -1; else *signgamp = 1; + if (q < 0x1p-120L) + return -__logl (q); z = q - p; if (z > 0.5L) { @@ -789,8 +791,6 @@ __ieee754_lgammal_r (long double x, int *signgamp) z = p - q; } z = q * __sinl (PIL * z); - if (z == 0.0L) - return (*signgamp * huge * huge); w = __ieee754_lgammal_r (q, &i); z = __logl (PIL / z) - w; return (z); diff --git a/libc/sysdeps/ieee754/ldbl-128/e_logl.c b/libc/sysdeps/ieee754/ldbl-128/e_logl.c index 395a76302..3d1034dd6 100644 --- a/libc/sysdeps/ieee754/ldbl-128/e_logl.c +++ b/libc/sysdeps/ieee754/ldbl-128/e_logl.c @@ -212,9 +212,8 @@ __ieee754_logl(long double x) } /* Extract exponent and reduce domain to 0.703125 <= u < 1.40625 */ - e = (int) (m >> 16) - (int) 0x3ffe; - m &= 0xffff; - u.parts32.w0 = m | 0x3ffe0000; + u.value = __frexpl (x, &e); + m = u.parts32.w0 & 0xffff; m |= 0x10000; /* Find lookup table index k from high order bits of the significand. */ if (m < 0x16800) diff --git a/libc/sysdeps/ieee754/ldbl-128/e_rem_pio2l.c b/libc/sysdeps/ieee754/ldbl-128/e_rem_pio2l.c index 644afc4b7..81c2dd894 100644 --- a/libc/sysdeps/ieee754/ldbl-128/e_rem_pio2l.c +++ b/libc/sysdeps/ieee754/ldbl-128/e_rem_pio2l.c @@ -1,5 +1,5 @@ /* Quad-precision floating point argument reduction. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jj@ultra.linux.cz> diff --git a/libc/sysdeps/ieee754/ldbl-128/gamma_productl.c b/libc/sysdeps/ieee754/ldbl-128/gamma_productl.c index 157dbab9f..2c8764967 100644 --- a/libc/sysdeps/ieee754/ldbl-128/gamma_productl.c +++ b/libc/sysdeps/ieee754/ldbl-128/gamma_productl.c @@ -1,5 +1,5 @@ /* Compute a product of X, X+1, ..., with an error estimate. - Copyright (C) 2013 Free Software Foundation, Inc. + Copyright (C) 2013-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128/ieee754.h b/libc/sysdeps/ieee754/ldbl-128/ieee754.h index 2324f0b30..d77b83592 100644 --- a/libc/sysdeps/ieee754/ldbl-128/ieee754.h +++ b/libc/sysdeps/ieee754/ldbl-128/ieee754.h @@ -1,4 +1,4 @@ -/* Copyright (C) 1992-2013 Free Software Foundation, Inc. +/* Copyright (C) 1992-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128/k_cosl.c b/libc/sysdeps/ieee754/ldbl-128/k_cosl.c index 997f7f115..0c76332bc 100644 --- a/libc/sysdeps/ieee754/ldbl-128/k_cosl.c +++ b/libc/sysdeps/ieee754/ldbl-128/k_cosl.c @@ -1,5 +1,5 @@ /* Quad-precision floating point cosine on <-pi/4,pi/4>. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jj@ultra.linux.cz> diff --git a/libc/sysdeps/ieee754/ldbl-128/k_sincosl.c b/libc/sysdeps/ieee754/ldbl-128/k_sincosl.c index 4f18deb6f..038718a29 100644 --- a/libc/sysdeps/ieee754/ldbl-128/k_sincosl.c +++ b/libc/sysdeps/ieee754/ldbl-128/k_sincosl.c @@ -1,5 +1,5 @@ /* Quad-precision floating point sine and cosine on <-pi/4,pi/4>. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jj@ultra.linux.cz> diff --git a/libc/sysdeps/ieee754/ldbl-128/k_sinl.c b/libc/sysdeps/ieee754/ldbl-128/k_sinl.c index e36a4024f..1c9c7c717 100644 --- a/libc/sysdeps/ieee754/ldbl-128/k_sinl.c +++ b/libc/sysdeps/ieee754/ldbl-128/k_sinl.c @@ -1,5 +1,5 @@ /* Quad-precision floating point sine on <-pi/4,pi/4>. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jj@ultra.linux.cz> diff --git a/libc/sysdeps/ieee754/ldbl-128/ldbl2mpn.c b/libc/sysdeps/ieee754/ldbl-128/ldbl2mpn.c index 64f98a59d..6d6ee2a55 100644 --- a/libc/sysdeps/ieee754/ldbl-128/ldbl2mpn.c +++ b/libc/sysdeps/ieee754/ldbl-128/ldbl2mpn.c @@ -1,4 +1,4 @@ -/* Copyright (C) 1995-2013 Free Software Foundation, Inc. +/* Copyright (C) 1995-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128/mpn2ldbl.c b/libc/sysdeps/ieee754/ldbl-128/mpn2ldbl.c index eabbdd229..6d17567fa 100644 --- a/libc/sysdeps/ieee754/ldbl-128/mpn2ldbl.c +++ b/libc/sysdeps/ieee754/ldbl-128/mpn2ldbl.c @@ -1,4 +1,4 @@ -/* Copyright (C) 1995-2013 Free Software Foundation, Inc. +/* Copyright (C) 1995-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128/printf_fphex.c b/libc/sysdeps/ieee754/ldbl-128/printf_fphex.c index e82228a53..5bd161510 100644 --- a/libc/sysdeps/ieee754/ldbl-128/printf_fphex.c +++ b/libc/sysdeps/ieee754/ldbl-128/printf_fphex.c @@ -1,6 +1,6 @@ /* Print floating point number in hexadecimal notation according to ISO C99. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128/s_erfl.c b/libc/sysdeps/ieee754/ldbl-128/s_erfl.c index d9d4195ad..ef65ed892 100644 --- a/libc/sysdeps/ieee754/ldbl-128/s_erfl.c +++ b/libc/sysdeps/ieee754/ldbl-128/s_erfl.c @@ -96,6 +96,7 @@ * erfc/erf(NaN) is NaN */ +#include <errno.h> #include <math.h> #include <math_private.h> @@ -918,14 +919,22 @@ __erfcl (long double x) r = __ieee754_expl (-z * z - 0.5625) * __ieee754_expl ((z - x) * (z + x) + p); if ((sign & 0x80000000) == 0) - return r / x; + { + long double ret = r / x; + if (ret == 0) + __set_errno (ERANGE); + return ret; + } else return two - r / x; } else { if ((sign & 0x80000000) == 0) - return tiny * tiny; + { + __set_errno (ERANGE); + return tiny * tiny; + } else return two - tiny; } diff --git a/libc/sysdeps/ieee754/ldbl-128/s_fma.c b/libc/sysdeps/ieee754/ldbl-128/s_fma.c index 8cc7db81c..e3cae61cf 100644 --- a/libc/sysdeps/ieee754/ldbl-128/s_fma.c +++ b/libc/sysdeps/ieee754/ldbl-128/s_fma.c @@ -1,5 +1,5 @@ /* Compute x * y + z as ternary operation. - Copyright (C) 2010-2013 Free Software Foundation, Inc. + Copyright (C) 2010-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jakub@redhat.com>, 2010. diff --git a/libc/sysdeps/ieee754/ldbl-128/s_fmal.c b/libc/sysdeps/ieee754/ldbl-128/s_fmal.c index 5ec0a0618..48b63ab51 100644 --- a/libc/sysdeps/ieee754/ldbl-128/s_fmal.c +++ b/libc/sysdeps/ieee754/ldbl-128/s_fmal.c @@ -1,5 +1,5 @@ /* Compute x * y + z as ternary operation. - Copyright (C) 2010-2013 Free Software Foundation, Inc. + Copyright (C) 2010-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jakub@redhat.com>, 2010. diff --git a/libc/sysdeps/ieee754/ldbl-128/s_fpclassifyl.c b/libc/sysdeps/ieee754/ldbl-128/s_fpclassifyl.c index dc57ba9ca..efe69f2d4 100644 --- a/libc/sysdeps/ieee754/ldbl-128/s_fpclassifyl.c +++ b/libc/sysdeps/ieee754/ldbl-128/s_fpclassifyl.c @@ -1,5 +1,5 @@ /* Return classification value corresponding to argument. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and Jakub Jelinek <jj@ultra.linux.cz>, 1999. diff --git a/libc/sysdeps/ieee754/ldbl-128/s_issignalingl.c b/libc/sysdeps/ieee754/ldbl-128/s_issignalingl.c index b69a66d0b..c3d77aff7 100644 --- a/libc/sysdeps/ieee754/ldbl-128/s_issignalingl.c +++ b/libc/sysdeps/ieee754/ldbl-128/s_issignalingl.c @@ -1,5 +1,5 @@ /* Test for signaling NaN. - Copyright (C) 2013 Free Software Foundation, Inc. + Copyright (C) 2013-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128/s_llrintl.c b/libc/sysdeps/ieee754/ldbl-128/s_llrintl.c index 9d56cf120..b3a2124f2 100644 --- a/libc/sysdeps/ieee754/ldbl-128/s_llrintl.c +++ b/libc/sysdeps/ieee754/ldbl-128/s_llrintl.c @@ -1,6 +1,6 @@ /* Round argument to nearest integral value according to current rounding direction. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and Jakub Jelinek <jj@ultra.linux.cz>, 1999. diff --git a/libc/sysdeps/ieee754/ldbl-128/s_llroundl.c b/libc/sysdeps/ieee754/ldbl-128/s_llroundl.c index c2136d990..8c2b48ea2 100644 --- a/libc/sysdeps/ieee754/ldbl-128/s_llroundl.c +++ b/libc/sysdeps/ieee754/ldbl-128/s_llroundl.c @@ -1,5 +1,5 @@ /* Round long double value to long long int. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and Jakub Jelinek <jj@ultra.linux.cz>, 1999. diff --git a/libc/sysdeps/ieee754/ldbl-128/s_lrintl.c b/libc/sysdeps/ieee754/ldbl-128/s_lrintl.c index 61c508f20..7dbab5cc1 100644 --- a/libc/sysdeps/ieee754/ldbl-128/s_lrintl.c +++ b/libc/sysdeps/ieee754/ldbl-128/s_lrintl.c @@ -1,6 +1,6 @@ /* Round argument to nearest integral value according to current rounding direction. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and Jakub Jelinek <jj@ultra.linux.cz>, 1999. diff --git a/libc/sysdeps/ieee754/ldbl-128/s_lroundl.c b/libc/sysdeps/ieee754/ldbl-128/s_lroundl.c index 68256ecbe..493592c4d 100644 --- a/libc/sysdeps/ieee754/ldbl-128/s_lroundl.c +++ b/libc/sysdeps/ieee754/ldbl-128/s_lroundl.c @@ -1,5 +1,5 @@ /* Round long double value to long int. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and Jakub Jelinek <jj@ultra.linux.cz>, 1999. diff --git a/libc/sysdeps/ieee754/ldbl-128/s_remquol.c b/libc/sysdeps/ieee754/ldbl-128/s_remquol.c index f916b64e4..d7b150385 100644 --- a/libc/sysdeps/ieee754/ldbl-128/s_remquol.c +++ b/libc/sysdeps/ieee754/ldbl-128/s_remquol.c @@ -1,5 +1,5 @@ /* Compute remainder and a congruent to the quotient. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and Jakub Jelinek <jj@ultra.linux.cz>, 1999. diff --git a/libc/sysdeps/ieee754/ldbl-128/s_roundl.c b/libc/sysdeps/ieee754/ldbl-128/s_roundl.c index 2445b76c8..4be1ad822 100644 --- a/libc/sysdeps/ieee754/ldbl-128/s_roundl.c +++ b/libc/sysdeps/ieee754/ldbl-128/s_roundl.c @@ -1,5 +1,5 @@ /* Round long double to integer away from zero. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and Jakub Jelinek <jj@ultra.linux.cz>, 1999. diff --git a/libc/sysdeps/ieee754/ldbl-128/s_signbitl.c b/libc/sysdeps/ieee754/ldbl-128/s_signbitl.c index e5fc78f95..434bc2ce6 100644 --- a/libc/sysdeps/ieee754/ldbl-128/s_signbitl.c +++ b/libc/sysdeps/ieee754/ldbl-128/s_signbitl.c @@ -1,5 +1,5 @@ /* Return nonzero value if number is negative. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-128/s_sincosl.c b/libc/sysdeps/ieee754/ldbl-128/s_sincosl.c index 260373bd3..591867acd 100644 --- a/libc/sysdeps/ieee754/ldbl-128/s_sincosl.c +++ b/libc/sysdeps/ieee754/ldbl-128/s_sincosl.c @@ -1,5 +1,5 @@ /* Compute sine and cosine of argument. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and Jakub Jelinek <jj@ultra.linux.cz>. diff --git a/libc/sysdeps/ieee754/ldbl-128/s_truncl.c b/libc/sysdeps/ieee754/ldbl-128/s_truncl.c index cd9f44802..2395b5a15 100644 --- a/libc/sysdeps/ieee754/ldbl-128/s_truncl.c +++ b/libc/sysdeps/ieee754/ldbl-128/s_truncl.c @@ -1,5 +1,5 @@ /* Truncate argument to nearest integral value not larger than the argument. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and Jakub Jelinek <jj@ultra.linux.cz>, 1999. diff --git a/libc/sysdeps/ieee754/ldbl-128/strtold_l.c b/libc/sysdeps/ieee754/ldbl-128/strtold_l.c index d3a1d1e86..b22cfe765 100644 --- a/libc/sysdeps/ieee754/ldbl-128/strtold_l.c +++ b/libc/sysdeps/ieee754/ldbl-128/strtold_l.c @@ -1,4 +1,4 @@ -/* Copyright (C) 1999-2013 Free Software Foundation, Inc. +/* Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128/t_expl.h b/libc/sysdeps/ieee754/ldbl-128/t_expl.h index 002fb8f17..6b94956ff 100644 --- a/libc/sysdeps/ieee754/ldbl-128/t_expl.h +++ b/libc/sysdeps/ieee754/ldbl-128/t_expl.h @@ -1,5 +1,5 @@ /* Accurate table for expl(). - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jj@ultra.linux.cz> diff --git a/libc/sysdeps/ieee754/ldbl-128/t_sincosl.c b/libc/sysdeps/ieee754/ldbl-128/t_sincosl.c index 95fe01ec9..b04321dc1 100644 --- a/libc/sysdeps/ieee754/ldbl-128/t_sincosl.c +++ b/libc/sysdeps/ieee754/ldbl-128/t_sincosl.c @@ -1,5 +1,5 @@ /* Quad-precision floating point sine and cosine tables. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jj@ultra.linux.cz> diff --git a/libc/sysdeps/ieee754/ldbl-128/x2y2m1l.c b/libc/sysdeps/ieee754/ldbl-128/x2y2m1l.c index 7bf9e297b..575bcd812 100644 --- a/libc/sysdeps/ieee754/ldbl-128/x2y2m1l.c +++ b/libc/sysdeps/ieee754/ldbl-128/x2y2m1l.c @@ -1,5 +1,5 @@ /* Compute x^2 + y^2 - 1, without large cancellation error. - Copyright (C) 2012-2013 Free Software Foundation, Inc. + Copyright (C) 2012-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/e_acoshl.c b/libc/sysdeps/ieee754/ldbl-128ibm/e_acoshl.c index 8a4a5bb7b..b0b33f752 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/e_acoshl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/e_acoshl.c @@ -29,7 +29,7 @@ static const long double one = 1.0L, -ln2 = 6.93147180559945286227e-01L; /* 0x3FE62E42, 0xFEFA39EF */ +ln2 = M_LN2l; long double __ieee754_acoshl(long double x) @@ -56,7 +56,7 @@ __ieee754_acoshl(long double x) return __ieee754_logl(2.0*x-one/(x+__ieee754_sqrtl(t-one))); } else { /* 1<x<2 */ t = x-one; - return __log1p(t+__ieee754_sqrtl(2.0*t+t*t)); + return __log1pl(t+__ieee754_sqrtl(2.0*t+t*t)); } } strong_alias (__ieee754_acoshl, __acoshl_finite) diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/e_exp10l.c b/libc/sysdeps/ieee754/ldbl-128ibm/e_exp10l.c index 49121ca31..5c68eb57b 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/e_exp10l.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/e_exp10l.c @@ -1,4 +1,4 @@ -/* Copyright (C) 2012-2013 Free Software Foundation, Inc. +/* Copyright (C) 2012-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/e_expl.c b/libc/sysdeps/ieee754/ldbl-128ibm/e_expl.c index f7c50bfd3..da71bb3b4 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/e_expl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/e_expl.c @@ -1,5 +1,5 @@ /* Quad-precision floating point e^x. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jj@ultra.linux.cz> Partly based on double-precision code @@ -134,18 +134,17 @@ static const long double C[] = { long double __ieee754_expl (long double x) { + long double result, x22; + union ibm_extended_long_double ex2_u, scale_u; + int unsafe; + /* Check for usual case. */ if (isless (x, himark) && isgreater (x, lomark)) { - int tval1, tval2, unsafe, n_i, exponent2; - long double x22, n, result, xl; - union ibm_extended_long_double ex2_u, scale_u; - fenv_t oldenv; - - feholdexcept (&oldenv); -#ifdef FE_TONEAREST - fesetround (FE_TONEAREST); -#endif + int tval1, tval2, n_i, exponent2; + long double n, xl; + + SET_RESTORE_ROUND (FE_TONEAREST); n = __roundl (x*M_1_LN2); x = x-n*M_LN2_0; @@ -201,11 +200,6 @@ __ieee754_expl (long double x) less than 4.8e-39. */ x22 = x + x*x*(P1+x*(P2+x*(P3+x*(P4+x*(P5+x*P6))))); - /* Return result. */ - fesetenv (&oldenv); - - result = x22 * ex2_u.ld + ex2_u.ld; - /* Now we can test whether the result is ultimate or if we are unsure. In the later case we should probably call a mpn based routine to give the ultimate result. @@ -235,10 +229,6 @@ __ieee754_expl (long double x) return __ieee754_expl_proc2 (origx); } */ - if (!unsafe) - return result; - else - return result * scale_u.ld; } /* Exceptional cases: */ else if (isless (x, himark)) @@ -253,5 +243,10 @@ __ieee754_expl (long double x) else /* Return x, if x is a NaN or Inf; or overflow, otherwise. */ return TWO1023*x; + + result = x22 * ex2_u.ld + ex2_u.ld; + if (!unsafe) + return result; + return result * scale_u.ld; } strong_alias (__ieee754_expl, __expl_finite) diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/e_gammal_r.c b/libc/sysdeps/ieee754/ldbl-128ibm/e_gammal_r.c index 84c13de9b..41a71dc86 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/e_gammal_r.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/e_gammal_r.c @@ -1,5 +1,5 @@ /* Implementation of gamma function according to ISO C. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and Jakub Jelinek <jj@ultra.linux.cz, 1999. diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/e_jnl.c b/libc/sysdeps/ieee754/ldbl-128ibm/e_jnl.c index 817977da5..6761a0d26 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/e_jnl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/e_jnl.c @@ -316,7 +316,7 @@ __ieee754_ynl (int n, long double x) if (x <= 0.0L) { if (x == 0.0L) - return -HUGE_VALL + x; + return ((n < 0 && (n & 1) != 0) ? 1.0L : -1.0L) / 0.0L; if (se & 0x80000000) return zero / (zero * x); } diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/e_logl.c b/libc/sysdeps/ieee754/ldbl-128ibm/e_logl.c index b7db2b978..58d6bc697 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/e_logl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/e_logl.c @@ -229,6 +229,14 @@ __ieee754_logl(long double x) /* Extract exponent and reduce domain to 0.703125 <= u < 1.40625 */ unsigned int w0; e = (int) (m >> 20) - (int) 0x3fe; + if (e == -1022) + { + x *= 0x1p106L; + xhi = ldbl_high (x); + EXTRACT_WORDS (hx, lx, xhi); + m = hx; + e = (int) (m >> 20) - (int) 0x3fe - 106; + } m &= 0xfffff; w0 = m | 0x3fe00000; m |= 0x100000; diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/e_rem_pio2l.c b/libc/sysdeps/ieee754/ldbl-128ibm/e_rem_pio2l.c index 36bc03226..c32490a61 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/e_rem_pio2l.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/e_rem_pio2l.c @@ -1,5 +1,5 @@ /* Quad-precision floating point argument reduction. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jj@ultra.linux.cz> diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/e_sqrtl.c b/libc/sysdeps/ieee754/ldbl-128ibm/e_sqrtl.c index 61feb367f..29e873565 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/e_sqrtl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/e_sqrtl.c @@ -1,7 +1,7 @@ /* * IBM Accurate Mathematical Library * written by International Business Machines Corp. - * Copyright (C) 2001-2013 Free Software Foundation, Inc. + * Copyright (C) 2001-2014 Free Software Foundation, Inc. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published by diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/gamma_productl.c b/libc/sysdeps/ieee754/ldbl-128ibm/gamma_productl.c index 7c6186d23..e94338592 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/gamma_productl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/gamma_productl.c @@ -1,5 +1,5 @@ /* Compute a product of X, X+1, ..., with an error estimate. - Copyright (C) 2013 Free Software Foundation, Inc. + Copyright (C) 2013-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/ieee754.h b/libc/sysdeps/ieee754/ldbl-128ibm/ieee754.h index 0c97a9920..52dbe8be9 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/ieee754.h +++ b/libc/sysdeps/ieee754/ldbl-128ibm/ieee754.h @@ -1,4 +1,4 @@ -/* Copyright (C) 1992-2013 Free Software Foundation, Inc. +/* Copyright (C) 1992-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/k_cosl.c b/libc/sysdeps/ieee754/ldbl-128ibm/k_cosl.c index 046f3b573..5409f2b94 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/k_cosl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/k_cosl.c @@ -1,5 +1,5 @@ /* Quad-precision floating point cosine on <-pi/4,pi/4>. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jj@ultra.linux.cz> diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/k_sincosl.c b/libc/sysdeps/ieee754/ldbl-128ibm/k_sincosl.c index 3ba9d7e90..14d135ec1 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/k_sincosl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/k_sincosl.c @@ -1,5 +1,5 @@ /* Quad-precision floating point sine and cosine on <-pi/4,pi/4>. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jj@ultra.linux.cz> diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/k_sinl.c b/libc/sysdeps/ieee754/ldbl-128ibm/k_sinl.c index b12ea134d..0e53a62a8 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/k_sinl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/k_sinl.c @@ -1,5 +1,5 @@ /* Quad-precision floating point sine on <-pi/4,pi/4>. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jj@ultra.linux.cz> diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/ldbl2mpn.c b/libc/sysdeps/ieee754/ldbl-128ibm/ldbl2mpn.c index e46fde74f..5cdfbe1be 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/ldbl2mpn.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/ldbl2mpn.c @@ -1,4 +1,4 @@ -/* Copyright (C) 1995-2013 Free Software Foundation, Inc. +/* Copyright (C) 1995-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/mpn2ldbl.c b/libc/sysdeps/ieee754/ldbl-128ibm/mpn2ldbl.c index c96852dfd..532a119b7 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/mpn2ldbl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/mpn2ldbl.c @@ -1,4 +1,4 @@ -/* Copyright (C) 1995-2013 Free Software Foundation, Inc. +/* Copyright (C) 1995-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/printf_fphex.c b/libc/sysdeps/ieee754/ldbl-128ibm/printf_fphex.c index e0ec422b0..fb4c9aca7 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/printf_fphex.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/printf_fphex.c @@ -1,5 +1,5 @@ /* Print floating point number in hexadecimal notation according to ISO C99. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_asinhl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_asinhl.c index 63c6edbb1..043b151ff 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_asinhl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_asinhl.c @@ -48,7 +48,7 @@ long double __asinhl(long double x) if(huge+x>one) return x; /* return x inexact except 0 */ } if(ix>0x41b0000000000000LL) { /* |x| > 2**28 */ - w = __ieee754_logl(fabs(x))+ln2; + w = __ieee754_logl(fabsl(x))+ln2; } else if (ix>0x4000000000000000LL) { /* 2**28 > |x| > 2.0 */ t = fabs(x); w = __ieee754_logl(2.0*t+one/(__ieee754_sqrtl(x*x+one)+t)); diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_ceill.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_ceill.c index c230d9e46..1389ccd68 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_ceill.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_ceill.c @@ -1,6 +1,6 @@ /* Ceil (round to +inf) long double floating-point values. IBM extended format long double version. - Copyright (C) 2006-2013 Free Software Foundation, Inc. + Copyright (C) 2006-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_cprojl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_cprojl.c index a344e9288..9edec7dbb 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_cprojl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_cprojl.c @@ -1,5 +1,5 @@ /* Compute projection of complex long double value to Riemann sphere. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_ctanhl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_ctanhl.c index 3cc0bcd2a..b2ac100d6 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_ctanhl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_ctanhl.c @@ -1,5 +1,5 @@ /* Complex hyperbole tangent for long double. IBM extended format version. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_ctanl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_ctanl.c index ce1ec8afb..c2c4c6f92 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_ctanl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_ctanl.c @@ -1,5 +1,5 @@ /* Complex tangent function for long double. IBM extended format version. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_erfl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_erfl.c index c861c65cc..95dc41584 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_erfl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_erfl.c @@ -101,6 +101,7 @@ * erfc/erf(NaN) is NaN */ +#include <errno.h> #include <math.h> #include <math_private.h> #include <math_ldbl_opt.h> @@ -940,14 +941,22 @@ __erfcl (long double x) r = __ieee754_expl (-z * z - 0.5625) * __ieee754_expl ((z - x) * (z + x) + p); if ((hx & 0x80000000) == 0) - return r / x; + { + long double ret = r / x; + if (ret == 0) + __set_errno (ERANGE); + return ret; + } else return two - r / x; } else { if ((hx & 0x80000000) == 0) - return tiny * tiny; + { + __set_errno (ERANGE); + return tiny * tiny; + } else return two - tiny; } diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_floorl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_floorl.c index 9f28d9810..50bc47750 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_floorl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_floorl.c @@ -1,6 +1,6 @@ /* Round to int long double floating-point values. IBM extended format long double version. - Copyright (C) 2006-2013 Free Software Foundation, Inc. + Copyright (C) 2006-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_fmal.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_fmal.c index 5ad93332e..46dedf83f 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_fmal.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_fmal.c @@ -1,5 +1,5 @@ /* Compute x * y + z as ternary operation. - Copyright (C) 2011-2013 Free Software Foundation, Inc. + Copyright (C) 2011-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by David Flaherty <flaherty@linux.vnet.ibm.com>. diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_fpclassifyl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_fpclassifyl.c index 90586e822..1fac09d88 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_fpclassifyl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_fpclassifyl.c @@ -1,5 +1,5 @@ /* Return classification value corresponding to argument. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and Jakub Jelinek <jj@ultra.linux.cz>, 1999. diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_issignalingl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_issignalingl.c index bdd58f8f2..ffc9bbe9a 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_issignalingl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_issignalingl.c @@ -1,5 +1,5 @@ /* Test for signaling NaN. - Copyright (C) 2013 Free Software Foundation, Inc. + Copyright (C) 2013-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_llrintl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_llrintl.c index 856034963..345f3905d 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_llrintl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_llrintl.c @@ -1,6 +1,6 @@ /* Round to long long int long double floating-point values. IBM extended format long double version. - Copyright (C) 2006-2013 Free Software Foundation, Inc. + Copyright (C) 2006-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or @@ -43,7 +43,7 @@ __llrintl (long double x) #endif ) { - save_round = fegetround (); + save_round = __fegetround (); if (__builtin_expect ((xh == -(double) (-__LONG_LONG_MAX__ - 1)), 0)) { diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_llroundl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_llroundl.c index c402b08d9..4a6e2d5f8 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_llroundl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_llroundl.c @@ -1,6 +1,6 @@ /* Round to long long int long double floating-point values. IBM extended format long double version. - Copyright (C) 2006-2013 Free Software Foundation, Inc. + Copyright (C) 2006-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_logbl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_logbl.c index e14028861..da8d71bde 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_logbl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_logbl.c @@ -44,5 +44,6 @@ __logbl (long double x) } return (long double) (rhx - 1023); } - +#ifndef __logbl long_double_symbol (libm, __logbl, logbl); +#endif diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_lrintl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_lrintl.c index 588098d09..56c3a81cd 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_lrintl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_lrintl.c @@ -1,6 +1,6 @@ /* Round to long int long double floating-point values. IBM extended format long double version. - Copyright (C) 2006-2013 Free Software Foundation, Inc. + Copyright (C) 2006-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or @@ -49,7 +49,7 @@ __lrintl (long double x) #endif ) { - save_round = fegetround (); + save_round = __fegetround (); #if __LONG_MAX__ == 2147483647 long long llhi = (long long) xh; diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_lroundl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_lroundl.c index c9f2e201a..107778df3 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_lroundl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_lroundl.c @@ -1,6 +1,6 @@ /* Round to long int long double floating-point values. IBM extended format long double version. - Copyright (C) 2006-2013 Free Software Foundation, Inc. + Copyright (C) 2006-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_nearbyintl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_nearbyintl.c index 92ced5218..4e997a68f 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_nearbyintl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_nearbyintl.c @@ -1,6 +1,6 @@ /* Round to int long double floating-point values without raising inexact. IBM extended format long double version. - Copyright (C) 2006-2013 Free Software Foundation, Inc. + Copyright (C) 2006-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_remquol.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_remquol.c index 195e108ca..aafffa529 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_remquol.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_remquol.c @@ -1,5 +1,5 @@ /* Compute remainder and a congruent to the quotient. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and Jakub Jelinek <jj@ultra.linux.cz>, 1999. diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_rintl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_rintl.c index 48dbe8569..20b75dd21 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_rintl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_rintl.c @@ -1,6 +1,6 @@ /* Round to int long double floating-point values. IBM extended format long double version. - Copyright (C) 2006-2013 Free Software Foundation, Inc. + Copyright (C) 2006-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or @@ -40,7 +40,7 @@ __rintl (long double x) __builtin_inf ()), 1)) { double orig_xh; - int save_round = fegetround (); + int save_round = __fegetround (); /* Long double arithmetic, including the canonicalisation below, only works in round-to-nearest mode. */ diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_roundl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_roundl.c index 8cc451b13..4dd0cc9ed 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_roundl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_roundl.c @@ -1,6 +1,6 @@ /* Round to int long double floating-point values. IBM extended format long double version. - Copyright (C) 2006-2013 Free Software Foundation, Inc. + Copyright (C) 2006-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_signbitl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_signbitl.c index aecb1fd79..6d2571106 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_signbitl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_signbitl.c @@ -1,5 +1,5 @@ /* Return nonzero value if number is negative. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_sincosl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_sincosl.c index a9e2f3d19..57971f9c9 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_sincosl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_sincosl.c @@ -1,5 +1,5 @@ /* Compute sine and cosine of argument. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and Jakub Jelinek <jj@ultra.linux.cz>. diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/s_truncl.c b/libc/sysdeps/ieee754/ldbl-128ibm/s_truncl.c index 2765b6313..490ad266d 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/s_truncl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/s_truncl.c @@ -1,6 +1,6 @@ /* Truncate (toward zero) long double floating-point values. IBM extended format long double version. - Copyright (C) 2006-2013 Free Software Foundation, Inc. + Copyright (C) 2006-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/strtold_l.c b/libc/sysdeps/ieee754/ldbl-128ibm/strtold_l.c index 93a80c5ee..0830a10e3 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/strtold_l.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/strtold_l.c @@ -1,4 +1,4 @@ -/* Copyright (C) 1999-2013 Free Software Foundation, Inc. +/* Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/t_sincosl.c b/libc/sysdeps/ieee754/ldbl-128ibm/t_sincosl.c index e4f54b05b..dcc36bda1 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/t_sincosl.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/t_sincosl.c @@ -1,5 +1,5 @@ /* Quad-precision floating point sine and cosine tables. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jj@ultra.linux.cz> diff --git a/libc/sysdeps/ieee754/ldbl-128ibm/x2y2m1l.c b/libc/sysdeps/ieee754/ldbl-128ibm/x2y2m1l.c index 06dcf02ff..d533339c9 100644 --- a/libc/sysdeps/ieee754/ldbl-128ibm/x2y2m1l.c +++ b/libc/sysdeps/ieee754/ldbl-128ibm/x2y2m1l.c @@ -1,5 +1,5 @@ /* Compute x^2 + y^2 - 1, without large cancellation error. - Copyright (C) 2012-2013 Free Software Foundation, Inc. + Copyright (C) 2012-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-64-128/strtold_l.c b/libc/sysdeps/ieee754/ldbl-64-128/strtold_l.c index e9b33f2d8..dfa79e99f 100644 --- a/libc/sysdeps/ieee754/ldbl-64-128/strtold_l.c +++ b/libc/sysdeps/ieee754/ldbl-64-128/strtold_l.c @@ -1,4 +1,4 @@ -/* Copyright (C) 1999-2013 Free Software Foundation, Inc. +/* Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-96/e_atan2l.c b/libc/sysdeps/ieee754/ldbl-96/e_atan2l.c deleted file mode 100644 index 209f29bbd..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/e_atan2l.c +++ /dev/null @@ -1,125 +0,0 @@ -/* e_atan2l.c -- long double version of e_atan2.c. - * Conversion to long double by Ulrich Drepper, - * Cygnus Support, drepper@cygnus.com. - */ - -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunPro, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ - -/* __ieee754_atan2l(y,x) - * Method : - * 1. Reduce y to positive by atan2(y,x)=-atan2(-y,x). - * 2. Reduce x to positive by (if x and y are unexceptional): - * ARG (x+iy) = arctan(y/x) ... if x > 0, - * ARG (x+iy) = pi - arctan[y/(-x)] ... if x < 0, - * - * Special cases: - * - * ATAN2((anything), NaN ) is NaN; - * ATAN2(NAN , (anything) ) is NaN; - * ATAN2(+-0, +(anything but NaN)) is +-0 ; - * ATAN2(+-0, -(anything but NaN)) is +-pi ; - * ATAN2(+-(anything but 0 and NaN), 0) is +-pi/2; - * ATAN2(+-(anything but INF and NaN), +INF) is +-0 ; - * ATAN2(+-(anything but INF and NaN), -INF) is +-pi; - * ATAN2(+-INF,+INF ) is +-pi/4 ; - * ATAN2(+-INF,-INF ) is +-3pi/4; - * ATAN2(+-INF, (anything but,0,NaN, and INF)) is +-pi/2; - * - * Constants: - * The hexadecimal values are the intended ones for the following - * constants. The decimal values may be used, provided that the - * compiler will convert from decimal to binary accurately enough - * to produce the hexadecimal values shown. - */ - -#include <math.h> -#include <math_private.h> - -static const long double -tiny = 1.0e-4900L, -zero = 0.0, -pi_o_4 = 7.85398163397448309628202E-01L, /* 0x3FFE, 0xC90FDAA2, 0x2168C235 */ -pi_o_2 = 1.5707963267948966192564E+00L, /* 0x3FFF, 0xC90FDAA2, 0x2168C235 */ -pi = 3.14159265358979323851281E+00L, /* 0x4000, 0xC90FDAA2, 0x2168C235 */ -pi_lo = -5.01655761266833202345176e-20L;/* 0xBFBE, 0xECE675D1, 0xFC8F8CBB */ - -long double -__ieee754_atan2l (long double y, long double x) -{ - long double z; - int32_t k,m,hx,hy,ix,iy; - u_int32_t sx,sy,lx,ly; - - GET_LDOUBLE_WORDS(sx,hx,lx,x); - ix = sx&0x7fff; - lx |= hx & 0x7fffffff; - GET_LDOUBLE_WORDS(sy,hy,ly,y); - iy = sy&0x7fff; - ly |= hy & 0x7fffffff; - if(((2*ix|((lx|-lx)>>31))>0xfffe)|| - ((2*iy|((ly|-ly)>>31))>0xfffe)) /* x or y is NaN */ - return x+y; - if(((sx-0x3fff)|lx)==0) return __atanl(y); /* x=1.0 */ - m = ((sy>>15)&1)|((sx>>14)&2); /* 2*sign(x)+sign(y) */ - - /* when y = 0 */ - if((iy|ly)==0) { - switch(m) { - case 0: - case 1: return y; /* atan(+-0,+anything)=+-0 */ - case 2: return pi+tiny;/* atan(+0,-anything) = pi */ - case 3: return -pi-tiny;/* atan(-0,-anything) =-pi */ - } - } - /* when x = 0 */ - if((ix|lx)==0) return (sy>=0x8000)? -pi_o_2-tiny: pi_o_2+tiny; - - /* when x is INF */ - if(ix==0x7fff) { - if(iy==0x7fff) { - switch(m) { - case 0: return pi_o_4+tiny;/* atan(+INF,+INF) */ - case 1: return -pi_o_4-tiny;/* atan(-INF,+INF) */ - case 2: return 3.0*pi_o_4+tiny;/*atan(+INF,-INF)*/ - case 3: return -3.0*pi_o_4-tiny;/*atan(-INF,-INF)*/ - } - } else { - switch(m) { - case 0: return zero ; /* atan(+...,+INF) */ - case 1: return -zero ; /* atan(-...,+INF) */ - case 2: return pi+tiny ; /* atan(+...,-INF) */ - case 3: return -pi-tiny ; /* atan(-...,-INF) */ - } - } - } - /* when y is INF */ - if(iy==0x7fff) return (sy>=0x8000)? -pi_o_2-tiny: pi_o_2+tiny; - - /* compute y/x */ - k = sy-sx; - if(k > 70) z=pi_o_2+0.5*pi_lo; /* |y/x| > 2**70 */ - else if(sx>=0x8000&&k<-70) z=0.0; /* |y|/x < -2**70 */ - else z=__atanl(fabsl(y/x)); /* safe to do y/x */ - switch (m) { - case 0: return z ; /* atan(+,+) */ - case 1: { - u_int32_t sz; - GET_LDOUBLE_EXP(sz,z); - SET_LDOUBLE_EXP(z,sz ^ 0x8000); - } - return z ; /* atan(-,+) */ - case 2: return pi-(z-pi_lo);/* atan(+,-) */ - default: /* case 3 */ - return (z-pi_lo)-pi;/* atan(-,-) */ - } -} -strong_alias (__ieee754_atan2l, __atan2l_finite) diff --git a/libc/sysdeps/ieee754/ldbl-96/e_gammal_r.c b/libc/sysdeps/ieee754/ldbl-96/e_gammal_r.c index 7cb3e8563..477c3a61d 100644 --- a/libc/sysdeps/ieee754/ldbl-96/e_gammal_r.c +++ b/libc/sysdeps/ieee754/ldbl-96/e_gammal_r.c @@ -1,5 +1,5 @@ /* Implementation of gamma function according to ISO C. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-96/e_hypotl.c b/libc/sysdeps/ieee754/ldbl-96/e_hypotl.c index 789548884..d3152f91e 100644 --- a/libc/sysdeps/ieee754/ldbl-96/e_hypotl.c +++ b/libc/sysdeps/ieee754/ldbl-96/e_hypotl.c @@ -89,6 +89,17 @@ long double __ieee754_hypotl(long double x, long double y) b *= t1; a *= t1; k -= 16382; + GET_LDOUBLE_EXP (ea, a); + GET_LDOUBLE_EXP (eb, b); + if (eb > ea) + { + t1 = a; + a = b; + b = t1; + j = ea; + ea = eb; + eb = j; + } } else { /* scale a and b by 2^9600 */ ea += 0x2580; /* a *= 2^9600 */ eb += 0x2580; /* b *= 2^9600 */ diff --git a/libc/sysdeps/ieee754/ldbl-96/e_ilogbl.c b/libc/sysdeps/ieee754/ldbl-96/e_ilogbl.c deleted file mode 100644 index 0c7d9d544..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/e_ilogbl.c +++ /dev/null @@ -1,59 +0,0 @@ -/* s_ilogbl.c -- long double version of s_ilogb.c. - * Conversion to long double by Ulrich Drepper, - * Cygnus Support, drepper@cygnus.com. - */ - -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunPro, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ - -#if defined(LIBM_SCCS) && !defined(lint) -static char rcsid[] = "$NetBSD: $"; -#endif - -/* ilogbl(long double x) - * return the binary exponent of non-zero x - * ilogbl(0) = FP_ILOGB0 - * ilogbl(NaN) = FP_ILOGBNAN (no signal is raised) - * ilogbl(+-Inf) = INT_MAX (no signal is raised) - */ - -#include <limits.h> -#include <math.h> -#include <math_private.h> - -int __ieee754_ilogbl (long double x) -{ - int32_t es,hx,lx,ix; - - GET_LDOUBLE_EXP(es,x); - es &= 0x7fff; - if(es==0) { - GET_LDOUBLE_WORDS(es,hx,lx,x); - if((hx|lx)==0) - return FP_ILOGB0; /* ilogbl(0) = FP_ILOGB0 */ - else /* subnormal x */ - if(hx==0) { - for (ix = -16415; lx>0; lx<<=1) ix -=1; - } else { - for (ix = -16383; hx>0; hx<<=1) ix -=1; - } - return ix; - } - else if (es<0x7fff) return es-0x3fff; - else if (FP_ILOGBNAN != INT_MAX) - { - GET_LDOUBLE_WORDS(es,hx,lx,x); - if (((hx & 0x7fffffff)|lx) == 0) - /* ISO C99 requires ilogbl(+-Inf) == INT_MAX. */ - return INT_MAX; - } - return FP_ILOGBNAN; -} diff --git a/libc/sysdeps/ieee754/ldbl-96/e_jnl.c b/libc/sysdeps/ieee754/ldbl-96/e_jnl.c index 58a9107f7..fa8e27efe 100644 --- a/libc/sysdeps/ieee754/ldbl-96/e_jnl.c +++ b/libc/sysdeps/ieee754/ldbl-96/e_jnl.c @@ -302,7 +302,8 @@ __ieee754_ynl (int n, long double x) if (__builtin_expect ((ix == 0x7fff) && ((i0 & 0x7fffffff) != 0), 0)) return x + x; if (__builtin_expect ((ix | i0 | i1) == 0, 0)) - return -HUGE_VALL + x; /* -inf and overflow exception. */ + /* -inf or inf and divide-by-zero exception. */ + return ((n < 0 && (n & 1) != 0) ? 1.0L : -1.0L) / 0.0L; if (__builtin_expect (se & 0x8000, 0)) return zero / (zero * x); sign = 1; diff --git a/libc/sysdeps/ieee754/ldbl-96/e_rem_pio2l.c b/libc/sysdeps/ieee754/ldbl-96/e_rem_pio2l.c index 18e6047b6..e18be6ee0 100644 --- a/libc/sysdeps/ieee754/ldbl-96/e_rem_pio2l.c +++ b/libc/sysdeps/ieee754/ldbl-96/e_rem_pio2l.c @@ -1,5 +1,5 @@ /* Extended-precision floating point argument reduction. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Based on quad-precision code by Jakub Jelinek <jj@ultra.linux.cz> diff --git a/libc/sysdeps/ieee754/ldbl-96/e_remainderl.c b/libc/sysdeps/ieee754/ldbl-96/e_remainderl.c deleted file mode 100644 index 290e483ae..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/e_remainderl.c +++ /dev/null @@ -1,72 +0,0 @@ -/* e_remainderl.c -- long double version of e_remainder.c. - * Conversion to long double by Ulrich Drepper, - * Cygnus Support, drepper@cygnus.com. - */ - -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunPro, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ - -/* __ieee754_remainderl(x,p) - * Return : - * returns x REM p = x - [x/p]*p as if in infinite - * precise arithmetic, where [x/p] is the (infinite bit) - * integer nearest x/p (in half way case choose the even one). - * Method : - * Based on fmod() return x-[x/p]chopped*p exactlp. - */ - -#include <math.h> -#include <math_private.h> - -static const long double zero = 0.0; - - -long double -__ieee754_remainderl(long double x, long double p) -{ - u_int32_t sx,sex,sep,x0,x1,p0,p1; - long double p_half; - - GET_LDOUBLE_WORDS(sex,x0,x1,x); - GET_LDOUBLE_WORDS(sep,p0,p1,p); - sx = sex&0x8000; - sep &= 0x7fff; - sex &= 0x7fff; - - /* purge off exception values */ - if((sep|p0|p1)==0) return (x*p)/(x*p); /* p = 0 */ - if((sex==0x7fff)|| /* x not finite */ - ((sep==0x7fff)&& /* p is NaN */ - ((p0|p1)!=0))) - return (x*p)/(x*p); - - - if (sep<0x7ffe) x = __ieee754_fmodl(x,p+p); /* now x < 2p */ - if (((sex-sep)|(x0-p0)|(x1-p1))==0) return zero*x; - x = fabsl(x); - p = fabsl(p); - if (sep<0x0002) { - if(x+x>p) { - x-=p; - if(x+x>=p) x -= p; - } - } else { - p_half = 0.5*p; - if(x>p_half) { - x-=p; - if(x>=p_half) x -= p; - } - } - GET_LDOUBLE_EXP(sex,x); - SET_LDOUBLE_EXP(x,sex^sx); - return x; -} -strong_alias (__ieee754_remainderl, __remainderl_finite) diff --git a/libc/sysdeps/ieee754/ldbl-96/gamma_product.c b/libc/sysdeps/ieee754/ldbl-96/gamma_product.c index d464e7084..214745624 100644 --- a/libc/sysdeps/ieee754/ldbl-96/gamma_product.c +++ b/libc/sysdeps/ieee754/ldbl-96/gamma_product.c @@ -1,5 +1,5 @@ /* Compute a product of X, X+1, ..., with an error estimate. - Copyright (C) 2013 Free Software Foundation, Inc. + Copyright (C) 2013-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-96/gamma_productl.c b/libc/sysdeps/ieee754/ldbl-96/gamma_productl.c index 157dbab9f..2c8764967 100644 --- a/libc/sysdeps/ieee754/ldbl-96/gamma_productl.c +++ b/libc/sysdeps/ieee754/ldbl-96/gamma_productl.c @@ -1,5 +1,5 @@ /* Compute a product of X, X+1, ..., with an error estimate. - Copyright (C) 2013 Free Software Foundation, Inc. + Copyright (C) 2013-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-96/k_cosl.c b/libc/sysdeps/ieee754/ldbl-96/k_cosl.c index 52ce076ca..daf7c060d 100644 --- a/libc/sysdeps/ieee754/ldbl-96/k_cosl.c +++ b/libc/sysdeps/ieee754/ldbl-96/k_cosl.c @@ -1,5 +1,5 @@ /* Extended-precision floating point cosine on <-pi/4,pi/4>. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Based on quad-precision cosine by Jakub Jelinek <jj@ultra.linux.cz> diff --git a/libc/sysdeps/ieee754/ldbl-96/k_sinl.c b/libc/sysdeps/ieee754/ldbl-96/k_sinl.c index 15581944a..f2d1c860e 100644 --- a/libc/sysdeps/ieee754/ldbl-96/k_sinl.c +++ b/libc/sysdeps/ieee754/ldbl-96/k_sinl.c @@ -1,5 +1,5 @@ /* Quad-precision floating point sine on <-pi/4,pi/4>. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Based on quad-precision sine by Jakub Jelinek <jj@ultra.linux.cz> diff --git a/libc/sysdeps/ieee754/ldbl-96/ldbl2mpn.c b/libc/sysdeps/ieee754/ldbl-96/ldbl2mpn.c index 3f85e0ae9..be8b93ffa 100644 --- a/libc/sysdeps/ieee754/ldbl-96/ldbl2mpn.c +++ b/libc/sysdeps/ieee754/ldbl-96/ldbl2mpn.c @@ -1,4 +1,4 @@ -/* Copyright (C) 1995-2013 Free Software Foundation, Inc. +/* Copyright (C) 1995-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-96/mpn2ldbl.c b/libc/sysdeps/ieee754/ldbl-96/mpn2ldbl.c index ab6f4e059..8fdb335fe 100644 --- a/libc/sysdeps/ieee754/ldbl-96/mpn2ldbl.c +++ b/libc/sysdeps/ieee754/ldbl-96/mpn2ldbl.c @@ -1,4 +1,4 @@ -/* Copyright (C) 1995-2013 Free Software Foundation, Inc. +/* Copyright (C) 1995-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-96/printf_fphex.c b/libc/sysdeps/ieee754/ldbl-96/printf_fphex.c index 715c93b50..798c78848 100644 --- a/libc/sysdeps/ieee754/ldbl-96/printf_fphex.c +++ b/libc/sysdeps/ieee754/ldbl-96/printf_fphex.c @@ -1,5 +1,5 @@ /* Print floating point number in hexadecimal notation according to ISO C99. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-96/s_cbrtl.c b/libc/sysdeps/ieee754/ldbl-96/s_cbrtl.c index e6b2f2540..49a689164 100644 --- a/libc/sysdeps/ieee754/ldbl-96/s_cbrtl.c +++ b/libc/sysdeps/ieee754/ldbl-96/s_cbrtl.c @@ -1,5 +1,5 @@ /* Compute cubic root of double value. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Dirk Alboth <dirka@uni-paderborn.de> and Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-96/s_ceill.c b/libc/sysdeps/ieee754/ldbl-96/s_ceill.c deleted file mode 100644 index aef8a32f6..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/s_ceill.c +++ /dev/null @@ -1,85 +0,0 @@ -/* s_ceill.c -- long double version of s_ceil.c. - * Conversion to long double by Ulrich Drepper, - * Cygnus Support, drepper@cygnus.com. - */ - -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunPro, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ - -#if defined(LIBM_SCCS) && !defined(lint) -static char rcsid[] = "$NetBSD: $"; -#endif - -/* - * ceill(x) - * Return x rounded toward -inf to integral value - * Method: - * Bit twiddling. - * Exception: - * Inexact flag raised if x not equal to ceil(x). - */ - -#include <math.h> -#include <math_private.h> - -static const long double huge = 1.0e4930; - -long double __ceill(long double x) -{ - int32_t i1,j0; - u_int32_t i,j,se,i0,sx; - GET_LDOUBLE_WORDS(se,i0,i1,x); - sx = (se>>15)&1; - j0 = (se&0x7fff)-0x3fff; - if(j0<31) { - if(j0<0) { /* raise inexact if x != 0 */ - if(huge+x>0.0) {/* return 0*sign(x) if |x|<1 */ - if(sx) {se=0x8000;i0=0;i1=0;} - else if((i0|i1)!=0) { se=0x3fff;i0=0;i1=0;} - } - } else { - i = (0x7fffffff)>>j0; - if(((i0&i)|i1)==0) return x; /* x is integral */ - if(huge+x>0.0) { /* raise inexact flag */ - if(sx==0) { - if (j0>0 && (i0+(0x80000000>>j0))>i0) - i0+=0x80000000>>j0; - else - { - i = 0x7fffffff; - ++se; - } - } - i0 &= (~i); i1=0; - } - } - } else if (j0>62) { - if(j0==0x4000) return x+x; /* inf or NaN */ - else return x; /* x is integral */ - } else { - i = ((u_int32_t)(0xffffffff))>>(j0-31); - if((i1&i)==0) return x; /* x is integral */ - if(huge+x>0.0) { /* raise inexact flag */ - if(sx==0) { - if(j0==31) i0+=1; - else { - j = i1 + (1<<(63-j0)); - if(j<i1) i0+=1; /* got a carry */ - i1 = j; - } - } - i1 &= (~i); - } - } - SET_LDOUBLE_WORDS(x,se,i0,i1); - return x; -} -weak_alias (__ceill, ceill) diff --git a/libc/sysdeps/ieee754/ldbl-96/s_erfl.c b/libc/sysdeps/ieee754/ldbl-96/s_erfl.c index 17d227852..47e4b9e90 100644 --- a/libc/sysdeps/ieee754/ldbl-96/s_erfl.c +++ b/libc/sysdeps/ieee754/ldbl-96/s_erfl.c @@ -104,6 +104,7 @@ */ +#include <errno.h> #include <math.h> #include <math_private.h> @@ -422,14 +423,22 @@ __erfcl (long double x) r = __ieee754_expl (-z * z - 0.5625) * __ieee754_expl ((z - x) * (z + x) + R / S); if ((se & 0x8000) == 0) - return r / x; + { + long double ret = r / x; + if (ret == 0) + __set_errno (ERANGE); + return ret; + } else return two - r / x; } else { if ((se & 0x8000) == 0) - return tiny * tiny; + { + __set_errno (ERANGE); + return tiny * tiny; + } else return two - tiny; } diff --git a/libc/sysdeps/ieee754/ldbl-96/s_fabsl.c b/libc/sysdeps/ieee754/ldbl-96/s_fabsl.c deleted file mode 100644 index fdc70e0dc..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/s_fabsl.c +++ /dev/null @@ -1,35 +0,0 @@ -/* s_fabsl.c -- long double version of s_fabs.c. - * Conversion to long double by Ulrich Drepper, - * Cygnus Support, drepper@cygnus.com. - */ - -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunPro, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ - -#if defined(LIBM_SCCS) && !defined(lint) -static char rcsid[] = "$NetBSD: $"; -#endif - -/* - * fabsl(x) returns the absolute value of x. - */ - -#include <math.h> -#include <math_private.h> - -long double __fabsl(long double x) -{ - u_int32_t exp; - GET_LDOUBLE_EXP(exp,x); - SET_LDOUBLE_EXP(x,exp&0x7fff); - return x; -} -weak_alias (__fabsl, fabsl) diff --git a/libc/sysdeps/ieee754/ldbl-96/s_finitel.c b/libc/sysdeps/ieee754/ldbl-96/s_finitel.c deleted file mode 100644 index fbf4cc691..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/s_finitel.c +++ /dev/null @@ -1,36 +0,0 @@ -/* s_finitel.c -- long double version of s_finite.c. - * Conversion to long double by Ulrich Drepper, - * Cygnus Support, drepper@cygnus.com. - */ - -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunPro, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ - -#if defined(LIBM_SCCS) && !defined(lint) -static char rcsid[] = "$NetBSD: $"; -#endif - -/* - * finitel(x) returns 1 is x is finite, else 0; - * no branching! - */ - -#include <math.h> -#include <math_private.h> - -int __finitel(long double x) -{ - int32_t exp; - GET_LDOUBLE_EXP(exp,x); - return (int)((u_int32_t)((exp&0x7fff)-0x7fff)>>31); -} -hidden_def (__finitel) -weak_alias (__finitel, finitel) diff --git a/libc/sysdeps/ieee754/ldbl-96/s_floorl.c b/libc/sysdeps/ieee754/ldbl-96/s_floorl.c deleted file mode 100644 index cad7935b3..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/s_floorl.c +++ /dev/null @@ -1,86 +0,0 @@ -/* s_floorl.c -- long double version of s_floor.c. - * Conversion to long double by Ulrich Drepper, - * Cygnus Support, drepper@cygnus.com. - */ - -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunPro, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ - -#if defined(LIBM_SCCS) && !defined(lint) -static char rcsid[] = "$NetBSD: $"; -#endif - -/* - * floorl(x) - * Return x rounded toward -inf to integral value - * Method: - * Bit twiddling. - * Exception: - * Inexact flag raised if x not equal to floor(x). - */ - -#include <math.h> -#include <math_private.h> - -static const long double huge = 1.0e4930; - -long double __floorl(long double x) -{ - int32_t i1,j0; - u_int32_t i,j,se,i0,sx; - GET_LDOUBLE_WORDS(se,i0,i1,x); - sx = (se>>15)&1; - j0 = (se&0x7fff)-0x3fff; - if(j0<31) { - if(j0<0) { /* raise inexact if x != 0 */ - if(huge+x>0.0) {/* return 0*sign(x) if |x|<1 */ - if(sx==0) {se=0;i0=i1=0;} - else if(((se&0x7fff)|i0|i1)!=0) - { se=0xbfff;i0=i1=0;} - } - } else { - i = (0x7fffffff)>>j0; - if(((i0&i)|i1)==0) return x; /* x is integral */ - if(huge+x>0.0) { /* raise inexact flag */ - if(sx) { - if (j0>0 && (i0+(0x80000000>>j0))>i0) - i0 += (0x80000000)>>j0; - else - { - i = 0x7fffffff; - ++se; - } - } - i0 &= (~i); i1=0; - } - } - } else if (j0>62) { - if(j0==0x4000) return x+x; /* inf or NaN */ - else return x; /* x is integral */ - } else { - i = ((u_int32_t)(0xffffffff))>>(j0-31); - if((i1&i)==0) return x; /* x is integral */ - if(huge+x>0.0) { /* raise inexact flag */ - if(sx) { - if(j0==31) i0+=1; - else { - j = i1+(1<<(63-j0)); - if(j<i1) i0 +=1 ; /* got a carry */ - i1=j; - } - } - i1 &= (~i); - } - } - SET_LDOUBLE_WORDS(x,se,i0,i1); - return x; -} -weak_alias (__floorl, floorl) diff --git a/libc/sysdeps/ieee754/ldbl-96/s_fma.c b/libc/sysdeps/ieee754/ldbl-96/s_fma.c index 88b902fcf..fde281104 100644 --- a/libc/sysdeps/ieee754/ldbl-96/s_fma.c +++ b/libc/sysdeps/ieee754/ldbl-96/s_fma.c @@ -1,5 +1,5 @@ /* Compute x * y + z as ternary operation. - Copyright (C) 2010-2013 Free Software Foundation, Inc. + Copyright (C) 2010-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jakub@redhat.com>, 2010. diff --git a/libc/sysdeps/ieee754/ldbl-96/s_fmal.c b/libc/sysdeps/ieee754/ldbl-96/s_fmal.c index a1cc0fd1f..056432135 100644 --- a/libc/sysdeps/ieee754/ldbl-96/s_fmal.c +++ b/libc/sysdeps/ieee754/ldbl-96/s_fmal.c @@ -1,5 +1,5 @@ /* Compute x * y + z as ternary operation. - Copyright (C) 2010-2013 Free Software Foundation, Inc. + Copyright (C) 2010-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jakub@redhat.com>, 2010. diff --git a/libc/sysdeps/ieee754/ldbl-96/s_fpclassifyl.c b/libc/sysdeps/ieee754/ldbl-96/s_fpclassifyl.c deleted file mode 100644 index d164f8633..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/s_fpclassifyl.c +++ /dev/null @@ -1,44 +0,0 @@ -/* Return classification value corresponding to argument. - Copyright (C) 1997-2013 Free Software Foundation, Inc. - This file is part of the GNU C Library. - Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. - Fixed by Andreas Schwab <schwab@issan.informatik.uni-dortmund.de>. - - The GNU C Library is free software; you can redistribute it and/or - modify it under the terms of the GNU Lesser General Public - License as published by the Free Software Foundation; either - version 2.1 of the License, or (at your option) any later version. - - The GNU C Library is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU - Lesser General Public License for more details. - - You should have received a copy of the GNU Lesser General Public - License along with the GNU C Library; if not, see - <http://www.gnu.org/licenses/>. */ - -#include <math.h> - -#include <math_private.h> - - -int -__fpclassifyl (long double x) -{ - u_int32_t ex, hx, lx, m; - int retval = FP_NORMAL; - - GET_LDOUBLE_WORDS (ex, hx, lx, x); - m = (hx & 0x7fffffff) | lx; - ex &= 0x7fff; - if ((ex | m) == 0) - retval = FP_ZERO; - else if (ex == 0 && (hx & 0x80000000) == 0) - retval = FP_SUBNORMAL; - else if (ex == 0x7fff) - retval = m != 0 ? FP_NAN : FP_INFINITE; - - return retval; -} -libm_hidden_def (__fpclassifyl) diff --git a/libc/sysdeps/ieee754/ldbl-96/s_isinfl.c b/libc/sysdeps/ieee754/ldbl-96/s_isinfl.c deleted file mode 100644 index 94639f00f..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/s_isinfl.c +++ /dev/null @@ -1,30 +0,0 @@ -/* - * Written by J.T. Conklin <jtc@netbsd.org>. - * Change for long double by Ulrich Drepper <drepper@cygnus.com>. - * Public domain. - */ - -#if defined(LIBM_SCCS) && !defined(lint) -static char rcsid[] = "$NetBSD: $"; -#endif - -/* - * isinfl(x) returns 1 if x is inf, -1 if x is -inf, else 0; - * no branching! - */ - -#include <math.h> -#include <math_private.h> - -int -__isinfl (long double x) -{ - int32_t se,hx,lx; - GET_LDOUBLE_WORDS(se,hx,lx,x); - lx |= (hx & 0x7fffffff) | ((se & 0x7fff) ^ 0x7fff); - lx |= -lx; - se &= 0x8000; - return ~(lx >> 31) & (1 - (se >> 14)); -} -hidden_def (__isinfl) -weak_alias (__isinfl, isinfl) diff --git a/libc/sysdeps/ieee754/ldbl-96/s_isnanl.c b/libc/sysdeps/ieee754/ldbl-96/s_isnanl.c deleted file mode 100644 index fd270fd84..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/s_isnanl.c +++ /dev/null @@ -1,40 +0,0 @@ -/* s_isnanl.c -- long double version of s_isnan.c. - * Conversion to long double by Ulrich Drepper, - * Cygnus Support, drepper@cygnus.com. - */ - -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunPro, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ - -#if defined(LIBM_SCCS) && !defined(lint) -static char rcsid[] = "$NetBSD: $"; -#endif - -/* - * isnanl(x) returns 1 is x is nan, else 0; - * no branching! - */ - -#include <math.h> -#include <math_private.h> - -int __isnanl(long double x) -{ - int32_t se,hx,lx; - GET_LDOUBLE_WORDS(se,hx,lx,x); - se = (se & 0x7fff) << 1; - lx |= hx & 0x7fffffff; - se |= (u_int32_t)(lx|(-lx))>>31; - se = 0xfffe - se; - return (int)(((u_int32_t)(se))>>31); -} -hidden_def (__isnanl) -weak_alias (__isnanl, isnanl) diff --git a/libc/sysdeps/ieee754/ldbl-96/s_issignalingl.c b/libc/sysdeps/ieee754/ldbl-96/s_issignalingl.c index 107aff151..bbb98735f 100644 --- a/libc/sysdeps/ieee754/ldbl-96/s_issignalingl.c +++ b/libc/sysdeps/ieee754/ldbl-96/s_issignalingl.c @@ -1,5 +1,5 @@ /* Test for signaling NaN. - Copyright (C) 2013 Free Software Foundation, Inc. + Copyright (C) 2013-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-96/s_llrintl.c b/libc/sysdeps/ieee754/ldbl-96/s_llrintl.c index 678ecd48b..0102f1824 100644 --- a/libc/sysdeps/ieee754/ldbl-96/s_llrintl.c +++ b/libc/sysdeps/ieee754/ldbl-96/s_llrintl.c @@ -1,6 +1,6 @@ /* Round argument to nearest integral value according to current rounding direction. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-96/s_llroundl.c b/libc/sysdeps/ieee754/ldbl-96/s_llroundl.c index 6074b962b..1559540a7 100644 --- a/libc/sysdeps/ieee754/ldbl-96/s_llroundl.c +++ b/libc/sysdeps/ieee754/ldbl-96/s_llroundl.c @@ -1,5 +1,5 @@ /* Round long double value to long long int. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-96/s_logbl.c b/libc/sysdeps/ieee754/ldbl-96/s_logbl.c deleted file mode 100644 index 4289be193..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/s_logbl.c +++ /dev/null @@ -1,51 +0,0 @@ -/* s_logbl.c -- long double version of s_logb.c. - * Conversion to long double by Ulrich Drepper, - * Cygnus Support, drepper@cygnus.com. - */ - -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunPro, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ - -/* - * long double logbl(x) - * IEEE 754 logb. Included to pass IEEE test suite. Not recommend. - * Use ilogb instead. - */ - -#include <math.h> -#include <math_private.h> - -long double -__logbl (long double x) -{ - int32_t es, lx, ix; - - GET_LDOUBLE_WORDS (es, ix, lx, x); - es &= 0x7fff; /* exponent */ - if ((es | ix | lx) == 0) - return -1.0 / fabs (x); - if (es == 0x7fff) - return x * x; - if (es == 0) /* IEEE 754 logb */ - { - /* POSIX specifies that denormal number is treated as - though it were normalized. */ - int ma; - if (ix == 0) - ma = __builtin_clz (lx) + 32; - else - ma = __builtin_clz (ix); - es -= ma - 1; - } - return (long double) (es - 16383); -} - -weak_alias (__logbl, logbl) diff --git a/libc/sysdeps/ieee754/ldbl-96/s_lrintl.c b/libc/sysdeps/ieee754/ldbl-96/s_lrintl.c index de1839673..a66893446 100644 --- a/libc/sysdeps/ieee754/ldbl-96/s_lrintl.c +++ b/libc/sysdeps/ieee754/ldbl-96/s_lrintl.c @@ -1,6 +1,6 @@ /* Round argument to nearest integral value according to current rounding direction. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-96/s_lroundl.c b/libc/sysdeps/ieee754/ldbl-96/s_lroundl.c index 790d5b763..a5c08d639 100644 --- a/libc/sysdeps/ieee754/ldbl-96/s_lroundl.c +++ b/libc/sysdeps/ieee754/ldbl-96/s_lroundl.c @@ -1,5 +1,5 @@ /* Round long double value to long int. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-96/s_nearbyintl.c b/libc/sysdeps/ieee754/ldbl-96/s_nearbyintl.c deleted file mode 100644 index c1d77f0c0..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/s_nearbyintl.c +++ /dev/null @@ -1,88 +0,0 @@ -/* s_rintl.c -- long double version of s_rint.c. - * Conversion to long double by Ulrich Drepper, - * Cygnus Support, drepper@cygnus.com. - */ -/* Adapted for use as nearbyint by Ulrich Drepper <drepper@cygnus.com>. */ - -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunPro, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ - -/* - * rintl(x) - * Return x rounded to integral value according to the prevailing - * rounding mode. - * Method: - * Using floating addition. - * Exception: - * Inexact flag raised if x not equal to rintl(x). - */ - -#include <fenv.h> -#include <math.h> -#include <math_private.h> - -static const long double -TWO63[2]={ - 9.223372036854775808000000e+18, /* 0x403E, 0x00000000, 0x00000000 */ - -9.223372036854775808000000e+18 /* 0xC03E, 0x00000000, 0x00000000 */ -}; - -long double __nearbyintl(long double x) -{ - fenv_t env; - int32_t se,j0,sx; - u_int32_t i,i0,i1; - long double w,t; - GET_LDOUBLE_WORDS(se,i0,i1,x); - sx = (se>>15)&1; - j0 = (se&0x7fff)-0x3fff; - if(j0<31) { - if(j0<0) { - if(((se&0x7fff)|i0|i1)==0) return x; - i1 |= i0; - i0 &= 0xe0000000; - i0 |= (i1|-i1)&0x80000000; - SET_LDOUBLE_MSW(x,i0); - feholdexcept (&env); - w = TWO63[sx]+x; - t = w-TWO63[sx]; - math_force_eval (t); - fesetenv (&env); - GET_LDOUBLE_EXP(i0,t); - SET_LDOUBLE_EXP(t,(i0&0x7fff)|(sx<<15)); - return t; - } else { - i = (0x7fffffff)>>j0; - if(((i0&i)|i1)==0) return x; /* x is integral */ - i>>=1; - if(((i0&i)|i1)!=0) { - if (j0==30) i1 = 0x40000000; else - i0 = (i0&(~i))|((0x20000000)>>j0); - } - } - } else if (j0>62) { - if(j0==0x4000) return x+x; /* inf or NaN */ - else return x; /* x is integral */ - } else { - i = ((u_int32_t)(0xffffffff))>>(j0-31); - if((i1&i)==0) return x; /* x is integral */ - i>>=1; - if((i1&i)!=0) i1 = (i1&(~i))|((0x40000000)>>(j0-31)); - } - SET_LDOUBLE_WORDS(x,se,i0,i1); - feholdexcept (&env); - w = TWO63[sx]+x; - t = w-TWO63[sx]; - math_force_eval (t); - fesetenv (&env); - return t; -} -weak_alias (__nearbyintl, nearbyintl) diff --git a/libc/sysdeps/ieee754/ldbl-96/s_nextafterl.c b/libc/sysdeps/ieee754/ldbl-96/s_nextafterl.c deleted file mode 100644 index 6859349b7..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/s_nextafterl.c +++ /dev/null @@ -1,96 +0,0 @@ -/* s_nextafterl.c -- long double version of s_nextafter.c. - * Conversion to long double by Ulrich Drepper, - * Cygnus Support, drepper@cygnus.com. - */ - -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunPro, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ - -#if defined(LIBM_SCCS) && !defined(lint) -static char rcsid[] = "$NetBSD: $"; -#endif - -/* IEEE functions - * nextafterl(x,y) - * return the next machine floating-point number of x in the - * direction toward y. - * Special cases: - */ - -#include <math.h> -#include <math_private.h> - -long double __nextafterl(long double x, long double y) -{ - int32_t hx,hy,ix,iy; - u_int32_t lx,ly,esx,esy; - - GET_LDOUBLE_WORDS(esx,hx,lx,x); - GET_LDOUBLE_WORDS(esy,hy,ly,y); - ix = esx&0x7fff; /* |x| */ - iy = esy&0x7fff; /* |y| */ - - if (((ix==0x7fff)&&((hx|lx)!=0)) || /* x is nan */ - ((iy==0x7fff)&&((hy|ly)!=0))) /* y is nan */ - return x+y; - if(x==y) return y; /* x=y, return y */ - if((ix|hx|lx)==0) { /* x == 0 */ - long double u; - SET_LDOUBLE_WORDS(x,esy&0x8000,0,1);/* return +-minsubnormal */ - u = math_opt_barrier (x); - u = u * u; - math_force_eval (u); /* raise underflow flag */ - return x; - } - if(esx<0x8000) { /* x > 0 */ - if(ix>iy||((ix==iy) && (hx>hy||((hx==hy)&&(lx>ly))))) { - /* x > y, x -= ulp */ - if(lx==0) { - if (hx==0) esx -= 1; - hx -= 1; - } - lx -= 1; - } else { /* x < y, x += ulp */ - lx += 1; - if(lx==0) { - hx += 1; - if (hx==0) - esx += 1; - } - } - } else { /* x < 0 */ - if(esy>=0||(ix>iy||((ix==iy)&&(hx>hy||((hx==hy)&&(lx>ly)))))){ - /* x < y, x -= ulp */ - if(lx==0) { - if (hx==0) esx -= 1; - hx -= 1; - } - lx -= 1; - } else { /* x > y, x += ulp */ - lx += 1; - if(lx==0) { - hx += 1; - if (hx==0) esx += 1; - } - } - } - esy = esx&0x7fff; - if(esy==0x7fff) return x+x; /* overflow */ - if(esy==0) { - long double u = x*x; /* underflow */ - math_force_eval (u); /* raise underflow flag */ - } - SET_LDOUBLE_WORDS(x,esx,hx,lx); - return x; -} -weak_alias (__nextafterl, nextafterl) -strong_alias (__nextafterl, __nexttowardl) -weak_alias (__nextafterl, nexttowardl) diff --git a/libc/sysdeps/ieee754/ldbl-96/s_remquol.c b/libc/sysdeps/ieee754/ldbl-96/s_remquol.c index c54682176..866935c6e 100644 --- a/libc/sysdeps/ieee754/ldbl-96/s_remquol.c +++ b/libc/sysdeps/ieee754/ldbl-96/s_remquol.c @@ -1,5 +1,5 @@ /* Compute remainder and a congruent to the quotient. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-96/s_rintl.c b/libc/sysdeps/ieee754/ldbl-96/s_rintl.c deleted file mode 100644 index b6f899d4e..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/s_rintl.c +++ /dev/null @@ -1,82 +0,0 @@ -/* s_rintl.c -- long double version of s_rint.c. - * Conversion to long double by Ulrich Drepper, - * Cygnus Support, drepper@cygnus.com. - */ - -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunPro, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ - -#if defined(LIBM_SCCS) && !defined(lint) -static char rcsid[] = "$NetBSD: $"; -#endif - -/* - * rintl(x) - * Return x rounded to integral value according to the prevailing - * rounding mode. - * Method: - * Using floating addition. - * Exception: - * Inexact flag raised if x not equal to rintl(x). - */ - -#include <math.h> -#include <math_private.h> - -static const long double -TWO63[2]={ - 9.223372036854775808000000e+18, /* 0x403E, 0x00000000, 0x00000000 */ - -9.223372036854775808000000e+18 /* 0xC03E, 0x00000000, 0x00000000 */ -}; - -long double __rintl(long double x) -{ - int32_t se,j0,sx; - u_int32_t i,i0,i1; - long double w,t; - GET_LDOUBLE_WORDS(se,i0,i1,x); - sx = (se>>15)&1; - j0 = (se&0x7fff)-0x3fff; - if(j0<31) { - if(j0<0) { - if(((se&0x7fff)|i0|i1)==0) return x; - i1 |= i0; - i0 &= 0xe0000000; - i0 |= (i1|-i1)&0x80000000; - SET_LDOUBLE_MSW(x,i0); - w = TWO63[sx]+x; - t = w-TWO63[sx]; - GET_LDOUBLE_EXP(i0,t); - SET_LDOUBLE_EXP(t,(i0&0x7fff)|(sx<<15)); - return t; - } else { - i = (0x7fffffff)>>j0; - if(((i0&i)|i1)==0) return x; /* x is integral */ - i>>=1; - if(((i0&i)|i1)!=0) { - if(j0==30) i1 = 0x40000000; else - i0 = (i0&(~i))|((0x20000000)>>j0); - } - } - } else if (j0>62) { - if(j0==0x4000) return x+x; /* inf or NaN */ - else return x; /* x is integral */ - } else { - i = ((u_int32_t)(0xffffffff))>>(j0-31); - if((i1&i)==0) return x; /* x is integral */ - i>>=1; - if((i1&i)!=0) i1 = (i1&(~i))|((0x40000000)>>(j0-31)); - } - SET_LDOUBLE_WORDS(x,se,i0,i1); - w = TWO63[sx]+x; - return w-TWO63[sx]; -} -weak_alias (__rintl, rintl) diff --git a/libc/sysdeps/ieee754/ldbl-96/s_roundl.c b/libc/sysdeps/ieee754/ldbl-96/s_roundl.c index 9e88bdda6..e16e9c6f6 100644 --- a/libc/sysdeps/ieee754/ldbl-96/s_roundl.c +++ b/libc/sysdeps/ieee754/ldbl-96/s_roundl.c @@ -1,5 +1,5 @@ /* Round long double to integer away from zero. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-96/s_scalbnl.c b/libc/sysdeps/ieee754/ldbl-96/s_scalbnl.c deleted file mode 100644 index 266a37b9c..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/s_scalbnl.c +++ /dev/null @@ -1,61 +0,0 @@ -/* s_scalbnl.c -- long double version of s_scalbn.c. - * Conversion to long double by Ulrich Drepper, - * Cygnus Support, drepper@cygnus.com. - */ - -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunPro, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ - -/* - * scalbnl (long double x, int n) - * scalbnl(x,n) returns x* 2**n computed by exponent - * manipulation rather than by actually performing an - * exponentiation or a multiplication. - */ - -#include <math.h> -#include <math_private.h> - -static const long double -two64 = 1.8446744073709551616e19L, -twom64 = 5.421010862427522170037e-20L, -huge = 1.0e+4900L, -tiny = 1.0e-4900L; - -long double -__scalbnl (long double x, int n) -{ - int32_t k,es,hx,lx; - GET_LDOUBLE_WORDS(es,hx,lx,x); - k = es&0x7fff; /* extract exponent */ - if (__builtin_expect(k==0, 0)) { /* 0 or subnormal x */ - if ((lx|(hx&0x7fffffff))==0) return x; /* +-0 */ - x *= two64; - GET_LDOUBLE_EXP(hx,x); - k = (hx&0x7fff) - 64; - } - if (__builtin_expect(k==0x7fff, 0)) return x+x; /* NaN or Inf */ - if (__builtin_expect(n< -50000, 0)) - return tiny*__copysignl(tiny,x); - if (__builtin_expect(n> 50000 || k+n > 0x7ffe, 0)) - return huge*__copysignl(huge,x); /* overflow */ - /* Now k and n are bounded we know that k = k+n does not - overflow. */ - k = k+n; - if (__builtin_expect(k > 0, 1)) /* normal result */ - {SET_LDOUBLE_EXP(x,(es&0x8000)|k); return x;} - if (k <= -64) - return tiny*__copysignl(tiny,x); /*underflow*/ - k += 64; /* subnormal result */ - SET_LDOUBLE_EXP(x,(es&0x8000)|k); - return x*twom64; -} -weak_alias (__scalbnl, scalbnl) diff --git a/libc/sysdeps/ieee754/ldbl-96/s_signbitl.c b/libc/sysdeps/ieee754/ldbl-96/s_signbitl.c index bc2ac303f..8256d9ddb 100644 --- a/libc/sysdeps/ieee754/ldbl-96/s_signbitl.c +++ b/libc/sysdeps/ieee754/ldbl-96/s_signbitl.c @@ -1,5 +1,5 @@ /* Return nonzero value if number is negative. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-96/s_sincosl.c b/libc/sysdeps/ieee754/ldbl-96/s_sincosl.c index 5a470ad3b..af8f899c2 100644 --- a/libc/sysdeps/ieee754/ldbl-96/s_sincosl.c +++ b/libc/sysdeps/ieee754/ldbl-96/s_sincosl.c @@ -1,5 +1,5 @@ /* Compute sine and cosine of argument. - Copyright (C) 1997-2013 Free Software Foundation, Inc. + Copyright (C) 1997-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. diff --git a/libc/sysdeps/ieee754/ldbl-96/s_truncl.c b/libc/sysdeps/ieee754/ldbl-96/s_truncl.c deleted file mode 100644 index cc58fb0b7..000000000 --- a/libc/sysdeps/ieee754/ldbl-96/s_truncl.c +++ /dev/null @@ -1,56 +0,0 @@ -/* Truncate argument to nearest integral value not larger than the argument. - Copyright (C) 1997-2013 Free Software Foundation, Inc. - This file is part of the GNU C Library. - Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. - - The GNU C Library is free software; you can redistribute it and/or - modify it under the terms of the GNU Lesser General Public - License as published by the Free Software Foundation; either - version 2.1 of the License, or (at your option) any later version. - - The GNU C Library is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU - Lesser General Public License for more details. - - You should have received a copy of the GNU Lesser General Public - License along with the GNU C Library; if not, see - <http://www.gnu.org/licenses/>. */ - -#include <math.h> - -#include <math_private.h> - - -long double -__truncl (long double x) -{ - int32_t i0, j0; - u_int32_t se, i1; - int sx; - - GET_LDOUBLE_WORDS (se, i0, i1, x); - sx = se & 0x8000; - j0 = (se & 0x7fff) - 0x3fff; - if (j0 < 31) - { - if (j0 < 0) - /* The magnitude of the number is < 1 so the result is +-0. */ - SET_LDOUBLE_WORDS (x, sx, 0, 0); - else - SET_LDOUBLE_WORDS (x, se, i0 & ~(0x7fffffff >> j0), 0); - } - else if (j0 > 63) - { - if (j0 == 0x4000) - /* x is inf or NaN. */ - return x + x; - } - else - { - SET_LDOUBLE_WORDS (x, se, i0, i1 & ~(0xffffffffu >> (j0 - 31))); - } - - return x; -} -weak_alias (__truncl, truncl) diff --git a/libc/sysdeps/ieee754/ldbl-96/strtold_l.c b/libc/sysdeps/ieee754/ldbl-96/strtold_l.c index dccf98c46..82c33e52d 100644 --- a/libc/sysdeps/ieee754/ldbl-96/strtold_l.c +++ b/libc/sysdeps/ieee754/ldbl-96/strtold_l.c @@ -1,4 +1,4 @@ -/* Copyright (C) 1999-2013 Free Software Foundation, Inc. +/* Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-96/t_sincosl.c b/libc/sysdeps/ieee754/ldbl-96/t_sincosl.c index 02a15009a..da8fcfa11 100644 --- a/libc/sysdeps/ieee754/ldbl-96/t_sincosl.c +++ b/libc/sysdeps/ieee754/ldbl-96/t_sincosl.c @@ -1,5 +1,5 @@ /* Extended-precision floating point sine and cosine tables. - Copyright (C) 1999-2013 Free Software Foundation, Inc. + Copyright (C) 1999-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Based on quad-precision tables by Jakub Jelinek <jj@ultra.linux.cz> diff --git a/libc/sysdeps/ieee754/ldbl-96/w_expl.c b/libc/sysdeps/ieee754/ldbl-96/w_expl.c index 8e8dd1864..af221184d 100644 --- a/libc/sysdeps/ieee754/ldbl-96/w_expl.c +++ b/libc/sysdeps/ieee754/ldbl-96/w_expl.c @@ -1,4 +1,4 @@ -/* Copyright (C) 2011-2013 Free Software Foundation, Inc. +/* Copyright (C) 2011-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@gmail.com>, 2011. diff --git a/libc/sysdeps/ieee754/ldbl-96/x2y2m1.c b/libc/sysdeps/ieee754/ldbl-96/x2y2m1.c index f5185b27d..6ae2a63cb 100644 --- a/libc/sysdeps/ieee754/ldbl-96/x2y2m1.c +++ b/libc/sysdeps/ieee754/ldbl-96/x2y2m1.c @@ -1,5 +1,5 @@ /* Compute x^2 + y^2 - 1, without large cancellation error. - Copyright (C) 2012-2013 Free Software Foundation, Inc. + Copyright (C) 2012-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-96/x2y2m1l.c b/libc/sysdeps/ieee754/ldbl-96/x2y2m1l.c index 7bf9e297b..575bcd812 100644 --- a/libc/sysdeps/ieee754/ldbl-96/x2y2m1l.c +++ b/libc/sysdeps/ieee754/ldbl-96/x2y2m1l.c @@ -1,5 +1,5 @@ /* Compute x^2 + y^2 - 1, without large cancellation error. - Copyright (C) 2012-2013 Free Software Foundation, Inc. + Copyright (C) 2012-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. The GNU C Library is free software; you can redistribute it and/or diff --git a/libc/sysdeps/ieee754/ldbl-opt/nldbl-compat.c b/libc/sysdeps/ieee754/ldbl-opt/nldbl-compat.c index 783846113..91900b861 100644 --- a/libc/sysdeps/ieee754/ldbl-opt/nldbl-compat.c +++ b/libc/sysdeps/ieee754/ldbl-opt/nldbl-compat.c @@ -1,5 +1,5 @@ /* *printf* family compatibility routines for IEEE double as long double - Copyright (C) 2006-2013 Free Software Foundation, Inc. + Copyright (C) 2006-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jakub@cygnus.com>, 2006. diff --git a/libc/sysdeps/ieee754/ldbl-opt/nldbl-compat.h b/libc/sysdeps/ieee754/ldbl-opt/nldbl-compat.h index 040853650..d1e09892f 100644 --- a/libc/sysdeps/ieee754/ldbl-opt/nldbl-compat.h +++ b/libc/sysdeps/ieee754/ldbl-opt/nldbl-compat.h @@ -1,5 +1,5 @@ /* Prototypes for compatibility double == long double entry points. - Copyright (C) 2006-2013 Free Software Foundation, Inc. + Copyright (C) 2006-2014 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Jakub Jelinek <jakub@cygnus.com>, 2006. |