Mon Mar 4 20:54:40 1996 Andreas Schwab <schwab@issan.informatik.uni-dortmund.de>

	* Makeconfig ($(common-objpfx)config.make): Depend on config.h.in.


Mon Mar  4 17:35:09 1996  Roland McGrath  <roland@charlie-brown.gnu.ai.mit.edu>

	* hurd/catch-signal.c (hurd_safe_memmove): New function.
	(hurd_safe_copyin, hurd_safe_copyout): New functions.
	* hurd/hurd/sigpreempt.h: Declare them.

Sun Mar  3 08:43:44 1996  Roland McGrath  <roland@charlie-brown.gnu.ai.mit.edu>

	Replace math code with fdlibm from Sun as modified for netbsd by
	JT Conklin and Ian Taylor, including x86 FPU support.
	* sysdeps/libm-ieee754, sysdeps/libm-i387: New directories.
	* math/math_private.h: New file.
	* sysdeps/i386/fpu/Implies: New file.
	* sysdeps/ieee754/Implies: New file.
	* math/machine/asm.h, math/machine/endian.h: New files.
	* math/Makefile, math/math.h: Rewritten.
	* mathcalls.h, math/mathcalls.h: New file, broken out of math.h.
	* math/finite.c: File removed.
	* sysdeps/generic/Makefile [$(subdir)=math]: Frobnication removed.

	* math/test-math.c: Include errno.h and string.h.

	* sysdeps/unix/bsd/dirstream.h: File removed.
	* sysdeps/unix/bsd/readdir.c: File removed.

7 files changed, 581 insertions, 282 deletions

diff --git a/math/Makefile b/math/Makefile
index c86d44a072..a2ceb06925 100644
--- a/math/Makefile
+++ b/math/Makefile
@@ -1,4 +1,6 @@
-# Copyright (C) 1991, 1992, 1994, 1995 Free Software Foundation, Inc.
+# Makefile for the math library.
+
+# Copyright (C) 1996 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
@@ -16,24 +18,58 @@
 # not, write to the Free Software Foundation, Inc., 675 Mass Ave,
 # Cambridge, MA 02139, USA.
 
-#
-#	Makefile for math.
-#
-subdir	:= math
+subdir		:= math
+
+# Installed header files.
+headers		:= math.h mathcalls.h __math.h huge_val.h nan.h
+
+# Internal header files.
+distribute	:= math_private.h machine/asm.h machine/endian.h
 
-headers	:= math.h __math.h huge_val.h nan.h
 
+# Build the -lm library.
 
-routines := acos asin atan cos sin tan cosh sinh tanh exp fabs ldexp	\
-	    log log10 floor sqrt fmod frexp pow atan2 ceil modf		\
-	    isinf isnan finite infnan copysign drem logb		\
-	    rint hypot cabs cbrt expm1 log1p acosh asinh atanh		\
-	    isinfl isnanl
-tests := # test-math
-install-lib := libm.a
+extra-libs	:= libm
+libm.so-version	:= .0
+libm-routines	:= e_acos e_acosf e_acosh e_acoshf e_asin e_asinf	\
+		   e_atan2 e_atan2f e_atanh e_atanhf e_cosh e_coshf	\
+		   e_exp e_expf e_fmod e_fmodf e_hypot e_hypotf e_j0	\
+		   e_j0f e_j1 e_j1f e_jn e_jnf e_lgamma_r e_lgammaf_r	\
+		   e_log e_log10 e_log10f e_logf e_pow e_powf		\
+		   e_rem_pio2 e_rem_pio2f e_remainder e_remainderf	\
+		   e_scalb e_scalbf e_sinh e_sinhf e_sqrt e_sqrtf k_cos	\
+		   k_cosf k_rem_pio2 k_rem_pio2f k_sin k_sinf		\
+		   k_standard k_tan k_tanf s_asinh s_asinhf s_atan	\
+		   s_atanf s_cbrt s_cbrtf s_ceil s_ceilf s_copysign	\
+		   s_copysignf s_cos s_cosf s_erf s_erff s_expm1	\
+		   s_expm1f s_fabs s_fabsf s_finite s_finitef s_floor	\
+		   s_floorf s_frexp s_frexpf s_ilogb s_ilogbf		\
+		   s_ldexp s_ldexpf s_lib_version s_log1p s_log1pf	\
+		   s_logb s_logbf s_matherr s_modf s_modff s_nextafter	\
+		   s_nextafterf s_rint s_rintf s_scalbn s_scalbnf	\
+		   s_signgam s_significand s_significandf s_sin s_sinf	\
+		   s_tan s_tanf s_tanh s_tanhf w_acos w_acosf w_acosh	\
+		   w_acoshf w_asin w_asinf w_atan2 w_atan2f w_atanh	\
+		   w_atanhf w_cabs w_cabsf w_cosh w_coshf w_drem	\
+		   w_dremf w_exp w_expf w_fmod w_fmodf w_gamma		\
+		   w_gamma_r w_gammaf w_gammaf_r w_hypot w_hypotf w_j0	\
+		   w_j0f w_j1 w_j1f w_jn w_jnf w_lgamma w_lgamma_r	\
+		   w_lgammaf w_lgammaf_r w_log w_log10 w_log10f w_logf	\
+		   w_pow w_powf w_remainder w_remainderf w_scalb	\
+		   w_scalbf w_sinh w_sinhf w_sqrt w_sqrtf
+
+# These functions are in libc instead of libm because __printf_fp
+# calls them, so any program using printf will need them linked in,
+# and we don't want to have to link every program with -lm.
+routines := s_isinf s_isinff s_isnan s_isnanf isinfl isnanl
 
 
 include ../Rules
 
-$(objpfx)libm.a: $(dep-dummy-lib); $(make-dummy-lib)
-lib: $(objpfx)libm.a
+
+# This file defines the default _LIB_VERSION variable that controls
+# the error return conventions for the math functions.
+CPPFLAGS-s_lib_version.c := -D_POSIX_MODE
+
+# The fdlibm code generates a lot of these warnings but is otherwise clean.
+override CFLAGS += -Wno-uninitialized -Wno-write-strings
diff --git a/math/finite.c b/math/finite.c
deleted file mode 100644
index 48f8f8011a..0000000000
--- a/math/finite.c
+++ /dev/null
@@ -1,31 +0,0 @@
-/* Copyright (C) 1991, 1992, 1995 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
-modify it under the terms of the GNU Library General Public License as
-published by the Free Software Foundation; either version 2 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
-Library General Public License for more details.
-
-You should have received a copy of the GNU Library General Public
-License along with the GNU C Library; see the file COPYING.LIB.  If
-not, write to the Free Software Foundation, Inc., 675 Mass Ave,
-Cambridge, MA 02139, USA.  */
-
-#include <ansidecl.h>
-#include <math.h>
-
-#undef	__finite
-
-/* Return nonzero if VALUE is finite and not NaN.  */
-int
-DEFUN(__finite, (value), double value)
-{
-  return !__isinf (value) && !__isnan (value);
-}
-
-weak_alias (__finite, finite)
diff --git a/math/machine/asm.h b/math/machine/asm.h
new file mode 100644
index 0000000000..548e248556
--- /dev/null
+++ b/math/machine/asm.h
@@ -0,0 +1,10 @@
+/* The libm assembly code wants to include <machine/asm.h> to define the
+   ENTRY macro.  We define assembly-related macros in sysdep.h and
+   asm-syntax.h.  */
+
+#include <sysdep.h>
+#include <asm-syntax.h>
+
+/* The libm assembly code uses this macro for RCSid strings.
+   We don't put RCSid strings into object files.  */
+#define RCSID(id) /* ignore them */
diff --git a/math/machine/endian.h b/math/machine/endian.h
new file mode 100644
index 0000000000..88804564bf
--- /dev/null
+++ b/math/machine/endian.h
@@ -0,0 +1,4 @@
+/* math_private.h wants to include <machine/endian.h>; we provide this
+   file so it can, leaving math_private.h unmodified from the original.  */
+
+#include <endian.h>
diff --git a/math/math.h b/math/math.h
index c30928de30..2693fcc9d0 100644
--- a/math/math.h
+++ b/math/math.h
@@ -1,4 +1,5 @@
-/* Copyright (C) 1991, 1992, 1993, 1995 Free Software Foundation, Inc.
+/* Declarations for math functions.
+Copyright (C) 1991, 92, 93, 95, 96 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
@@ -30,7 +31,8 @@ __BEGIN_DECLS
 #define	__need_Emath
 #include <errno.h>
 
-/* Get machine-dependent HUGE_VAL value (returned on overflow).  */
+/* Get machine-dependent HUGE_VAL value (returned on overflow).
+   On all IEEE754 machines, this is +Infinity.  */
 #include <huge_val.h>
 
 /* Get machine-dependent NAN value (returned for some domain errors).  */
@@ -39,260 +41,84 @@ __BEGIN_DECLS
 #endif
 
 
-/* Trigonometric functions.  */
-
-/* Arc cosine of X.  */
-extern double acos __P ((double __x)) __attribute__ ((__const__));
-/* Arc sine of X.  */
-extern double asin __P ((double __x)) __attribute__ ((__const__));
-/* Arc tangent of X.  */
-extern double atan __P ((double __x)) __attribute__ ((__const__));
-/* Arc tangent of Y/X.  */
-extern double atan2 __P ((double __y, double __x)) __attribute__ ((__const__));
-
-/* Cosine of X.  */
-extern double cos __P ((double __x)) __attribute__ ((__const__));
-/* Sine of X.  */
-extern double sin __P ((double __x)) __attribute__ ((__const__));
-/* Tangent of X.  */
-extern double tan __P ((double __x)) __attribute__ ((__const__));
-
-
-/* Hyperbolic functions.  */
-
-/* Hyperbolic cosine of X.  */
-extern double cosh __P ((double __x)) __attribute__ ((__const__));
-/* Hyperbolic sine of X.  */
-extern double sinh __P ((double __x)) __attribute__ ((__const__));
-/* Hyperbolic tangent of X.  */
-extern double tanh __P ((double __x)) __attribute__ ((__const__));
-
-#ifdef	__USE_MISC
-/* Hyperbolic arc cosine of X.  */
-extern double acosh __P ((double __x)) __attribute__ ((__const__));
-/* Hyperbolic arc sine of X.  */
-extern double asinh __P ((double __x)) __attribute__ ((__const__));
-/* Hyperbolic arc tangent of X.  */
-extern double atanh __P ((double __x)) __attribute__ ((__const__));
+/* The file <mathcalls.h> contains the prototypes for all the actual
+   math functions.  These macros are used for those prototypes, so
+   we can easily declare each function as both `name' and `__name',
+   and can declare the float versions `namef' and `__namef'.  */
+
+#define __MATHCALL(function,suffix, args)	\
+  __MATHDECL (_Mdouble_, function,suffix, args)
+#define __MATHDECL(type, function,suffix, args) \
+  __MATHDECL_1(type, function,suffix, args); \
+  __MATHDECL_1(type, __##function,suffix, args)
+#define __MATHDECL_1(type, function,suffix, args) \
+  extern type __MATH_PRECNAME(function,suffix) args
+
+#define _Mdouble_ 		double
+#define __MATH_PRECNAME(name,r)	name##r
+#include <mathcalls.h>
+#undef	_Mdouble_
+#undef	__MATH_PRECNAME
+
+#ifdef __USE_MISC
+/* Include the file of declarations again, this type using `float'
+   instead of `double' and appending f to each function name.  */
+
+#define _Mdouble_ 		float
+#define __MATH_PRECNAME(name,r)	name##f##r
+#include <mathcalls.h>
+#undef	_Mdouble_
+#undef	__MATH_PRECNAME
 #endif
 
-/* Exponential and logarithmic functions.  */
-
-/* Exponentional function of X.  */
-extern double exp __P ((double __x)) __attribute__ ((__const__));
-
-/* Break VALUE into a normalized fraction and an integral power of 2.  */
-extern double frexp __P ((double __value, int *__exp));
-
-/* X times (two to the EXP power).  */
-extern double ldexp __P ((double __x, int __exp)) __attribute__ ((__const__));
-
-/* Natural logarithm of X.  */
-extern double log __P ((double __x)) __attribute__ ((__const__));
-
-/* Base-ten logarithm of X.  */
-extern double log10 __P ((double __x)) __attribute__ ((__const__));
 
 #ifdef	__USE_MISC
-/* Return exp(X) - 1.  */
-extern double __expm1 __P ((double __x)) __attribute__ ((__const__));
-extern double expm1 __P ((double __x)) __attribute__ ((__const__));
-
-/* Return log(1 + X).  */
-extern double log1p __P ((double __x)) __attribute__ ((__const__));
-#endif
-
-/* Break VALUE into integral and fractional parts.  */
-extern double modf __P ((double __value, double *__iptr));
-
+/* Support for various different standard error handling behaviors.  */
 
-/* Power functions.  */
+typedef enum { _IEEE_ = -1, _SVID_, _XOPEN_, _POSIX_ } _LIB_VERSION_TYPE;
 
-/* Return X to the Y power.  */
-extern double pow __P ((double __x, double __y)) __attribute__ ((__const__));
-
-/* Return the square root of X.  */
-extern double sqrt __P ((double __x)) __attribute__ ((__const__));
-
-#ifdef	__USE_MISC
-/* Return the cube root of X.  */
-extern double cbrt __P ((double __x)) __attribute__ ((__const__));
+/* This variable can be changed at run-time to any of the values above to
+   affect floating point error handling behavior (it may also be necessary
+   to change the hardware FPU exception settings).  */
+extern _LIB_VERSION_TYPE _LIB_VERSION;
 #endif
 
 
-/* Nearest integer, absolute value, and remainder functions.  */
-
-/* Smallest integral value not less than X.  */
-extern double ceil __P ((double __x)) __attribute__ ((__const__));
-
-/* Absolute value of X.  */
-extern double fabs __P ((double __x)) __attribute__ ((__const__));
-
-/* Largest integer not greater than X.  */
-extern double floor __P ((double __x)) __attribute__ ((__const__));
-
-/* Floating-point modulo remainder of X/Y.  */
-extern double fmod __P ((double __x, double __y)) __attribute__ ((__const__));
-
+#ifdef __USE_SVID
+/* In SVID error handling, `matherr' is called with this description
+   of the exceptional condition.  */
+struct exception
+  {
+    int type;
+    char *name;
+    double arg1;
+    double arg2;
+    double retval;
+  };
 
-/* Return 0 if VALUE is finite or NaN, +1 if it
-   is +Infinity, -1 if it is -Infinity.  */
-extern int __isinf __P ((double __value)) __attribute__ ((__const__));
+extern int matherr __P ((struct exception *));
 
-/* Return nonzero if VALUE is not a number.  */
-extern int __isnan __P ((double __value)) __attribute__ ((__const__));
+#define X_TLOSS		1.41484755040568800000e+16
 
-/* Return nonzero if VALUE is finite and not NaN.  */
-extern int __finite __P ((double __value)) __attribute__ ((__const__));
-#ifdef	__OPTIMIZE__
-#define	__finite(value)	(!__isinf(value))
-#endif
-
-/* Deal with an infinite or NaN result.
-   If ERROR is ERANGE, result is +Inf;
-   if ERROR is - ERANGE, result is -Inf;
-   otherwise result is NaN.
-   This will set `errno' to either ERANGE or EDOM,
-   and may return an infinity or NaN, or may do something else.  */
-extern double __infnan __P ((int __error));
+/* Types of exceptions in the `type' field.  */
+#define	DOMAIN		1
+#define	SING		2
+#define	OVERFLOW	3
+#define	UNDERFLOW	4
+#define	TLOSS		5
+#define	PLOSS		6
 
-/* Return X with its signed changed to Y's.  */
-extern double __copysign __P ((double __x, double __y))
-     __attribute__ ((__const__));
+/* SVID mode specifies returning this large value instead of infinity.  */
+#define HUGE		FLT_MAX
+#include <float.h>		/* Defines FLT_MAX.  */
 
-/* Return X times (2 to the Nth power).  */
-extern double __scalb __P ((double __x, int __n))
-     __attribute__ ((__const__));
-
-#ifdef	__OPTIMIZE__
-#define	__scalb(x, n)	ldexp ((x), (n))
 #endif
 
-/* Return the remainder of X/Y.  */
-extern double __drem __P ((double __x, double __y))
-     __attribute__ ((__const__));
-
-/* Return the base 2 signed integral exponent of X.  */
-extern double __logb __P ((double __x)) __attribute__ ((__const__));
-
-#ifdef	__USE_MISC
-
-/* Return the integer nearest X in the direction of the
-   prevailing rounding mode.  */
-extern double __rint __P ((double __x)) __attribute__ ((__const__));
-extern double rint __P ((double __x)) __attribute__ ((__const__));
-
-/* Return `sqrt(X*X + Y*Y)'.  */
-extern double hypot __P ((double __x, double __y)) __attribute__ ((__const__));
-
-struct __cabs_complex
-{
-  double __x, __y;
-};
-
-/* Return `sqrt(X*X + Y*Y)'.  */
-extern double cabs __P ((struct __cabs_complex)) __attribute__ ((__const__));
-
-extern int isinf __P ((double __value)) __attribute__ ((__const__));
-extern int isnan __P ((double __value)) __attribute__ ((__const__));
-extern int finite __P ((double __value)) __attribute__ ((__const__));
-extern double infnan __P ((int __error)) __attribute__ ((__const__));
-extern double copysign __P ((double __x, double __y))
-     __attribute__ ((__const__));
-extern double scalb __P ((double __x, int __n)) __attribute__ ((__const__));
-extern double drem __P ((double __x, double __y)) __attribute__ ((__const__));
-extern double logb __P ((double __x)) __attribute__ ((__const__));
-
-#ifdef	__OPTIMIZE__
-#define	isinf(value)	__isinf(value)
-#define	isnan(value)	__isnan(value)
-#define	infnan(error)	__infnan(error)
-#define	finite(value)	__finite(value)
-#define	copysign(x, y)	__copysign((x), (y))
-#define	scalb(x, n)	__scalb((x), (n))
-#define	drem(x, y)	__drem((x), (y))
-#define	logb(x)		__logb(x)
-#endif /* Optimizing.  */
-
-#endif /* Use misc.  */
-
-
-#if 0
-/* The "Future Library Directions" section of the
-   ANSI Standard reserves these as `float' and
-   `long double' versions of the above functions.  */
-
-extern float acosf __P ((float __x)) __attribute__ ((__const__));
-extern float asinf __P ((float __x)) __attribute__ ((__const__));
-extern float atanf __P ((float __x)) __attribute__ ((__const__));
-extern float atan2f __P ((float __y, float __x)) __attribute__ ((__const__));
-extern float cosf __P ((float __x)) __attribute__ ((__const__));
-extern float sinf __P ((float __x)) __attribute__ ((__const__));
-extern float tanf __P ((float __x)) __attribute__ ((__const__));
-extern float coshf __P ((float __x)) __attribute__ ((__const__));
-extern float sinhf __P ((float __x)) __attribute__ ((__const__));
-extern float tanhf __P ((float __x)) __attribute__ ((__const__));
-extern float expf __P ((float __x)) __attribute__ ((__const__));
-extern float frexpf __P ((float __value, int *__exp));
-extern float ldexpf __P ((float __x, int __exp)) __attribute__ ((__const__));
-extern float logf __P ((float __x)) __attribute__ ((__const__));
-extern float log10f __P ((float __x)) __attribute__ ((__const__));
-extern float modff __P ((float __value, float *__iptr));
-extern float powf __P ((float __x, float __y)) __attribute__ ((__const__));
-extern float sqrtf __P ((float __x)) __attribute__ ((__const__));
-extern float ceilf __P ((float __x)) __attribute__ ((__const__));
-extern float fabsf __P ((float __x)) __attribute__ ((__const__));
-extern float floorf __P ((float __x)) __attribute__ ((__const__));
-extern float fmodf __P ((float __x, float __y)) __attribute__ ((__const__));
-
-extern __long_double_t acosl __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t asinl __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t atanl __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t atan2l __P ((__long_double_t __y, __long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t cosl __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t sinl __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t tanl __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t coshl __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t sinhl __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t tanhl __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t expl __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t frexpl __P ((__long_double_t __value, int *__exp));
-extern __long_double_t ldexpl __P ((__long_double_t __x, int __exp))
-     __attribute__ ((__const__));
-extern __long_double_t logl __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t log10l __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t modfl __P ((__long_double_t __value,
-				   __long_double_t *__ip));
-extern __long_double_t powl __P ((__long_double_t __x, __long_double_t __y))
-     __attribute__ ((__const__));
-extern __long_double_t sqrtl __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t ceill __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t fabsl __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t floorl __P ((__long_double_t __x))
-     __attribute__ ((__const__));
-extern __long_double_t fmodl __P ((__long_double_t __x, __long_double_t __y))
-     __attribute__ ((__const__));
-#endif /* 0 */
 
 /* Get machine-dependent inline versions (if there are any).  */
 #include <__math.h>
 
+
 __END_DECLS
 
 
diff --git a/math/math_private.h b/math/math_private.h
new file mode 100644
index 0000000000..95f65a293b
--- /dev/null
+++ b/math/math_private.h
@@ -0,0 +1,222 @@
+/*
+ * ====================================================
+ * 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.
+ * ====================================================
+ */
+
+/*
+ * from: @(#)fdlibm.h 5.1 93/09/24
+ * $Id$
+ */
+
+#ifndef _MATH_PRIVATE_H_
+#define _MATH_PRIVATE_H_
+
+#include <machine/endian.h>
+#include <sys/types.h>
+
+/* The original fdlibm code used statements like:
+	n0 = ((*(int*)&one)>>29)^1;		* index of high word *
+	ix0 = *(n0+(int*)&x);			* high word of x *
+	ix1 = *((1-n0)+(int*)&x);		* low word of x *
+   to dig two 32 bit words out of the 64 bit IEEE floating point
+   value.  That is non-ANSI, and, moreover, the gcc instruction
+   scheduler gets it wrong.  We instead use the following macros.
+   Unlike the original code, we determine the endianness at compile
+   time, not at run time; I don't see much benefit to selecting
+   endianness at run time.  */
+
+/* A union which permits us to convert between a double and two 32 bit
+   ints.  */
+
+#if BYTE_ORDER == BIG_ENDIAN
+
+typedef union
+{
+  double value;
+  struct
+  {
+    u_int32_t msw;
+    u_int32_t lsw;
+  } parts;
+} ieee_double_shape_type;
+
+#endif
+
+#if BYTE_ORDER == LITTLE_ENDIAN
+
+typedef union
+{
+  double value;
+  struct
+  {
+    u_int32_t lsw;
+    u_int32_t msw;
+  } parts;
+} ieee_double_shape_type;
+
+#endif
+
+/* Get two 32 bit ints from a double.  */
+
+#define EXTRACT_WORDS(ix0,ix1,d)				\
+do {								\
+  ieee_double_shape_type ew_u;					\
+  ew_u.value = (d);						\
+  (ix0) = ew_u.parts.msw;					\
+  (ix1) = ew_u.parts.lsw;					\
+} while (0)
+
+/* Get the more significant 32 bit int from a double.  */
+
+#define GET_HIGH_WORD(i,d)					\
+do {								\
+  ieee_double_shape_type gh_u;					\
+  gh_u.value = (d);						\
+  (i) = gh_u.parts.msw;						\
+} while (0)
+
+/* Get the less significant 32 bit int from a double.  */
+
+#define GET_LOW_WORD(i,d)					\
+do {								\
+  ieee_double_shape_type gl_u;					\
+  gl_u.value = (d);						\
+  (i) = gl_u.parts.lsw;						\
+} while (0)
+
+/* Set a double from two 32 bit ints.  */
+
+#define INSERT_WORDS(d,ix0,ix1)					\
+do {								\
+  ieee_double_shape_type iw_u;					\
+  iw_u.parts.msw = (ix0);					\
+  iw_u.parts.lsw = (ix1);					\
+  (d) = iw_u.value;						\
+} while (0)
+
+/* Set the more significant 32 bits of a double from an int.  */
+
+#define SET_HIGH_WORD(d,v)					\
+do {								\
+  ieee_double_shape_type sh_u;					\
+  sh_u.value = (d);						\
+  sh_u.parts.msw = (v);						\
+  (d) = sh_u.value;						\
+} while (0)
+
+/* Set the less significant 32 bits of a double from an int.  */
+
+#define SET_LOW_WORD(d,v)					\
+do {								\
+  ieee_double_shape_type sl_u;					\
+  sl_u.value = (d);						\
+  sl_u.parts.lsw = (v);						\
+  (d) = sl_u.value;						\
+} while (0)
+
+/* A union which permits us to convert between a float and a 32 bit
+   int.  */
+
+typedef union
+{
+  float value;
+  u_int32_t word;
+} ieee_float_shape_type;
+
+/* Get a 32 bit int from a float.  */
+
+#define GET_FLOAT_WORD(i,d)					\
+do {								\
+  ieee_float_shape_type gf_u;					\
+  gf_u.value = (d);						\
+  (i) = gf_u.word;						\
+} while (0)
+
+/* Set a float from a 32 bit int.  */
+
+#define SET_FLOAT_WORD(d,i)					\
+do {								\
+  ieee_float_shape_type sf_u;					\
+  sf_u.word = (i);						\
+  (d) = sf_u.value;						\
+} while (0)
+
+/* ieee style elementary functions */
+extern double __ieee754_sqrt __P((double));
+extern double __ieee754_acos __P((double));
+extern double __ieee754_acosh __P((double));
+extern double __ieee754_log __P((double));
+extern double __ieee754_atanh __P((double));
+extern double __ieee754_asin __P((double));
+extern double __ieee754_atan2 __P((double,double));
+extern double __ieee754_exp __P((double));
+extern double __ieee754_cosh __P((double));
+extern double __ieee754_fmod __P((double,double));
+extern double __ieee754_pow __P((double,double));
+extern double __ieee754_lgamma_r __P((double,int *));
+extern double __ieee754_gamma_r __P((double,int *));
+extern double __ieee754_lgamma __P((double));
+extern double __ieee754_gamma __P((double));
+extern double __ieee754_log10 __P((double));
+extern double __ieee754_sinh __P((double));
+extern double __ieee754_hypot __P((double,double));
+extern double __ieee754_j0 __P((double));
+extern double __ieee754_j1 __P((double));
+extern double __ieee754_y0 __P((double));
+extern double __ieee754_y1 __P((double));
+extern double __ieee754_jn __P((int,double));
+extern double __ieee754_yn __P((int,double));
+extern double __ieee754_remainder __P((double,double));
+extern int    __ieee754_rem_pio2 __P((double,double*));
+extern double __ieee754_scalb __P((double,double));
+
+/* fdlibm kernel function */
+extern double __kernel_standard __P((double,double,int));
+extern double __kernel_sin __P((double,double,int));
+extern double __kernel_cos __P((double,double));
+extern double __kernel_tan __P((double,double,int));
+extern int    __kernel_rem_pio2 __P((double*,double*,int,int,int,const int*));
+
+
+/* ieee style elementary float functions */
+extern float __ieee754_sqrtf __P((float));
+extern float __ieee754_acosf __P((float));
+extern float __ieee754_acoshf __P((float));
+extern float __ieee754_logf __P((float));
+extern float __ieee754_atanhf __P((float));
+extern float __ieee754_asinf __P((float));
+extern float __ieee754_atan2f __P((float,float));
+extern float __ieee754_expf __P((float));
+extern float __ieee754_coshf __P((float));
+extern float __ieee754_fmodf __P((float,float));
+extern float __ieee754_powf __P((float,float));
+extern float __ieee754_lgammaf_r __P((float,int *));
+extern float __ieee754_gammaf_r __P((float,int *));
+extern float __ieee754_lgammaf __P((float));
+extern float __ieee754_gammaf __P((float));
+extern float __ieee754_log10f __P((float));
+extern float __ieee754_sinhf __P((float));
+extern float __ieee754_hypotf __P((float,float));
+extern float __ieee754_j0f __P((float));
+extern float __ieee754_j1f __P((float));
+extern float __ieee754_y0f __P((float));
+extern float __ieee754_y1f __P((float));
+extern float __ieee754_jnf __P((int,float));
+extern float __ieee754_ynf __P((int,float));
+extern float __ieee754_remainderf __P((float,float));
+extern int   __ieee754_rem_pio2f __P((float,float*));
+extern float __ieee754_scalbf __P((float,float));
+
+/* float versions of fdlibm kernel functions */
+extern float __kernel_sinf __P((float,float,int));
+extern float __kernel_cosf __P((float,float));
+extern float __kernel_tanf __P((float,float,int));
+extern int   __kernel_rem_pio2f __P((float*,float*,int,int,int,const int*));
+
+#endif /* _MATH_PRIVATE_H_ */
diff --git a/math/mathcalls.h b/math/mathcalls.h
new file mode 100644
index 0000000000..656352ec2c
--- /dev/null
+++ b/math/mathcalls.h
@@ -0,0 +1,232 @@
+/* Prototype declarations for math functions; helper file for <math.h>.
+Copyright (C) 1996 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
+modify it under the terms of the GNU Library General Public License as
+published by the Free Software Foundation; either version 2 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
+Library General Public License for more details.
+
+You should have received a copy of the GNU Library General Public
+License along with the GNU C Library; see the file COPYING.LIB.  If
+not, write to the Free Software Foundation, Inc., 675 Mass Ave,
+Cambridge, MA 02139, USA.  */
+
+/* NOTE: Because of the special way this file is used by <math.h>, this
+   file must NOT be protected from multiple inclusion as header files
+   usually are.
+
+   This file provides prototype declarations for the math functions.
+   Most functions are declared using the macro:
+
+   __MATHCALL (NAME,[_r], (ARGS...));
+
+   This means there is a function `NAME' returning `double' and a function
+   `NAMEf' returning `float'.  Each place `_Mdouble_' appears in the
+   prototype, that is actually `double' in the prototype for `NAME' and
+   `float' in the prototype for `NAMEf'.  Reentrant variant functions are
+   called `NAME_r' and `NAMEf_r'.
+
+   Functions returning other types like `int' are declared using the macro:
+
+   __MATHDECL (TYPE, NAME,[_r], (ARGS...));
+
+   This is just like __MATHCALL but for a function returning `TYPE'
+   instead of `_Mdouble_'.  In all of these cases, there is still
+   both a `NAME' and a `NAMEf' that takes `float' arguments.  */
+
+#ifndef _MATH_H
+ #error "Never include mathcalls.h directly; include <math.h> instead."
+#endif
+
+
+/* Trigonometric functions.  */
+
+/* Arc cosine of X.  */
+__MATHCALL (acos,, (_Mdouble_ __x));
+/* Arc sine of X.  */
+__MATHCALL (asin,, (_Mdouble_ __x));
+/* Arc tangent of X.  */
+__MATHCALL (atan,, (_Mdouble_ __x));
+/* Arc tangent of Y/X.  */
+__MATHCALL (atan2,, (_Mdouble_ __y, _Mdouble_ __x));
+
+/* Cosine of X.  */
+__MATHCALL (cos,, (_Mdouble_ __x));
+/* Sine of X.  */
+__MATHCALL (sin,, (_Mdouble_ __x));
+/* Tangent of X.  */
+__MATHCALL (tan,, (_Mdouble_ __x));
+
+
+/* Hyperbolic functions.  */
+
+/* Hyperbolic cosine of X.  */
+__MATHCALL (cosh,, (_Mdouble_ __x));
+/* Hyperbolic sine of X.  */
+__MATHCALL (sinh,, (_Mdouble_ __x));
+/* Hyperbolic tangent of X.  */
+__MATHCALL (tanh,, (_Mdouble_ __x));
+
+#ifdef	__USE_MISC
+/* Hyperbolic arc cosine of X.  */
+__MATHCALL (acosh,, (_Mdouble_ __x));
+/* Hyperbolic arc sine of X.  */
+__MATHCALL (asinh,, (_Mdouble_ __x));
+/* Hyperbolic arc tangent of X.  */
+__MATHCALL (atanh,, (_Mdouble_ __x));
+#endif
+
+/* Exponential and logarithmic functions.  */
+
+/* Exponentional function of X.  */
+__MATHCALL (exp,, (_Mdouble_ __x));
+
+/* Break VALUE into a normalized fraction and an integral power of 2.  */
+__MATHCALL (frexp,, (_Mdouble_ __value, int *__exp));
+
+/* X times (two to the EXP power).  */
+__MATHCALL (ldexp,, (_Mdouble_ __x, int __exp));
+
+/* Natural logarithm of X.  */
+__MATHCALL (log,, (_Mdouble_ __x));
+
+/* Base-ten logarithm of X.  */
+__MATHCALL (log10,, (_Mdouble_ __x));
+
+#ifdef	__USE_MISC
+/* Return exp(X) - 1.  */
+__MATHCALL (expm1,, (_Mdouble_ __x));
+
+/* Return log(1 + X).  */
+__MATHCALL (log1p,, (_Mdouble_ __x));
+#endif
+
+/* Break VALUE into integral and fractional parts.  */
+__MATHCALL (modf,, (_Mdouble_ __value, _Mdouble_ *__iptr));
+
+
+/* Power functions.  */
+
+/* Return X to the Y power.  */
+__MATHCALL (pow,, (_Mdouble_ __x, _Mdouble_ __y));
+
+/* Return the square root of X.  */
+__MATHCALL (sqrt,, (_Mdouble_ __x));
+
+#ifdef	__USE_MISC
+/* Return the cube root of X.  */
+__MATHCALL (cbrt,, (_Mdouble_ __x));
+#endif
+
+
+/* Nearest integer, absolute value, and remainder functions.  */
+
+/* Smallest integral value not less than X.  */
+__MATHCALL (ceil,, (_Mdouble_ __x));
+
+/* Absolute value of X.  */
+__MATHCALL (fabs,, (_Mdouble_ __x));
+
+/* Largest integer not greater than X.  */
+__MATHCALL (floor,, (_Mdouble_ __x));
+
+/* Floating-point modulo remainder of X/Y.  */
+__MATHCALL (fmod,, (_Mdouble_ __x, _Mdouble_ __y));
+
+
+#ifdef __USE_MISC
+
+/* Return 0 if VALUE is finite or NaN, +1 if it
+   is +Infinity, -1 if it is -Infinity.  */
+__MATHDECL (int, isinf,, (_Mdouble_ __value));
+
+/* Return nonzero if VALUE is not a number.  */
+__MATHDECL (int, isnan,, (_Mdouble_ __value));
+
+/* Return nonzero if VALUE is finite and not NaN.  */
+__MATHDECL (int, finite,, (_Mdouble_ __value));
+
+/* Deal with an infinite or NaN result.
+   If ERROR is ERANGE, result is +Inf;
+   if ERROR is - ERANGE, result is -Inf;
+   otherwise result is NaN.
+   This will set `errno' to either ERANGE or EDOM,
+   and may return an infinity or NaN, or may do something else.  */
+__MATHCALL (infnan,, (int __error));
+
+/* Return X with its signed changed to Y's.  */
+__MATHCALL (copysign,, (_Mdouble_ __x, _Mdouble_ __y));
+
+/* Return X times (2 to the Nth power).  */
+__MATHCALL (scalb,, (_Mdouble_ __x, _Mdouble_ __n));
+
+/* Return X times (2 to the Nth power).  */
+__MATHCALL (scalbn,, (_Mdouble_ __x, int __n));
+
+/* Return the remainder of X/Y.  */
+__MATHCALL (drem,, (_Mdouble_ __x, _Mdouble_ __y));
+
+/* Return the base 2 signed integral exponent of X.  */
+__MATHCALL (logb,, (_Mdouble_ __x));
+
+/* Return the integer nearest X in the direction of the
+   prevailing rounding mode.  */
+__MATHCALL (rint,, (_Mdouble_ __x));
+
+/* Return `sqrt(X*X + Y*Y)'.  */
+__MATHCALL (hypot,, (_Mdouble_ __x, _Mdouble_ __y));
+
+struct __MATH_PRECNAME(__cabs_complex,)
+{
+  _Mdouble_ x, y;
+};
+
+/* Return `sqrt(X*X + Y*Y)'.  */
+__MATHCALL (cabs,, (struct __MATH_PRECNAME(__cabs_complex,)));
+
+
+/* Return X + epsilon if X < Y, X - epsilon if X > Y.  */
+__MATHCALL (nextafter,, (_Mdouble_ __x, _Mdouble_ __y));
+
+/* Return the remainder of integer divison X / Y with infinite precision.  */
+__MATHCALL (remainder,, (_Mdouble_ __x, _Mdouble_ __y));
+
+/* Return the binary exponent of X, which must be nonzero.  */
+__MATHDECL (int, ilogb,, (_Mdouble_ __x));
+
+/* Return the fractional part of X after dividing out `ilogb (X)'.  */
+__MATHCALL (significand,, (_Mdouble_ __x));
+
+
+
+/* Error, gamma, and Bessel functions.  */
+__MATHCALL (erf,, (_Mdouble_));
+__MATHCALL (erfc,, (_Mdouble_));
+__MATHCALL (gamma,, (_Mdouble_));
+__MATHCALL (j0,, (_Mdouble_));
+__MATHCALL (j1,, (_Mdouble_));
+__MATHCALL (jn,, (int, _Mdouble_));
+__MATHCALL (lgamma,, (_Mdouble_));
+__MATHCALL (y0,, (_Mdouble_));
+__MATHCALL (y1,, (_Mdouble_));
+__MATHCALL (yn,, (int, _Mdouble_));
+
+/* This variable is used by `gamma' and `lgamma'.  */
+extern int signgam;
+
+#ifdef __USE_REENTRANT
+
+/* Reentrant versions of gamma and lgamma.  Those functions use the global
+   variable `signgam'.  The reentrant versions instead take a pointer and
+   store the value through it.  */
+__MATHCALL (gamma,_r, (_Mdouble_, int *));
+__MATHCALL (lgamma,_r, (_Mdouble_, int *));
+#endif
+
+#endif /* Use misc.  */