[Openmcl-cvs-notifications] r10966 - /trunk/source/lisp-kernel/windows-calls.c

[gb at clozure.com]

Author: gb
Date: Fri Oct  3 19:53:46 2008
New Revision: 10966

Log:
Provide win64 versions of acosh[f],asinh[f],atanh[f], using public-domain
code from mingw32.

Modified:
    trunk/source/lisp-kernel/windows-calls.c

Modified: trunk/source/lisp-kernel/windows-calls.c
=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=
=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=
=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=
=3D=3D=3D
--- trunk/source/lisp-kernel/windows-calls.c (original)
+++ trunk/source/lisp-kernel/windows-calls.c Fri Oct  3 19:53:46 2008
@@ -455,6 +455,161 @@
 {
   return gettimeofday(tp, tzp);
 }
+
+#ifdef WIN_64
+
+/* Make sure that the lisp calls these functions, when they do something */
+/* This code is taken from the 32-bit mingw library and is in the
+   public domain */
+double
+acosh(double x)
+{
+  if (isnan (x)) =

+    return x;
+
+  if (x < 1.0)
+    {
+      errno =3D EDOM;
+      return nan("");
+    }
+
+  if (x > 0x1p32)
+    /*  Avoid overflow (and unnecessary calculation when
+        sqrt (x * x - 1) =3D=3D x). GCC optimizes by replacing
+        the long double M_LN2 const with a fldln2 insn.  */ =

+    return log (x) + 6.9314718055994530941723E-1L;
+
+  /* Since  x >=3D 1, the arg to log will always be greater than
+     the fyl2xp1 limit (approx 0.29) so just use logl. */ =

+  return log (x + sqrt((x + 1.0) * (x - 1.0)));
+}
+
+float
+acoshf(float x)
+{
+  if (isnan (x)) =

+    return x;
+  if (x < 1.0f)
+    {
+      errno =3D EDOM;
+      return nan("");
+    }
+
+ if (x > 0x1p32f)
+    /*  Avoid overflow (and unnecessary calculation when
+        sqrt (x * x - 1) =3D=3D x). GCC optimizes by replacing
+        the long double M_LN2 const with a fldln2 insn.  */ =

+    return log (x) + 6.9314718055994530941723E-1L;
+
+  /* Since  x >=3D 1, the arg to log will always be greater than
+     the fyl2xp1 limit (approx 0.29) so just use logl. */ =

+  return log (x + sqrt((x + 1.0) * (x - 1.0)));
+}
+
+double
+asinh(double x)
+{
+  double z;
+  if (!isfinite (x))
+    return x;
+  z =3D fabs (x);
+
+  /* Avoid setting FPU underflow exception flag in x * x. */
+#if 0
+  if ( z < 0x1p-32)
+    return x;
+#endif
+
+  /* Use log1p to avoid cancellation with small x. Put
+     x * x in denom, so overflow is harmless. =

+     asinh(x) =3D log1p (x + sqrt (x * x + 1.0) - 1.0)
+              =3D log1p (x + x * x / (sqrt (x * x + 1.0) + 1.0))  */
+
+  z =3D log1p (z + z * z / (sqrt (z * z + 1.0) + 1.0));
+
+  return ( x > 0.0 ? z : -z);
+}
+
+float
+asinhf(float x)
+{
+  float z;
+  if (!isfinite (x))
+    return x;
+  z =3D fabsf (x);
+
+  /* Avoid setting FPU underflow exception flag in x * x. */
+#if 0
+  if ( z < 0x1p-32)
+    return x;
+#endif
+
+
+  /* Use log1p to avoid cancellation with small x. Put
+     x * x in denom, so overflow is harmless. =

+     asinh(x) =3D log1p (x + sqrt (x * x + 1.0) - 1.0)
+              =3D log1p (x + x * x / (sqrt (x * x + 1.0) + 1.0))  */
+
+  z =3D log1p (z + z * z / (sqrt (z * z + 1.0) + 1.0));
+
+  return ( x > 0.0 ? z : -z);
+}
+
+double
+atanh(double x)
+{
+  double z;
+  if (isnan (x))
+    return x;
+  z =3D fabs (x);
+  if (z =3D=3D 1.0)
+    {
+      errno  =3D ERANGE;
+      return (x > 0 ? INFINITY : -INFINITY);
+    }
+  if (z > 1.0)
+    {
+      errno =3D EDOM;
+      return nan("");
+    }
+  /* Rearrange formula to avoid precision loss for small x.
+
+  atanh(x) =3D 0.5 * log ((1.0 + x)/(1.0 - x))
+	   =3D 0.5 * log1p ((1.0 + x)/(1.0 - x) - 1.0)
+           =3D 0.5 * log1p ((1.0 + x - 1.0 + x) /(1.0 - x)) =

+           =3D 0.5 * log1p ((2.0 * x ) / (1.0 - x))  */
+  z =3D 0.5 * log1p ((z + z) / (1.0 - z));
+  return x >=3D 0 ? z : -z;
+}
+
+float
+atanhf(float x)
+{
+  float z;
+  if (isnan (x))
+    return x;
+  z =3D fabsf (x);
+  if (z =3D=3D 1.0)
+    {
+      errno  =3D ERANGE;
+      return (x > 0 ? INFINITY : -INFINITY);
+    }
+  if ( z > 1.0)
+    {
+      errno =3D EDOM;
+      return nanf("");
+    }
+  /* Rearrange formula to avoid precision loss for small x.
+
+  atanh(x) =3D 0.5 * log ((1.0 + x)/(1.0 - x))
+	   =3D 0.5 * log1p ((1.0 + x)/(1.0 - x) - 1.0)
+           =3D 0.5 * log1p ((1.0 + x - 1.0 + x) /(1.0 - x)) =

+           =3D 0.5 * log1p ((2.0 * x ) / (1.0 - x))  */
+  z =3D 0.5 * log1p ((z + z) / (1.0 - z));
+  return x >=3D 0 ? z : -z;
+}
+
+#endif
 =

 typedef struct {
   char *name;