about summary refs log tree commit diff
path: root/sysdeps/powerpc/powerpc64/fpu
diff options
context:
space:
mode:
authorJakub Jelinek <jakub@redhat.com>2007-07-12 18:26:36 +0000
committerJakub Jelinek <jakub@redhat.com>2007-07-12 18:26:36 +0000
commit0ecb606cb6cf65de1d9fc8a919bceb4be476c602 (patch)
tree2ea1f8305970753e4a657acb2ccc15ca3eec8e2c /sysdeps/powerpc/powerpc64/fpu
parent7d58530341304d403a6626d7f7a1913165fe2f32 (diff)
downloadglibc-0ecb606cb6cf65de1d9fc8a919bceb4be476c602.tar.gz
glibc-0ecb606cb6cf65de1d9fc8a919bceb4be476c602.tar.xz
glibc-0ecb606cb6cf65de1d9fc8a919bceb4be476c602.zip
2.5-18.1
Diffstat (limited to 'sysdeps/powerpc/powerpc64/fpu')
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/e_sqrt.c29
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/e_sqrtf.c29
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_ceil.S19
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_ceilf.S17
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_ceill.S133
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_copysign.S20
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_copysignl.S51
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_fabs.S5
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_fabsl.S36
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_fdim.c5
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_floor.S15
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_floorf.S13
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_floorl.S134
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_fmax.S5
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_fmin.S5
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_isnan.c7
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_llrint.S7
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_llround.S7
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_llroundf.S6
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_nearbyintl.S114
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_rint.S19
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_rintf.S17
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_round.S21
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_roundf.S23
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_roundl.S133
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_trunc.S21
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_truncf.S21
-rw-r--r--sysdeps/powerpc/powerpc64/fpu/s_truncl.S121
28 files changed, 940 insertions, 93 deletions
diff --git a/sysdeps/powerpc/powerpc64/fpu/e_sqrt.c b/sysdeps/powerpc/powerpc64/fpu/e_sqrt.c
new file mode 100644
index 0000000000..0a229cbe27
--- /dev/null
+++ b/sysdeps/powerpc/powerpc64/fpu/e_sqrt.c
@@ -0,0 +1,29 @@
+/* Double-precision floating point square root.
+   Copyright (C) 1997, 2002, 2003, 2004, 2006 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 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.  */
+
+#include <math.h>
+#include <math_private.h>
+
+double
+__ieee754_sqrt (double x)
+{
+  double z;
+  __asm __volatile ("fsqrt %0,%1" : "=f" (z) : "f" (x));
+  return z;
+}
diff --git a/sysdeps/powerpc/powerpc64/fpu/e_sqrtf.c b/sysdeps/powerpc/powerpc64/fpu/e_sqrtf.c
new file mode 100644
index 0000000000..0f17a64a8a
--- /dev/null
+++ b/sysdeps/powerpc/powerpc64/fpu/e_sqrtf.c
@@ -0,0 +1,29 @@
+/* Single-precision floating point square root.
+   Copyright (C) 1997, 2003, 2004, 2006 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 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.  */
+
+#include <math.h>
+#include <math_private.h>
+
+float
+__ieee754_sqrtf (float x)
+{
+  double z;
+  __asm ("fsqrts %0,%1" : "=f" (z) : "f" (x));
+  return z;
+}
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_ceil.S b/sysdeps/powerpc/powerpc64/fpu/s_ceil.S
index a1bfaa70c2..02b70940ee 100644
--- a/sysdeps/powerpc/powerpc64/fpu/s_ceil.S
+++ b/sysdeps/powerpc/powerpc64/fpu/s_ceil.S
@@ -1,5 +1,5 @@
 /* ceil function.  PowerPC64 version.
-   Copyright (C) 2004 Free Software Foundation, Inc.
+   Copyright (C) 2004, 2006 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
@@ -18,15 +18,14 @@
    02111-1307 USA.  */
 
 #include <sysdep.h>
+#include <math_ldbl_opt.h>
 
 	.section	".toc","aw"
 .LC0:	/* 2**52 */
 	.tc FD_43300000_0[TC],0x4330000000000000
-.LC1:	/* -0.0 */
-	.tc FD_80000000_0[TC],0x8000000000000000
 	.section	".text"
 
-ENTRY (__ceil)
+EALIGN (__ceil, 4, 0)
 	CALL_MCOUNT 0
 	mffs	fp11		/* Save current FPU rounding mode.  */
 	lfd	fp13,.LC0@toc(2)
@@ -39,17 +38,18 @@ ENTRY (__ceil)
 	ble-	cr6,.L4
 	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
 	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
-.L9:	
+	fabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = 0.0; */
 	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
 	blr
 .L4:
 	bge-	cr6,.L9		/* if (x < 0.0)  */
 	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
 	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
-	fcmpu	cr5,fp1,fp12	/* if (x > 0.0)  */
+	fnabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = -0.0; */
+.L9:
 	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
-	bnelr+	cr5
-	lfd	fp1,.LC1@toc(2)	/* x must be -0.0 for the 0.0 case.  */	
 	blr
 	END (__ceil)
 
@@ -59,3 +59,6 @@ weak_alias (__ceil, ceil)
 weak_alias (__ceil, ceill)
 strong_alias (__ceil, __ceill)
 #endif
+#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_0)
+compat_symbol (libm, __ceil, ceill, GLIBC_2_0)
+#endif
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_ceilf.S b/sysdeps/powerpc/powerpc64/fpu/s_ceilf.S
index 42eb274389..1ccd133b66 100644
--- a/sysdeps/powerpc/powerpc64/fpu/s_ceilf.S
+++ b/sysdeps/powerpc/powerpc64/fpu/s_ceilf.S
@@ -21,15 +21,13 @@
 
 	.section	".toc","aw"
 .LC0:	/* 2**23 */
-	.tc FD_41600000_0[TC],0x4160000000000000
-.LC1:	/* -0.0 */
-	.tc FD_80000000_0[TC],0x8000000000000000
+	.tc FD_4b000000_0[TC],0x4b00000000000000
 	.section	".text"
 
-ENTRY (__ceilf)
+EALIGN (__ceilf, 4, 0)
 	CALL_MCOUNT 0
 	mffs	fp11		/* Save current FPU rounding mode.  */
-	lfd	fp13,.LC0@toc(2)
+	lfs	fp13,.LC0@toc(2)
 	fabs	fp0,fp1
 	fsubs	fp12,fp13,fp13	/* generate 0.0  */
 	fcmpu	cr7,fp0,fp13	/* if (fabs(x) > TWO23)  */
@@ -39,17 +37,18 @@ ENTRY (__ceilf)
 	ble-	cr6,.L4
 	fadds	fp1,fp1,fp13	/* x+= TWO23;  */
 	fsubs	fp1,fp1,fp13	/* x-= TWO23;  */
-.L9:	
+	fabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = 0.0; */
 	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
 	blr
 .L4:
 	bge-	cr6,.L9		/* if (x < 0.0)  */
 	fsubs	fp1,fp1,fp13	/* x-= TWO23;  */
 	fadds	fp1,fp1,fp13	/* x+= TWO23;  */
-	fcmpu	cr5,fp1,fp12	/* if (x > 0.0)  */
+	fnabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = -0.0; */
+.L9:
 	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
-	bnelr+	cr5
-	lfd	fp1,.LC1@toc(2)	/* x must be -0.0 for the 0.0 case.  */	
 	blr
 	END (__ceilf)
 
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_ceill.S b/sysdeps/powerpc/powerpc64/fpu/s_ceill.S
new file mode 100644
index 0000000000..a8f8a0afc5
--- /dev/null
+++ b/sysdeps/powerpc/powerpc64/fpu/s_ceill.S
@@ -0,0 +1,133 @@
+/* s_ceill.S IBM extended format long double version.
+   Copyright (C) 2004, 2006 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 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.  */
+
+#include <sysdep.h>
+#include <math_ldbl_opt.h>
+
+	.section	".toc","aw"
+.LC0:	/* 2**52 */
+	.tc FD_43300000_0[TC],0x4330000000000000
+
+	.section	".text"
+
+/* long double [fp1,fp2] ceill (long double x [fp1,fp2])
+   IEEE 1003.1 ceil function.
+
+   PowerPC64 long double uses the IBM extended format which is
+   represented two 64-floating point double values. The values are
+   non-overlapping giving an effective precision of 106 bits. The first
+   double contains the high order bits of mantisa and is always ceiled
+   to represent a normal ceiling of long double to double. Since the
+   long double value is sum of the high and low values, the low double
+   normally has the opposite sign to compensate for the this ceiling.
+
+   For long double there are two cases:
+   1) |x| < 2**52, all the integer bits are in the high double.
+      ceil the high double and set the low double to -0.0.
+   2) |x| >= 2**52, ceiling involves both doubles.
+      See the comment before lable .L2 for details.
+   */
+
+ENTRY (__ceill)
+	mffs	fp11		/* Save current FPU rounding mode.  */
+	lfd	fp13,.LC0@toc(2)
+	fabs	fp0,fp1
+	fabs	fp9,fp2
+	fsub	fp12,fp13,fp13	/* generate 0.0  */
+	fcmpu	cr7,fp0,fp13	/* if (fabs(x) > TWO52)  */
+	fcmpu	cr6,fp1,fp12	/* if (x > 0.0)  */
+	bnl-	cr7,.L2
+	mtfsfi	7,2		/* Set rounding mode toward +inf.  */
+	fneg	fp2,fp12
+	ble-	cr6,.L1
+	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
+	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
+	fabs	fp1,fp1		/* if (x == 0.0)  */
+.L0:
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
+	blr			/* x = 0.0; */
+.L1:
+	bge-	cr6,.L0		/* if (x < 0.0)  */
+	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
+	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
+	fcmpu	cr5,fp1,fp12	/* if (x > 0.0)  */
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
+	fnabs	fp1,fp1		/* if (x == 0.0)  */
+	blr			/* x = -0.0; */
+	
+/* The high double is > TWO52 so we need to round the low double and
+   perhaps the high double.  In this case we have to round the low
+   double and handle any adjustment to the high double that may be
+   caused by rounding (up).  This is complicated by the fact that the
+   high double may already be rounded and the low double may have the
+   opposite sign to compensate.This gets a bit tricky so we use the
+   following algorithm:
+
+   tau = floor(x_high/TWO52);
+   x0 = x_high - tau;
+   x1 = x_low + tau;
+   r1 = rint(x1);
+   y_high = x0 + r1;
+   y_low = x0 - y_high + r1;
+   return y;  */
+.L2:
+	fcmpu	cr7,fp9,fp13	/* if (|x_low| > TWO52)  */
+	fcmpu	cr0,fp9,fp12	/* || (|x_low| == 0.0)  */
+	fcmpu	cr5,fp2,fp12	/* if (x_low > 0.0)  */
+	bgelr-	cr7		/*   return x;	*/
+	beqlr-  cr0
+	mtfsfi	7,2		/* Set rounding mode toward +inf.  */
+	fdiv	fp8,fp1,fp13	/* x_high/TWO52  */
+	
+	bng-	cr6,.L6		/* if (x > 0.0)  */
+	fctidz	fp0,fp8
+	fcfid	fp8,fp0		/* tau = floor(x_high/TWO52);  */
+	bng	cr5,.L4		/* if (x_low > 0.0)  */
+	fmr	fp3,fp1
+	fmr	fp4,fp2
+	b	.L5
+.L4:				/* if (x_low < 0.0)  */
+	fsub	fp3,fp1,fp8	/* x0 = x_high - tau;  */
+	fadd	fp4,fp2,fp8	/* x1 = x_low + tau;  */
+.L5:
+	fadd	fp5,fp4,fp13	/* r1 = r1 + TWO52;  */
+	fsub	fp5,fp5,fp13	/* r1 = r1 - TWO52;  */
+	b	.L9
+.L6:				/* if (x < 0.0)  */
+	fctidz	fp0,fp8
+	fcfid	fp8,fp0		/* tau = floor(x_high/TWO52);  */	
+	bnl	cr5,.L7		/* if (x_low < 0.0)  */
+	fmr	fp3,fp1
+	fmr	fp4,fp2
+	b	.L8
+.L7:				/* if (x_low > 0.0)  */
+	fsub	fp3,fp1,fp8	/* x0 = x_high - tau;  */
+	fadd	fp4,fp2,fp8	/* x1 = x_low + tau;  */
+.L8:
+	fsub	fp5,fp4,fp13	/* r1-= TWO52;  */
+	fadd	fp5,fp5,fp13	/* r1+= TWO52;  */
+.L9:
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
+	fadd	fp1,fp3,fp5	/* y_high = x0 + r1;  */
+	fsub	fp2,fp3,fp1	/* y_low = x0 - y_high + r1;  */
+	fadd	fp2,fp2,fp5
+	blr
+END (__ceill)
+
+long_double_symbol (libm, __ceill, ceill)
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_copysign.S b/sysdeps/powerpc/powerpc64/fpu/s_copysign.S
index a43ed12cf0..38171e31d7 100644
--- a/sysdeps/powerpc/powerpc64/fpu/s_copysign.S
+++ b/sysdeps/powerpc/powerpc64/fpu/s_copysign.S
@@ -1,5 +1,5 @@
 /* Copy a sign bit between floating-point values.  PowerPC64 version.
-   Copyright (C) 1997, 1999, 2000, 2002 Free Software Foundation, Inc.
+   Copyright (C) 1997, 1999, 2000, 2002, 2006 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
@@ -21,6 +21,7 @@
    when it's coded in C.  */
 
 #include <sysdep.h>
+#include <math_ldbl_opt.h>
 
 ENTRY(__copysign)
 	CALL_MCOUNT 0
@@ -28,7 +29,11 @@ ENTRY(__copysign)
    copysign(x,y) returns a value with the magnitude of x and
    with the sign bit of y.  */
 	stdu	r1,-48(r1)
+	cfi_adjust_cfa_offset (48)
 	stfd	fp2,24(r1)
+	nop
+	nop
+	nop
 	ld	r3,24(r1)
 	cmpdi   r3,0
 	addi    r1,r1,48
@@ -39,13 +44,20 @@ L(0):	fnabs   fp1,fp1
 	blr
 	END (__copysign)
 
-weak_alias(__copysign,copysign)
+weak_alias (__copysign,copysign)
 
 /* It turns out that it's safe to use this code even for single-precision.  */
-weak_alias(__copysign,copysignf)
+weak_alias (__copysign,copysignf)
 strong_alias(__copysign,__copysignf)
 
 #ifdef NO_LONG_DOUBLE
-weak_alias(__copysign,copysignl)
+weak_alias (__copysign,copysignl)
 strong_alias(__copysign,__copysignl)
 #endif
+#ifdef IS_IN_libm
+# if LONG_DOUBLE_COMPAT(libm, GLIBC_2_0)
+compat_symbol (libm, __copysign, copysignl, GLIBC_2_0)
+# endif
+#elif LONG_DOUBLE_COMPAT(libc, GLIBC_2_0)
+compat_symbol (libc, __copysign, copysignl, GLIBC_2_0)
+#endif
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_copysignl.S b/sysdeps/powerpc/powerpc64/fpu/s_copysignl.S
new file mode 100644
index 0000000000..b2c62eacba
--- /dev/null
+++ b/sysdeps/powerpc/powerpc64/fpu/s_copysignl.S
@@ -0,0 +1,51 @@
+/* Copy a sign bit between floating-point values.
+   IBM extended format long double version.
+   Copyright (C) 2004, 2006 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 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.  */
+
+#include <sysdep.h>
+#include <math_ldbl_opt.h>
+
+ENTRY(__copysignl)
+/* long double [f1,f2] copysign (long double [f1,f2] x, long double [f3,f4] y);
+   copysign(x,y) returns a value with the magnitude of x and
+   with the sign bit of y.  */
+	stfd	fp3,-16(r1)
+	ld	r3,-16(r1)
+	cmpdi   r3,0
+	blt     L(0)
+	fmr	fp0,fp1
+	fabs	fp1,fp1
+	fcmpu	cr1,fp0,fp1
+	beqlr	cr1
+	fneg	fp2,fp2
+	blr
+L(0):
+	fmr	fp0,fp1
+	fnabs   fp1,fp1
+	fcmpu	cr1,fp0,fp1
+	beqlr	cr1
+	fneg	fp2,fp2
+	blr
+END (__copysignl)
+
+#ifdef IS_IN_libm
+long_double_symbol (libm, __copysignl, copysignl)
+#else
+long_double_symbol (libc, __copysignl, copysignl)
+#endif
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_fabs.S b/sysdeps/powerpc/powerpc64/fpu/s_fabs.S
new file mode 100644
index 0000000000..53d21301ee
--- /dev/null
+++ b/sysdeps/powerpc/powerpc64/fpu/s_fabs.S
@@ -0,0 +1,5 @@
+#include <math_ldbl_opt.h>
+#include <sysdeps/powerpc/fpu/s_fabs.S>
+#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_0)
+compat_symbol (libm, __fabs, fabsl, GLIBC_2_0)
+#endif
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_fabsl.S b/sysdeps/powerpc/powerpc64/fpu/s_fabsl.S
new file mode 100644
index 0000000000..3655e5b2f3
--- /dev/null
+++ b/sysdeps/powerpc/powerpc64/fpu/s_fabsl.S
@@ -0,0 +1,36 @@
+/* Copy a sign bit between floating-point values.
+   IBM extended format long double version.
+   Copyright (C) 2004, 2006 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 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.  */
+
+#include <sysdep.h>
+#include <math_ldbl_opt.h>
+
+ENTRY(__fabsl)
+/* long double [f1,f2] fabs (long double [f1,f2] x);
+   fabs(x,y) returns a value with the magnitude of x and
+   with the sign bit of y.  */
+	fmr	fp0,fp1
+	fabs	fp1,fp1
+	fcmpu	cr1,fp0,fp1
+	beqlr	cr1
+	fneg	fp2,fp2
+	blr
+END (__fabsl)
+
+long_double_symbol (libm, __fabsl, fabsl)
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_fdim.c b/sysdeps/powerpc/powerpc64/fpu/s_fdim.c
new file mode 100644
index 0000000000..e34b51ee54
--- /dev/null
+++ b/sysdeps/powerpc/powerpc64/fpu/s_fdim.c
@@ -0,0 +1,5 @@
+#include <math_ldbl_opt.h>
+#include <sysdeps/powerpc/fpu/s_fdim.c>
+#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_1)
+compat_symbol (libm, __fdim, fdiml, GLIBC_2_1);
+#endif
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_floor.S b/sysdeps/powerpc/powerpc64/fpu/s_floor.S
index 80cbdc5709..65a2848b67 100644
--- a/sysdeps/powerpc/powerpc64/fpu/s_floor.S
+++ b/sysdeps/powerpc/powerpc64/fpu/s_floor.S
@@ -1,5 +1,5 @@
 /* Floor function.  PowerPC64 version.
-   Copyright (C) 2004 Free Software Foundation, Inc.
+   Copyright (C) 2004, 2006 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
@@ -18,13 +18,14 @@
    02111-1307 USA.  */
 
 #include <sysdep.h>
+#include <math_ldbl_opt.h>
 
 	.section	".toc","aw"
 .LC0:	/* 2**52 */
 	.tc FD_43300000_0[TC],0x4330000000000000
 	.section	".text"
 
-ENTRY (__floor)
+EALIGN (__floor, 4, 0)
 	CALL_MCOUNT 0
 	mffs	fp11		/* Save current FPU rounding mode.  */
 	lfd	fp13,.LC0@toc(2)
@@ -37,15 +38,16 @@ ENTRY (__floor)
 	ble-	cr6,.L4
 	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
 	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
-	fcmpu	cr5,fp1,fp12	/* if (x > 0.0)  */
+	fabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = 0.0; */
 	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
-	bnelr+	cr5
-	fmr	fp1,fp12	/* x must be +0.0 for the 0.0 case.  */
 	blr
 .L4:
 	bge-	cr6,.L9		/* if (x < 0.0)  */
 	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
 	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
+	fnabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = -0.0; */
 .L9:
 	mtfsf	0x01,fp11	/* restore previous rounding mode.  */	
 	blr
@@ -57,3 +59,6 @@ weak_alias (__floor, floor)
 weak_alias (__floor, floorl)
 strong_alias (__floor, __floorl)
 #endif
+#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_0)
+compat_symbol (libm, __floor, floorl, GLIBC_2_0)
+#endif
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_floorf.S b/sysdeps/powerpc/powerpc64/fpu/s_floorf.S
index 20cbb15ebd..bcdbf7823d 100644
--- a/sysdeps/powerpc/powerpc64/fpu/s_floorf.S
+++ b/sysdeps/powerpc/powerpc64/fpu/s_floorf.S
@@ -21,13 +21,13 @@
 
 	.section	".toc","aw"
 .LC0:	/* 2**23 */
-	.tc FD_41600000_0[TC],0x4160000000000000
+	.tc FD_4b000000_0[TC],0x4b00000000000000
 	.section	".text"
 
-ENTRY (__floorf)
+EALIGN (__floorf, 4, 0)
 	CALL_MCOUNT 0
 	mffs	fp11		/* Save current FPU rounding mode.  */
-	lfd	fp13,.LC0@toc(2)
+	lfs	fp13,.LC0@toc(2)
 	fabs	fp0,fp1
 	fsubs	fp12,fp13,fp13	/* generate 0.0  */
 	fcmpu	cr7,fp0,fp13	/* if (fabs(x) > TWO23)  */
@@ -37,15 +37,16 @@ ENTRY (__floorf)
 	ble-	cr6,.L4
 	fadds	fp1,fp1,fp13	/* x+= TWO23;  */
 	fsubs	fp1,fp1,fp13	/* x-= TWO23;  */
-	fcmpu	cr5,fp1,fp12	/* if (x > 0.0)  */
+	fabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = 0.0; */
 	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
-	bnelr+	cr5
-	fmr	fp1,fp12	/* x must be +0.0 for the 0.0 case.  */
 	blr
 .L4:
 	bge-	cr6,.L9		/* if (x < 0.0)  */
 	fsubs	fp1,fp1,fp13	/* x-= TWO23;  */
 	fadds	fp1,fp1,fp13	/* x+= TWO23;  */
+	fnabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = -0.0; */
 .L9:
 	mtfsf	0x01,fp11	/* restore previous rounding mode.  */	
 	blr
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_floorl.S b/sysdeps/powerpc/powerpc64/fpu/s_floorl.S
new file mode 100644
index 0000000000..01b3c2101d
--- /dev/null
+++ b/sysdeps/powerpc/powerpc64/fpu/s_floorl.S
@@ -0,0 +1,134 @@
+/* long double floor function.
+   IBM extended format long double version.
+   Copyright (C) 2004, 2006 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 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.  */
+
+#include <sysdep.h>
+#include <math_ldbl_opt.h>
+
+	.section	".toc","aw"
+.LC0:	/* 2**52 */
+	.tc FD_43300000_0[TC],0x4330000000000000
+
+	.section	".text"
+/* long double [fp1,fp2] floorl (long double x [fp1,fp2])
+   IEEE 1003.1 floor function.
+
+   PowerPC64 long double uses the IBM extended format which is
+   represented two 64-floating point double values. The values are
+   non-overlapping giving an effective precision of 106 bits. The first
+   double contains the high order bits of mantisa and is always rounded
+   to represent a normal rounding of long double to double. Since the
+   long double value is sum of the high and low values, the low double
+   normally has the opposite sign to compensate for the this rounding.
+
+   For long double there are two cases:
+   1) |x| < 2**52, all the integer bits are in the high double.
+      floor the high double and set the low double to -0.0.
+   2) |x| >= 2**52, Rounding involves both doubles.
+      See the comment before lable .L2 for details.
+   */
+
+ENTRY (__floorl)
+	mffs	fp11		/* Save current FPU rounding mode.  */
+	lfd	fp13,.LC0@toc(2)
+	fabs	fp0,fp1
+	fabs	fp9,fp2
+	fsub	fp12,fp13,fp13	/* generate 0.0  */
+	fcmpu	cr7,fp0,fp13	/* if (fabs(x) > TWO52)  */
+	fcmpu	cr6,fp1,fp12	/* if (x > 0.0)  */
+	bnl-	cr7,.L2
+	mtfsfi	7,3		/* Set rounding mode toward -inf.  */
+	fneg	fp2,fp12	/* set low double to -0.0.  */
+	ble-	cr6,.L0
+	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
+	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
+	fcmpu	cr5,fp1,fp12	/* if (x > 0.0)  */
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
+	bnelr+	cr5
+	fmr	fp1,fp12	/* x must be +0.0 for the 0.0 case.  */
+	blr
+.L0:
+	bge-	cr6,.L1		/* if (x < 0.0)  */
+	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
+	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
+.L1:
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
+	blr
+	
+
+/* The high double is > TWO52 so we need to round the low double and
+   perhaps the high double.  In this case we have to round the low
+   double and handle any adjustment to the high double that may be
+   caused by rounding (up).  This is complicated by the fact that the
+   high double may already be rounded and the low double may have the
+   opposite sign to compensate.This gets a bit tricky so we use the
+   following algorithm:
+
+   tau = floor(x_high/TWO52);
+   x0 = x_high - tau;
+   x1 = x_low + tau;
+   r1 = rint(x1);
+   y_high = x0 + r1;
+   y_low = x0 - y_high + r1;
+   return y;  */
+.L2:
+	fcmpu	cr7,fp9,fp13	/* if (|x_low| > TWO52)  */
+	fcmpu	cr0,fp9,fp12	/* || (|x_low| == 0.0)  */
+	fcmpu	cr5,fp2,fp12	/* if (x_low > 0.0)  */
+	bgelr-	cr7		/*   return x;	*/
+	beqlr-  cr0
+	mtfsfi	7,3		/* Set rounding mode toward -inf.  */
+	fdiv	fp8,fp1,fp13	/* x_high/TWO52  */
+	
+	bng-	cr6,.L6		/* if (x > 0.0)  */
+	fctidz	fp0,fp8
+	fcfid	fp8,fp0		/* tau = floor(x_high/TWO52);  */
+	bng	cr5,.L4		/* if (x_low > 0.0)  */
+	fmr	fp3,fp1
+	fmr	fp4,fp2
+	b	.L5
+.L4:				/* if (x_low < 0.0)  */
+	fsub	fp3,fp1,fp8	/* x0 = x_high - tau;  */
+	fadd	fp4,fp2,fp8	/* x1 = x_low + tau;  */
+.L5:
+	fadd	fp5,fp4,fp13	/* r1 = r1 + TWO52;  */
+	fsub	fp5,fp5,fp13	/* r1 = r1 - TWO52;  */
+	b	.L9
+.L6:				/* if (x < 0.0)  */
+	fctidz	fp0,fp8
+	fcfid	fp8,fp0		/* tau = floor(x_high/TWO52);  */	
+	bnl	cr5,.L7		/* if (x_low < 0.0)  */
+	fmr	fp3,fp1
+	fmr	fp4,fp2
+	b	.L8
+.L7:				/* if (x_low > 0.0)  */
+	fsub	fp3,fp1,fp8	/* x0 = x_high - tau;  */
+	fadd	fp4,fp2,fp8	/* x1 = x_low + tau;  */
+.L8:
+	fsub	fp5,fp4,fp13	/* r1-= TWO52;  */
+	fadd	fp5,fp5,fp13	/* r1+= TWO52;  */
+.L9:
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
+	fadd	fp1,fp3,fp5	/* y_high = x0 + r1;  */
+	fsub	fp2,fp3,fp1	/* y_low = x0 - y_high + r1;  */
+	fadd	fp2,fp2,fp5
+	blr
+END (__floorl)
+
+long_double_symbol (libm, __floorl, floorl)
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_fmax.S b/sysdeps/powerpc/powerpc64/fpu/s_fmax.S
new file mode 100644
index 0000000000..69735761ab
--- /dev/null
+++ b/sysdeps/powerpc/powerpc64/fpu/s_fmax.S
@@ -0,0 +1,5 @@
+#include <math_ldbl_opt.h>
+#include <sysdeps/powerpc/fpu/s_fmax.S>
+#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_1)
+compat_symbol (libm, __fmax, fmaxl, GLIBC_2_1)
+#endif
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_fmin.S b/sysdeps/powerpc/powerpc64/fpu/s_fmin.S
new file mode 100644
index 0000000000..6d4a0a946c
--- /dev/null
+++ b/sysdeps/powerpc/powerpc64/fpu/s_fmin.S
@@ -0,0 +1,5 @@
+#include <math_ldbl_opt.h>
+#include <sysdeps/powerpc/fpu/s_fmin.S>
+#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_1)
+compat_symbol (libm, __fmin, fminl, GLIBC_2_1)
+#endif
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_isnan.c b/sysdeps/powerpc/powerpc64/fpu/s_isnan.c
new file mode 100644
index 0000000000..397717ba9c
--- /dev/null
+++ b/sysdeps/powerpc/powerpc64/fpu/s_isnan.c
@@ -0,0 +1,7 @@
+#include <sysdeps/powerpc/fpu/s_isnan.c>
+#ifndef IS_IN_libm
+# if LONG_DOUBLE_COMPAT(libc, GLIBC_2_0)
+compat_symbol (libc, __isnan, __isnanl, GLIBC_2_0);
+compat_symbol (libc, isnan, isnanl, GLIBC_2_0);
+# endif
+#endif
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_llrint.S b/sysdeps/powerpc/powerpc64/fpu/s_llrint.S
index 610b561f25..ff0ba948a5 100644
--- a/sysdeps/powerpc/powerpc64/fpu/s_llrint.S
+++ b/sysdeps/powerpc/powerpc64/fpu/s_llrint.S
@@ -1,5 +1,5 @@
 /* Round double to long int.  PowerPC64 version.
-   Copyright (C) 2004 Free Software Foundation, Inc.
+   Copyright (C) 2004, 2006 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
@@ -18,6 +18,7 @@
    02111-1307 USA.  */
 
 #include <sysdep.h>
+#include <math_ldbl_opt.h>
 
 /* long long int[r3] __llrint (double x[fp1])  */
 ENTRY (__llrint)	
@@ -41,3 +42,7 @@ weak_alias (__llrint, llrintl)
 strong_alias (__lrint, __lrintl)
 weak_alias (__lrint, lrintl)
 #endif
+#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_1)
+compat_symbol (libm, __llrint, llrintl, GLIBC_2_1)
+compat_symbol (libm, __lrint, lrintl, GLIBC_2_1)
+#endif
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_llround.S b/sysdeps/powerpc/powerpc64/fpu/s_llround.S
index a3dcd4c33d..d023b8f2c0 100644
--- a/sysdeps/powerpc/powerpc64/fpu/s_llround.S
+++ b/sysdeps/powerpc/powerpc64/fpu/s_llround.S
@@ -1,5 +1,5 @@
 /* llround function.  PowerPC64 version.
-   Copyright (C) 2004 Free Software Foundation, Inc.
+   Copyright (C) 2004, 2006 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
@@ -18,6 +18,7 @@
    02111-1307 USA.  */
 
 #include <sysdep.h>
+#include <math_ldbl_opt.h>
 
 	.section	".toc","aw"
 .LC0:	/* -0.0 */
@@ -66,3 +67,7 @@ strong_alias (__llround, __llroundl)
 weak_alias (__lround, lroundl)
 strong_alias (__lround, __lroundl)
 #endif
+#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_1)
+compat_symbol (libm, __llround, llroundl, GLIBC_2_1)
+compat_symbol (libm, __lround, lroundl, GLIBC_2_1)
+#endif
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_llroundf.S b/sysdeps/powerpc/powerpc64/fpu/s_llroundf.S
index b5ca43bf20..bbbd05492e 100644
--- a/sysdeps/powerpc/powerpc64/fpu/s_llroundf.S
+++ b/sysdeps/powerpc/powerpc64/fpu/s_llroundf.S
@@ -1,5 +1,5 @@
 /* llroundf function.  PowerPC64 version.
-   Copyright (C) 2004 Free Software Foundation, Inc.
+   Copyright (C) 2004, 2006 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
@@ -28,8 +28,8 @@
 	
 /* long long [r3] llroundf (float x [fp1])
    IEEE 1003.1 llroundf function.  IEEE specifies "roundf to the nearest 
-   integer value, roundfing halfway cases away from zero, regardless of
-   the current roundfing mode."  However PowerPC Architecture defines
+   integer value, rounding halfway cases away from zero, regardless of
+   the current rounding mode."  However PowerPC Architecture defines
    "roundf to Nearest" as "Choose the best approximation. In case of a 
    tie, choose the one that is even (least significant bit o).". 
    So we can't use the PowerPC "round to Nearest" mode. Instead we set
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_nearbyintl.S b/sysdeps/powerpc/powerpc64/fpu/s_nearbyintl.S
new file mode 100644
index 0000000000..0d0eb36f98
--- /dev/null
+++ b/sysdeps/powerpc/powerpc64/fpu/s_nearbyintl.S
@@ -0,0 +1,114 @@
+/* nearbyint long double.
+   IBM extended format long double version.
+   Copyright (C) 2004, 2006 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 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.  */
+
+#include <sysdep.h>
+#include <math_ldbl_opt.h>
+
+	.section	".toc","aw"
+.LC0:	/* 2**52 */
+	.tc FD_43300000_0[TC],0x4330000000000000
+	.section	".text"
+
+/* long double [fp1,fp2] nearbyintl (long double x [fp1,fp2])
+   IEEE 1003.1 nearbyintl function.  nearbyintl is simular to the rintl
+   but does raise the "inexact" exception.  This implementation is
+   based on rintl but explicitly maskes the inexact exception on entry
+   and clears any pending inexact before restoring the exception mask
+   on exit.
+
+   PowerPC64 long double uses the IBM extended format which is
+   represented two 64-floating point double values. The values are
+   non-overlapping giving an effective precision of 106 bits. The first
+   double contains the high order bits of mantisa and is always rounded
+   to represent a normal rounding of long double to double. Since the
+   long double value is sum of the high and low values, the low double
+   normally has the opposite sign to compensate for the this rounding.
+
+   For long double there are two cases:
+   1) |x| < 2**52, all the integer bits are in the high double.
+      floor the high double and set the low double to -0.0.
+   2) |x| >= 2**52, Rounding involves both doubles.
+      See the comment before lable .L2 for details.
+   */
+ENTRY (__nearbyintl)
+	mffs	fp11		/* Save current FPSCR.  */
+	lfd	fp13,.LC0@toc(2)
+	fabs	fp0,fp1
+	mtfsb0  28		/* Disable "inexact" exceptions.  */
+	fsub	fp12,fp13,fp13	/* generate 0.0  */
+	fabs	fp9,fp2
+	fcmpu	cr7,fp0,fp13	/* if (fabs(x) > TWO52)  */
+	fcmpu	cr6,fp1,fp12	/* if (x > 0.0)  */
+	bnl-	cr7,.L2
+	fmr	fp2,fp12
+	bng-	cr6,.L4
+	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
+	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
+	b	.L9
+.L4:
+	bnl-	cr6,.L9		/* if (x < 0.0)  */
+	fsub	fp1,fp13,fp1	/* x = TWO52 - x;  */
+	fsub	fp0,fp1,fp13	/* x = - (x - TWO52);  */
+	fneg	fp1,fp0
+.L9:
+	mtfsb0	6		/* Clear any pending "inexact" exceptions.  */
+	mtfsf	0x01,fp11	/* restore exception mask.  */
+	blr
+
+/* The high double is > TWO52 so we need to round the low double and
+   perhaps the high double.  This gets a bit tricky so we use the
+   following algorithm:
+
+   tau = floor(x_high/TWO52);
+   x0 = x_high - tau;
+   x1 = x_low + tau;
+   r1 = nearbyint(x1);
+   y_high = x0 + r1;
+   y_low = r1 - tau;
+   return y;  */
+.L2:
+	fcmpu	cr7,fp9,fp13	/* if (|x_low| > TWO52)  */
+	fcmpu	cr0,fp9,fp12	/* || (|x_low| == 0.0)  */
+	bge-	cr7,.L9		/*   return x;	*/
+	beq-  cr0,.L9
+	fdiv	fp8,fp1,fp13	/* x_high/TWO52  */
+	fctidz	fp0,fp8
+	fcfid	fp8,fp0		/* tau = floor(x_high/TWO52);  */
+	fsub	fp3,fp1,fp8	/* x0 = x_high - tau;  */
+	fadd	fp4,fp2,fp8	/* x1 = x_low + tau;  */
+
+	fcmpu	cr6,fp4,fp12	/* if (x1 > 0.0)  */
+	bng-	cr6,.L8
+	fadd	fp5,fp4,fp13	/* r1 = x1 + TWO52;  */
+	fsub	fp5,fp5,fp13	/* r1 = r1 - TWO52;  */
+	b	.L6
+.L8:
+	fmr	fp5,fp4
+	bge-	cr6,.L6		/* if (x1 < 0.0)  */
+	fsub	fp5,fp13,fp4	/* r1 = TWO52 - x1;  */
+	fsub	fp0,fp5,fp13	/* r1 = - (r1 - TWO52);  */
+	fneg	fp5,fp0
+.L6:
+	fadd	fp1,fp3,fp5	/* y_high = x0 + r1;  */
+	fsub	fp2,fp5,fp8	/* y_low = r1 - tau;  */
+	b	.L9
+END (__nearbyintl)
+
+long_double_symbol (libm, __nearbyintl, nearbyintl)
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_rint.S b/sysdeps/powerpc/powerpc64/fpu/s_rint.S
index 79e807269d..b4fbc0b51b 100644
--- a/sysdeps/powerpc/powerpc64/fpu/s_rint.S
+++ b/sysdeps/powerpc/powerpc64/fpu/s_rint.S
@@ -1,5 +1,5 @@
 /* Round to int floating-point values.  PowerPC64 version.
-   Copyright (C) 2004 Free Software Foundation, Inc.
+   Copyright (C) 2004, 2006 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
@@ -21,13 +21,14 @@
    when it's coded in C.  */
 
 #include <sysdep.h>
+#include <math_ldbl_opt.h>
 
 	.section	".toc","aw"
 .LC0:	/* 2**52 */
 	.tc FD_43300000_0[TC],0x4330000000000000
 	.section	".text"
 
-ENTRY (__rint)
+EALIGN (__rint, 4, 0)
 	CALL_MCOUNT 0
 	lfd	fp13,.LC0@toc(2)
 	fabs	fp0,fp1
@@ -38,13 +39,14 @@ ENTRY (__rint)
 	bng-	cr6,.L4
 	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
 	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
-	blr
+	fabs	fp1,fp1		/* if (x == 0.0)  */
+	blr			/* x = 0.0; */
 .L4:
 	bnllr-	cr6		/* if (x < 0.0)  */
-	fsub	fp1,fp13,fp1	/* x = TWO52 - x;  */
-	fsub	fp0,fp1,fp13	/* x = - (x - TWO52);  */
-	fneg	fp1,fp0
-	blr
+	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
+	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
+	fnabs	fp1,fp1		/* if (x == 0.0)  */
+	blr			/* x = -0.0; */
 	END (__rint)
 
 weak_alias (__rint, rint)
@@ -53,3 +55,6 @@ weak_alias (__rint, rint)
 weak_alias (__rint, rintl)
 strong_alias (__rint, __rintl)
 #endif
+#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_0)
+compat_symbol (libm, __rint, rintl, GLIBC_2_0)
+#endif
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_rintf.S b/sysdeps/powerpc/powerpc64/fpu/s_rintf.S
index eb34dd5e77..e4fa9ba2e6 100644
--- a/sysdeps/powerpc/powerpc64/fpu/s_rintf.S
+++ b/sysdeps/powerpc/powerpc64/fpu/s_rintf.S
@@ -21,12 +21,12 @@
 
 	.section	".toc","aw"
 .LC0:	/* 2**23 */
-	.tc FD_41600000_0[TC],0x4160000000000000
+	.tc FD_4b000000_0[TC],0x4b00000000000000
 	.section	".text"
 
-ENTRY (__rintf)
+EALIGN (__rintf, 4, 0)
 	CALL_MCOUNT 0
-	lfd	fp13,.LC0@toc(2)
+	lfs	fp13,.LC0@toc(2)
 	fabs	fp0,fp1
 	fsubs	fp12,fp13,fp13	/* generate 0.0  */
 	fcmpu	cr7,fp0,fp13	/* if (fabs(x) > TWO23)  */
@@ -35,13 +35,14 @@ ENTRY (__rintf)
 	bng-	cr6,.L4
 	fadds	fp1,fp1,fp13	/* x+= TWO23;  */
 	fsubs	fp1,fp1,fp13	/* x-= TWO23;  */
-	blr
+	fabs	fp1,fp1		/* if (x == 0.0)  */
+	blr			/* x = 0.0; */
 .L4:
 	bnllr-	cr6		/* if (x < 0.0)  */
-	fsubs	fp1,fp13,fp1	/* x = TWO23 - x;  */
-	fsubs	fp0,fp1,fp13	/* x = - (x - TWO23);  */
-	fneg	fp1,fp0
-	blr
+	fsubs	fp1,fp1,fp13	/* x-= TWO23;  */
+	fadds	fp1,fp1,fp13	/* x+= TWO23;  */
+	fnabs	fp1,fp1		/* if (x == 0.0)  */
+	blr			/* x = -0.0; */
 	END (__rintf)
 
 weak_alias (__rintf, rintf)
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_round.S b/sysdeps/powerpc/powerpc64/fpu/s_round.S
index c0b6d46fea..15afca1543 100644
--- a/sysdeps/powerpc/powerpc64/fpu/s_round.S
+++ b/sysdeps/powerpc/powerpc64/fpu/s_round.S
@@ -1,5 +1,5 @@
 /* round function.  PowerPC64 version.
-   Copyright (C) 2004 Free Software Foundation, Inc.
+   Copyright (C) 2004, 2006 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
@@ -18,14 +18,13 @@
    02111-1307 USA.  */
 
 #include <sysdep.h>
+#include <math_ldbl_opt.h>
 
 	.section	".toc","aw"
 .LC0:	/* 2**52 */
 	.tc FD_43300000_0[TC],0x4330000000000000
 .LC1:	/* 0.5 */
 	.tc FD_3fe00000_0[TC],0x3fe0000000000000
-.LC2:	/* -0.0 */
-	.tc FD_80000000_0[TC],0x8000000000000000
 	.section	".text"
 	
 /* double [fp1] round (double x [fp1])
@@ -38,7 +37,7 @@
    "Round toward Zero" mode and round by adding +-0.5 before rounding
    to the integer value.  */
 
-ENTRY (__round)
+EALIGN (__round, 4, 0)
 	CALL_MCOUNT 0
 	mffs	fp11		/* Save current FPU rounding mode.  */
 	lfd	fp13,.LC0@toc(2)
@@ -53,7 +52,8 @@ ENTRY (__round)
 	fadd	fp1,fp1,fp10	/* x+= 0.5;  */
 	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
 	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
-.L9:	
+	fabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = 0.0; */
 	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
 	blr
 .L4:
@@ -61,10 +61,10 @@ ENTRY (__round)
 	bge-	cr6,.L9		/* if (x < 0.0)  */
 	fsub	fp1,fp9,fp13	/* x-= TWO52;  */
 	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
-	fcmpu	cr5,fp1,fp12	/* if (x > 0.0)  */
-	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
-	bnelr+	cr5
-	lfd	fp1,.LC2@toc(2)	/* x must be -0.0 for the 0.0 case.  */	
+	fnabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = -0.0; */
+.L9:
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */	
 	blr
 	END (__round)
 
@@ -74,3 +74,6 @@ weak_alias (__round, round)
 weak_alias (__round, roundl)
 strong_alias (__round, __roundl)
 #endif
+#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_1)
+compat_symbol (libm, __round, roundl, GLIBC_2_1)
+#endif
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_roundf.S b/sysdeps/powerpc/powerpc64/fpu/s_roundf.S
index 23ee4c052b..d2e29fdb8f 100644
--- a/sysdeps/powerpc/powerpc64/fpu/s_roundf.S
+++ b/sysdeps/powerpc/powerpc64/fpu/s_roundf.S
@@ -21,11 +21,9 @@
 
 	.section	".toc","aw"
 .LC0:	/* 2**23 */
-	.tc FD_41600000_0[TC],0x4160000000000000
+	.tc FD_4b000000_0[TC],0x4b00000000000000
 .LC1:	/* 0.5 */
-	.tc FD_3fe00000_0[TC],0x3fe0000000000000
-.LC2:	/* -0.0 */
-	.tc FD_80000000_0[TC],0x8000000000000000
+	.tc FD_3f000000_0[TC],0x3f00000000000000
 	.section	".text"
 	
 /* float [fp1] roundf  (float x [fp1])
@@ -38,22 +36,23 @@
    "Round toward Zero" mode and round by adding +-0.5 before rounding
    to the integer value.  */
 
-ENTRY (__roundf )
+EALIGN (__roundf, 4, 0)
 	CALL_MCOUNT 0
 	mffs	fp11		/* Save current FPU rounding mode.  */
-	lfd	fp13,.LC0@toc(2)
+	lfs	fp13,.LC0@toc(2)
 	fabs	fp0,fp1
 	fsubs	fp12,fp13,fp13	/* generate 0.0  */
 	fcmpu	cr7,fp0,fp13	/* if (fabs(x) > TWO23)  */
 	fcmpu	cr6,fp1,fp12	/* if (x > 0.0)  */
 	bnllr-	cr7
 	mtfsfi	7,1		/* Set rounding mode toward 0.  */
-	lfd	fp10,.LC1@toc(2)
+	lfs	fp10,.LC1@toc(2)
 	ble-	cr6,.L4
 	fadds	fp1,fp1,fp10	/* x+= 0.5;  */
 	fadds	fp1,fp1,fp13	/* x+= TWO23;  */
 	fsubs	fp1,fp1,fp13	/* x-= TWO23;  */
-.L9:	
+	fabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = 0.0; */
 	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
 	blr
 .L4:
@@ -61,10 +60,10 @@ ENTRY (__roundf )
 	bge-	cr6,.L9		/* if (x < 0.0)  */
 	fsubs	fp1,fp9,fp13	/* x-= TWO23;  */
 	fadds	fp1,fp1,fp13	/* x+= TWO23;  */
-	fcmpu	cr5,fp1,fp12	/* if (x > 0.0)  */
-	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
-	bnelr+	cr5
-	lfd	fp1,.LC2@toc(2)	/* x must be -0.0 for the 0.0 case.  */	
+	fnabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = -0.0; */
+.L9:
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */	
 	blr
 	END (__roundf)
 
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_roundl.S b/sysdeps/powerpc/powerpc64/fpu/s_roundl.S
new file mode 100644
index 0000000000..20da828a82
--- /dev/null
+++ b/sysdeps/powerpc/powerpc64/fpu/s_roundl.S
@@ -0,0 +1,133 @@
+/* long double round function.
+   IBM extended format long double version.
+   Copyright (C) 2004, 2006 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 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.  */
+
+#include <sysdep.h>
+#include <math_ldbl_opt.h>
+
+	.section	".toc","aw"
+.LC0:	/* 2**52 */
+	.tc FD_43300000_0[TC],0x4330000000000000
+.LC1:	/* 0.5 */
+	.tc FD_3fe00000_0[TC],0x3fe0000000000000
+	.section	".text"
+
+/* long double [fp1,fp2] roundl (long double x [fp1,fp2])
+   IEEE 1003.1 round function.  IEEE specifies "round to the nearest
+   integer value, rounding halfway cases away from zero, regardless of
+   the current rounding mode."  However PowerPC Architecture defines
+   "Round to Nearest" as "Choose the best approximation. In case of a
+   tie, choose the one that is even (least significant bit o).".
+   So we can't use the PowerPC "Round to Nearest" mode. Instead we set
+   "Round toward Zero" mode and round by adding +-0.5 before rounding
+   to the integer value.  */
+
+ENTRY (__roundl)
+	mffs	fp11		/* Save current FPU rounding mode.  */
+	lfd	fp13,.LC0@toc(2)
+	fabs	fp0,fp1
+	fabs	fp9,fp2
+	fsub	fp12,fp13,fp13	/* generate 0.0  */
+	fcmpu	cr7,fp0,fp13	/* if (fabs(x) > TWO52)  */
+	fcmpu	cr6,fp1,fp12	/* if (x > 0.0)  */
+	bnl-	cr7,.L2
+	mtfsfi	7,1		/* Set rounding mode toward 0.  */
+	lfd	fp10,.LC1@toc(2)
+	ble-	cr6,.L1
+	fneg	fp2,fp12
+	fadd	fp1,fp1,fp10	/* x+= 0.5;  */
+	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
+	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
+	fabs	fp1,fp1		/* if (x == 0.0) x = 0.0;  */
+.L0:
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
+	blr
+.L1:
+	fsub	fp9,fp1,fp10	/* x-= 0.5;  */
+	fneg	fp2,fp12
+	bge-	cr6,.L0		/* if (x < 0.0)  */
+	fsub	fp1,fp9,fp13	/* x-= TWO52;  */
+	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
+	fnabs	fp1,fp1		/* if (x == 0.0) x = -0.0;  */
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
+	blr
+
+/* The high double is > TWO52 so we need to round the low double and
+   perhaps the high double.  In this case we have to round the low
+   double and handle any adjustment to the high double that may be
+   caused by rounding (up).  This is complicated by the fact that the
+   high double may already be rounded and the low double may have the
+   opposite sign to compensate.This gets a bit tricky so we use the
+   following algorithm:
+
+   tau = floor(x_high/TWO52);
+   x0 = x_high - tau;
+   x1 = x_low + tau;
+   r1 = rint(x1);
+   y_high = x0 + r1;
+   y_low = x0 - y_high + r1;
+   return y;  */
+.L2:
+	fcmpu	cr7,fp9,fp13	/* if (|x_low| > TWO52)  */
+	fcmpu	cr0,fp9,fp12	/* || (|x_low| == 0.0)  */
+	fcmpu	cr5,fp2,fp12	/* if (x_low > 0.0)  */
+	lfd	fp10,.LC1@toc(2)
+	bgelr-	cr7		/*   return x;	*/
+	beqlr-  cr0
+	mtfsfi	7,1		/* Set rounding mode toward 0.  */
+	fdiv	fp8,fp1,fp13	/* x_high/TWO52  */
+	
+	bng-	cr6,.L6		/* if (x > 0.0)  */
+	fctidz	fp0,fp8
+	fcfid	fp8,fp0		/* tau = floor(x_high/TWO52);  */
+	bng	cr5,.L4		/* if (x_low > 0.0)  */
+	fmr	fp3,fp1
+	fmr	fp4,fp2
+	b	.L5
+.L4:				/* if (x_low < 0.0)  */
+	fsub	fp3,fp1,fp8	/* x0 = x_high - tau;  */
+	fadd	fp4,fp2,fp8	/* x1 = x_low + tau;  */
+.L5:
+	fadd	fp5,fp4,fp10	/* r1 = x1 + 0.5;  */
+	fadd	fp5,fp5,fp13	/* r1 = r1 + TWO52;  */
+	fsub	fp5,fp5,fp13	/* r1 = r1 - TWO52;  */
+	b	.L9
+.L6:				/* if (x < 0.0)  */
+	fctidz	fp0,fp8
+	fcfid	fp8,fp0		/* tau = floor(x_high/TWO52);  */	
+	bnl	cr5,.L7		/* if (x_low < 0.0)  */
+	fmr	fp3,fp1
+	fmr	fp4,fp2
+	b	.L8
+.L7:				/* if (x_low > 0.0)  */
+	fsub	fp3,fp1,fp8	/* x0 = x_high - tau;  */
+	fadd	fp4,fp2,fp8	/* x1 = x_low + tau;  */
+.L8:
+	fsub	fp5,fp4,fp10	/* r1 = x1 - 0.5;  */
+	fsub	fp5,fp5,fp13	/* r1-= TWO52;  */
+	fadd	fp5,fp5,fp13	/* r1+= TWO52;  */
+.L9:
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
+	fadd	fp1,fp3,fp5	/* y_high = x0 + r1;  */
+	fsub	fp2,fp3,fp1	/* y_low = x0 - y_high + r1;  */
+	fadd	fp2,fp2,fp5
+	blr
+END (__roundl)
+
+long_double_symbol (libm, __roundl, roundl)
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_trunc.S b/sysdeps/powerpc/powerpc64/fpu/s_trunc.S
index 3ddd298525..086ed0025e 100644
--- a/sysdeps/powerpc/powerpc64/fpu/s_trunc.S
+++ b/sysdeps/powerpc/powerpc64/fpu/s_trunc.S
@@ -1,5 +1,5 @@
 /* trunc function.  PowerPC64 version.
-   Copyright (C) 2004 Free Software Foundation, Inc.
+   Copyright (C) 2004, 2006 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
@@ -18,12 +18,11 @@
    02111-1307 USA.  */
 
 #include <sysdep.h>
+#include <math_ldbl_opt.h>
 
 	.section	".toc","aw"
 .LC0:	/* 2**52 */
 	.tc FD_43300000_0[TC],0x4330000000000000
-.LC2:	/* -0.0 */
-	.tc FD_80000000_0[TC],0x8000000000000000
 	.section	".text"
 	
 /* double [fp1] trunc (double x [fp1])
@@ -33,7 +32,7 @@
    We set "round toward Zero" mode and trunc by adding +-2**52 then
    subtracting +-2**52.  */
 
-ENTRY (__trunc)
+EALIGN (__trunc, 4, 0)
 	CALL_MCOUNT 0
 	mffs	fp11		/* Save current FPU rounding mode.  */
 	lfd	fp13,.LC0@toc(2)
@@ -46,17 +45,18 @@ ENTRY (__trunc)
 	ble-	cr6,.L4
 	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
 	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
-.L9:	
-	mtfsf	0x01,fp11	/* restore previous truncing mode.  */
+	fabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = 0.0; */
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
 	blr
 .L4:
 	bge-	cr6,.L9		/* if (x < 0.0)  */
 	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
 	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
-	fcmpu	cr5,fp1,fp12	/* if (x > 0.0)  */
+	fnabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = -0.0; */
+.L9:
 	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
-	bnelr+	cr5
-	lfd	fp1,.LC2@toc(2)	/* x must be -0.0 for the 0.0 case.  */	
 	blr
 	END (__trunc)
 
@@ -66,3 +66,6 @@ weak_alias (__trunc, trunc)
 weak_alias (__trunc, truncl)
 strong_alias (__trunc, __truncl)
 #endif
+#if LONG_DOUBLE_COMPAT(libm, GLIBC_2_1)
+compat_symbol (libm, __trunc, truncl, GLIBC_2_1)
+#endif
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_truncf.S b/sysdeps/powerpc/powerpc64/fpu/s_truncf.S
index b38b722a6f..1456e7f434 100644
--- a/sysdeps/powerpc/powerpc64/fpu/s_truncf.S
+++ b/sysdeps/powerpc/powerpc64/fpu/s_truncf.S
@@ -1,5 +1,5 @@
 /* truncf function.  PowerPC64 version.
-   Copyright (C) 2004 Free Software Foundation, Inc.
+   Copyright (C) 2004, 2006 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
@@ -21,9 +21,7 @@
 
 	.section	".toc","aw"
 .LC0:	/* 2**23 */
-	.tc FD_41600000_0[TC],0x4160000000000000
-.LC2:	/* -0.0 */
-	.tc FD_80000000_0[TC],0x8000000000000000
+	.tc FD_4b000000_0[TC],0x4b00000000000000
 	.section	".text"
 	
 /* float [fp1] truncf (float x [fp1])
@@ -33,10 +31,10 @@
    We set "round toward Zero" mode and trunc by adding +-2**23 then
    subtracting +-2**23.  */
 
-ENTRY (__truncf)
+EALIGN (__truncf, 4, 0)
 	CALL_MCOUNT 0
 	mffs	fp11		/* Save current FPU rounding mode.  */
-	lfd	fp13,.LC0@toc(2)
+	lfs	fp13,.LC0@toc(2)
 	fabs	fp0,fp1
 	fsubs	fp12,fp13,fp13	/* generate 0.0  */
 	fcmpu	cr7,fp0,fp13	/* if (fabs(x) > TWO23)  */
@@ -46,17 +44,18 @@ ENTRY (__truncf)
 	ble-	cr6,.L4
 	fadds	fp1,fp1,fp13	/* x+= TWO23;  */
 	fsubs	fp1,fp1,fp13	/* x-= TWO23;  */
-.L9:	
-	mtfsf	0x01,fp11	/* restore previous truncing mode.  */
+	fabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = 0.0; */
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
 	blr
 .L4:
 	bge-	cr6,.L9		/* if (x < 0.0)  */
 	fsubs	fp1,fp1,fp13	/* x-= TWO23;  */
 	fadds	fp1,fp1,fp13	/* x+= TWO23;  */
-	fcmpu	cr5,fp1,fp12	/* if (x > 0.0)  */
+	fnabs	fp1,fp1		/* if (x == 0.0)  */
+				/* x = -0.0; */
+.L9:
 	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
-	bnelr+	cr5
-	lfd	fp1,.LC2@toc(2)	/* x must be -0.0 for the 0.0 case.  */	
 	blr
 	END (__truncf)
 
diff --git a/sysdeps/powerpc/powerpc64/fpu/s_truncl.S b/sysdeps/powerpc/powerpc64/fpu/s_truncl.S
new file mode 100644
index 0000000000..1864fc14b7
--- /dev/null
+++ b/sysdeps/powerpc/powerpc64/fpu/s_truncl.S
@@ -0,0 +1,121 @@
+/* long double trunc function.
+   IBM extended format long double version.
+   Copyright (C) 2004, 2006 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 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.  */
+
+#include <sysdep.h>
+#include <math_ldbl_opt.h>
+
+	.section	".toc","aw"
+.LC0:	/* 2**52 */
+	.tc FD_43300000_0[TC],0x4330000000000000
+.LC1:	/* 0.5 */
+	.tc FD_3fe00000_0[TC],0x3fe0000000000000
+	.section	".text"
+
+/* long double [fp1,fp2] truncl (long double x [fp1,fp2])  */
+
+ENTRY (__truncl)
+	mffs	fp11		/* Save current FPU rounding mode.  */
+	lfd	fp13,.LC0@toc(2)
+	fabs	fp0,fp1
+	fabs	fp9,fp2
+	fsub	fp12,fp13,fp13	/* generate 0.0  */
+	fcmpu	cr7,fp0,fp13	/* if (fabs(x) > TWO52)  */
+	fcmpu	cr6,fp1,fp12	/* if (x > 0.0)  */
+	bnl-	cr7,.L2
+	mtfsfi	7,1		/* Set rounding mode toward 0.  */
+	ble-	cr6,.L1
+	fneg	fp2,fp12
+	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
+	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
+	fabs	fp1,fp1		/* if (x == 0.0) x = 0.0;  */
+.L0:
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
+	blr
+.L1:
+	fneg	fp2,fp12
+	bge-	cr6,.L0		/* if (x < 0.0)  */
+	fsub	fp1,fp1,fp13	/* x-= TWO52;  */
+	fadd	fp1,fp1,fp13	/* x+= TWO52;  */
+	fnabs	fp1,fp1		/* if (x == 0.0) x = -0.0;  */
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
+	blr
+
+/* The high double is > TWO52 so we need to round the low double and
+   perhaps the high double.  In this case we have to round the low
+   double and handle any adjustment to the high double that may be
+   caused by rounding (up).  This is complicated by the fact that the
+   high double may already be rounded and the low double may have the
+   opposite sign to compensate.This gets a bit tricky so we use the
+   following algorithm:
+
+   tau = floor(x_high/TWO52);
+   x0 = x_high - tau;
+   x1 = x_low + tau;
+   r1 = rint(x1);
+   y_high = x0 + r1;
+   y_low = x0 - y_high + r1;
+   return y;  */
+.L2:
+	fcmpu	cr7,fp9,fp13	/* if (|x_low| > TWO52)  */
+	fcmpu	cr0,fp9,fp12	/* || (|x_low| == 0.0)  */
+	fcmpu	cr5,fp2,fp12	/* if (x_low > 0.0)  */
+	bgelr-	cr7		/*   return x;	*/
+	beqlr-  cr0
+	mtfsfi	7,1		/* Set rounding mode toward 0.  */
+	fdiv	fp8,fp1,fp13	/* x_high/TWO52  */
+	
+	bng-	cr6,.L6		/* if (x > 0.0)  */
+	fctidz	fp0,fp8
+	fcfid	fp8,fp0		/* tau = floor(x_high/TWO52);  */
+	fadd	fp8,fp8,fp8	/* tau++; Make tau even  */
+	bng	cr5,.L4		/* if (x_low > 0.0)  */
+	fmr	fp3,fp1
+	fmr	fp4,fp2
+	b	.L5
+.L4:				/* if (x_low < 0.0)  */
+	fsub	fp3,fp1,fp8	/* x0 = x_high - tau;  */
+	fadd	fp4,fp2,fp8	/* x1 = x_low + tau;  */
+.L5:
+	fadd	fp5,fp4,fp13	/* r1 = r1 + TWO52;  */
+	fsub	fp5,fp5,fp13	/* r1 = r1 - TWO52;  */
+	b	.L9
+.L6:				/* if (x < 0.0)  */
+	fctidz	fp0,fp8
+	fcfid	fp8,fp0		/* tau = floor(x_high/TWO52);  */
+	fadd	fp8,fp8,fp8	/* tau++; Make tau even  */	
+	bnl	cr5,.L7		/* if (x_low < 0.0)  */
+	fmr	fp3,fp1
+	fmr	fp4,fp2
+	b	.L8
+.L7:				/* if (x_low > 0.0)  */
+	fsub	fp3,fp1,fp8	/* x0 = x_high - tau;  */
+	fadd	fp4,fp2,fp8	/* x1 = x_low + tau;  */
+.L8:
+	fsub	fp5,fp4,fp13	/* r1-= TWO52;  */
+	fadd	fp5,fp5,fp13	/* r1+= TWO52;  */
+.L9:
+	mtfsf	0x01,fp11	/* restore previous rounding mode.  */
+	fadd	fp1,fp3,fp5	/* y_high = x0 + r1;  */
+	fsub	fp2,fp3,fp1	/* y_low = x0 - y_high + r1;  */
+	fadd	fp2,fp2,fp5
+	blr
+END (__truncl)
+
+long_double_symbol (libm, __truncl, truncl)