- ldexp implementation has issues. switch to the one from libm

- remove frexp in hppa64, cloned from hppa
- move generic ieee754 implementations of modf and ldexp to gen
ok kettenis@, "looks good" millert@

28 files changed, 23 insertions, 3111 deletions

diff --git a/lib/libc/arch/alpha/gen/Makefile.inc b/lib/libc/arch/alpha/gen/Makefile.inc
index 2f82f39d175..3e4eac62983 100644
--- a/lib/libc/arch/alpha/gen/Makefile.inc
+++ b/lib/libc/arch/alpha/gen/Makefile.inc
@@ -1,7 +1,7 @@
-#	$OpenBSD: Makefile.inc,v 1.8 2008/07/24 09:31:06 martynas Exp $
+#	$OpenBSD: Makefile.inc,v 1.9 2009/04/19 16:42:05 martynas Exp $
 #	$NetBSD: Makefile.inc,v 1.3 1995/04/29 05:09:14 cgd Exp $
 
-SRCS+=	_setjmp.S fabs.S infinity.c ldexp.c modf.c nan.c setjmp.S
+SRCS+=	_setjmp.S fabs.S infinity.c nan.c setjmp.S
 SRCS+=	flt_rounds.c fpgetmask.c fpgetround.c fpgetsticky.c fpsetmask.c \
 	fpsetround.c fpsetsticky.c
 SRCS+=	sigsetjmp.S
diff --git a/lib/libc/arch/alpha/gen/ldexp.c b/lib/libc/arch/alpha/gen/ldexp.c
deleted file mode 100644
index 66cf9ccd311..00000000000
--- a/lib/libc/arch/alpha/gen/ldexp.c
+++ /dev/null
@@ -1,138 +0,0 @@
-/*	$OpenBSD: ldexp.c,v 1.6 2008/12/10 01:15:01 martynas Exp $	*/
-/*	$NetBSD: ldexp.c,v 1.1 1995/02/10 17:50:24 cgd Exp $	*/
-
-/*
- * Copyright (c) 1994, 1995 Carnegie-Mellon University.
- * All rights reserved.
- *
- * Author: Chris G. Demetriou
- * 
- * Permission to use, copy, modify and distribute this software and
- * its documentation is hereby granted, provided that both the copyright
- * notice and this permission notice appear in all copies of the
- * software, derivative works or modified versions, and any portions
- * thereof, and that both notices appear in supporting documentation.
- * 
- * CARNEGIE MELLON ALLOWS FREE USE OF THIS SOFTWARE IN ITS "AS IS" 
- * CONDITION.  CARNEGIE MELLON DISCLAIMS ANY LIABILITY OF ANY KIND 
- * FOR ANY DAMAGES WHATSOEVER RESULTING FROM THE USE OF THIS SOFTWARE.
- * 
- * Carnegie Mellon requests users of this software to return to
- *
- *  Software Distribution Coordinator  or  Software.Distribution@CS.CMU.EDU
- *  School of Computer Science
- *  Carnegie Mellon University
- *  Pittsburgh PA 15213-3890
- *
- * any improvements or extensions that they make and grant Carnegie the
- * rights to redistribute these changes.
- */
-
-#include <sys/types.h>
-#include <sys/cdefs.h>
-#include <machine/ieee.h>
-#include <errno.h>
-#include <math.h>
-
-/*
- * double ldexp(double val, int exp)
- * returns: val * (2**exp)
- */
-double
-ldexp(val, exp)
-	double val;
-	int exp;
-{
-	register int oldexp, newexp, mulexp;
-	union doub {
-		double v;
-		struct ieee_double s;
-	} u, mul;
-
-	/*
-	 * If input is zero, or no change, just return input.
-	 * Likewise, if input is Inf or NaN, just return it.
-	 */
-	u.v = val;
-	oldexp = u.s.dbl_exp;
-	if (val == 0 || exp == 0 || oldexp == DBL_EXP_INFNAN)
-		return (val);
-
-	/*
-	 * Compute new exponent and check for over/under flow.
-	 * Underflow, unfortunately, could mean switching to denormal.
-	 * If result out of range, set ERANGE and return 0 if too small
-	 * or Inf if too big, with the same sign as the input value.
-	 */
-	newexp = oldexp + exp;
-	if (newexp >= DBL_EXP_INFNAN) {
-		/* u.s.dbl_sign = val < 0; -- already set */
-		u.s.dbl_exp = DBL_EXP_INFNAN;
-		u.s.dbl_frach = u.s.dbl_fracl = 0;
-		errno = ERANGE;
-		return (u.v);		/* Inf */
-	}
-	if (newexp <= 0) {
-		/*
-		 * The output number is either a denormal or underflows
-		 * (see comments in machine/ieee.h).
-		 */
-		if (newexp <= -DBL_FRACBITS) {
-			/* u.s.dbl_sign = val < 0; -- already set */
-			u.s.dbl_exp = 0;
-			u.s.dbl_frach = u.s.dbl_fracl = 0;
-			errno = ERANGE;
-			return (u.v);		/* zero */
-		}
-		/*
-		 * We are going to produce a denorm.  Our `exp' argument
-		 * might be as small as -2097, and we cannot compute
-		 * 2^-2097, so we may have to do this as many as three
-		 * steps (not just two, as for positive `exp's below).
-		 */
-		mul.v = 0;
-		while (exp <= -DBL_EXP_BIAS) {
-			mul.s.dbl_exp = 1;
-			val *= mul.v;
-			exp += DBL_EXP_BIAS - 1;
-		}
-		mul.s.dbl_exp = exp + DBL_EXP_BIAS;
-		val *= mul.v;
-		return (val);
-	}
-
-	/*
-	 * Newexp is positive.
-	 *
-	 * If oldexp is zero, we are starting with a denorm, and simply
-	 * adjusting the exponent will produce bogus answers.  We need
-	 * to fix that first.
-	 */
-	if (oldexp == 0) {
-		/*
-		 * Multiply by 2^mulexp to make the number normalizable.
-		 * We cannot multiply by more than 2^1023, but `exp'
-		 * argument might be as large as 2046.  A single
-		 * adjustment, however, will normalize the number even
-		 * for huge `exp's, and then we can use exponent
-		 * arithmetic just as for normal `double's.
-		 */
-		mulexp = exp <= DBL_EXP_BIAS ? exp : DBL_EXP_BIAS;
-		mul.v = 0;
-		mul.s.dbl_exp = mulexp + DBL_EXP_BIAS;
-		val *= mul.v;
-		if (mulexp == exp)
-			return (val);
-		u.v = val;
-		newexp -= mulexp;
-	}
-
-	/*
-	 * Both oldexp and newexp are positive; just replace the
-	 * old exponent with the new one.
-	 */
-	u.s.dbl_exp = newexp;
-	return (u.v);
-}
-
-__weak_alias(ldexpl, ldexp);
diff --git a/lib/libc/arch/alpha/gen/modf.c b/lib/libc/arch/alpha/gen/modf.c
deleted file mode 100644
index 2315ef19034..00000000000
--- a/lib/libc/arch/alpha/gen/modf.c
+++ /dev/null
@@ -1,105 +0,0 @@
-/*	$OpenBSD: modf.c,v 1.5 2005/08/07 16:40:13 espie Exp $	*/
-/*	$NetBSD: modf.c,v 1.1 1995/02/10 17:50:25 cgd Exp $	*/
-
-/*
- * Copyright (c) 1994, 1995 Carnegie-Mellon University.
- * All rights reserved.
- *
- * Author: Chris G. Demetriou
- * 
- * Permission to use, copy, modify and distribute this software and
- * its documentation is hereby granted, provided that both the copyright
- * notice and this permission notice appear in all copies of the
- * software, derivative works or modified versions, and any portions
- * thereof, and that both notices appear in supporting documentation.
- * 
- * CARNEGIE MELLON ALLOWS FREE USE OF THIS SOFTWARE IN ITS "AS IS" 
- * CONDITION.  CARNEGIE MELLON DISCLAIMS ANY LIABILITY OF ANY KIND 
- * FOR ANY DAMAGES WHATSOEVER RESULTING FROM THE USE OF THIS SOFTWARE.
- * 
- * Carnegie Mellon requests users of this software to return to
- *
- *  Software Distribution Coordinator  or  Software.Distribution@CS.CMU.EDU
- *  School of Computer Science
- *  Carnegie Mellon University
- *  Pittsburgh PA 15213-3890
- *
- * any improvements or extensions that they make and grant Carnegie the
- * rights to redistribute these changes.
- */
-
-#include <sys/types.h>
-#include <machine/ieee.h>
-#include <errno.h>
-#include <math.h>
-
-/*
- * double modf(double val, double *iptr)
- * returns: f and i such that |f| < 1.0, (f + i) = val, and
- *	sign(f) == sign(i) == sign(val).
- *
- * Beware signedness when doing subtraction, and also operand size!
- */
-double
-modf(val, iptr)
-	double val, *iptr;
-{
-	union doub {
-		double v;
-		struct ieee_double s;
-	} u, v;
-	u_int64_t frac;
-
-	/*
-	 * If input is Inf or NaN, return it and leave i alone.
-	 */
-	u.v = val;
-	if (u.s.dbl_exp == DBL_EXP_INFNAN)
-		return (u.v);
-
-	/*
-	 * If input can't have a fractional part, return
-	 * (appropriately signed) zero, and make i be the input.
-	 */
-	if ((int)u.s.dbl_exp - DBL_EXP_BIAS > DBL_FRACBITS - 1) {
-		*iptr = u.v;
-		v.v = 0.0;
-		v.s.dbl_sign = u.s.dbl_sign;
-		return (v.v);
-	}
-
-	/*
-	 * If |input| < 1.0, return it, and set i to the appropriately
-	 * signed zero.
-	 */
-	if (u.s.dbl_exp < DBL_EXP_BIAS) {
-		v.v = 0.0;
-		v.s.dbl_sign = u.s.dbl_sign;
-		*iptr = v.v;
-		return (u.v);
-	}
-
-	/*
-	 * There can be a fractional part of the input.
-	 * If you look at the math involved for a few seconds, it's
-	 * plain to see that the integral part is the input, with the
-	 * low (DBL_FRACBITS - (exponent - DBL_EXP_BIAS)) bits zeroed,
-	 * the fractional part is the part with the rest of the
-	 * bits zeroed.  Just zeroing the high bits to get the
-	 * fractional part would yield a fraction in need of
-	 * normalization.  Therefore, we take the easy way out, and
-	 * just use subtraction to get the fractional part.
-	 */
-	v.v = u.v;
-	/* Zero the low bits of the fraction, the sleazy way. */
-	frac = ((u_int64_t)v.s.dbl_frach << 32) + v.s.dbl_fracl;
-	frac >>= DBL_FRACBITS - (u.s.dbl_exp - DBL_EXP_BIAS);
-	frac <<= DBL_FRACBITS - (u.s.dbl_exp - DBL_EXP_BIAS);
-	v.s.dbl_fracl = frac & 0xffffffff;
-	v.s.dbl_frach = frac >> 32;
-	*iptr = v.v;
-
-	u.v -= v.v;
-	u.s.dbl_sign = v.s.dbl_sign;
-	return (u.v);
-}
diff --git a/lib/libc/arch/amd64/gen/Makefile.inc b/lib/libc/arch/amd64/gen/Makefile.inc
index 2d8fa9e3c0c..77e7de0706f 100644
--- a/lib/libc/arch/amd64/gen/Makefile.inc
+++ b/lib/libc/arch/amd64/gen/Makefile.inc
@@ -1,6 +1,6 @@
-#	$OpenBSD: Makefile.inc,v 1.6 2008/12/09 19:52:33 martynas Exp $
+#	$OpenBSD: Makefile.inc,v 1.7 2009/04/19 16:42:05 martynas Exp $
 
-SRCS+=	_setjmp.S fabs.S infinity.c ldexp.c modf.S nan.c setjmp.S \
+SRCS+=	_setjmp.S fabs.S infinity.c modf.S nan.c setjmp.S \
 	sigsetjmp.S
 SRCS+=	fpclassifyl.c isfinitel.c isinfl.c isnanl.c isnormall.c signbitl.c
 SRCS+=	flt_rounds.S fpgetmask.S fpgetround.S fpgetsticky.S fpsetmask.S \
diff --git a/lib/libc/arch/arm/gen/Makefile.inc b/lib/libc/arch/arm/gen/Makefile.inc
index 9877978c8dd..7c5b71185eb 100644
--- a/lib/libc/arch/arm/gen/Makefile.inc
+++ b/lib/libc/arch/arm/gen/Makefile.inc
@@ -1,9 +1,9 @@
-# $OpenBSD: Makefile.inc,v 1.9 2008/07/24 09:31:06 martynas Exp $
+# $OpenBSD: Makefile.inc,v 1.10 2009/04/19 16:42:05 martynas Exp $
 # $NetBSD: Makefile.inc,v 1.6 2003/08/01 17:03:47 lukem Exp $
 
 SRCS+=	byte_swap_2.S byte_swap_4.S divsi3.S fabs.c flt_rounds.c infinity.c \
 	nan.c
-SRCS+=	setjmp.S _setjmp.S sigsetjmp.S modf.c ldexp.c
+SRCS+=	setjmp.S _setjmp.S sigsetjmp.S
 
 SRCS+=	alloca.S
 LSRCS+= alloca.c
diff --git a/lib/libc/arch/arm/gen/ldexp.c b/lib/libc/arch/arm/gen/ldexp.c
deleted file mode 100644
index 90e528024da..00000000000
--- a/lib/libc/arch/arm/gen/ldexp.c
+++ /dev/null
@@ -1,146 +0,0 @@
-/*	$OpenBSD: ldexp.c,v 1.4 2008/12/09 20:32:06 martynas Exp $	*/
-/*	$NetBSD: ldexp.c,v 1.2 2001/11/08 22:45:45 bjh21 Exp $	*/
-
-/*-
- * Copyright (c) 1999 The NetBSD Foundation, Inc.
- * All rights reserved.
- *
- * This code is derived from software contributed to The NetBSD Foundation
- * by Charles M. Hannum.
- *
- * Redistribution and use in source and binary forms, with or without
- * modification, are permitted provided that the following conditions
- * are met:
- * 1. Redistributions of source code must retain the above copyright
- *    notice, this list of conditions and the following disclaimer.
- * 2. Redistributions in binary form must reproduce the above copyright
- *    notice, this list of conditions and the following disclaimer in the
- *    documentation and/or other materials provided with the distribution.
- *
- * THIS SOFTWARE IS PROVIDED BY THE NETBSD FOUNDATION, INC. AND CONTRIBUTORS
- * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
- * TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
- * PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR CONTRIBUTORS
- * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
- * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
- * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
- * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
- * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
- * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
- * POSSIBILITY OF SUCH DAMAGE.
- */
-
-#include <sys/cdefs.h>
-#include <sys/types.h>
-#include <machine/ieee.h>
-#include <errno.h>
-#include <math.h>
-
-/*
- * Multiply the given value by 2^exponent.
- */
-double
-ldexp(val, expo)
-	double val;
-	int expo;
-{
-	register int oldexp, newexp;
-	union {
-		double v;
-		struct ieee_double s;
-	} u, mul;
-
-	u.v = val;
-	oldexp = u.s.dbl_exp;
-
-	/*
-	 * If input is zero, Inf or NaN, just return it.
-	 */
-	if (u.v == 0.0 || oldexp == DBL_EXP_INFNAN)
-		return (val);
-
-	if (oldexp == 0) {
-		/*
-		 * u.v is denormal.  We must adjust it so that the exponent
-		 * arithmetic below will work.
-		 */
-		if (expo <= DBL_EXP_BIAS) {
-			/*
-			 * Optimization: if the scaling can be done in a single
-			 * multiply, or underflows, just do it now.
-			 */
-			if (expo <= -DBL_FRACBITS) {
-				errno = ERANGE;
-				return (0.0);
-			}
-			mul.v = 0.0;
-			mul.s.dbl_exp = expo + DBL_EXP_BIAS;
-			u.v *= mul.v;
-			if (u.v == 0.0) {
-				errno = ERANGE;
-				return (0.0);
-			}
-			return (u.v);
-		} else {
-			/*
-			 * We know that expo is very large, and therefore the
-			 * result cannot be denormal (though it may be Inf).
-			 * Shift u.v by just enough to make it normal.
-			 */
-			mul.v = 0.0;
-			mul.s.dbl_exp = DBL_FRACBITS + DBL_EXP_BIAS;
-			u.v *= mul.v;
-			expo -= DBL_FRACBITS;
-			oldexp = u.s.dbl_exp;
-		}
-	}
-
-	/*
-	 * u.v is now normalized and oldexp has been adjusted if necessary.
-	 * Calculate the new exponent and check for underflow and overflow.
-	 */
-	newexp = oldexp + expo;
-
-	if (newexp <= 0) {
-		/*
-		 * The output number is either denormal or underflows (see
-		 * comments in machine/ieee.h).
-		 */
-		if (newexp <= -DBL_FRACBITS) {
-			errno = ERANGE;
-			return (0.0);
-		}
-		/*
-		 * Denormalize the result.  We do this with a multiply. If expo
-		 * is very large, it won't fit in a double, so we have to
-		 * adjust the exponent first.  This is safe because we know
-		 * that u.v is normal at this point.
-		 */
-		if (expo <= -DBL_EXP_BIAS) {
-			u.s.dbl_exp = 1;
-			expo += oldexp - 1;
-		}
-		mul.v = 0.0;
-		mul.s.dbl_exp = expo + DBL_EXP_BIAS;
-		u.v *= mul.v;
-		return (u.v);
-	} else if (newexp >= DBL_EXP_INFNAN) {
-		/*
-		 * The result overflowed; return +/-Inf.
-		 */
-		u.s.dbl_exp = DBL_EXP_INFNAN;
-		u.s.dbl_frach = 0;
-		u.s.dbl_fracl = 0;
-		errno = ERANGE;
-		return (u.v);
-	} else {
-		/*
-		 * The result is normal; just replace the old exponent with the
-		 * new one.
-		 */
-		u.s.dbl_exp = newexp;
-		return (u.v);
-	}
-}
-
-__weak_alias(ldexpl, ldexp);
diff --git a/lib/libc/arch/arm/gen/modf.c b/lib/libc/arch/arm/gen/modf.c
deleted file mode 100644
index e675c8a1921..00000000000
--- a/lib/libc/arch/arm/gen/modf.c
+++ /dev/null
@@ -1,105 +0,0 @@
-/*	$OpenBSD: modf.c,v 1.2 2005/08/07 16:40:14 espie Exp $	*/
-/*	$NetBSD: modf.c,v 1.1 1995/02/10 17:50:25 cgd Exp $	*/
-
-/*
- * Copyright (c) 1994, 1995 Carnegie-Mellon University.
- * All rights reserved.
- *
- * Author: Chris G. Demetriou
- * 
- * Permission to use, copy, modify and distribute this software and
- * its documentation is hereby granted, provided that both the copyright
- * notice and this permission notice appear in all copies of the
- * software, derivative works or modified versions, and any portions
- * thereof, and that both notices appear in supporting documentation.
- * 
- * CARNEGIE MELLON ALLOWS FREE USE OF THIS SOFTWARE IN ITS "AS IS" 
- * CONDITION.  CARNEGIE MELLON DISCLAIMS ANY LIABILITY OF ANY KIND 
- * FOR ANY DAMAGES WHATSOEVER RESULTING FROM THE USE OF THIS SOFTWARE.
- * 
- * Carnegie Mellon requests users of this software to return to
- *
- *  Software Distribution Coordinator  or  Software.Distribution@CS.CMU.EDU
- *  School of Computer Science
- *  Carnegie Mellon University
- *  Pittsburgh PA 15213-3890
- *
- * any improvements or extensions that they make and grant Carnegie the
- * rights to redistribute these changes.
- */
-
-#include <sys/types.h>
-#include <machine/ieee.h>
-#include <errno.h>
-#include <math.h>
-
-/*
- * double modf(double val, double *iptr)
- * returns: f and i such that |f| < 1.0, (f + i) = val, and
- *	sign(f) == sign(i) == sign(val).
- *
- * Beware signedness when doing subtraction, and also operand size!
- */
-double
-modf(val, iptr)
-	double val, *iptr;
-{
-	union doub {
-		double v;
-		struct ieee_double s;
-	} u, v;
-	u_int64_t frac;
-
-	/*
-	 * If input is Inf or NaN, return it and leave i alone.
-	 */
-	u.v = val;
-	if (u.s.dbl_exp == DBL_EXP_INFNAN)
-		return (u.v);
-
-	/*
-	 * If input can't have a fractional part, return
-	 * (appropriately signed) zero, and make i be the input.
-	 */
-	if ((int)u.s.dbl_exp - DBL_EXP_BIAS > DBL_FRACBITS - 1) {
-		*iptr = u.v;
-		v.v = 0.0;
-		v.s.dbl_sign = u.s.dbl_sign;
-		return (v.v);
-	}
-
-	/*
-	 * If |input| < 1.0, return it, and set i to the appropriately
-	 * signed zero.
-	 */
-	if (u.s.dbl_exp < DBL_EXP_BIAS) {
-		v.v = 0.0;
-		v.s.dbl_sign = u.s.dbl_sign;
-		*iptr = v.v;
-		return (u.v);
-	}
-
-	/*
-	 * There can be a fractional part of the input.
-	 * If you look at the math involved for a few seconds, it's
-	 * plain to see that the integral part is the input, with the
-	 * low (DBL_FRACBITS - (exponent - DBL_EXP_BIAS)) bits zeroed,
-	 * the fractional part is the part with the rest of the
-	 * bits zeroed.  Just zeroing the high bits to get the
-	 * fractional part would yield a fraction in need of
-	 * normalization.  Therefore, we take the easy way out, and
-	 * just use subtraction to get the fractional part.
-	 */
-	v.v = u.v;
-	/* Zero the low bits of the fraction, the sleazy way. */
-	frac = ((u_int64_t)v.s.dbl_frach << 32) + v.s.dbl_fracl;
-	frac >>= DBL_FRACBITS - (u.s.dbl_exp - DBL_EXP_BIAS);
-	frac <<= DBL_FRACBITS - (u.s.dbl_exp - DBL_EXP_BIAS);
-	v.s.dbl_fracl = frac & 0xffffffff;
-	v.s.dbl_frach = frac >> 32;
-	*iptr = v.v;
-
-	u.v -= v.v;
-	u.s.dbl_sign = v.s.dbl_sign;
-	return (u.v);
-}
diff --git a/lib/libc/arch/hppa/gen/Makefile.inc b/lib/libc/arch/hppa/gen/Makefile.inc
index dbe7648d951..28670ad29a4 100644
--- a/lib/libc/arch/hppa/gen/Makefile.inc
+++ b/lib/libc/arch/hppa/gen/Makefile.inc
@@ -1,10 +1,9 @@
-#	$OpenBSD: Makefile.inc,v 1.9 2008/07/24 09:31:06 martynas Exp $
+#	$OpenBSD: Makefile.inc,v 1.10 2009/04/19 16:42:05 martynas Exp $
 
 SRCS+=	setjmp.S
-SRCS+=	fabs.c ldexp.c
+SRCS+=	fabs.c
 SRCS+=	infinity.c nan.c setjmp.S
 SRCS+=	flt_rounds.c fpgetmask.c fpgetround.c fpgetsticky.c fpsetmask.c \
 	fpsetround.c fpsetsticky.c
-SRCS+=	modf.c
 
 SRCS+=	alloca.c
diff --git a/lib/libc/arch/hppa/gen/ldexp.c b/lib/libc/arch/hppa/gen/ldexp.c
deleted file mode 100644
index de38f23e1ca..00000000000
--- a/lib/libc/arch/hppa/gen/ldexp.c
+++ /dev/null
@@ -1,142 +0,0 @@
-/*	$OpenBSD: ldexp.c,v 1.5 2008/12/10 00:59:07 deraadt Exp $	*/
-/*
- * Copyright (c) 1992, 1993
- *	The Regents of the University of California.  All rights reserved.
- *
- * This software was developed by the Computer Systems Engineering group
- * at Lawrence Berkeley Laboratory under DARPA contract BG 91-66 and
- * contributed to Berkeley.
- *
- * Redistribution and use in source and binary forms, with or without
- * modification, are permitted provided that the following conditions
- * are met:
- * 1. Redistributions of source code must retain the above copyright
- *    notice, this list of conditions and the following disclaimer.
- * 2. Redistributions in binary form must reproduce the above copyright
- *    notice, this list of conditions and the following disclaimer in the
- *    documentation and/or other materials provided with the distribution.
- * 3. Neither the name of the University nor the names of its contributors
- *    may be used to endorse or promote products derived from this software
- *    without specific prior written permission.
- *
- * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
- * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
- * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
- * ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
- * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
- * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
- * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
- * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
- * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
- * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
- * SUCH DAMAGE.
- *
- * from: Header: ldexp.c,v 1.1 91/07/07 04:28:19 torek Exp
- */
-
-#include <sys/types.h>
-#include <machine/ieee.h>
-#include <errno.h>
-
-/*
- * double ldexp(double val, int exp)
- * returns: val * (2**exp)
- */
-double
-ldexp(val, exp)
-	double val;
-	int exp;
-{
-	register int oldexp, newexp, mulexp;
-	union doub {
-		double v;
-		struct ieee_double s;
-	} u, mul;
-
-	/*
-	 * If input is zero, or no change, just return input.
-	 * Likewise, if input is Inf or NaN, just return it.
-	 */
-	u.v = val;
-	oldexp = u.s.dbl_exp;
-	if (val == 0 || exp == 0 || oldexp == DBL_EXP_INFNAN)
-		return (val);
-
-	/*
-	 * Compute new exponent and check for over/under flow.
-	 * Underflow, unfortunately, could mean switching to denormal.
-	 * If result out of range, set ERANGE and return 0 if too small
-	 * or Inf if too big, with the same sign as the input value.
-	 */
-	newexp = oldexp + exp;
-	if (newexp >= DBL_EXP_INFNAN) {
-		/* u.s.dbl_sign = val < 0; -- already set */
-		u.s.dbl_exp = DBL_EXP_INFNAN;
-		u.s.dbl_frach = u.s.dbl_fracl = 0;
-		errno = ERANGE;
-		return (u.v);		/* Inf */
-	}
-	if (newexp <= 0) {
-		/*
-		 * The output number is either a denormal or underflows
-		 * (see comments in machine/ieee.h).
-		 */
-		if (newexp <= -DBL_FRACBITS) {
-			/* u.s.dbl_sign = val < 0; -- already set */
-			u.s.dbl_exp = 0;
-			u.s.dbl_frach = u.s.dbl_fracl = 0;
-			errno = ERANGE;
-			return (u.v);		/* zero */
-		}
-		/*
-		 * We are going to produce a denorm.  Our `exp' argument
-		 * might be as small as -2097, and we cannot compute
-		 * 2^-2097, so we may have to do this as many as three
-		 * steps (not just two, as for positive `exp's below).
-		 */
-		mul.v = 0;
-		while (exp <= -DBL_EXP_BIAS) {
-			mul.s.dbl_exp = 1;
-			val *= mul.v;
-			exp += DBL_EXP_BIAS - 1;
-		}
-		mul.s.dbl_exp = exp + DBL_EXP_BIAS;
-		val *= mul.v;
-		return (val);
-	}
-
-	/*
-	 * Newexp is positive.
-	 *
-	 * If oldexp is zero, we are starting with a denorm, and simply
-	 * adjusting the exponent will produce bogus answers.  We need
-	 * to fix that first.
-	 */
-	if (oldexp == 0) {
-		/*
-		 * Multiply by 2^mulexp to make the number normalizable.
-		 * We cannot multiply by more than 2^1023, but `exp'
-		 * argument might be as large as 2046.  A single
-		 * adjustment, however, will normalize the number even
-		 * for huge `exp's, and then we can use exponent
-		 * arithmetic just as for normal `double's.
-		 */
-		mulexp = exp <= DBL_EXP_BIAS ? exp : DBL_EXP_BIAS;
-		mul.v = 0;
-		mul.s.dbl_exp = mulexp + DBL_EXP_BIAS;
-		val *= mul.v;
-		if (mulexp == exp)
-			return (val);
-		u.v = val;
-		newexp -= mulexp;
-	}
-
-	/*
-	 * Both oldexp and newexp are positive; just replace the
-	 * old exponent with the new one.
-	 */
-	u.s.dbl_exp = newexp;
-	return (u.v);
-}
-
-__weak_alias(ldexpl, ldexp);
diff --git a/lib/libc/arch/hppa/gen/modf.c b/lib/libc/arch/hppa/gen/modf.c
deleted file mode 100644
index 424bf36ddf5..00000000000
--- a/lib/libc/arch/hppa/gen/modf.c
+++ /dev/null
@@ -1,305 +0,0 @@
-/*	$OpenBSD: modf.c,v 1.2 2005/08/07 16:40:14 espie Exp $	*/
-/* @(#)s_modf.c 5.1 93/09/24 */
-/*
- * ====================================================
- * 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.
- * ====================================================
- */
-
-/*
- * modf(double x, double *iptr) 
- * return fraction part of x, and return x's integral part in *iptr.
- * Method:
- *	Bit twiddling.
- *
- * Exception:
- *	No exception.
- */
-
-#include "math.h"
-
-/*
- * ====================================================
- * 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
- */
-
-#ifndef _MATH_PRIVATE_H_
-#define _MATH_PRIVATE_H_
-
-#include <sys/types.h>
-#include <machine/endian.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.  */
-
-/*
- * The arm32 port is little endian except for the FP word order which is
- * big endian.
- */
-
-#if (BYTE_ORDER == BIG_ENDIAN) || defined(arm32)
-
-typedef union 
-{
-  double value;
-  struct 
-  {
-    u_int32_t msw;
-    u_int32_t lsw;
-  } parts;
-} ieee_double_shape_type;
-
-#endif
-
-#if (BYTE_ORDER == LITTLE_ENDIAN) && !defined(arm32)
-
-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(double);			
-extern double __ieee754_acos(double);			
-extern double __ieee754_acosh(double);			
-extern double __ieee754_log(double);			
-extern double __ieee754_atanh(double);			
-extern double __ieee754_asin(double);			
-extern double __ieee754_atan2(double,double);			
-extern double __ieee754_exp(double);
-extern double __ieee754_cosh(double);
-extern double __ieee754_fmod(double,double);
-extern double __ieee754_pow(double,double);
-extern double __ieee754_lgamma_r(double,int *);
-extern double __ieee754_gamma_r(double,int *);
-extern double __ieee754_lgamma(double);
-extern double __ieee754_gamma(double);
-extern double __ieee754_log10(double);
-extern double __ieee754_sinh(double);
-extern double __ieee754_hypot(double,double);
-extern double __ieee754_j0(double);
-extern double __ieee754_j1(double);
-extern double __ieee754_y0(double);
-extern double __ieee754_y1(double);
-extern double __ieee754_jn(int,double);
-extern double __ieee754_yn(int,double);
-extern double __ieee754_remainder(double,double);
-extern int    __ieee754_rem_pio2(double,double*);
-extern double __ieee754_scalb(double,double);
-
-/* fdlibm kernel function */
-extern double __kernel_standard(double,double,int);	
-extern double __kernel_sin(double,double,int);
-extern double __kernel_cos(double,double);
-extern double __kernel_tan(double,double,int);
-extern int    __kernel_rem_pio2(double*,double*,int,int,int,const int*);
-
-
-/* ieee style elementary float functions */
-extern float __ieee754_sqrtf(float);			
-extern float __ieee754_acosf(float);			
-extern float __ieee754_acoshf(float);			
-extern float __ieee754_logf(float);			
-extern float __ieee754_atanhf(float);			
-extern float __ieee754_asinf(float);			
-extern float __ieee754_atan2f(float,float);			
-extern float __ieee754_expf(float);
-extern float __ieee754_coshf(float);
-extern float __ieee754_fmodf(float,float);
-extern float __ieee754_powf(float,float);
-extern float __ieee754_lgammaf_r(float,int *);
-extern float __ieee754_gammaf_r(float,int *);
-extern float __ieee754_lgammaf(float);
-extern float __ieee754_gammaf(float);
-extern float __ieee754_log10f(float);
-extern float __ieee754_sinhf(float);
-extern float __ieee754_hypotf(float,float);
-extern float __ieee754_j0f(float);
-extern float __ieee754_j1f(float);
-extern float __ieee754_y0f(float);
-extern float __ieee754_y1f(float);
-extern float __ieee754_jnf(int,float);
-extern float __ieee754_ynf(int,float);
-extern float __ieee754_remainderf(float,float);
-extern int   __ieee754_rem_pio2f(float,float*);
-extern float __ieee754_scalbf(float,float);
-
-/* float versions of fdlibm kernel functions */
-extern float __kernel_sinf(float,float,int);
-extern float __kernel_cosf(float,float);
-extern float __kernel_tanf(float,float,int);
-extern int   __kernel_rem_pio2f(float*,float*,int,int,int,const int*);
-
-#endif /* _MATH_PRIVATE_H_ */
-#ifdef __STDC__
-static const double one = 1.0;
-#else
-static double one = 1.0;
-#endif
-
-#ifdef __STDC__
-	double modf(double x, double *iptr)
-#else
-	double modf(x, iptr)
-	double x,*iptr;
-#endif
-{
-	int32_t i0,i1,j0;
-	u_int32_t i;
-	EXTRACT_WORDS(i0,i1,x);
-	j0 = ((i0>>20)&0x7ff)-0x3ff;	/* exponent of x */
-	if(j0<20) {			/* integer part in high x */
-	    if(j0<0) {			/* |x|<1 */
-		INSERT_WORDS(*iptr,i0&0x80000000,0);	/* *iptr = +-0 */
-		return x;
-	    } else {
-		i = (0x000fffff)>>j0;
-		if(((i0&i)|i1)==0) {		/* x is integral */
-		    u_int32_t high;
-		    *iptr = x;
-		    GET_HIGH_WORD(high,x);
-		    INSERT_WORDS(x,high&0x80000000,0);	/* return +-0 */
-		    return x;
-		} else {
-		    INSERT_WORDS(*iptr,i0&(~i),0);
-		    return x - *iptr;
-		}
-	    }
-	} else if (j0>51) {		/* no fraction part */
-	    u_int32_t high;
-	    *iptr = x*one;
-	    GET_HIGH_WORD(high,x);
-	    INSERT_WORDS(x,high&0x80000000,0);	/* return +-0 */
-	    return x;
-	} else {			/* fraction part in low x */
-	    i = ((u_int32_t)(0xffffffff))>>(j0-20);
-	    if((i1&i)==0) { 		/* x is integral */
-		u_int32_t high;
-		*iptr = x;
-		GET_HIGH_WORD(high,x);
-		INSERT_WORDS(x,high&0x80000000,0);	/* return +-0 */
-		return x;
-	    } else {
-		INSERT_WORDS(*iptr,i0,i1&(~i));
-		return x - *iptr;
-	    }
-	}
-}
diff --git a/lib/libc/arch/hppa64/gen/Makefile.inc b/lib/libc/arch/hppa64/gen/Makefile.inc
index da20439174e..7dd2ab91212 100644
--- a/lib/libc/arch/hppa64/gen/Makefile.inc
+++ b/lib/libc/arch/hppa64/gen/Makefile.inc
@@ -1,11 +1,10 @@
-#	$OpenBSD: Makefile.inc,v 1.5 2008/12/09 19:52:33 martynas Exp $
+#	$OpenBSD: Makefile.inc,v 1.6 2009/04/19 16:42:05 martynas Exp $
 
 SRCS+=	setjmp.S
-SRCS+=	fabs.c frexp.c ldexp.c
+SRCS+=	fabs.c
 SRCS+=	infinity.c nan.c setjmp.S
 SRCS+=	flt_rounds.c fpgetmask.c fpgetround.c fpgetsticky.c fpsetmask.c \
 	fpsetround.c fpsetsticky.c
 SRCS+=	fpclassifyl.c isfinitel.c isinfl.c isnanl.c isnormall.c signbitl.c
-SRCS+=	modf.c
 
 SRCS+=	alloca.c
diff --git a/lib/libc/arch/hppa64/gen/frexp.c b/lib/libc/arch/hppa64/gen/frexp.c
deleted file mode 100644
index d2bdcaf59d5..00000000000
--- a/lib/libc/arch/hppa64/gen/frexp.c
+++ /dev/null
@@ -1,70 +0,0 @@
-/*	$OpenBSD: frexp.c,v 1.2 2005/08/07 16:40:14 espie Exp $	*/
-
-/*
- * Copyright (c) 1992, 1993
- *	The Regents of the University of California.  All rights reserved.
- *
- * This software was developed by the Computer Systems Engineering group
- * at Lawrence Berkeley Laboratory under DARPA contract BG 91-66 and
- * contributed to Berkeley.
- *
- * Redistribution and use in source and binary forms, with or without
- * modification, are permitted provided that the following conditions
- * are met:
- * 1. Redistributions of source code must retain the above copyright
- *    notice, this list of conditions and the following disclaimer.
- * 2. Redistributions in binary form must reproduce the above copyright
- *    notice, this list of conditions and the following disclaimer in the
- *    documentation and/or other materials provided with the distribution.
- * 3. Neither the name of the University nor the names of its contributors
- *    may be used to endorse or promote products derived from this software
- *    without specific prior written permission.
- *
- * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
- * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
- * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
- * ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
- * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
- * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
- * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
- * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
- * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
- * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
- * SUCH DAMAGE.
- */
-
-#include <sys/types.h>
-#include <machine/ieee.h>
-
-/*
- * Split the given value into a fraction in the range [0.5, 1.0) and
- * an exponent, such that frac * (2^exp) == value.  If value is 0,
- * return 0.
- */
-double
-frexp(value, eptr)
-	double value;
-	int *eptr;
-{
-	union {
-		double v;
-		struct ieee_double s;
-	} u;
-
-	if (value) {
-		/*
-		 * Fractions in [0.5..1.0) have an exponent of 2^-1.
-		 * Leave Inf and NaN alone, however.
-		 * WHAT ABOUT DENORMS?
-		 */
-		u.v = value;
-		if (u.s.dbl_exp != DBL_EXP_INFNAN) {
-			*eptr = u.s.dbl_exp - (DBL_EXP_BIAS - 1);
-			u.s.dbl_exp = DBL_EXP_BIAS - 1;
-		}
-		return (u.v);
-	} else {
-		*eptr = 0;
-		return ((double)0);
-	}
-}
diff --git a/lib/libc/arch/hppa64/gen/ldexp.c b/lib/libc/arch/hppa64/gen/ldexp.c
deleted file mode 100644
index 35c425edd50..00000000000
--- a/lib/libc/arch/hppa64/gen/ldexp.c
+++ /dev/null
@@ -1,140 +0,0 @@
-/*	$OpenBSD: ldexp.c,v 1.2 2005/08/07 16:40:14 espie Exp $	*/
-/*
- * Copyright (c) 1992, 1993
- *	The Regents of the University of California.  All rights reserved.
- *
- * This software was developed by the Computer Systems Engineering group
- * at Lawrence Berkeley Laboratory under DARPA contract BG 91-66 and
- * contributed to Berkeley.
- *
- * Redistribution and use in source and binary forms, with or without
- * modification, are permitted provided that the following conditions
- * are met:
- * 1. Redistributions of source code must retain the above copyright
- *    notice, this list of conditions and the following disclaimer.
- * 2. Redistributions in binary form must reproduce the above copyright
- *    notice, this list of conditions and the following disclaimer in the
- *    documentation and/or other materials provided with the distribution.
- * 3. Neither the name of the University nor the names of its contributors
- *    may be used to endorse or promote products derived from this software
- *    without specific prior written permission.
- *
- * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
- * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
- * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
- * ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
- * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
- * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
- * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
- * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
- * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
- * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
- * SUCH DAMAGE.
- *
- * from: Header: ldexp.c,v 1.1 91/07/07 04:28:19 torek Exp
- */
-
-#include <sys/types.h>
-#include <machine/ieee.h>
-#include <errno.h>
-
-/*
- * double ldexp(double val, int exp)
- * returns: val * (2**exp)
- */
-double
-ldexp(val, exp)
-	double val;
-	int exp;
-{
-	register int oldexp, newexp, mulexp;
-	union doub {
-		double v;
-		struct ieee_double s;
-	} u, mul;
-
-	/*
-	 * If input is zero, or no change, just return input.
-	 * Likewise, if input is Inf or NaN, just return it.
-	 */
-	u.v = val;
-	oldexp = u.s.dbl_exp;
-	if (val == 0 || exp == 0 || oldexp == DBL_EXP_INFNAN)
-		return (val);
-
-	/*
-	 * Compute new exponent and check for over/under flow.
-	 * Underflow, unfortunately, could mean switching to denormal.
-	 * If result out of range, set ERANGE and return 0 if too small
-	 * or Inf if too big, with the same sign as the input value.
-	 */
-	newexp = oldexp + exp;
-	if (newexp >= DBL_EXP_INFNAN) {
-		/* u.s.dbl_sign = val < 0; -- already set */
-		u.s.dbl_exp = DBL_EXP_INFNAN;
-		u.s.dbl_frach = u.s.dbl_fracl = 0;
-		errno = ERANGE;
-		return (u.v);		/* Inf */
-	}
-	if (newexp <= 0) {
-		/*
-		 * The output number is either a denormal or underflows
-		 * (see comments in machine/ieee.h).
-		 */
-		if (newexp <= -DBL_FRACBITS) {
-			/* u.s.dbl_sign = val < 0; -- already set */
-			u.s.dbl_exp = 0;
-			u.s.dbl_frach = u.s.dbl_fracl = 0;
-			errno = ERANGE;
-			return (u.v);		/* zero */
-		}
-		/*
-		 * We are going to produce a denorm.  Our `exp' argument
-		 * might be as small as -2097, and we cannot compute
-		 * 2^-2097, so we may have to do this as many as three
-		 * steps (not just two, as for positive `exp's below).
-		 */
-		mul.v = 0;
-		while (exp <= -DBL_EXP_BIAS) {
-			mul.s.dbl_exp = 1;
-			val *= mul.v;
-			exp += DBL_EXP_BIAS - 1;
-		}
-		mul.s.dbl_exp = exp + DBL_EXP_BIAS;
-		val *= mul.v;
-		return (val);
-	}
-
-	/*
-	 * Newexp is positive.
-	 *
-	 * If oldexp is zero, we are starting with a denorm, and simply
-	 * adjusting the exponent will produce bogus answers.  We need
-	 * to fix that first.
-	 */
-	if (oldexp == 0) {
-		/*
-		 * Multiply by 2^mulexp to make the number normalizable.
-		 * We cannot multiply by more than 2^1023, but `exp'
-		 * argument might be as large as 2046.  A single
-		 * adjustment, however, will normalize the number even
-		 * for huge `exp's, and then we can use exponent
-		 * arithmetic just as for normal `double's.
-		 */
-		mulexp = exp <= DBL_EXP_BIAS ? exp : DBL_EXP_BIAS;
-		mul.v = 0;
-		mul.s.dbl_exp = mulexp + DBL_EXP_BIAS;
-		val *= mul.v;
-		if (mulexp == exp)
-			return (val);
-		u.v = val;
-		newexp -= mulexp;
-	}
-
-	/*
-	 * Both oldexp and newexp are positive; just replace the
-	 * old exponent with the new one.
-	 */
-	u.s.dbl_exp = newexp;
-	return (u.v);
-}
diff --git a/lib/libc/arch/hppa64/gen/modf.c b/lib/libc/arch/hppa64/gen/modf.c
deleted file mode 100644
index 424bf36ddf5..00000000000
--- a/lib/libc/arch/hppa64/gen/modf.c
+++ /dev/null
@@ -1,305 +0,0 @@
-/*	$OpenBSD: modf.c,v 1.2 2005/08/07 16:40:14 espie Exp $	*/
-/* @(#)s_modf.c 5.1 93/09/24 */
-/*
- * ====================================================
- * 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.
- * ====================================================
- */
-
-/*
- * modf(double x, double *iptr) 
- * return fraction part of x, and return x's integral part in *iptr.
- * Method:
- *	Bit twiddling.
- *
- * Exception:
- *	No exception.
- */
-
-#include "math.h"
-
-/*
- * ====================================================
- * 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
- */
-
-#ifndef _MATH_PRIVATE_H_
-#define _MATH_PRIVATE_H_
-
-#include <sys/types.h>
-#include <machine/endian.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.  */
-
-/*
- * The arm32 port is little endian except for the FP word order which is
- * big endian.
- */
-
-#if (BYTE_ORDER == BIG_ENDIAN) || defined(arm32)
-
-typedef union 
-{
-  double value;
-  struct 
-  {
-    u_int32_t msw;
-    u_int32_t lsw;
-  } parts;
-} ieee_double_shape_type;
-
-#endif
-
-#if (BYTE_ORDER == LITTLE_ENDIAN) && !defined(arm32)
-
-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(double);			
-extern double __ieee754_acos(double);			
-extern double __ieee754_acosh(double);			
-extern double __ieee754_log(double);			
-extern double __ieee754_atanh(double);			
-extern double __ieee754_asin(double);			
-extern double __ieee754_atan2(double,double);			
-extern double __ieee754_exp(double);
-extern double __ieee754_cosh(double);
-extern double __ieee754_fmod(double,double);
-extern double __ieee754_pow(double,double);
-extern double __ieee754_lgamma_r(double,int *);
-extern double __ieee754_gamma_r(double,int *);
-extern double __ieee754_lgamma(double);
-extern double __ieee754_gamma(double);
-extern double __ieee754_log10(double);
-extern double __ieee754_sinh(double);
-extern double __ieee754_hypot(double,double);
-extern double __ieee754_j0(double);
-extern double __ieee754_j1(double);
-extern double __ieee754_y0(double);
-extern double __ieee754_y1(double);
-extern double __ieee754_jn(int,double);
-extern double __ieee754_yn(int,double);
-extern double __ieee754_remainder(double,double);
-extern int    __ieee754_rem_pio2(double,double*);
-extern double __ieee754_scalb(double,double);
-
-/* fdlibm kernel function */
-extern double __kernel_standard(double,double,int);	
-extern double __kernel_sin(double,double,int);
-extern double __kernel_cos(double,double);
-extern double __kernel_tan(double,double,int);
-extern int    __kernel_rem_pio2(double*,double*,int,int,int,const int*);
-
-
-/* ieee style elementary float functions */
-extern float __ieee754_sqrtf(float);			
-extern float __ieee754_acosf(float);			
-extern float __ieee754_acoshf(float);			
-extern float __ieee754_logf(float);			
-extern float __ieee754_atanhf(float);			
-extern float __ieee754_asinf(float);			
-extern float __ieee754_atan2f(float,float);			
-extern float __ieee754_expf(float);
-extern float __ieee754_coshf(float);
-extern float __ieee754_fmodf(float,float);
-extern float __ieee754_powf(float,float);
-extern float __ieee754_lgammaf_r(float,int *);
-extern float __ieee754_gammaf_r(float,int *);
-extern float __ieee754_lgammaf(float);
-extern float __ieee754_gammaf(float);
-extern float __ieee754_log10f(float);
-extern float __ieee754_sinhf(float);
-extern float __ieee754_hypotf(float,float);
-extern float __ieee754_j0f(float);
-extern float __ieee754_j1f(float);
-extern float __ieee754_y0f(float);
-extern float __ieee754_y1f(float);
-extern float __ieee754_jnf(int,float);
-extern float __ieee754_ynf(int,float);
-extern float __ieee754_remainderf(float,float);
-extern int   __ieee754_rem_pio2f(float,float*);
-extern float __ieee754_scalbf(float,float);
-
-/* float versions of fdlibm kernel functions */
-extern float __kernel_sinf(float,float,int);
-extern float __kernel_cosf(float,float);
-extern float __kernel_tanf(float,float,int);
-extern int   __kernel_rem_pio2f(float*,float*,int,int,int,const int*);
-
-#endif /* _MATH_PRIVATE_H_ */
-#ifdef __STDC__
-static const double one = 1.0;
-#else
-static double one = 1.0;
-#endif
-
-#ifdef __STDC__
-	double modf(double x, double *iptr)
-#else
-	double modf(x, iptr)
-	double x,*iptr;
-#endif
-{
-	int32_t i0,i1,j0;
-	u_int32_t i;
-	EXTRACT_WORDS(i0,i1,x);
-	j0 = ((i0>>20)&0x7ff)-0x3ff;	/* exponent of x */
-	if(j0<20) {			/* integer part in high x */
-	    if(j0<0) {			/* |x|<1 */
-		INSERT_WORDS(*iptr,i0&0x80000000,0);	/* *iptr = +-0 */
-		return x;
-	    } else {
-		i = (0x000fffff)>>j0;
-		if(((i0&i)|i1)==0) {		/* x is integral */
-		    u_int32_t high;
-		    *iptr = x;
-		    GET_HIGH_WORD(high,x);
-		    INSERT_WORDS(x,high&0x80000000,0);	/* return +-0 */
-		    return x;
-		} else {
-		    INSERT_WORDS(*iptr,i0&(~i),0);
-		    return x - *iptr;
-		}
-	    }
-	} else if (j0>51) {		/* no fraction part */
-	    u_int32_t high;
-	    *iptr = x*one;
-	    GET_HIGH_WORD(high,x);
-	    INSERT_WORDS(x,high&0x80000000,0);	/* return +-0 */
-	    return x;
-	} else {			/* fraction part in low x */
-	    i = ((u_int32_t)(0xffffffff))>>(j0-20);
-	    if((i1&i)==0) { 		/* x is integral */
-		u_int32_t high;
-		*iptr = x;
-		GET_HIGH_WORD(high,x);
-		INSERT_WORDS(x,high&0x80000000,0);	/* return +-0 */
-		return x;
-	    } else {
-		INSERT_WORDS(*iptr,i0,i1&(~i));
-		return x - *iptr;
-	    }
-	}
-}
diff --git a/lib/libc/arch/i386/gen/Makefile.inc b/lib/libc/arch/i386/gen/Makefile.inc
index 9e4fdf7184d..69ecd98a6ad 100644
--- a/lib/libc/arch/i386/gen/Makefile.inc
+++ b/lib/libc/arch/i386/gen/Makefile.inc
@@ -1,6 +1,6 @@
-#	$OpenBSD: Makefile.inc,v 1.8 2008/12/09 19:52:33 martynas Exp $
+#	$OpenBSD: Makefile.inc,v 1.9 2009/04/19 16:42:05 martynas Exp $
 
-SRCS+=	_setjmp.S alloca.S fabs.S infinity.c ldexp.c \
+SRCS+=	_setjmp.S alloca.S fabs.S infinity.c \
 	modf.S nan.c setjmp.S sigsetjmp.S
 SRCS+=	fpclassifyl.c isfinitel.c isinfl.c isnanl.c isnormall.c signbitl.c
 SRCS+=	flt_rounds.S fpgetmask.S fpgetround.S fpgetsticky.S fpsetmask.S \
diff --git a/lib/libc/arch/m88k/gen/Makefile.inc b/lib/libc/arch/m88k/gen/Makefile.inc
index 255c6918e7a..475a3faa124 100644
--- a/lib/libc/arch/m88k/gen/Makefile.inc
+++ b/lib/libc/arch/m88k/gen/Makefile.inc
@@ -1,16 +1,15 @@
-#	$OpenBSD: Makefile.inc,v 1.8 2008/12/09 19:52:33 martynas Exp $
+#	$OpenBSD: Makefile.inc,v 1.9 2009/04/19 16:42:05 martynas Exp $
 #	$NetBSD: Makefile.inc,v 1.3 1995/04/10 21:09:06 jtc Exp $
 
-#SRCS+=	_setjmp.S fabs.S infinity.c ldexp.c modf.S nan.c
+#SRCS+=	_setjmp.S fabs.S infinity.c modf.S nan.c
 #SRCS+=	flt_rounds.c fpgetmask.c fpgetround.c fpgetsticky.c fpsetmask.c \
 #	fpsetround.c fpsetsticky.c
 #SRCS+=	fixunsdfsi.S mul.S umul.S saveregs.S setjmp.S sigsetjmp.S
 
-SRCS+=	_setjmp.S fabs.S frexp.c infinity.c ldexp.c nan.c
+SRCS+=	_setjmp.S fabs.S infinity.c nan.c
 SRCS+=	flt_rounds.c fpgetmask.c fpgetround.c fpgetsticky.c fpsetmask.c \
 	fpsetround.c fpsetsticky.c
 SRCS+=	fpclassifyl.c isfinitel.c isinfl.c isnanl.c isnormall.c signbitl.c
 SRCS+=	setjmp.S sigsetjmp.S
-SRCS+=	modf.c
 
 SRCS+=	alloca.c
diff --git a/lib/libc/arch/m88k/gen/ldexp.c b/lib/libc/arch/m88k/gen/ldexp.c
deleted file mode 100644
index c216a416497..00000000000
--- a/lib/libc/arch/m88k/gen/ldexp.c
+++ /dev/null
@@ -1,140 +0,0 @@
-/*	$OpenBSD: ldexp.c,v 1.5 2005/08/07 16:40:14 espie Exp $	*/
-/*
- * Copyright (c) 1992, 1993
- *	The Regents of the University of California.  All rights reserved.
- *
- * This software was developed by the Computer Systems Engineering group
- * at Lawrence Berkeley Laboratory under DARPA contract BG 91-66 and
- * contributed to Berkeley.
- *
- * Redistribution and use in source and binary forms, with or without
- * modification, are permitted provided that the following conditions
- * are met:
- * 1. Redistributions of source code must retain the above copyright
- *    notice, this list of conditions and the following disclaimer.
- * 2. Redistributions in binary form must reproduce the above copyright
- *    notice, this list of conditions and the following disclaimer in the
- *    documentation and/or other materials provided with the distribution.
- * 3. Neither the name of the University nor the names of its contributors
- *    may be used to endorse or promote products derived from this software
- *    without specific prior written permission.
- *
- * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
- * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
- * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
- * ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
- * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
- * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
- * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
- * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
- * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
- * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
- * SUCH DAMAGE.
- *
- * from: Header: ldexp.c,v 1.1 91/07/07 04:28:19 torek Exp
- */
-
-#include <sys/types.h>
-#include <machine/ieee.h>
-#include <errno.h>
-
-/*
- * double ldexp(double val, int exp)
- * returns: val * (2**exp)
- */
-double
-ldexp(val, exp)
-	double val;
-	int exp;
-{
-	register int oldexp, newexp, mulexp;
-	union doub {
-		double v;
-		struct ieee_double s;
-	} u, mul;
-
-	/*
-	 * If input is zero, or no change, just return input.
-	 * Likewise, if input is Inf or NaN, just return it.
-	 */
-	u.v = val;
-	oldexp = u.s.dbl_exp;
-	if (val == 0 || exp == 0 || oldexp == DBL_EXP_INFNAN)
-		return (val);
-
-	/*
-	 * Compute new exponent and check for over/under flow.
-	 * Underflow, unfortunately, could mean switching to denormal.
-	 * If result out of range, set ERANGE and return 0 if too small
-	 * or Inf if too big, with the same sign as the input value.
-	 */
-	newexp = oldexp + exp;
-	if (newexp >= DBL_EXP_INFNAN) {
-		/* u.s.dbl_sign = val < 0; -- already set */
-		u.s.dbl_exp = DBL_EXP_INFNAN;
-		u.s.dbl_frach = u.s.dbl_fracl = 0;
-		errno = ERANGE;
-		return (u.v);		/* Inf */
-	}
-	if (newexp <= 0) {
-		/*
-		 * The output number is either a denormal or underflows
-		 * (see comments in machine/ieee.h).
-		 */
-		if (newexp <= -DBL_FRACBITS) {
-			/* u.s.dbl_sign = val < 0; -- already set */
-			u.s.dbl_exp = 0;
-			u.s.dbl_frach = u.s.dbl_fracl = 0;
-			errno = ERANGE;
-			return (u.v);		/* zero */
-		}
-		/*
-		 * We are going to produce a denorm.  Our `exp' argument
-		 * might be as small as -2097, and we cannot compute
-		 * 2^-2097, so we may have to do this as many as three
-		 * steps (not just two, as for positive `exp's below).
-		 */
-		mul.v = 0;
-		while (exp <= -DBL_EXP_BIAS) {
-			mul.s.dbl_exp = 1;
-			val *= mul.v;
-			exp += DBL_EXP_BIAS - 1;
-		}
-		mul.s.dbl_exp = exp + DBL_EXP_BIAS;
-		val *= mul.v;
-		return (val);
-	}
-
-	/*
-	 * Newexp is positive.
-	 *
-	 * If oldexp is zero, we are starting with a denorm, and simply
-	 * adjusting the exponent will produce bogus answers.  We need
-	 * to fix that first.
-	 */
-	if (oldexp == 0) {
-		/*
-		 * Multiply by 2^mulexp to make the number normalizable.
-		 * We cannot multiply by more than 2^1023, but `exp'
-		 * argument might be as large as 2046.  A single
-		 * adjustment, however, will normalize the number even
-		 * for huge `exp's, and then we can use exponent
-		 * arithmetic just as for normal `double's.
-		 */
-		mulexp = exp <= DBL_EXP_BIAS ? exp : DBL_EXP_BIAS;
-		mul.v = 0;
-		mul.s.dbl_exp = mulexp + DBL_EXP_BIAS;
-		val *= mul.v;
-		if (mulexp == exp)
-			return (val);
-		u.v = val;
-		newexp -= mulexp;
-	}
-
-	/*
-	 * Both oldexp and newexp are positive; just replace the
-	 * old exponent with the new one.
-	 */
-	u.s.dbl_exp = newexp;
-	return (u.v);
-}
diff --git a/lib/libc/arch/m88k/gen/modf.c b/lib/libc/arch/m88k/gen/modf.c
deleted file mode 100644
index 102451e870f..00000000000
--- a/lib/libc/arch/m88k/gen/modf.c
+++ /dev/null
@@ -1,305 +0,0 @@
-/*	$OpenBSD: modf.c,v 1.7 2005/08/07 16:40:14 espie Exp $	*/
-/* @(#)s_modf.c 5.1 93/09/24 */
-/*
- * ====================================================
- * 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.
- * ====================================================
- */
-
-/*
- * modf(double x, double *iptr) 
- * return fraction part of x, and return x's integral part in *iptr.
- * Method:
- *	Bit twiddling.
- *
- * Exception:
- *	No exception.
- */
-
-#include "math.h"
-
-/*
- * ====================================================
- * 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
- */
-
-#ifndef _MATH_PRIVATE_H_
-#define _MATH_PRIVATE_H_
-
-#include <sys/types.h>
-#include <machine/endian.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.  */
-
-/*
- * The arm32 port is little endian except for the FP word order which is
- * big endian.
- */
-
-#if (BYTE_ORDER == BIG_ENDIAN) || defined(arm32)
-
-typedef union 
-{
-  double value;
-  struct 
-  {
-    u_int32_t msw;
-    u_int32_t lsw;
-  } parts;
-} ieee_double_shape_type;
-
-#endif
-
-#if (BYTE_ORDER == LITTLE_ENDIAN) && !defined(arm32)
-
-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(double);			
-extern double __ieee754_acos(double);			
-extern double __ieee754_acosh(double);			
-extern double __ieee754_log(double);			
-extern double __ieee754_atanh(double);			
-extern double __ieee754_asin(double);			
-extern double __ieee754_atan2(double,double);			
-extern double __ieee754_exp(double);
-extern double __ieee754_cosh(double);
-extern double __ieee754_fmod(double,double);
-extern double __ieee754_pow(double,double);
-extern double __ieee754_lgamma_r(double,int *);
-extern double __ieee754_gamma_r(double,int *);
-extern double __ieee754_lgamma(double);
-extern double __ieee754_gamma(double);
-extern double __ieee754_log10(double);
-extern double __ieee754_sinh(double);
-extern double __ieee754_hypot(double,double);
-extern double __ieee754_j0(double);
-extern double __ieee754_j1(double);
-extern double __ieee754_y0(double);
-extern double __ieee754_y1(double);
-extern double __ieee754_jn(int,double);
-extern double __ieee754_yn(int,double);
-extern double __ieee754_remainder(double,double);
-extern int    __ieee754_rem_pio2(double,double*);
-extern double __ieee754_scalb(double,double);
-
-/* fdlibm kernel function */
-extern double __kernel_standard(double,double,int);	
-extern double __kernel_sin(double,double,int);
-extern double __kernel_cos(double,double);
-extern double __kernel_tan(double,double,int);
-extern int    __kernel_rem_pio2(double*,double*,int,int,int,const int*);
-
-
-/* ieee style elementary float functions */
-extern float __ieee754_sqrtf(float);			
-extern float __ieee754_acosf(float);			
-extern float __ieee754_acoshf(float);			
-extern float __ieee754_logf(float);			
-extern float __ieee754_atanhf(float);			
-extern float __ieee754_asinf(float);			
-extern float __ieee754_atan2f(float,float);			
-extern float __ieee754_expf(float);
-extern float __ieee754_coshf(float);
-extern float __ieee754_fmodf(float,float);
-extern float __ieee754_powf(float,float);
-extern float __ieee754_lgammaf_r(float,int *);
-extern float __ieee754_gammaf_r(float,int *);
-extern float __ieee754_lgammaf(float);
-extern float __ieee754_gammaf(float);
-extern float __ieee754_log10f(float);
-extern float __ieee754_sinhf(float);
-extern float __ieee754_hypotf(float,float);
-extern float __ieee754_j0f(float);
-extern float __ieee754_j1f(float);
-extern float __ieee754_y0f(float);
-extern float __ieee754_y1f(float);
-extern float __ieee754_jnf(int,float);
-extern float __ieee754_ynf(int,float);
-extern float __ieee754_remainderf(float,float);
-extern int   __ieee754_rem_pio2f(float,float*);
-extern float __ieee754_scalbf(float,float);
-
-/* float versions of fdlibm kernel functions */
-extern float __kernel_sinf(float,float,int);
-extern float __kernel_cosf(float,float);
-extern float __kernel_tanf(float,float,int);
-extern int   __kernel_rem_pio2f(float*,float*,int,int,int,const int*);
-
-#endif /* _MATH_PRIVATE_H_ */
-#ifdef __STDC__
-static const double one = 1.0;
-#else
-static double one = 1.0;
-#endif
-
-#ifdef __STDC__
-	double modf(double x, double *iptr)
-#else
-	double modf(x, iptr)
-	double x,*iptr;
-#endif
-{
-	int32_t i0,i1,j0;
-	u_int32_t i;
-	EXTRACT_WORDS(i0,i1,x);
-	j0 = ((i0>>20)&0x7ff)-0x3ff;	/* exponent of x */
-	if(j0<20) {			/* integer part in high x */
-	    if(j0<0) {			/* |x|<1 */
-		INSERT_WORDS(*iptr,i0&0x80000000,0);	/* *iptr = +-0 */
-		return x;
-	    } else {
-		i = (0x000fffff)>>j0;
-		if(((i0&i)|i1)==0) {		/* x is integral */
-		    u_int32_t high;
-		    *iptr = x;
-		    GET_HIGH_WORD(high,x);
-		    INSERT_WORDS(x,high&0x80000000,0);	/* return +-0 */
-		    return x;
-		} else {
-		    INSERT_WORDS(*iptr,i0&(~i),0);
-		    return x - *iptr;
-		}
-	    }
-	} else if (j0>51) {		/* no fraction part */
-	    u_int32_t high;
-	    *iptr = x*one;
-	    GET_HIGH_WORD(high,x);
-	    INSERT_WORDS(x,high&0x80000000,0);	/* return +-0 */
-	    return x;
-	} else {			/* fraction part in low x */
-	    i = ((u_int32_t)(0xffffffff))>>(j0-20);
-	    if((i1&i)==0) { 		/* x is integral */
-		u_int32_t high;
-		*iptr = x;
-		GET_HIGH_WORD(high,x);
-		INSERT_WORDS(x,high&0x80000000,0);	/* return +-0 */
-		return x;
-	    } else {
-		INSERT_WORDS(*iptr,i0,i1&(~i));
-		return x - *iptr;
-	    }
-	}
-}
diff --git a/lib/libc/arch/powerpc/gen/Makefile.inc b/lib/libc/arch/powerpc/gen/Makefile.inc
index 183826e016c..c00548b79d9 100644
--- a/lib/libc/arch/powerpc/gen/Makefile.inc
+++ b/lib/libc/arch/powerpc/gen/Makefile.inc
@@ -1,5 +1,5 @@
-SRCS+= infinity.c setjmp.S sigsetjmp.S flt_rounds.c modf.c nan.c
-SRCS+= ldexp.c fabs.c
+SRCS+= infinity.c setjmp.S sigsetjmp.S flt_rounds.c nan.c
+SRCS+= fabs.c
 SRCS+= fpgetmask.c fpsetmask.c
 SRCS+= fpgetround.c fpsetround.c
 SRCS+= fpgetsticky.c fpsetsticky.c
diff --git a/lib/libc/arch/powerpc/gen/ldexp.c b/lib/libc/arch/powerpc/gen/ldexp.c
deleted file mode 100644
index 0a7bf70558a..00000000000
--- a/lib/libc/arch/powerpc/gen/ldexp.c
+++ /dev/null
@@ -1,143 +0,0 @@
-/*	$OpenBSD: ldexp.c,v 1.6 2008/12/10 01:15:02 martynas Exp $	*/
-/*
- * Copyright (c) 1992, 1993
- *	The Regents of the University of California.  All rights reserved.
- *
- * This software was developed by the Computer Systems Engineering group
- * at Lawrence Berkeley Laboratory under DARPA contract BG 91-66 and
- * contributed to Berkeley.
- *
- * Redistribution and use in source and binary forms, with or without
- * modification, are permitted provided that the following conditions
- * are met:
- * 1. Redistributions of source code must retain the above copyright
- *    notice, this list of conditions and the following disclaimer.
- * 2. Redistributions in binary form must reproduce the above copyright
- *    notice, this list of conditions and the following disclaimer in the
- *    documentation and/or other materials provided with the distribution.
- * 3. Neither the name of the University nor the names of its contributors
- *    may be used to endorse or promote products derived from this software
- *    without specific prior written permission.
- *
- * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
- * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
- * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
- * ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
- * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
- * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
- * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
- * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
- * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
- * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
- * SUCH DAMAGE.
- *
- * from: Header: ldexp.c,v 1.1 91/07/07 04:28:19 torek Exp
- */
-
-#include <sys/types.h>
-#include <sys/cdefs.h>
-#include <machine/ieee.h>
-#include <errno.h>
-
-/*
- * double ldexp(double val, int exp)
- * returns: val * (2**exp)
- */
-double
-ldexp(val, exp)
-	double val;
-	int exp;
-{
-	register int oldexp, newexp, mulexp;
-	union doub {
-		double v;
-		struct ieee_double s;
-	} u, mul;
-
-	/*
-	 * If input is zero, or no change, just return input.
-	 * Likewise, if input is Inf or NaN, just return it.
-	 */
-	u.v = val;
-	oldexp = u.s.dbl_exp;
-	if (val == 0 || exp == 0 || oldexp == DBL_EXP_INFNAN)
-		return (val);
-
-	/*
-	 * Compute new exponent and check for over/under flow.
-	 * Underflow, unfortunately, could mean switching to denormal.
-	 * If result out of range, set ERANGE and return 0 if too small
-	 * or Inf if too big, with the same sign as the input value.
-	 */
-	newexp = oldexp + exp;
-	if (newexp >= DBL_EXP_INFNAN) {
-		/* u.s.dbl_sign = val < 0; -- already set */
-		u.s.dbl_exp = DBL_EXP_INFNAN;
-		u.s.dbl_frach = u.s.dbl_fracl = 0;
-		errno = ERANGE;
-		return (u.v);		/* Inf */
-	}
-	if (newexp <= 0) {
-		/*
-		 * The output number is either a denormal or underflows
-		 * (see comments in machine/ieee.h).
-		 */
-		if (newexp <= -DBL_FRACBITS) {
-			/* u.s.dbl_sign = val < 0; -- already set */
-			u.s.dbl_exp = 0;
-			u.s.dbl_frach = u.s.dbl_fracl = 0;
-			errno = ERANGE;
-			return (u.v);		/* zero */
-		}
-		/*
-		 * We are going to produce a denorm.  Our `exp' argument
-		 * might be as small as -2097, and we cannot compute
-		 * 2^-2097, so we may have to do this as many as three
-		 * steps (not just two, as for positive `exp's below).
-		 */
-		mul.v = 0;
-		while (exp <= -DBL_EXP_BIAS) {
-			mul.s.dbl_exp = 1;
-			val *= mul.v;
-			exp += DBL_EXP_BIAS - 1;
-		}
-		mul.s.dbl_exp = exp + DBL_EXP_BIAS;
-		val *= mul.v;
-		return (val);
-	}
-
-	/*
-	 * Newexp is positive.
-	 *
-	 * If oldexp is zero, we are starting with a denorm, and simply
-	 * adjusting the exponent will produce bogus answers.  We need
-	 * to fix that first.
-	 */
-	if (oldexp == 0) {
-		/*
-		 * Multiply by 2^mulexp to make the number normalizable.
-		 * We cannot multiply by more than 2^1023, but `exp'
-		 * argument might be as large as 2046.  A single
-		 * adjustment, however, will normalize the number even
-		 * for huge `exp's, and then we can use exponent
-		 * arithmetic just as for normal `double's.
-		 */
-		mulexp = exp <= DBL_EXP_BIAS ? exp : DBL_EXP_BIAS;
-		mul.v = 0;
-		mul.s.dbl_exp = mulexp + DBL_EXP_BIAS;
-		val *= mul.v;
-		if (mulexp == exp)
-			return (val);
-		u.v = val;
-		newexp -= mulexp;
-	}
-
-	/*
-	 * Both oldexp and newexp are positive; just replace the
-	 * old exponent with the new one.
-	 */
-	u.s.dbl_exp = newexp;
-	return (u.v);
-}
-
-__weak_alias(ldexpl, ldexp);
diff --git a/lib/libc/arch/powerpc/gen/modf.c b/lib/libc/arch/powerpc/gen/modf.c
deleted file mode 100644
index e0e748864d3..00000000000
--- a/lib/libc/arch/powerpc/gen/modf.c
+++ /dev/null
@@ -1,305 +0,0 @@
-/*	$OpenBSD: modf.c,v 1.8 2005/08/07 16:40:15 espie Exp $	*/
-/* @(#)s_modf.c 5.1 93/09/24 */
-/*
- * ====================================================
- * 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.
- * ====================================================
- */
-
-/*
- * modf(double x, double *iptr) 
- * return fraction part of x, and return x's integral part in *iptr.
- * Method:
- *	Bit twiddling.
- *
- * Exception:
- *	No exception.
- */
-
-#include "math.h"
-
-/*
- * ====================================================
- * 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
- */
-
-#ifndef _MATH_PRIVATE_H_
-#define _MATH_PRIVATE_H_
-
-#include <sys/types.h>
-#include <machine/endian.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.  */
-
-/*
- * The arm32 port is little endian except for the FP word order which is
- * big endian.
- */
-
-#if (BYTE_ORDER == BIG_ENDIAN) || defined(arm32)
-
-typedef union 
-{
-  double value;
-  struct 
-  {
-    u_int32_t msw;
-    u_int32_t lsw;
-  } parts;
-} ieee_double_shape_type;
-
-#endif
-
-#if (BYTE_ORDER == LITTLE_ENDIAN) && !defined(arm32)
-
-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(double);			
-extern double __ieee754_acos(double);			
-extern double __ieee754_acosh(double);			
-extern double __ieee754_log(double);			
-extern double __ieee754_atanh(double);			
-extern double __ieee754_asin(double);			
-extern double __ieee754_atan2(double,double);			
-extern double __ieee754_exp(double);
-extern double __ieee754_cosh(double);
-extern double __ieee754_fmod(double,double);
-extern double __ieee754_pow(double,double);
-extern double __ieee754_lgamma_r(double,int *);
-extern double __ieee754_gamma_r(double,int *);
-extern double __ieee754_lgamma(double);
-extern double __ieee754_gamma(double);
-extern double __ieee754_log10(double);
-extern double __ieee754_sinh(double);
-extern double __ieee754_hypot(double,double);
-extern double __ieee754_j0(double);
-extern double __ieee754_j1(double);
-extern double __ieee754_y0(double);
-extern double __ieee754_y1(double);
-extern double __ieee754_jn(int,double);
-extern double __ieee754_yn(int,double);
-extern double __ieee754_remainder(double,double);
-extern int    __ieee754_rem_pio2(double,double*);
-extern double __ieee754_scalb(double,double);
-
-/* fdlibm kernel function */
-extern double __kernel_standard(double,double,int);	
-extern double __kernel_sin(double,double,int);
-extern double __kernel_cos(double,double);
-extern double __kernel_tan(double,double,int);
-extern int    __kernel_rem_pio2(double*,double*,int,int,int,const int*);
-
-
-/* ieee style elementary float functions */
-extern float __ieee754_sqrtf(float);			
-extern float __ieee754_acosf(float);			
-extern float __ieee754_acoshf(float);			
-extern float __ieee754_logf(float);			
-extern float __ieee754_atanhf(float);			
-extern float __ieee754_asinf(float);			
-extern float __ieee754_atan2f(float,float);			
-extern float __ieee754_expf(float);
-extern float __ieee754_coshf(float);
-extern float __ieee754_fmodf(float,float);
-extern float __ieee754_powf(float,float);
-extern float __ieee754_lgammaf_r(float,int *);
-extern float __ieee754_gammaf_r(float,int *);
-extern float __ieee754_lgammaf(float);
-extern float __ieee754_gammaf(float);
-extern float __ieee754_log10f(float);
-extern float __ieee754_sinhf(float);
-extern float __ieee754_hypotf(float,float);
-extern float __ieee754_j0f(float);
-extern float __ieee754_j1f(float);
-extern float __ieee754_y0f(float);
-extern float __ieee754_y1f(float);
-extern float __ieee754_jnf(int,float);
-extern float __ieee754_ynf(int,float);
-extern float __ieee754_remainderf(float,float);
-extern int   __ieee754_rem_pio2f(float,float*);
-extern float __ieee754_scalbf(float,float);
-
-/* float versions of fdlibm kernel functions */
-extern float __kernel_sinf(float,float,int);
-extern float __kernel_cosf(float,float);
-extern float __kernel_tanf(float,float,int);
-extern int   __kernel_rem_pio2f(float*,float*,int,int,int,const int*);
-
-#endif /* _MATH_PRIVATE_H_ */
-#ifdef __STDC__
-static const double one = 1.0;
-#else
-static double one = 1.0;
-#endif
-
-#ifdef __STDC__
-	double modf(double x, double *iptr)
-#else
-	double modf(x, iptr)
-	double x,*iptr;
-#endif
-{
-	int32_t i0,i1,j0;
-	u_int32_t i;
-	EXTRACT_WORDS(i0,i1,x);
-	j0 = ((i0>>20)&0x7ff)-0x3ff;	/* exponent of x */
-	if(j0<20) {			/* integer part in high x */
-	    if(j0<0) {			/* |x|<1 */
-		INSERT_WORDS(*iptr,i0&0x80000000,0);	/* *iptr = +-0 */
-		return x;
-	    } else {
-		i = (0x000fffff)>>j0;
-		if(((i0&i)|i1)==0) {		/* x is integral */
-		    u_int32_t high;
-		    *iptr = x;
-		    GET_HIGH_WORD(high,x);
-		    INSERT_WORDS(x,high&0x80000000,0);	/* return +-0 */
-		    return x;
-		} else {
-		    INSERT_WORDS(*iptr,i0&(~i),0);
-		    return x - *iptr;
-		}
-	    }
-	} else if (j0>51) {		/* no fraction part */
-	    u_int32_t high;
-	    *iptr = x*one;
-	    GET_HIGH_WORD(high,x);
-	    INSERT_WORDS(x,high&0x80000000,0);	/* return +-0 */
-	    return x;
-	} else {			/* fraction part in low x */
-	    i = ((u_int32_t)(0xffffffff))>>(j0-20);
-	    if((i1&i)==0) { 		/* x is integral */
-		u_int32_t high;
-		*iptr = x;
-		GET_HIGH_WORD(high,x);
-		INSERT_WORDS(x,high&0x80000000,0);	/* return +-0 */
-		return x;
-	    } else {
-		INSERT_WORDS(*iptr,i0,i1&(~i));
-		return x - *iptr;
-	    }
-	}
-}
diff --git a/lib/libc/arch/sh/gen/Makefile.inc b/lib/libc/arch/sh/gen/Makefile.inc
index 393c2fa35c1..386a3a47b4b 100644
--- a/lib/libc/arch/sh/gen/Makefile.inc
+++ b/lib/libc/arch/sh/gen/Makefile.inc
@@ -1,7 +1,5 @@
-#	$OpenBSD: Makefile.inc,v 1.6 2008/07/24 09:31:06 martynas Exp $
-
-SRCS+= flt_rounds.c infinity.c nan.c setjmp.S _setjmp.S sigsetjmp.S \
-	modf.c ldexp.c
+#	$OpenBSD: Makefile.inc,v 1.7 2009/04/19 16:42:06 martynas Exp $
 
+SRCS+= flt_rounds.c infinity.c nan.c setjmp.S _setjmp.S sigsetjmp.S
 SRCS+=	fabs.c fpgetmask.c fpgetround.c fpgetsticky.c \
 	fpsetmask.c fpsetround.c fpsetsticky.c
diff --git a/lib/libc/arch/sh/gen/ldexp.c b/lib/libc/arch/sh/gen/ldexp.c
deleted file mode 100644
index 574a3aae0b7..00000000000
--- a/lib/libc/arch/sh/gen/ldexp.c
+++ /dev/null
@@ -1,143 +0,0 @@
-/*	$OpenBSD: ldexp.c,v 1.3 2008/12/10 01:15:02 martynas Exp $	*/
-/*
- * Copyright (c) 1992, 1993
- *	The Regents of the University of California.  All rights reserved.
- *
- * This software was developed by the Computer Systems Engineering group
- * at Lawrence Berkeley Laboratory under DARPA contract BG 91-66 and
- * contributed to Berkeley.
- *
- * Redistribution and use in source and binary forms, with or without
- * modification, are permitted provided that the following conditions
- * are met:
- * 1. Redistributions of source code must retain the above copyright
- *    notice, this list of conditions and the following disclaimer.
- * 2. Redistributions in binary form must reproduce the above copyright
- *    notice, this list of conditions and the following disclaimer in the
- *    documentation and/or other materials provided with the distribution.
- * 3. Neither the name of the University nor the names of its contributors
- *    may be used to endorse or promote products derived from this software
- *    without specific prior written permission.
- *
- * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
- * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
- * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
- * ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
- * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
- * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
- * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
- * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
- * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
- * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
- * SUCH DAMAGE.
- *
- * from: Header: ldexp.c,v 1.1 91/07/07 04:28:19 torek Exp
- */
-
-#include <sys/types.h>
-#include <sys/cdefs.h>
-#include <machine/ieee.h>
-#include <errno.h>
-
-/*
- * double ldexp(double val, int exp)
- * returns: val * (2**exp)
- */
-double
-ldexp(val, exp)
-	double val;
-	int exp;
-{
-	register int oldexp, newexp, mulexp;
-	union doub {
-		double v;
-		struct ieee_double s;
-	} u, mul;
-
-	/*
-	 * If input is zero, or no change, just return input.
-	 * Likewise, if input is Inf or NaN, just return it.
-	 */
-	u.v = val;
-	oldexp = u.s.dbl_exp;
-	if (val == 0 || exp == 0 || oldexp == DBL_EXP_INFNAN)
-		return (val);
-
-	/*
-	 * Compute new exponent and check for over/under flow.
-	 * Underflow, unfortunately, could mean switching to denormal.
-	 * If result out of range, set ERANGE and return 0 if too small
-	 * or Inf if too big, with the same sign as the input value.
-	 */
-	newexp = oldexp + exp;
-	if (newexp >= DBL_EXP_INFNAN) {
-		/* u.s.dbl_sign = val < 0; -- already set */
-		u.s.dbl_exp = DBL_EXP_INFNAN;
-		u.s.dbl_frach = u.s.dbl_fracl = 0;
-		errno = ERANGE;
-		return (u.v);		/* Inf */
-	}
-	if (newexp <= 0) {
-		/*
-		 * The output number is either a denormal or underflows
-		 * (see comments in machine/ieee.h).
-		 */
-		if (newexp <= -DBL_FRACBITS) {
-			/* u.s.dbl_sign = val < 0; -- already set */
-			u.s.dbl_exp = 0;
-			u.s.dbl_frach = u.s.dbl_fracl = 0;
-			errno = ERANGE;
-			return (u.v);		/* zero */
-		}
-		/*
-		 * We are going to produce a denorm.  Our `exp' argument
-		 * might be as small as -2097, and we cannot compute
-		 * 2^-2097, so we may have to do this as many as three
-		 * steps (not just two, as for positive `exp's below).
-		 */
-		mul.v = 0;
-		while (exp <= -DBL_EXP_BIAS) {
-			mul.s.dbl_exp = 1;
-			val *= mul.v;
-			exp += DBL_EXP_BIAS - 1;
-		}
-		mul.s.dbl_exp = exp + DBL_EXP_BIAS;
-		val *= mul.v;
-		return (val);
-	}
-
-	/*
-	 * Newexp is positive.
-	 *
-	 * If oldexp is zero, we are starting with a denorm, and simply
-	 * adjusting the exponent will produce bogus answers.  We need
-	 * to fix that first.
-	 */
-	if (oldexp == 0) {
-		/*
-		 * Multiply by 2^mulexp to make the number normalizable.
-		 * We cannot multiply by more than 2^1023, but `exp'
-		 * argument might be as large as 2046.  A single
-		 * adjustment, however, will normalize the number even
-		 * for huge `exp's, and then we can use exponent
-		 * arithmetic just as for normal `double's.
-		 */
-		mulexp = exp <= DBL_EXP_BIAS ? exp : DBL_EXP_BIAS;
-		mul.v = 0;
-		mul.s.dbl_exp = mulexp + DBL_EXP_BIAS;
-		val *= mul.v;
-		if (mulexp == exp)
-			return (val);
-		u.v = val;
-		newexp -= mulexp;
-	}
-
-	/*
-	 * Both oldexp and newexp are positive; just replace the
-	 * old exponent with the new one.
-	 */
-	u.s.dbl_exp = newexp;
-	return (u.v);
-}
-
-__weak_alias(ldexpl, ldexp);
diff --git a/lib/libc/arch/sh/gen/modf.c b/lib/libc/arch/sh/gen/modf.c
deleted file mode 100644
index da646541432..00000000000
--- a/lib/libc/arch/sh/gen/modf.c
+++ /dev/null
@@ -1,305 +0,0 @@
-/*	$OpenBSD: modf.c,v 1.1 2006/10/10 22:07:10 miod Exp $	*/
-/* @(#)s_modf.c 5.1 93/09/24 */
-/*
- * ====================================================
- * 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.
- * ====================================================
- */
-
-/*
- * modf(double x, double *iptr) 
- * return fraction part of x, and return x's integral part in *iptr.
- * Method:
- *	Bit twiddling.
- *
- * Exception:
- *	No exception.
- */
-
-#include "math.h"
-
-/*
- * ====================================================
- * 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
- */
-
-#ifndef _MATH_PRIVATE_H_
-#define _MATH_PRIVATE_H_
-
-#include <sys/types.h>
-#include <machine/endian.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.  */
-
-/*
- * The arm32 port is little endian except for the FP word order which is
- * big endian.
- */
-
-#if (BYTE_ORDER == BIG_ENDIAN) || defined(arm32)
-
-typedef union 
-{
-  double value;
-  struct 
-  {
-    u_int32_t msw;
-    u_int32_t lsw;
-  } parts;
-} ieee_double_shape_type;
-
-#endif
-
-#if (BYTE_ORDER == LITTLE_ENDIAN) && !defined(arm32)
-
-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(double);			
-extern double __ieee754_acos(double);			
-extern double __ieee754_acosh(double);			
-extern double __ieee754_log(double);			
-extern double __ieee754_atanh(double);			
-extern double __ieee754_asin(double);			
-extern double __ieee754_atan2(double,double);			
-extern double __ieee754_exp(double);
-extern double __ieee754_cosh(double);
-extern double __ieee754_fmod(double,double);
-extern double __ieee754_pow(double,double);
-extern double __ieee754_lgamma_r(double,int *);
-extern double __ieee754_gamma_r(double,int *);
-extern double __ieee754_lgamma(double);
-extern double __ieee754_gamma(double);
-extern double __ieee754_log10(double);
-extern double __ieee754_sinh(double);
-extern double __ieee754_hypot(double,double);
-extern double __ieee754_j0(double);
-extern double __ieee754_j1(double);
-extern double __ieee754_y0(double);
-extern double __ieee754_y1(double);
-extern double __ieee754_jn(int,double);
-extern double __ieee754_yn(int,double);
-extern double __ieee754_remainder(double,double);
-extern int    __ieee754_rem_pio2(double,double*);
-extern double __ieee754_scalb(double,double);
-
-/* fdlibm kernel function */
-extern double __kernel_standard(double,double,int);	
-extern double __kernel_sin(double,double,int);
-extern double __kernel_cos(double,double);
-extern double __kernel_tan(double,double,int);
-extern int    __kernel_rem_pio2(double*,double*,int,int,int,const int*);
-
-
-/* ieee style elementary float functions */
-extern float __ieee754_sqrtf(float);			
-extern float __ieee754_acosf(float);			
-extern float __ieee754_acoshf(float);			
-extern float __ieee754_logf(float);			
-extern float __ieee754_atanhf(float);			
-extern float __ieee754_asinf(float);			
-extern float __ieee754_atan2f(float,float);			
-extern float __ieee754_expf(float);
-extern float __ieee754_coshf(float);
-extern float __ieee754_fmodf(float,float);
-extern float __ieee754_powf(float,float);
-extern float __ieee754_lgammaf_r(float,int *);
-extern float __ieee754_gammaf_r(float,int *);
-extern float __ieee754_lgammaf(float);
-extern float __ieee754_gammaf(float);
-extern float __ieee754_log10f(float);
-extern float __ieee754_sinhf(float);
-extern float __ieee754_hypotf(float,float);
-extern float __ieee754_j0f(float);
-extern float __ieee754_j1f(float);
-extern float __ieee754_y0f(float);
-extern float __ieee754_y1f(float);
-extern float __ieee754_jnf(int,float);
-extern float __ieee754_ynf(int,float);
-extern float __ieee754_remainderf(float,float);
-extern int   __ieee754_rem_pio2f(float,float*);
-extern float __ieee754_scalbf(float,float);
-
-/* float versions of fdlibm kernel functions */
-extern float __kernel_sinf(float,float,int);
-extern float __kernel_cosf(float,float);
-extern float __kernel_tanf(float,float,int);
-extern int   __kernel_rem_pio2f(float*,float*,int,int,int,const int*);
-
-#endif /* _MATH_PRIVATE_H_ */
-#ifdef __STDC__
-static const double one = 1.0;
-#else
-static double one = 1.0;
-#endif
-
-#ifdef __STDC__
-	double modf(double x, double *iptr)
-#else
-	double modf(x, iptr)
-	double x,*iptr;
-#endif
-{
-	int32_t i0,i1,j0;
-	u_int32_t i;
-	EXTRACT_WORDS(i0,i1,x);
-	j0 = ((i0>>20)&0x7ff)-0x3ff;	/* exponent of x */
-	if(j0<20) {			/* integer part in high x */
-	    if(j0<0) {			/* |x|<1 */
-		INSERT_WORDS(*iptr,i0&0x80000000,0);	/* *iptr = +-0 */
-		return x;
-	    } else {
-		i = (0x000fffff)>>j0;
-		if(((i0&i)|i1)==0) {		/* x is integral */
-		    u_int32_t high;
-		    *iptr = x;
-		    GET_HIGH_WORD(high,x);
-		    INSERT_WORDS(x,high&0x80000000,0);	/* return +-0 */
-		    return x;
-		} else {
-		    INSERT_WORDS(*iptr,i0&(~i),0);
-		    return x - *iptr;
-		}
-	    }
-	} else if (j0>51) {		/* no fraction part */
-	    u_int32_t high;
-	    *iptr = x*one;
-	    GET_HIGH_WORD(high,x);
-	    INSERT_WORDS(x,high&0x80000000,0);	/* return +-0 */
-	    return x;
-	} else {			/* fraction part in low x */
-	    i = ((u_int32_t)(0xffffffff))>>(j0-20);
-	    if((i1&i)==0) { 		/* x is integral */
-		u_int32_t high;
-		*iptr = x;
-		GET_HIGH_WORD(high,x);
-		INSERT_WORDS(x,high&0x80000000,0);	/* return +-0 */
-		return x;
-	    } else {
-		INSERT_WORDS(*iptr,i0,i1&(~i));
-		return x - *iptr;
-	    }
-	}
-}
diff --git a/lib/libc/arch/sparc/gen/Makefile.inc b/lib/libc/arch/sparc/gen/Makefile.inc
index a2475620099..993c7121866 100644
--- a/lib/libc/arch/sparc/gen/Makefile.inc
+++ b/lib/libc/arch/sparc/gen/Makefile.inc
@@ -1,6 +1,6 @@
-#	$OpenBSD: Makefile.inc,v 1.8 2008/12/09 19:52:34 martynas Exp $
+#	$OpenBSD: Makefile.inc,v 1.9 2009/04/19 16:42:06 martynas Exp $
 
-SRCS+=	_setjmp.S fabs.S infinity.c ldexp.c modf.S nan.c
+SRCS+=	_setjmp.S fabs.S infinity.c modf.S nan.c
 SRCS+=	flt_rounds.c fpgetmask.c fpgetround.c fpgetsticky.c fpsetmask.c \
 	fpsetround.c fpsetsticky.c
 SRCS+=	fixunsdfsi.S mul.S umul.S saveregs.S setjmp.S sigsetjmp.S
diff --git a/lib/libc/arch/sparc/gen/ldexp.c b/lib/libc/arch/sparc/gen/ldexp.c
deleted file mode 100644
index a12b322d122..00000000000
--- a/lib/libc/arch/sparc/gen/ldexp.c
+++ /dev/null
@@ -1,141 +0,0 @@
-/*	$OpenBSD: ldexp.c,v 1.6 2008/12/10 01:15:02 martynas Exp $ */
-/*
- * Copyright (c) 1992, 1993
- *	The Regents of the University of California.  All rights reserved.
- *
- * This software was developed by the Computer Systems Engineering group
- * at Lawrence Berkeley Laboratory under DARPA contract BG 91-66 and
- * contributed to Berkeley.
- *
- * Redistribution and use in source and binary forms, with or without
- * modification, are permitted provided that the following conditions
- * are met:
- * 1. Redistributions of source code must retain the above copyright
- *    notice, this list of conditions and the following disclaimer.
- * 2. Redistributions in binary form must reproduce the above copyright
- *    notice, this list of conditions and the following disclaimer in the
- *    documentation and/or other materials provided with the distribution.
- * 3. Neither the name of the University nor the names of its contributors
- *    may be used to endorse or promote products derived from this software
- *    without specific prior written permission.
- *
- * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
- * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
- * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
- * ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
- * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
- * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
- * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
- * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
- * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
- * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
- * SUCH DAMAGE.
- */
-
-#include <sys/types.h>
-#include <sys/cdefs.h>
-#include <machine/ieee.h>
-#include <errno.h>
-
-/*
- * double ldexp(double val, int exp)
- * returns: val * (2**exp)
- */
-double
-ldexp(val, exp)
-	double val;
-	int exp;
-{
-	register int oldexp, newexp, mulexp;
-	union doub {
-		double v;
-		struct ieee_double s;
-	} u, mul;
-
-	/*
-	 * If input is zero, or no change, just return input.
-	 * Likewise, if input is Inf or NaN, just return it.
-	 */
-	u.v = val;
-	oldexp = u.s.dbl_exp;
-	if (val == 0 || exp == 0 || oldexp == DBL_EXP_INFNAN)
-		return (val);
-
-	/*
-	 * Compute new exponent and check for over/under flow.
-	 * Underflow, unfortunately, could mean switching to denormal.
-	 * If result out of range, set ERANGE and return 0 if too small
-	 * or Inf if too big, with the same sign as the input value.
-	 */
-	newexp = oldexp + exp;
-	if (newexp >= DBL_EXP_INFNAN) {
-		/* u.s.dbl_sign = val < 0; -- already set */
-		u.s.dbl_exp = DBL_EXP_INFNAN;
-		u.s.dbl_frach = u.s.dbl_fracl = 0;
-		errno = ERANGE;
-		return (u.v);		/* Inf */
-	}
-	if (newexp <= 0) {
-		/*
-		 * The output number is either a denormal or underflows
-		 * (see comments in machine/ieee.h).
-		 */
-		if (newexp <= -DBL_FRACBITS) {
-			/* u.s.dbl_sign = val < 0; -- already set */
-			u.s.dbl_exp = 0;
-			u.s.dbl_frach = u.s.dbl_fracl = 0;
-			errno = ERANGE;
-			return (u.v);		/* zero */
-		}
-		/*
-		 * We are going to produce a denorm.  Our `exp' argument
-		 * might be as small as -2097, and we cannot compute
-		 * 2^-2097, so we may have to do this as many as three
-		 * steps (not just two, as for positive `exp's below).
-		 */
-		mul.v = 0;
-		while (exp <= -DBL_EXP_BIAS) {
-			mul.s.dbl_exp = 1;
-			val *= mul.v;
-			exp += DBL_EXP_BIAS - 1;
-		}
-		mul.s.dbl_exp = exp + DBL_EXP_BIAS;
-		val *= mul.v;
-		return (val);
-	}
-
-	/*
-	 * Newexp is positive.
-	 *
-	 * If oldexp is zero, we are starting with a denorm, and simply
-	 * adjusting the exponent will produce bogus answers.  We need
-	 * to fix that first.
-	 */
-	if (oldexp == 0) {
-		/*
-		 * Multiply by 2^mulexp to make the number normalizable.
-		 * We cannot multiply by more than 2^1023, but `exp'
-		 * argument might be as large as 2046.  A single
-		 * adjustment, however, will normalize the number even
-		 * for huge `exp's, and then we can use exponent
-		 * arithmetic just as for normal `double's.
-		 */
-		mulexp = exp <= DBL_EXP_BIAS ? exp : DBL_EXP_BIAS;
-		mul.v = 0;
-		mul.s.dbl_exp = mulexp + DBL_EXP_BIAS;
-		val *= mul.v;
-		if (mulexp == exp)
-			return (val);
-		u.v = val;
-		newexp -= mulexp;
-	}
-
-	/*
-	 * Both oldexp and newexp are positive; just replace the
-	 * old exponent with the new one.
-	 */
-	u.s.dbl_exp = newexp;
-	return (u.v);
-}
-
-__weak_alias(ldexpl, ldexp);
diff --git a/lib/libc/arch/sparc64/gen/Makefile.inc b/lib/libc/arch/sparc64/gen/Makefile.inc
index 5ae927d7383..95bc28135ac 100644
--- a/lib/libc/arch/sparc64/gen/Makefile.inc
+++ b/lib/libc/arch/sparc64/gen/Makefile.inc
@@ -1,6 +1,6 @@
-#	$OpenBSD: Makefile.inc,v 1.7 2008/12/09 19:52:34 martynas Exp $
+#	$OpenBSD: Makefile.inc,v 1.8 2009/04/19 16:42:06 martynas Exp $
 
-SRCS+=	_setjmp.S fabs.S infinity.c ldexp.c modf.S nan.c
+SRCS+=	_setjmp.S fabs.S infinity.c modf.S nan.c
 SRCS+=	flt_rounds.c fpgetmask.c fpgetround.c fpgetsticky.c fpsetmask.c \
 	fpsetround.c fpsetsticky.c
 SRCS+=	fpclassifyl.c isfinitel.c isinfl.c isnanl.c isnormall.c signbitl.c
diff --git a/lib/libc/arch/sparc64/gen/ldexp.c b/lib/libc/arch/sparc64/gen/ldexp.c
deleted file mode 100644
index aa7ea600734..00000000000
--- a/lib/libc/arch/sparc64/gen/ldexp.c
+++ /dev/null
@@ -1,145 +0,0 @@
-/*	$OpenBSD: ldexp.c,v 1.2 2008/06/26 05:42:05 ray Exp $	*/
-/*	$NetBSD: ldexp.c,v 1.8 1999/08/30 18:28:26 mycroft Exp $	*/
-
-/*-
- * Copyright (c) 1999 The NetBSD Foundation, Inc.
- * All rights reserved.
- *
- * This code is derived from software contributed to The NetBSD Foundation
- * by Charles M. Hannum.
- *
- * Redistribution and use in source and binary forms, with or without
- * modification, are permitted provided that the following conditions
- * are met:
- * 1. Redistributions of source code must retain the above copyright
- *    notice, this list of conditions and the following disclaimer.
- * 2. Redistributions in binary form must reproduce the above copyright
- *    notice, this list of conditions and the following disclaimer in the
- *    documentation and/or other materials provided with the distribution.
- *
- * THIS SOFTWARE IS PROVIDED BY THE NETBSD FOUNDATION, INC. AND CONTRIBUTORS
- * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
- * TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
- * PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR CONTRIBUTORS
- * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
- * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
- * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
- * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
- * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
- * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
- * POSSIBILITY OF SUCH DAMAGE.
- */
-
-#include <sys/cdefs.h>
-
-#include <sys/types.h>
-#include <machine/ieee.h>
-#include <errno.h>
-#include <math.h>
-
-/*
- * Multiply the given value by 2^exp.
- */
-double
-ldexp(val, exp)
-	double val;
-	int exp;
-{
-	register int oldexp, newexp;
-	union {
-		double v;
-		struct ieee_double s;
-	} u, mul;
-
-	u.v = val;
-	oldexp = u.s.dbl_exp;
-
-	/*
-	 * If input is zero, Inf or NaN, just return it.
-	 */
-	if (u.v == 0.0 || oldexp == DBL_EXP_INFNAN)
-		return (val);
-
-	if (oldexp == 0) {
-		/*
-		 * u.v is denormal.  We must adjust it so that the exponent
-		 * arithmetic below will work.
-		 */
-		if (exp <= DBL_EXP_BIAS) {
-			/*
-			 * Optimization: if the scaling can be done in a single
-			 * multiply, or underflows, just do it now.
-			 */
-			if (exp <= -DBL_FRACBITS) {
-				errno = ERANGE;
-				return (0.0);
-			}
-			mul.v = 0.0;
-			mul.s.dbl_exp = exp + DBL_EXP_BIAS;
-			u.v *= mul.v;
-			if (u.v == 0.0) {
-				errno = ERANGE;
-				return (0.0);
-			}
-			return (u.v);
-		} else {
-			/*
-			 * We know that exp is very large, and therefore the
-			 * result cannot be denormal (though it may be Inf).
-			 * Shift u.v by just enough to make it normal.
-			 */
-			mul.v = 0.0;
-			mul.s.dbl_exp = DBL_FRACBITS + DBL_EXP_BIAS;
-			u.v *= mul.v;
-			exp -= DBL_FRACBITS;
-			oldexp = u.s.dbl_exp;
-		}
-	}
-
-	/*
-	 * u.v is now normalized and oldexp has been adjusted if necessary.
-	 * Calculate the new exponent and check for underflow and overflow.
-	 */
-	newexp = oldexp + exp;
-
-	if (newexp <= 0) {
-		/*
-		 * The output number is either denormal or underflows (see
-		 * comments in machine/ieee.h).
-		 */
-		if (newexp <= -DBL_FRACBITS) {
-			errno = ERANGE;
-			return (0.0);
-		}
-		/*
-		 * Denormalize the result.  We do this with a multiply.  If exp
-		 * is very large, it won't fit in a double, so we have to
-		 * adjust the exponent first.  This is safe because we know
-		 * that u.v is normal at this point.
-		 */
-		if (exp <= -DBL_EXP_BIAS) {
-			u.s.dbl_exp = 1;
-			exp += oldexp - 1;
-		}
-		mul.v = 0.0;
-		mul.s.dbl_exp = exp + DBL_EXP_BIAS;
-		u.v *= mul.v;
-		return (u.v);
-	} else if (newexp >= DBL_EXP_INFNAN) {
-		/*
-		 * The result overflowed; return +/-Inf.
-		 */
-		u.s.dbl_exp = DBL_EXP_INFNAN;
-		u.s.dbl_frach = 0;
-		u.s.dbl_fracl = 0;
-		errno = ERANGE;
-		return (u.v);
-	} else {
-		/*
-		 * The result is normal; just replace the old exponent with the
-		 * new one.
-		 */
-		u.s.dbl_exp = newexp;
-		return (u.v);
-	}
-}