Author: kargl
Date: Mon Sep 1 18:57:13 2014
New Revision: 270947
URL: http://svnweb.freebsd.org/changeset/base/270947
Log:
According to the ISO C standard, lgamma(-integer) returns
inf and raises the divided-by-zero exception. Compilers
constant fold one/zero to inf but do not raise the exception.
Introduce a volatile vzero to prevent the constant folding.
Move the declaration of zero into the main declaration block.
While here, fix a nearby disordering of 'lx,ix'
Discussed with: bde
Modified:
head/lib/msun/src/e_lgamma_r.c
head/lib/msun/src/e_lgammaf_r.c
Modified: head/lib/msun/src/e_lgamma_r.c
==============================================================================
--- head/lib/msun/src/e_lgamma_r.c Mon Sep 1 18:56:22 2014
(r270946)
+++ head/lib/msun/src/e_lgamma_r.c Mon Sep 1 18:57:13 2014
(r270947)
@@ -86,7 +86,10 @@ __FBSDID("$FreeBSD$");
#include "math.h"
#include "math_private.h"
+static const volatile double vzero = 0;
+
static const double
+zero= 0.00000000000000000000e+00,
two52= 4.50359962737049600000e+15, /* 0x43300000, 0x00000000 */
half= 5.00000000000000000000e-01, /* 0x3FE00000, 0x00000000 */
one = 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
@@ -154,8 +157,6 @@ w4 = -5.95187557450339963135e-04, /* 0x
w5 = 8.36339918996282139126e-04, /* 0x3F4B67BA, 0x4CDAD5D1 */
w6 = -1.63092934096575273989e-03; /* 0xBF5AB89D, 0x0B9E43E4 */
-static const double zero= 0.00000000000000000000e+00;
-
/*
* Compute sin(pi*x) without actually doing the pi*x multiplication.
* sin_pi(x) is only called for x < 0 and |x| < 2**(p-1) where p is
@@ -204,7 +205,7 @@ __ieee754_lgamma_r(double x, int *signga
{
double t,y,z,nadj,p,p1,p2,p3,q,r,w;
int32_t hx;
- int i,lx,ix;
+ int i,ix,lx;
EXTRACT_WORDS(hx,lx,x);
@@ -212,7 +213,7 @@ __ieee754_lgamma_r(double x, int *signga
*signgamp = 1;
ix = hx&0x7fffffff;
if(ix>=0x7ff00000) return x*x;
- if((ix|lx)==0) return one/zero;
+ if((ix|lx)==0) return one/vzero;
if(ix<0x3b900000) { /* |x|<2**-70, return -log(|x|) */
if(hx<0) {
*signgamp = -1;
@@ -221,9 +222,9 @@ __ieee754_lgamma_r(double x, int *signga
}
if(hx<0) {
if(ix>=0x43300000) /* |x|>=2**52, must be -integer */
- return one/zero;
+ return one/vzero;
t = sin_pi(x);
- if(t==zero) return one/zero; /* -integer */
+ if(t==zero) return one/vzero; /* -integer */
nadj = __ieee754_log(pi/fabs(t*x));
if(t<zero) *signgamp = -1;
x = -x;
Modified: head/lib/msun/src/e_lgammaf_r.c
==============================================================================
--- head/lib/msun/src/e_lgammaf_r.c Mon Sep 1 18:56:22 2014
(r270946)
+++ head/lib/msun/src/e_lgammaf_r.c Mon Sep 1 18:57:13 2014
(r270947)
@@ -19,7 +19,10 @@ __FBSDID("$FreeBSD$");
#include "math.h"
#include "math_private.h"
+static const volatile float vzero = 0;
+
static const float
+zero= 0.0000000000e+00,
two23= 8.3886080000e+06, /* 0x4b000000 */
half= 5.0000000000e-01, /* 0x3f000000 */
one = 1.0000000000e+00, /* 0x3f800000 */
@@ -87,8 +90,6 @@ w4 = -5.9518753551e-04, /* 0xba1c065c *
w5 = 8.3633989561e-04, /* 0x3a5b3dd2 */
w6 = -1.6309292987e-03; /* 0xbad5c4e8 */
-static const float zero= 0.0000000000e+00;
-
static float
sin_pif(float x)
{
@@ -140,7 +141,7 @@ __ieee754_lgammaf_r(float x, int *signga
*signgamp = 1;
ix = hx&0x7fffffff;
if(ix>=0x7f800000) return x*x;
- if(ix==0) return one/zero;
+ if(ix==0) return one/vzero;
if(ix<0x35000000) { /* |x|<2**-21, return -log(|x|) */
if(hx<0) {
*signgamp = -1;
@@ -149,9 +150,9 @@ __ieee754_lgammaf_r(float x, int *signga
}
if(hx<0) {
if(ix>=0x4b000000) /* |x|>=2**23, must be -integer */
- return one/zero;
+ return one/vzero;
t = sin_pif(x);
- if(t==zero) return one/zero; /* -integer */
+ if(t==zero) return one/vzero; /* -integer */
nadj = __ieee754_logf(pi/fabsf(t*x));
if(t<zero) *signgamp = -1;
x = -x;
_______________________________________________
svn-src-all@freebsd.org mailing list
http://lists.freebsd.org/mailman/listinfo/svn-src-all
To unsubscribe, send any mail to "svn-src-all-unsubscr...@freebsd.org"
