Optimize libm

libm is now somewhat integrated with gcc's -ffinite-math-only option
and lots of the wrapper functions have been optimized.

1 files changed, 62 insertions, 50 deletions

diff --git a/sysdeps/ieee754/flt-32/e_atanhf.c b/sysdeps/ieee754/flt-32/e_atanhf.c
index f26a15bbc1..ddd18ab300 100644
--- a/sysdeps/ieee754/flt-32/e_atanhf.c
+++ b/sysdeps/ieee754/flt-32/e_atanhf.c
@@ -1,58 +1,70 @@
-/* e_atanhf.c -- float version of e_atanh.c.
- * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
- */
+/* Copyright (C) 2011 Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+   Contributed by Ulrich Drepper <drepper@gmail.com>, 2011.
 
-/*
- * ====================================================
- * 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.
- * ====================================================
- */
+   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, write to the Free
+   Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
+   02111-1307 USA.  */
+
+
+/* __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))
 
-#if defined(LIBM_SCCS) && !defined(lint)
-static char rcsid[] = "$NetBSD: e_atanhf.c,v 1.4 1995/05/10 20:44:56 jtc Exp $";
-#endif
+   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 <inttypes.h>
 #include "math.h"
 #include "math_private.h"
 
-#ifdef __STDC__
-static const float one = 1.0, huge = 1e30;
-#else
-static float one = 1.0, huge = 1e30;
-#endif
-
-#ifdef __STDC__
-static const float zero = 0.0;
-#else
-static float zero = 0.0;
-#endif
-
-#ifdef __STDC__
-	float __ieee754_atanhf(float x)
-#else
-	float __ieee754_atanhf(x)
-	float x;
-#endif
+static const float huge = 1e30;
+
+float
+__ieee754_atanhf (float x)
 {
-	float t;
-	int32_t hx,ix;
-	GET_FLOAT_WORD(hx,x);
-	ix = hx&0x7fffffff;
-	if (ix>0x3f800000) 		/* |x|>1 */
-	    return (x-x)/(x-x);
-	if(ix==0x3f800000) 
-	    return x/zero;
-	if(ix<0x31800000&&(huge+x)>zero) return x;	/* x<2**-28 */
-	SET_FLOAT_WORD(x,ix);
-	if(ix<0x3f000000) {		/* x < 0.5 */
-	    t = x+x;
-	    t = (float)0.5*__log1pf(t+t*x/(one-x));
-	} else 
-	    t = (float)0.5*__log1pf((x+x)/(one-x));
-	if(hx>=0) return t; else return -t;
+  float xa = fabsf (x);
+  float t;
+  if (xa < 0.5f)
+    {
+      if (__builtin_expect (xa < 0x1.0p-28f, 0) && (huge + x) > 0.0f)
+	return x;
+
+      t = xa + xa;
+      t = 0.5f * __log1pf (t + t * xa / (1.0f - xa));
+    }
+  else if (__builtin_expect (xa < 1.0f, 1))
+    t = 0.5f * __log1pf ((xa + xa) / (1.0f - xa));
+  else
+    {
+      if (xa > 1.0f)
+	return (x - x) / (x - x);
+
+      return x / 0.0f;
+    }
+
+  return __copysignf (t, x);
 }
+strong_alias (__ieee754_atanhf, __atanhf_finite)