[cmucl-commit] [git] CMU Common Lisp branch master updated. snapshot-2014-06-68-g60d71b6
Raymond Toy
rtoy at common-lisp.net
Sat Aug 2 22:37:44 UTC 2014
This is an automated email from the git hooks/post-receive script. It was
generated because a ref change was pushed to the repository containing
the project "CMU Common Lisp".
The branch, master has been updated
via 60d71b68ef6af3e5fb3e045f7ce2a42a5ff23ac1 (commit)
via 1b876280ab1774cb36fb63906c58cd5af1e5d230 (commit)
from 21aad0ecc9f96297deedf9362c19fd877910592a (commit)
Those revisions listed above that are new to this repository have
not appeared on any other notification email; so we list those
revisions in full, below.
- Log -----------------------------------------------------------------
commit 60d71b68ef6af3e5fb3e045f7ce2a42a5ff23ac1
Author: Raymond Toy <toy.raymond at gmail.com>
Date: Sat Aug 2 15:30:51 2014 -0700
Fix up the inverse hyperbolics.
* fdlibm.h:
* Declare fdlibm_log1p and __ieee754_log
* e_acosh.c:
* Use unions to get the parts of a double.
* Use fdlibm_log1p instead of log1p.
* e_atanh.c:
* Use unions to get the parts of a double.
* Use fdlibm_log1p instead of log1p.
* s_asinh.c:
* Rename from asinh to fdlibm_asinh.
* Use unions to get the parts of a double.
* Use fdlibm_log1p instead of log1p.
diff --git a/src/lisp/e_acosh.c b/src/lisp/e_acosh.c
index 6839554..fb248d5 100644
--- a/src/lisp/e_acosh.c
+++ b/src/lisp/e_acosh.c
@@ -45,7 +45,10 @@ ln2 = 6.93147180559945286227e-01; /* 0x3FE62E42, 0xFEFA39EF */
{
double t;
int hx;
- hx = __HI(x);
+ union { int i[2]; double d; } ux;
+
+ ux.d = x;
+ hx = ux.i[HIWORD];
if(hx<0x3ff00000) { /* x < 1 */
return (x-x)/(x-x);
} else if(hx >=0x41b00000) { /* x > 2**28 */
@@ -53,13 +56,13 @@ ln2 = 6.93147180559945286227e-01; /* 0x3FE62E42, 0xFEFA39EF */
return x+x;
} else
return __ieee754_log(x)+ln2; /* acosh(huge)=log(2x) */
- } else if(((hx-0x3ff00000)|__LO(x))==0) {
+ } else if(((hx-0x3ff00000)|ux.i[LOWORD])==0) {
return 0.0; /* acosh(1) = 0 */
} else if (hx > 0x40000000) { /* 2**28 > x > 2 */
t=x*x;
return __ieee754_log(2.0*x-one/(x+sqrt(t-one)));
} else { /* 1<x<2 */
t = x-one;
- return log1p(t+sqrt(2.0*t+t*t));
+ return fdlibm_log1p(t+sqrt(2.0*t+t*t));
}
}
diff --git a/src/lisp/e_atanh.c b/src/lisp/e_atanh.c
index 302a89f..3588436 100644
--- a/src/lisp/e_atanh.c
+++ b/src/lisp/e_atanh.c
@@ -50,19 +50,24 @@ static double zero = 0.0;
double t;
int hx,ix;
unsigned lx;
- hx = __HI(x); /* high word */
- lx = __LO(x); /* low word */
+ union { int i[2]; double d; } ux;
+
+ ux.d = x;
+ hx = ux.i[HIWORD]; /* high word */
+ lx = ux.i[LOWORD]; /* low word */
ix = hx&0x7fffffff;
if ((ix|((lx|(-lx))>>31))>0x3ff00000) /* |x|>1 */
return (x-x)/(x-x);
if(ix==0x3ff00000)
return x/zero;
if(ix<0x3e300000&&(huge+x)>zero) return x; /* x<2**-28 */
- __HI(x) = ix; /* x <- |x| */
+ ux.d = x;
+ ux.i[HIWORD] = ix; /* x <- |x| */
+ x = ux.d;
if(ix<0x3fe00000) { /* x < 0.5 */
t = x+x;
- t = 0.5*log1p(t+t*x/(one-x));
+ t = 0.5*fdlibm_log1p(t+t*x/(one-x));
} else
- t = 0.5*log1p((x+x)/(one-x));
+ t = 0.5*fdlibm_log1p((x+x)/(one-x));
if(hx>=0) return t; else return -t;
}
diff --git a/src/lisp/fdlibm.h b/src/lisp/fdlibm.h
index 325d7d4..72f57f8 100644
--- a/src/lisp/fdlibm.h
+++ b/src/lisp/fdlibm.h
@@ -52,7 +52,9 @@ extern double fdlibm_sin(double x);
extern double fdlibm_cos(double x);
extern double fdlibm_tan(double x);
extern double fdlibm_expm1(double x);
+extern double fdlibm_log1p(double x);
extern double __ieee754_exp(double x);
+extern double __ieee754_log(double x);
#endif
diff --git a/src/lisp/s_asinh.c b/src/lisp/s_asinh.c
index c809467..6cae81f 100644
--- a/src/lisp/s_asinh.c
+++ b/src/lisp/s_asinh.c
@@ -34,15 +34,18 @@ ln2 = 6.93147180559945286227e-01, /* 0x3FE62E42, 0xFEFA39EF */
huge= 1.00000000000000000000e+300;
#ifdef __STDC__
- double asinh(double x)
+ double fdlibm_asinh(double x)
#else
- double asinh(x)
+ double fdlibm_asinh(x)
double x;
#endif
{
double t,w;
int hx,ix;
- hx = __HI(x);
+ union { int i[2]; double d; } ux;
+
+ ux.d = x;
+ hx = ux.i[HIWORD];
ix = hx&0x7fffffff;
if(ix>=0x7ff00000) return x+x; /* x is inf or NaN */
if(ix< 0x3e300000) { /* |x|<2**-28 */
@@ -55,7 +58,7 @@ huge= 1.00000000000000000000e+300;
w = __ieee754_log(2.0*t+one/(sqrt(x*x+one)+t));
} else { /* 2.0 > |x| > 2**-28 */
t = x*x;
- w =log1p(fabs(x)+t/(one+sqrt(one+t)));
+ w =fdlibm_log1p(fabs(x)+t/(one+sqrt(one+t)));
}
if(hx>0) return w; else return -w;
}
commit 1b876280ab1774cb36fb63906c58cd5af1e5d230
Author: Raymond Toy <toy.raymond at gmail.com>
Date: Sat Aug 2 15:17:58 2014 -0700
Import inverse hyperbolic functions from fdlibm, as is.
diff --git a/src/lisp/e_acosh.c b/src/lisp/e_acosh.c
new file mode 100644
index 0000000..6839554
--- /dev/null
+++ b/src/lisp/e_acosh.c
@@ -0,0 +1,65 @@
+
+/* @(#)e_acosh.c 1.3 95/01/18 */
+/*
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ *
+ */
+
+/* __ieee754_acosh(x)
+ * Method :
+ * Based on
+ * acosh(x) = log [ x + sqrt(x*x-1) ]
+ * we have
+ * acosh(x) := log(x)+ln2, if x is large; else
+ * acosh(x) := log(2x-1/(sqrt(x*x-1)+x)) if x>2; else
+ * acosh(x) := log1p(t+sqrt(2.0*t+t*t)); where t=x-1.
+ *
+ * Special cases:
+ * acosh(x) is NaN with signal if x<1.
+ * acosh(NaN) is NaN without signal.
+ */
+
+#include "fdlibm.h"
+
+#ifdef __STDC__
+static const double
+#else
+static double
+#endif
+one = 1.0,
+ln2 = 6.93147180559945286227e-01; /* 0x3FE62E42, 0xFEFA39EF */
+
+#ifdef __STDC__
+ double __ieee754_acosh(double x)
+#else
+ double __ieee754_acosh(x)
+ double x;
+#endif
+{
+ double t;
+ int hx;
+ hx = __HI(x);
+ if(hx<0x3ff00000) { /* x < 1 */
+ return (x-x)/(x-x);
+ } else if(hx >=0x41b00000) { /* x > 2**28 */
+ if(hx >=0x7ff00000) { /* x is inf of NaN */
+ return x+x;
+ } else
+ return __ieee754_log(x)+ln2; /* acosh(huge)=log(2x) */
+ } else if(((hx-0x3ff00000)|__LO(x))==0) {
+ return 0.0; /* acosh(1) = 0 */
+ } else if (hx > 0x40000000) { /* 2**28 > x > 2 */
+ t=x*x;
+ return __ieee754_log(2.0*x-one/(x+sqrt(t-one)));
+ } else { /* 1<x<2 */
+ t = x-one;
+ return log1p(t+sqrt(2.0*t+t*t));
+ }
+}
diff --git a/src/lisp/e_atanh.c b/src/lisp/e_atanh.c
new file mode 100644
index 0000000..302a89f
--- /dev/null
+++ b/src/lisp/e_atanh.c
@@ -0,0 +1,68 @@
+
+/* @(#)e_atanh.c 1.3 95/01/18 */
+/*
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ *
+ */
+
+/* __ieee754_atanh(x)
+ * Method :
+ * 1.Reduced x to positive by atanh(-x) = -atanh(x)
+ * 2.For x>=0.5
+ * 1 2x x
+ * atanh(x) = --- * log(1 + -------) = 0.5 * log1p(2 * --------)
+ * 2 1 - x 1 - x
+ *
+ * For x<0.5
+ * atanh(x) = 0.5*log1p(2x+2x*x/(1-x))
+ *
+ * Special cases:
+ * atanh(x) is NaN if |x| > 1 with signal;
+ * atanh(NaN) is that NaN with no signal;
+ * atanh(+-1) is +-INF with signal.
+ *
+ */
+
+#include "fdlibm.h"
+
+#ifdef __STDC__
+static const double one = 1.0, huge = 1e300;
+#else
+static double one = 1.0, huge = 1e300;
+#endif
+
+static double zero = 0.0;
+
+#ifdef __STDC__
+ double __ieee754_atanh(double x)
+#else
+ double __ieee754_atanh(x)
+ double x;
+#endif
+{
+ double t;
+ int hx,ix;
+ unsigned lx;
+ hx = __HI(x); /* high word */
+ lx = __LO(x); /* low word */
+ ix = hx&0x7fffffff;
+ if ((ix|((lx|(-lx))>>31))>0x3ff00000) /* |x|>1 */
+ return (x-x)/(x-x);
+ if(ix==0x3ff00000)
+ return x/zero;
+ if(ix<0x3e300000&&(huge+x)>zero) return x; /* x<2**-28 */
+ __HI(x) = ix; /* x <- |x| */
+ if(ix<0x3fe00000) { /* x < 0.5 */
+ t = x+x;
+ t = 0.5*log1p(t+t*x/(one-x));
+ } else
+ t = 0.5*log1p((x+x)/(one-x));
+ if(hx>=0) return t; else return -t;
+}
diff --git a/src/lisp/s_asinh.c b/src/lisp/s_asinh.c
new file mode 100644
index 0000000..c809467
--- /dev/null
+++ b/src/lisp/s_asinh.c
@@ -0,0 +1,61 @@
+
+/* @(#)s_asinh.c 1.3 95/01/18 */
+/*
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+/* asinh(x)
+ * Method :
+ * Based on
+ * asinh(x) = sign(x) * log [ |x| + sqrt(x*x+1) ]
+ * we have
+ * asinh(x) := x if 1+x*x=1,
+ * := sign(x)*(log(x)+ln2)) for large |x|, else
+ * := sign(x)*log(2|x|+1/(|x|+sqrt(x*x+1))) if|x|>2, else
+ * := sign(x)*log1p(|x| + x^2/(1 + sqrt(1+x^2)))
+ */
+
+#include "fdlibm.h"
+
+#ifdef __STDC__
+static const double
+#else
+static double
+#endif
+one = 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
+ln2 = 6.93147180559945286227e-01, /* 0x3FE62E42, 0xFEFA39EF */
+huge= 1.00000000000000000000e+300;
+
+#ifdef __STDC__
+ double asinh(double x)
+#else
+ double asinh(x)
+ double x;
+#endif
+{
+ double t,w;
+ int hx,ix;
+ hx = __HI(x);
+ ix = hx&0x7fffffff;
+ if(ix>=0x7ff00000) return x+x; /* x is inf or NaN */
+ if(ix< 0x3e300000) { /* |x|<2**-28 */
+ if(huge+x>one) return x; /* return x inexact except 0 */
+ }
+ if(ix>0x41b00000) { /* |x| > 2**28 */
+ w = __ieee754_log(fabs(x))+ln2;
+ } else if (ix>0x40000000) { /* 2**28 > |x| > 2.0 */
+ t = fabs(x);
+ w = __ieee754_log(2.0*t+one/(sqrt(x*x+one)+t));
+ } else { /* 2.0 > |x| > 2**-28 */
+ t = x*x;
+ w =log1p(fabs(x)+t/(one+sqrt(one+t)));
+ }
+ if(hx>0) return w; else return -w;
+}
-----------------------------------------------------------------------
Summary of changes:
src/lisp/e_acosh.c | 68 ++++++++++++++++++++++++++++++++++++++++++++++++
src/lisp/e_atanh.c | 73 ++++++++++++++++++++++++++++++++++++++++++++++++++++
src/lisp/fdlibm.h | 2 ++
src/lisp/s_asinh.c | 64 +++++++++++++++++++++++++++++++++++++++++++++
4 files changed, 207 insertions(+)
create mode 100644 src/lisp/e_acosh.c
create mode 100644 src/lisp/e_atanh.c
create mode 100644 src/lisp/s_asinh.c
hooks/post-receive
--
CMU Common Lisp
More information about the cmucl-commit
mailing list