¸

1998-08-23  Ricahrd Henderson  <rth@cygnus.com>

	* sysdeps/alpha/elf/crtbegin.S: Fix .prologue; no pv used.
	* sysdeps/alpha/elf/crtend.S: Likewise.
	* sysdeps/alpha/elf/start.S: Likewise.

	* sysdeps/unix/sysv/linux/alpha/rt_sigaction.S: Fix .prologue;
	non-standard pv usage.

	* sysdeps/unix/sysv/linux/alpha/brk.S: Use jmp macro for relaxation.
	* sysdeps/unix/sysv/linux/alpha/getitimer.S: Likewise.
	* sysdeps/unix/sysv/linux/alpha/getrusage.S: Likewise.
	* sysdeps/unix/sysv/linux/alpha/gettimeofday.S: Likewise.
	* sysdeps/unix/sysv/linux/alpha/ieee_get_fp_control.S: Likewise.
	* sysdeps/unix/sysv/linux/alpha/ieee_set_fp_control.S: Likewise.
	* sysdeps/unix/sysv/linux/alpha/select.S: Likewise.
	* sysdeps/unix/sysv/linux/alpha/setitimer.S: Likewise.
	* sysdeps/unix/sysv/linux/alpha/settimeofday.S: Likewise.
	* sysdeps/unix/sysv/linux/alpha/utimes.S: Likewise.
	* sysdeps/unix/sysv/linux/alpha/wait4.S: Likewise.

	* sysdeps/alpha/fpu/e_sqrt.c: Use the asm version when the input is
	a finite non-denormal, deferring to the full IEEE version otherwise.

	* sysdeps/alpha/fpu/bits/mathinline.h (__floorf, __floor):
	Early out for -0.  Optimize for !_IEEE_FP_INEXACT.
	* sysdeps/alpha/fpu/s_floor.c: New.
	* sysdeps/alpha/fpu/s_floorf.c: New.
	* sysdeps/alpha/fpu/s_ceil.c: New.
	* sysdeps/alpha/fpu/s_ceilf.c: New.

9 files changed, 289 insertions, 183 deletions

diff --git a/sysdeps/alpha/elf/crtbegin.S b/sysdeps/alpha/elf/crtbegin.S
index e6147c2ec3..25ddaafd6e 100644
--- a/sysdeps/alpha/elf/crtbegin.S
+++ b/sysdeps/alpha/elf/crtbegin.S
@@ -74,7 +74,7 @@ __do_global_dtors_aux:
 	stq     $9,8($sp)
 	stq     $26,0($sp)
 	.mask   (1<<26)|(1<<9), -16
-	.prologue 1
+	.prologue 0
 
 	lda     $9,__DTOR_LIST__
 	br      1f
diff --git a/sysdeps/alpha/elf/crtend.S b/sysdeps/alpha/elf/crtend.S
index 7f51d81da1..198aba8a5b 100644
--- a/sysdeps/alpha/elf/crtend.S
+++ b/sysdeps/alpha/elf/crtend.S
@@ -74,7 +74,7 @@ __do_global_ctors_aux:
 	stq     $9,8($sp)
 	stq     $26,0($sp)
 	.mask   (1<<26)|(1<<9), -16
-	.prologue 1
+	.prologue 0
 
 	lda     $9,__CTOR_END__
 	br      1f
diff --git a/sysdeps/alpha/elf/start.S b/sysdeps/alpha/elf/start.S
index ecb8174152..b5e5df102a 100644
--- a/sysdeps/alpha/elf/start.S
+++ b/sysdeps/alpha/elf/start.S
@@ -30,7 +30,7 @@ _start:
 	br	gp, 1f
 1:	ldgp	gp, 0(gp)
 	subq	sp, 16, sp
-	.prologue 1
+	.prologue 0
 
   /* Load address of the user's main function.  */
 	lda	a0, main
diff --git a/sysdeps/alpha/fpu/bits/mathinline.h b/sysdeps/alpha/fpu/bits/mathinline.h
index 492d9f18a7..681ea70cf9 100644
--- a/sysdeps/alpha/fpu/bits/mathinline.h
+++ b/sysdeps/alpha/fpu/bits/mathinline.h
@@ -75,7 +75,7 @@ __inline_copysign(copysign, double)
 #undef __MATH_INLINE_copysign
 
 
-#if defined __GNUC__ && (__GNUC__ > 2 || __GNUC__ == 2 && __GNUC_MINOR__ >= 8)
+#if defined __GNUC__ && (__GNUC__ > 2 || (__GNUC__ == 2 && __GNUC_MINOR__ >= 8))
 __MATH_INLINE float __fabsf (float __x) { return __builtin_fabsf (__x); }
 __MATH_INLINE float fabsf (float __x) { return __builtin_fabsf (__x); }
 __MATH_INLINE double __fabs (double __x) { return __builtin_fabs (__x); }
@@ -106,7 +106,8 @@ __inline_fabs(fabs, double)
 __MATH_INLINE float
 __floorf (float __x)
 {
-  if (fabsf (__x) < 16777216.0f)  /* 1 << FLT_MANT_DIG */
+  /* Check not zero since floor(-0) == -0.  */
+  if (__x != 0 && fabsf (__x) < 16777216.0f)  /* 1 << FLT_MANT_DIG */
     {
       /* Note that Alpha S_Floating is stored in registers in a
 	 restricted T_Floating format, so we don't even need to
@@ -116,10 +117,13 @@ __floorf (float __x)
       float __tmp1, __tmp2;
 
       __asm ("cvtst/s %3,%2\n\t"
+#ifdef _IEEE_FP_INEXACT
 	     "cvttq/svim %2,%1\n\t"
-	     "cvtqt/suim %1,%0\n\t"
-	     "trapb"
-	     : "=&f"(__x), "=&f"(__tmp1), "=&f"(__tmp2)
+#else
+	     "cvttq/svm %2,%1\n\t"
+#endif
+	     "cvtqt/m %1,%0\n\t"
+	     : "=f"(__x), "=&f"(__tmp1), "=&f"(__tmp2)
 	     : "f"(__x));
     }
   return __x;
@@ -128,13 +132,17 @@ __floorf (float __x)
 __MATH_INLINE double
 __floor (double __x)
 {
-  if (fabs (__x) < 9007199254740992.0)  /* 1 << DBL_MANT_DIG */
+  if (__x != 0 && fabs (__x) < 9007199254740992.0)  /* 1 << DBL_MANT_DIG */
     {
       double __tmp1;
-      __asm ("cvttq/svim %2,%1\n\t"
-	     "cvtqt/suim %1,%0\n\t"
-	     "trapb"
-	     : "=&f"(__x), "=&f"(__tmp1)
+      __asm (
+#ifdef _IEEE_FP_INEXACT
+	     "cvttq/svim %2,%1\n\t"
+#else
+	     "cvttq/svm %2,%1\n\t"
+#endif
+	     "cvtqt/m %1,%0\n\t"
+	     : "=f"(__x), "=&f"(__tmp1)
 	     : "f"(__x));
     }
   return __x;
diff --git a/sysdeps/alpha/fpu/e_sqrt.c b/sysdeps/alpha/fpu/e_sqrt.c
index 58de39f392..7b4e596664 100644
--- a/sysdeps/alpha/fpu/e_sqrt.c
+++ b/sysdeps/alpha/fpu/e_sqrt.c
@@ -1,4 +1,4 @@
-/* Copyright (C) 1996, 1997 Free Software Foundation, Inc.
+/* Copyright (C) 1996, 1997, 1998 Free Software Foundation, Inc.
    Contributed by David Mosberger (davidm@cs.arizona.edu).
 
    This file is part of the GNU C Library.
@@ -18,16 +18,15 @@
    write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
    Boston, MA 02111-1307, USA.  */
 
-/*
- * We have three versions, depending on how exact we need the results.
- */
-
-#if defined(_IEEE_FP) && defined(_IEEE_FP_INEXACT)
 
-/* Most demanding: go to the original source.  */
-#include <libm-ieee754/e_sqrt.c>
+#if !defined(_IEEE_FP_INEXACT)
 
-#else
+/*
+ * This version is much faster than generic sqrt implementation, but
+ * it doesn't handle the inexact flag.  It doesn't handle exceptional
+ * values either, but will defer to the full ieee754_sqrt routine which
+ * can.
+ */
 
 /* Careful with rearranging this without consulting the assembly below.  */
 const static struct sqrt_data_struct {
@@ -54,112 +53,6 @@ const static struct sqrt_data_struct {
 	0x1527f,0x1334a,0x11051,0xe951, 0xbe01, 0x8e0d, 0x5924, 0x1edd }
 };
 
-#ifdef _IEEE_FP
-/*
- * This version is much faster than the standard one included above,
- * but it doesn't maintain the inexact flag.
- */
-
-#define lobits(x) (((unsigned int *)&x)[0])
-#define hibits(x) (((unsigned int *)&x)[1])
-
-static inline double initial_guess(double x, unsigned int k,
-	const struct sqrt_data_struct * const ptr)
-{
-	double ret = 0.0;
-
-	k = 0x5fe80000 - (k >> 1);
-	k = k - ptr->T2[63&(k>>14)];
-	hibits(ret) = k;
-	return ret;
-}
-
-/* up = nextafter(1,+Inf), dn = nextafter(1,-Inf) */
-
-#define __half			(ptr->half)
-#define __one_and_a_half	(ptr->one_and_a_half)
-#define __two_to_minus_30	(ptr->two_to_minus_30)
-#define __one			(ptr->one)
-#define __up			(ptr->up)
-#define __dn			(ptr->dn)
-#define __Nan			(ptr->nan)
-
-#define Double(x) (*(double *)&x)
-
-/* Multiply with chopping rounding.. */
-#define choppedmul(a,b,c) \
-  __asm__("multc %1,%2,%0":"=&f" (c):"f" (a), "f" (b))
-
-double
-__ieee754_sqrt(double x)
-{
-  const struct sqrt_data_struct * const ptr = &sqrt_data;
-  unsigned long k, bits;
-  double y, z, zp, zn;
-  double dn, up, low, high;
-  double half, one_and_a_half, one, two_to_minus_30;
-
-  *(double *)&bits = x;
-  k = bits;
-
-  /* Negative or NaN or Inf */
-  if ((k >> 52) >= 0x7ff)
-    goto special;
-  y = initial_guess(x, k >> 32, ptr);
-  half = Double(__half);
-  one_and_a_half = Double(__one_and_a_half);
-  y = y*(one_and_a_half - half*x*y*y);
-  dn = Double(__dn);
-  two_to_minus_30 = Double(__two_to_minus_30);
-  y = y*((one_and_a_half - two_to_minus_30) - half*x*y*y);
-  up = Double(__up);
-  z = x*y;
-  one = Double(__one);
-  z = z + half*z*(one-z*y);
-
-  choppedmul(z,dn,zp);
-  choppedmul(z,up,zn);
-
-  choppedmul(z,zp,low);
-  low = low - x;
-  choppedmul(z,zn,high);
-  high = high - x;
-
-  /* I can't get gcc to use fcmov's.. */
-  __asm__("fcmovge %2,%3,%0"
-	  :"=f" (z)
-	  :"0" (z), "f" (low), "f" (zp));
-  __asm__("fcmovlt %2,%3,%0"
-	  :"=f" (z)
-	  :"0" (z), "f" (high), "f" (zn));
-  return z;	/* Argh! gcc jumps to end here */
-
-special:
-  /* throw away sign bit */
-  k <<= 1;
-  /* -0 */
-  if (!k)
-    return x;
-  /* special? */
-  if ((k >> 53) == 0x7ff) {
-    /* NaN? */
-    if (k << 11)
-      return x;
-    /* sqrt(+Inf) = +Inf */
-    if (x > 0)
-      return x;
-  }
-
-  x = Double(__Nan);
-  return x;
-}
-
-#else
-/*
- * This version is much faster than generic sqrt implementation, but
- * it doesn't handle exceptional values or the inexact flag.
- */
-
 asm ("\
   /* Define offsets into the structure defined in C above.  */
 	$DN = 0*8
@@ -174,7 +67,7 @@ asm ("\
 	$Y = 8
 
 	.text
-	.align	3
+	.align	5
 	.globl	__ieee754_sqrt
 	.ent	__ieee754_sqrt
 __ieee754_sqrt:
@@ -187,72 +80,86 @@ __ieee754_sqrt:
 #endif
 "	.prologue 1
 
-	stt	$f16, $K($sp)
-	lda	$4, sqrt_data			# load base address into t3
-	fblt	$f16, $negative
-
-  /* Compute initial guess.  */
+	.align 4
+	stt	$f16, $K($sp)		# e0    :
+	mult	$f31, $f31, $f31	# .. fm :
+	lda	$4, sqrt_data		# e0    :
+	fblt	$f16, $fixup		# .. fa :
 
-	.align 3
-
-	ldah	$2, 0x5fe8			# e0    :
-	ldq	$3, $K($sp)			# .. e1 :
-	ldt	$f12, $HALF($4)			# e0    :
+	ldah	$2, 0x5fe8		# e0    :
+	ldq	$3, $K($sp)		# .. e1 :
+	ldt	$f12, $HALF($4)		# e0    :
 	ldt	$f18, $ALMOST_THREE_HALF($4)	# .. e1 :
-	srl	$3, 33, $1			# e0    :
-	mult	$f16, $f12, $f11		# .. fm : $f11 = x * 0.5
-	subl	$2, $1, $2			# e0    :
-	addt	$f12, $f12, $f17		# .. fa : $f17 = 1.0
-	srl	$2, 12, $1			# e0    :
-	and	$1, 0xfc, $1			# .. e1 :
-	addq	$1, $4, $1			# e0    :
-	ldl	$1, $T2($1)			# .. e1 :
-	addt	$f12, $f17, $f15		# fa    : $f15 = 1.5
-	subl	$2, $1, $2			# .. e1 :
-	sll	$2, 32, $2			# e0    :
-	ldt	$f14, $DN($4)			# .. e1 :
-	stq	$2, $Y($sp)			# e0    :
-	nop					# .. e1 : avoid pipe flash
-	nop					# e0    :
-	ldt	$f13, $Y($sp)			# .. e1 :
 
-	mult/su	$f11, $f13, $f10	# fm    : $f10 = (x * 0.5) * y
-	mult	$f10, $f13, $f10	# fm    : $f10 = ((x * 0.5) * y) * y
-	subt	$f15, $f10, $f1		# fa    : $f1 = (1.5 - 0.5*x*y*y)
-	mult	$f13, $f1, $f13         # fm    : yp = y*(1.5 - 0.5*x*y*y)
- 	mult/su	$f11, $f13, $f1		# fm    : $f11 = x * 0.5 * yp
-	mult	$f1, $f13, $f11		# fm    : $f11 = (x * 0.5 * yp) * yp
-	subt	$f18, $f11, $f1		# fa    : $f1= (1.5-2^-30) - 0.5*x*yp*yp
-	mult	$f13, $f1, $f13		# fm    : ypp = $f13 = yp*$f1
-	subt	$f15, $f12, $f1		# fa    : $f1 = (1.5 - 0.5)
-	ldt	$f15, $UP($4)		# .. e1 :
-	mult/su	$f16, $f13, $f10	# fm    : z = $f10 = x * ypp
-	mult	$f10, $f13, $f11	# fm    : $f11 = z*ypp
+	sll	$3, 52, $5		# e0    :
+	lda	$6, 0x7fd		# .. e1 :
+	fnop				# .. fa :
+	fnop				# .. fm :
+
+	subq	$5, 1, $5		# e1    :
+	srl	$3, 33, $1		# .. e0 :
+	cmpule	$5, $6, $5		# e0    :
+	beq	$5, $fixup		# .. e1 :
+
+	mult	$f16, $f12, $f11	# fm    : $f11 = x * 0.5
+	subl	$2, $1, $2		# .. e0 :
+	addt	$f12, $f12, $f17	# .. fa : $f17 = 1.0
+	srl	$2, 12, $1		# e0    :
+
+	and	$1, 0xfc, $1		# e0    :
+	addq	$1, $4, $1		# e1    :
+	ldl	$1, $T2($1)		# e0    :
+	addt	$f12, $f17, $f15	# .. fa : $f15 = 1.5
+
+	subl	$2, $1, $2		# e0    :
+	ldt	$f14, $DN($4)		# .. e1 :
+	sll	$2, 32, $2		# e0    :
+	stq	$2, $Y($sp)		# e0    :
+
+	ldt	$f13, $Y($sp)		# e0    :
+	mult/su	$f11, $f13, $f10	# fm   2: $f10 = (x * 0.5) * y
+	mult	$f10, $f13, $f10	# fm   4: $f10 = ((x * 0.5) * y) * y
+	subt	$f15, $f10, $f1		# fa   4: $f1 = (1.5 - 0.5*x*y*y)
+
+	mult	$f13, $f1, $f13         # fm   4: yp = y*(1.5 - 0.5*x*y*y)
+ 	mult/su	$f11, $f13, $f1		# fm   4: $f11 = x * 0.5 * yp
+	mult	$f1, $f13, $f11		# fm   4: $f11 = (x * 0.5 * yp) * yp
+	subt	$f18, $f11, $f1		# fa   4: $f1= (1.5-2^-30) - 0.5*x*yp*yp
+
+	mult	$f13, $f1, $f13		# fm   4: ypp = $f13 = yp*$f1
+	subt	$f15, $f12, $f1		# .. fa : $f1 = (1.5 - 0.5)
+	ldt	$f15, $UP($4)		# .. e0 :
+	mult/su	$f16, $f13, $f10	# fm   4: z = $f10 = x * ypp
+
+	mult	$f10, $f13, $f11	# fm   4: $f11 = z*ypp
 	mult	$f10, $f12, $f12	# fm    : $f12 = z*0.5
-	subt	$f1, $f11, $f1		# .. fa : $f1 = 1 - z*ypp
-	mult	$f12, $f1, $f12		# fm    : $f12 = z*0.5*(1 - z*ypp)
-	addt	$f10, $f12, $f0		# fa    : zp=res=$f0= z + z*0.5*(1 - z*ypp)
+	subt	$f1, $f11, $f1		# fa   4: $f1 = 1 - z*ypp
+	mult	$f12, $f1, $f12		# fm   4: $f12 = z*0.5*(1 - z*ypp)
 
-	mult/c	$f0, $f14, $f12		# fm    : zmi = zp * DN
+	addt	$f10, $f12, $f0		# fa   4: zp=res= z + z*0.5*(1 - z*ypp)
+	mult/c	$f0, $f14, $f12		# fm   4: zmi = zp * DN
 	mult/c	$f0, $f15, $f11		# fm    : zpl = zp * UP
 	mult/c	$f0, $f12, $f1		# fm    : $f1 = zp * zmi
-	mult/c	$f0, $f11, $f15		# fm    : $f15 = zp * zpl
 
-	subt/su	$f1, $f16, $f13		# fa    : y1 = zp*zmi - x
-	subt/su	$f15, $f16, $f14	# fa    : y2 = zp*zpl - x
-
-	fcmovge	$f13, $f12, $f0		# res = (y1 >= 0) ? zmi : res
-	fcmovlt	$f14, $f11, $f0		# res = (y2 <  0) ? zpl : res
+	mult/c	$f0, $f11, $f15		# fm    : $f15 = zp * zpl
+	subt/su	$f1, $f16, $f13		# .. fa : y1 = zp*zmi - x
+	subt/su	$f15, $f16, $f14	# fa   4: y2 = zp*zpl - x
+	fcmovge	$f13, $f12, $f0		# fa   3: res = (y1 >= 0) ? zmi : res
 
-	addq	$sp, 16, $sp		# e0    :
+	fcmovlt	$f14, $f11, $f0		# fa   4: res = (y2 <  0) ? zpl : res
+	addq	$sp, 16, $sp		# .. e0 :
 	ret				# .. e1 :
 
-$negative:
-	ldt	$f0, $NAN($4)
+	.align 4
+$fixup:
 	addq	$sp, 16, $sp
-	ret
+	br	"ASM_ALPHA_NG_SYMBOL_PREFIX"__full_ieee754_sqrt..ng
 
 	.end	__ieee754_sqrt");
 
-#endif /* _IEEE_FP */
-#endif /* _IEEE_FP && _IEEE_FP_INEXACT */
+static double __full_ieee754_sqrt(double) __attribute__((unused));
+#define __ieee754_sqrt __full_ieee754_sqrt
+
+#endif /* _IEEE_FP_INEXACT */
+
+#include <sysdeps/libm-ieee754/e_sqrt.c>
diff --git a/sysdeps/alpha/fpu/s_ceil.c b/sysdeps/alpha/fpu/s_ceil.c
new file mode 100644
index 0000000000..23491db5ac
--- /dev/null
+++ b/sysdeps/alpha/fpu/s_ceil.c
@@ -0,0 +1,58 @@
+/* Copyright (C) 1998 Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+   Contributed by Richard Henderson.
+
+   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., 59 Temple Place - Suite 330,
+   Boston, MA 02111-1307, USA.  */
+
+#include <math.h>
+
+double
+__ceil (double x)
+{
+  if (x != 0 && fabs (x) < 9007199254740992.0)  /* 1 << DBL_MANT_DIG */
+    {
+      double tmp1;
+      unsigned long fpcr0, fpcr1;
+      unsigned long pinf = 3UL << 58;
+
+      /* Set round to +inf.  */
+      __asm __volatile("excb; mf_fpcr %0" : "=f"(fpcr0));
+      __asm __volatile("mt_fpcr %0; excb" : : "f"(fpcr0 | pinf));
+
+      /* Calculate!  */
+#ifdef _IEEE_FP_INEXACT
+      __asm("cvttq/svid %2,%1\n\tcvtqt/suid %1,%0"
+	    : "=f"(x), "=&f"(tmp1)
+	    : "f"(x));
+#else
+      __asm("cvttq/svd %2,%1\n\tcvtqt/d %1,%0"
+	    : "=f"(x), "=&f"(tmp1)
+	    : "f"(x));
+#endif
+
+      /* Reset rounding mode, while retaining new exception bits.  */
+      __asm __volatile("excb; mf_fpcr %0" : "=f"(fpcr1));
+      fpcr0 = (fpcr0 & pinf) | (fpcr1 & ~pinf);
+      __asm __volatile("mt_fpcr %0; excb" : : "f"(fpcr0));
+    }
+  return x;
+}
+
+weak_alias (__ceil, ceil)
+#ifdef NO_LONG_DOUBLE
+strong_alias (__ceil, __ceill)
+weak_alias (__ceil, ceill)
+#endif
diff --git a/sysdeps/alpha/fpu/s_ceilf.c b/sysdeps/alpha/fpu/s_ceilf.c
new file mode 100644
index 0000000000..3defaeb01e
--- /dev/null
+++ b/sysdeps/alpha/fpu/s_ceilf.c
@@ -0,0 +1,59 @@
+/* Copyright (C) 1998 Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+   Contributed by Richard Henderson.
+
+   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., 59 Temple Place - Suite 330,
+   Boston, MA 02111-1307, USA.  */
+
+#include <math.h>
+
+float
+__ceilf (float x)
+{
+  if (x != 0 && fabsf (x) < 16777216.0f)  /* 1 << FLT_MANT_DIG */
+    {
+      float tmp1, tmp2;
+      unsigned long fpcr0, fpcr1;
+      unsigned long pinf = 3UL << 58;
+
+      /* Set round to +inf.  */
+      __asm __volatile("excb; mf_fpcr %0" : "=f"(fpcr0));
+      __asm __volatile("mt_fpcr %0; excb" : : "f"(fpcr0 | pinf));
+
+      /* Calculate!  
+         Note that Alpha S_Floating is stored in registers in a
+         restricted T_Floating format, so we don't even need to
+         convert back to S_Floating in the end.  The initial
+         conversion to T_Floating is needed to handle denormals.  */
+
+#ifdef _IEEE_FP_INEXACT
+      __asm("cvtst/s %3,%2\n\tcvttq/svid %2,%1\n\tcvtqt/suid %1,%0"
+	    : "=f"(x), "=&f"(tmp1), "=&f"(tmp2)
+	    : "f"(x));
+#else
+      __asm("cvtst/s %3,%2\n\tcvttq/svd %2,%1\n\tcvtqt/d %1,%0"
+	    : "=f"(x), "=&f"(tmp1), "=&f"(tmp2)
+	    : "f"(x));
+#endif
+
+      /* Reset rounding mode, while retaining new exception bits.  */
+      __asm __volatile("excb; mf_fpcr %0" : "=f"(fpcr1));
+      fpcr0 = (fpcr0 & pinf) | (fpcr1 & ~pinf);
+      __asm __volatile("mt_fpcr %0; excb" : : "f"(fpcr0));
+    }
+  return x;
+}
+
+weak_alias (__ceilf, ceilf)
diff --git a/sysdeps/alpha/fpu/s_floor.c b/sysdeps/alpha/fpu/s_floor.c
new file mode 100644
index 0000000000..7b647923e2
--- /dev/null
+++ b/sysdeps/alpha/fpu/s_floor.c
@@ -0,0 +1,39 @@
+/* Copyright (C) 1998 Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+   Contributed by Richard Henderson.
+
+   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., 59 Temple Place - Suite 330,
+   Boston, MA 02111-1307, USA.  */
+
+#ifndef __USE_EXTERN_INLINES
+#define __USE_EXTERN_INLINES
+#endif
+#define __floor __i_floor
+
+#include <math.h>
+
+#undef __floor
+
+double
+__floor (double x)
+{
+  return __i_floor(x);
+}
+
+weak_alias (__floor, floor)
+#ifdef NO_LONG_DOUBLE
+strong_alias (__floor, __floorl)
+weak_alias (__floor, floorl)
+#endif
diff --git a/sysdeps/alpha/fpu/s_floorf.c b/sysdeps/alpha/fpu/s_floorf.c
new file mode 100644
index 0000000000..d25643d19e
--- /dev/null
+++ b/sysdeps/alpha/fpu/s_floorf.c
@@ -0,0 +1,35 @@
+/* Copyright (C) 1998 Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+   Contributed by Richard Henderson.
+
+   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., 59 Temple Place - Suite 330,
+   Boston, MA 02111-1307, USA.  */
+
+#ifndef __USE_EXTERN_INLINES
+#define __USE_EXTERN_INLINES
+#endif
+#define __floorf __i_floorf
+
+#include <math.h>
+
+#undef __floorf
+
+float
+__floorf (float x)
+{
+  return __i_floorf(x);
+}
+
+weak_alias (__floorf, floorf)