Add and document ilogbl(), a long double version of ilogb().

4 files changed, 65 insertions, 4 deletions

diff --git a/lib/msun/Makefile b/lib/msun/Makefile
index 496e740..712f0c2 100644
--- a/lib/msun/Makefile
+++ b/lib/msun/Makefile
@@ -85,7 +85,7 @@ COMMON_SRCS= b_exp.c b_log.c b_tgamma.c \
 	s_expm1.c s_expm1f.c s_fabsf.c s_fdim.c s_finite.c s_finitef.c \
 	s_floor.c s_floorf.c s_fmax.c s_fmaxf.c s_fmaxl.c s_fmin.c \
 	s_fminf.c s_fminl.c s_frexp.c s_frexpf.c s_ilogb.c s_ilogbf.c \
-	s_isfinite.c s_isnan.c s_isnormal.c s_ldexpf.c \
+	s_ilogbl.c s_isfinite.c s_isnan.c s_isnormal.c s_ldexpf.c \
 	s_lib_version.c s_log1p.c \
 	s_log1pf.c s_logb.c s_logbf.c s_matherr.c s_modff.c \
 	s_nearbyint.c s_nextafter.c s_nextafterf.c \
@@ -167,7 +167,8 @@ MLINKS+=fmod.3 fmodf.3
 MLINKS+=hypot.3 cabs.3 hypot.3 cabsf.3 hypot.3 hypotf.3
 MLINKS+=ieee.3 copysign.3 ieee.3 copysignf.3 ieee.3 copysignl.3 \
 	ieee.3 finite.3 ieee.3 finitef.3 \
-	ieee.3 ilogb.3 ieee.3 ilogbf.3 ieee.3 nextafter.3 ieee.3 nextafterf.3 \
+	ieee.3 ilogb.3 ieee.3 ilogbf.3 ieee.3 ilogbl.3 \
+	ieee.3 nextafter.3 ieee.3 nextafterf.3 \
 	ieee.3 remainder.3 ieee.3 remainderf.3 \
 	ieee.3 scalbln.3 ieee.3 scalblnf.3 ieee.3 scalbn.3 ieee.3 scalbnf.3
 MLINKS+=ieee_test.3 logb.3 ieee_test.3 logbf.3
diff --git a/lib/msun/man/ieee.3 b/lib/msun/man/ieee.3
index 73a39f6..527bdbc 100644
--- a/lib/msun/man/ieee.3
+++ b/lib/msun/man/ieee.3
@@ -43,6 +43,7 @@
 .Nm finitef ,
 .Nm ilogb ,
 .Nm ilogbf ,
+.Nm ilogbl ,
 .Nm nextafter ,
 .Nm nextafterf ,
 .Nm remainder ,
@@ -70,6 +71,8 @@
 .Fn ilogb "double x"
 .Ft int
 .Fn ilogbf "float x"
+.Ft int
+.Fn ilogbl "long double x"
 .Ft double
 .Fn nextafter "double x" "double y"
 .Ft float
@@ -116,9 +119,10 @@ zero is returned
 .Fa x
 is \*(Na).
 .Pp
-.Fn ilogb
-and
+.Fn ilogb ,
 .Fn ilogbf
+and
+.Fn ilogbl
 return
 .Fa x Ns 's exponent,
 in integer format.
diff --git a/lib/msun/src/math.h b/lib/msun/src/math.h
index bea5f2b..7d85545 100644
--- a/lib/msun/src/math.h
+++ b/lib/msun/src/math.h
@@ -447,7 +447,9 @@ long double	fminl(long double, long double) __pure2;
 long double	fmodl(long double, long double);
 long double	frexpl(long double	value, int *);
 long double	hypotl(long double, long double);
+#endif
 int		ilogbl(long double);
+#if 0
 long double	ldexpl(long double, int);
 long double	lgammal(long double);
 long long	llrintl(long double);
diff --git a/lib/msun/src/s_ilogbl.c b/lib/msun/src/s_ilogbl.c
new file mode 100644
index 0000000..a6a4708
--- /dev/null
+++ b/lib/msun/src/s_ilogbl.c
@@ -0,0 +1,54 @@
+/*
+ * From: @(#)s_ilogb.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.
+ * ====================================================
+ */
+
+#ifndef lint
+static char rcsid[] = "$FreeBSD$";
+#endif
+
+#include <float.h>
+#include <limits.h>
+#include <math.h>
+
+#include "fpmath.h"
+
+int
+ilogbl(long double x)
+{
+	union IEEEl2bits u;
+	unsigned long m;
+	int b;
+
+	u.e = x;
+	if (u.bits.exp == 0) {
+		if ((u.bits.manl | u.bits.manh) == 0)
+			return (FP_ILOGB0);
+		/* denormalized */
+		if (u.bits.manh == 0) {
+			m = 1lu << (LDBL_MANL_SIZE - 1);
+			for (b = LDBL_MANH_SIZE; !(u.bits.manl & m); m >>= 1)
+				b++;
+		} else {
+			m = 1lu << (LDBL_MANH_SIZE - 1);
+			for (b = 0; !(u.bits.manh & m); m >>= 1)
+				b++;
+		}
+#ifdef LDBL_IMPLICIT_NBIT
+		b++;
+#endif
+		return (LDBL_MIN_EXP - b - 1);
+	} else if (u.bits.exp < (LDBL_MAX_EXP << 1) - 1)
+		return (u.bits.exp - LDBL_MAX_EXP + 1);
+	else if (u.bits.manl != 0 || u.bits.manh != 0)
+		return (FP_ILOGBNAN);
+	else
+		return (INT_MAX);
+}