Add docs and tests for fp classificators

2 files changed, 534 insertions, 304 deletions

diff --git a/source/tanya/math/fp.d b/source/tanya/math/fp.d
deleted file mode 100644
index 7206c22..0000000
--- a/source/tanya/math/fp.d
+++ /dev/null
@@ -1,304 +0,0 @@
-/* This Source Code Form is subject to the terms of the Mozilla Public
- * License, v. 2.0. If a copy of the MPL was not distributed with this
- * file, You can obtain one at http://mozilla.org/MPL/2.0/. */
-
-/**
- * Copyright: Eugene Wissner 2017.
- * License: $(LINK2 https://www.mozilla.org/en-US/MPL/2.0/,
- *                  Mozilla Public License, v. 2.0).
- * Authors: $(LINK2 mailto:info@caraus.de, Eugene Wissner)
- * Source: $(LINK2 https://github.com/caraus-ecms/tanya/blob/master/source/tanya/math/fp.d,
- *                 tanya/math/fp.d)
- */
-module tanya.math.fp;
-
-import tanya.math.nbtheory;
-
-/**
- * Floating-point number precisions according to IEEE-754.
- */
-enum IEEEPrecision : ubyte
-{
-    /// Single precision: 64-bit.
-    single = 4,
-
-    /// Single precision: 64-bit.
-    double_ = 8,
-
-    /// Extended precision: 80-bit.
-    extended = 10,
-}
-
-/**
- * Tests the precision of floating-point type $(D_PARAM F).
- *
- * For $(D_KEYWORD float), $(D_PSYMBOL ieeePrecision) always evaluates to
- * $(D_INLINECODE IEEEPrecision.single); for $(D_KEYWORD double) - to
- * $(D_INLINECODE IEEEPrecision.double). It returns different values only
- * for $(D_KEYWORD real), since $(D_KEYWORD real) is a platform-dependent type.
- *
- * If $(D_PARAM F) is a $(D_KEYWORD real) and the target platform isn't
- * currently supported, static assertion error will be raised (you can use
- * $(D_INLINECODE is(typeof(ieeePrecision!F))) for testing the platform support
- * without a compilation error).
- *
- * Params:
- *  F = Type to be tested.
- *
- * Returns: Precision according to IEEE-754.
- *
- * See_Also: $(D_PSYMBOL IEEEPrecision).
- */
-template ieeePrecision(F)
-if (isFloatingPoint!F)
-{
-    static if (F.sizeof == float.sizeof)
-    {
-        enum IEEEPrecision ieeePrecision = IEEEPrecision.single;
-    }
-    else static if (F.sizeof == double.sizeof)
-    {
-        enum IEEEPrecision ieeePrecision = IEEEPrecision.double_;
-    }
-    else version (X86)
-    {
-        enum IEEEPrecision ieeePrecision = IEEEPrecision.extended;
-    }
-    else version (X86_64)
-    {
-        enum IEEEPrecision ieeePrecision = IEEEPrecision.extended;
-    }
-    else
-    {
-        static assert(false, "Unsupported IEEE 754 floating point precision");
-    }
-}
-
-private union FloatBits(F)
-{
-    F floating;
-    static if (ieeePrecision!F == IEEEPrecision.single)
-    {
-        uint integral;
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.double_)
-    {
-        ulong integral;
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.extended)
-    {
-        struct // Little-endian.
-        {
-            ulong mantissa;
-            ushort exp;
-        }
-    }
-    else
-    {
-        static assert(false, "Unsupported IEEE 754 floating point precision");
-    }
-}
-
-enum FloatingPointClass : ubyte
-{
-    nan,
-    zero,
-    infinite,
-    subnormal,
-    normal,
-}
-
-FloatingPointClass classify(F)(F x)
-if (isFloatingPoint!F)
-{
-    if (x == 0)
-    {
-        return FloatingPointClass.zero;
-    }
-    FloatBits!F bits;
-    bits.floating = abs(x);
-
-    static if (ieeePrecision!F == IEEEPrecision.single)
-    {
-        if (bits.integral > 0x7f800000)
-        {
-            return FloatingPointClass.nan;
-        }
-        else if (bits.integral == 0x7f800000)
-        {
-            return FloatingPointClass.infinite;
-        }
-        else if (bits.integral < 0x800000)
-        {
-            return FloatingPointClass.subnormal;
-        }
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.double_)
-    {
-        if (bits.integral > 0x7ff0000000000000)
-        {
-            return FloatingPointClass.nan;
-        }
-        else if (bits.integral == 0x7ff0000000000000)
-        {
-            return FloatingPointClass.infinite;
-        }
-        else if (bits.integral < 0x10000000000000)
-        {
-            return FloatingPointClass.subnormal;
-        }
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.extended)
-    {
-        if (bits.exp == 0x7fff)
-        {
-            if ((bits.mantissa & 0x7fffffffffffffff) == 0)
-            {
-                return FloatingPointClass.infinite;
-            }
-            else
-            {
-                return FloatingPointClass.nan;
-            }
-        }
-        else if (bits.exp == 0)
-        {
-            return FloatingPointClass.subnormal;
-        }
-        else if (bits.mantissa < 0x8000000000000000) // "Unnormal".
-        {
-            return FloatingPointClass.nan;
-        }
-    }
-
-    return FloatingPointClass.normal;
-}
-
-bool isFinite(F)(F x)
-if (isFloatingPoint!F)
-{
-    FloatBits!F bits;
-    static if (ieeePrecision!F == IEEEPrecision.single)
-    {
-        bits.floating = x;
-        bits.integral &= 0x7f800000;
-        return bits.integral != 0x7f800000;
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.double_)
-    {
-        bits.floating = x;
-        bits.integral &= 0x7ff0000000000000;
-        return bits.integral != 0x7ff0000000000000;
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.extended)
-    {
-        bits.floating = abs(x);
-        return (bits.exp != 0x7fff) && (bits.mantissa >= 0x8000000000000000);
-    }
-}
-
-bool isNaN(F)(F x)
-if (isFloatingPoint!F)
-{
-    FloatBits!F bits;
-    bits.floating = abs(x);
-
-    static if (ieeePrecision!F == IEEEPrecision.single)
-    {
-        return bits.integral > 0x7f800000;
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.double_)
-    {
-        return bits.integral > 0x7ff0000000000000;
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.extended)
-    {
-        if ((bits.exp == 0x7fff && (bits.mantissa & 0x7fffffffffffffff) != 0)
-         || ((bits.exp != 0) && (bits.mantissa < 0x8000000000000000)))
-        {
-            return true;
-        }
-        return false;
-    }
-}
-
-bool isInfinity(F)(F x)
-if (isFloatingPoint!F)
-{
-    FloatBits!F bits;
-    bits.floating = abs(x);
-    static if (ieeePrecision!F == IEEEPrecision.single)
-    {
-        return bits.integral == 0x7f800000;
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.double_)
-    {
-        return bits.integral == 0x7ff0000000000000;
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.extended)
-    {
-        return (bits.exp == 0x7fff)
-            && ((bits.mantissa & 0x7fffffffffffffff) == 0);
-    }
-}
-
-bool isSubnormal(F)(F x)
-if (isFloatingPoint!F)
-{
-    FloatBits!F bits;
-    bits.floating = abs(x);
-    static if (ieeePrecision!F == IEEEPrecision.single)
-    {
-        return bits.integral < 0x800000;
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.double_)
-    {
-        return bits.integral < 0x10000000000000;
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.extended)
-    {
-        return bits.exp == 0;
-    }
-}
-
-bool isNormal(F)(F x)
-if (isFloatingPoint!F)
-{
-    static if (ieeePrecision!F == IEEEPrecision.single)
-    {
-        FloatBits!F bits;
-        bits.floating = x;
-        bits.integral &= 0x7f800000;
-        return bits.integral != 0 && bits.integral != 0x7f800000;
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.double_)
-    {
-        FloatBits!F bits;
-        bits.floating = x;
-        bits.integral &= 0x7ff0000000000000;
-        return bits.integral != 0 && bits.integral != 0x7ff0000000000000;
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.extended)
-    {
-        return classify(x) == FloatingPointClass.normal;
-    }
-}
-
-bool signBit(F)(F x)
-if (isFloatingPoint!F)
-{
-    FloatBits!F bits;
-    bits.floating = x;
-    static if (ieeePrecision!F == IEEEPrecision.single)
-    {
-        return (bits.integral & (1 << 31)) != 0;
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.double_)
-    {
-        return (bits.integral & (1 << 63)) != 0;
-    }
-    else static if (ieeePrecision!F == IEEEPrecision.extended)
-    {
-        return (bits.exp & (1 << 15)) != 0;
-    }
-}
diff --git a/source/tanya/math/package.d b/source/tanya/math/package.d
index 86f1158..5d7c37d 100644
--- a/source/tanya/math/package.d
+++ b/source/tanya/math/package.d
@@ -19,6 +19,540 @@ public import tanya.math.nbtheory;
 public import tanya.math.random;
 import tanya.meta.trait;
 
+/// Floating-point number precisions according to IEEE-754.
+enum IEEEPrecision : ubyte
+{
+    single = 4, /// Single precision: 64-bit.
+    double_ = 8, /// Single precision: 64-bit.
+    doubleExtended = 10, /// Double extended precision: 80-bit.
+}
+
+/**
+ * Tests the precision of floating-point type $(D_PARAM F).
+ *
+ * For $(D_KEYWORD float), $(D_PSYMBOL ieeePrecision) always evaluates to
+ * $(D_INLINECODE IEEEPrecision.single); for $(D_KEYWORD double) - to
+ * $(D_INLINECODE IEEEPrecision.double). It returns different values only
+ * for $(D_KEYWORD real), since $(D_KEYWORD real) is a platform-dependent type.
+ *
+ * If $(D_PARAM F) is a $(D_KEYWORD real) and the target platform isn't
+ * currently supported, static assertion error will be raised (you can use
+ * $(D_INLINECODE is(typeof(ieeePrecision!F))) for testing the platform support
+ * without a compilation error).
+ *
+ * Params:
+ *  F = Type to be tested.
+ *
+ * Returns: Precision according to IEEE-754.
+ *
+ * See_Also: $(D_PSYMBOL IEEEPrecision).
+ */
+template ieeePrecision(F)
+if (isFloatingPoint!F)
+{
+    static if (F.sizeof == float.sizeof)
+    {
+        enum IEEEPrecision ieeePrecision = IEEEPrecision.single;
+    }
+    else static if (F.sizeof == double.sizeof)
+    {
+        enum IEEEPrecision ieeePrecision = IEEEPrecision.double_;
+    }
+    else version (X86)
+    {
+        enum IEEEPrecision ieeePrecision = IEEEPrecision.doubleExtended;
+    }
+    else version (X86_64)
+    {
+        enum IEEEPrecision ieeePrecision = IEEEPrecision.doubleExtended;
+    }
+    else
+    {
+        static assert(false, "Unsupported IEEE 754 floating point precision");
+    }
+}
+
+///
+pure nothrow @safe @nogc unittest
+{
+    static assert(ieeePrecision!float == IEEEPrecision.single);
+    static assert(ieeePrecision!double == IEEEPrecision.double_);
+}
+
+private union FloatBits(F)
+{
+    F floating;
+    static if (ieeePrecision!F == IEEEPrecision.single)
+    {
+        uint integral;
+        enum uint expMask = 0x7f800000;
+    }
+    else static if (ieeePrecision!F == IEEEPrecision.double_)
+    {
+        ulong integral;
+        enum ulong expMask = 0x7ff0000000000000;
+    }
+    else static if (ieeePrecision!F == IEEEPrecision.doubleExtended)
+    {
+        struct // Little-endian.
+        {
+            ulong mantissa;
+            ushort exp;
+        }
+        enum ulong mantissaMask = 0x7fffffffffffffff;
+        enum uint expMask = 0x7fff;
+    }
+    else
+    {
+        static assert(false, "Unsupported IEEE 754 floating point precision");
+    }
+}
+
+/**
+ * Floating-point number classifications.
+ */
+enum FloatingPointClass : ubyte
+{
+    /**
+     * Not a Number.
+     *
+     * See_Also: $(D_PSYMBOL isNaN).
+     */
+    nan,
+
+    /// Zero.
+    zero,
+
+    /**
+     * Infinity.
+     *
+     * See_Also: $(D_PSYMBOL isInfinity).
+     */
+    infinite,
+
+    /**
+     * Denormalized number.
+     *
+     * See_Also: $(D_PSYMBOL isSubnormal).
+     */
+    subnormal,
+
+    /**
+     * Normalized number.
+     *
+     * See_Also: $(D_PSYMBOL isNormal).
+     */
+    normal,
+}
+
+/**
+ * Returns whether $(D_PARAM x) is a NaN, zero, infinity, subnormal or
+ * normalized number.
+ *
+ * This function doesn't distinguish between negative and positive infinity,
+ * negative and positive NaN or negative and positive zero.
+ *
+ * Params:
+ *  F = Type of the floating point number.
+ *  x = Floating point number.
+ *
+ * Returns: Classification of $(D_PARAM x).
+ */
+FloatingPointClass classify(F)(F x)
+if (isFloatingPoint!F)
+{
+    if (x == 0)
+    {
+        return FloatingPointClass.zero;
+    }
+    FloatBits!F bits;
+    bits.floating = abs(x);
+
+    static if (ieeePrecision!F == IEEEPrecision.single)
+    {
+        if (bits.integral > bits.expMask)
+        {
+            return FloatingPointClass.nan;
+        }
+        else if (bits.integral == bits.expMask)
+        {
+            return FloatingPointClass.infinite;
+        }
+        else if (bits.integral < (1 << 23))
+        {
+            return FloatingPointClass.subnormal;
+        }
+    }
+    else static if (ieeePrecision!F == IEEEPrecision.double_)
+    {
+        if (bits.integral > bits.expMask)
+        {
+            return FloatingPointClass.nan;
+        }
+        else if (bits.integral == bits.expMask)
+        {
+            return FloatingPointClass.infinite;
+        }
+        else if (bits.integral < (1L << 52))
+        {
+            return FloatingPointClass.subnormal;
+        }
+    }
+    else static if (ieeePrecision!F == IEEEPrecision.doubleExtended)
+    {
+        if (bits.exp == bits.expMask)
+        {
+            if ((bits.mantissa & bits.mantissaMask) == 0)
+            {
+                return FloatingPointClass.infinite;
+            }
+            else
+            {
+                return FloatingPointClass.nan;
+            }
+        }
+        else if (bits.exp == 0)
+        {
+            return FloatingPointClass.subnormal;
+        }
+        else if (bits.mantissa < (1L << 63)) // "Unnormal".
+        {
+            return FloatingPointClass.nan;
+        }
+    }
+
+    return FloatingPointClass.normal;
+}
+
+///
+pure nothrow @safe @nogc unittest
+{
+    assert(classify(0.0) == FloatingPointClass.zero);
+    assert(classify(double.nan) == FloatingPointClass.nan);
+    assert(classify(double.infinity) == FloatingPointClass.infinite);
+    assert(classify(-double.infinity) == FloatingPointClass.infinite);
+    assert(classify(1.4) == FloatingPointClass.normal);
+    assert(classify(1.11254e-307 / 10) == FloatingPointClass.subnormal);
+
+    assert(classify(0.0f) == FloatingPointClass.zero);
+    assert(classify(float.nan) == FloatingPointClass.nan);
+    assert(classify(float.infinity) == FloatingPointClass.infinite);
+    assert(classify(-float.infinity) == FloatingPointClass.infinite);
+    assert(classify(0.3) == FloatingPointClass.normal);
+    assert(classify(5.87747e-38f / 10) == FloatingPointClass.subnormal);
+
+    assert(classify(0.0L) == FloatingPointClass.zero);
+    assert(classify(real.nan) == FloatingPointClass.nan);
+    assert(classify(real.infinity) == FloatingPointClass.infinite);
+    assert(classify(-real.infinity) == FloatingPointClass.infinite);
+}
+
+private pure nothrow @nogc @safe unittest
+{
+    static if (ieeePrecision!float == IEEEPrecision.doubleExtended)
+    {
+        assert(classify(1.68105e-10) == FloatingPointClass.normal);
+        assert(classify(1.68105e-4932L) == FloatingPointClass.subnormal);
+
+        // Emulate unnormals, because they aren't generated anymore since i386
+        FloatBits!real unnormal;
+        unnormal.exp = 0x123;
+        unnormal.mantissa = 0x1;
+        assert(classify(unnormal) == FloatingPointClass.subnormal);
+    }
+}
+
+/**
+ * Determines whether $(D_PARAM x) is a finite number.
+ *
+ * Params:
+ *  F = Type of the floating point number.
+ *  x = Floating point number.
+ *
+ * Returns: $(D_KEYWORD true) if $(D_PARAM x) is a finite number,
+ *          $(D_KEYWORD false) otherwise.
+ *
+ * See_Also: $(D_PSYMBOL isInfinity).
+ */
+bool isFinite(F)(F x)
+if (isFloatingPoint!F)
+{
+    FloatBits!F bits;
+    static if (ieeePrecision!F == IEEEPrecision.single
+            || ieeePrecision!F == IEEEPrecision.double_)
+    {
+        bits.floating = x;
+        bits.integral &= bits.expMask;
+        return bits.integral != bits.expMask;
+    }
+    else static if (ieeePrecision!F == IEEEPrecision.doubleExtended)
+    {
+        bits.floating = abs(x);
+        return (bits.exp != bits.expMask)
+            && (bits.exp == 0 || bits.mantissa >= (1L << 63));
+    }
+}
+
+///
+pure nothrow @safe @nogc unittest
+{
+    assert(!isFinite(float.infinity));
+    assert(!isFinite(-double.infinity));
+    assert(isFinite(0.0));
+    assert(!isFinite(float.nan));
+    assert(isFinite(5.87747e-38f / 10));
+    assert(isFinite(1.11254e-307 / 10));
+    assert(isFinite(0.5));
+}
+
+/**
+ * Determines whether $(D_PARAM x) is $(B n)ot $(B a) $(B n)umber (NaN).
+ *
+ * Params:
+ *  F = Type of the floating point number.
+ *  x = Floating point number.
+ *
+ * Returns: $(D_KEYWORD true) if $(D_PARAM x) is not a number,
+ *          $(D_KEYWORD false) otherwise.
+ */
+bool isNaN(F)(F x)
+if (isFloatingPoint!F)
+{
+    FloatBits!F bits;
+    bits.floating = abs(x);
+
+    static if (ieeePrecision!F == IEEEPrecision.single
+            || ieeePrecision!F == IEEEPrecision.double_)
+    {
+        return bits.integral > bits.expMask;
+    }
+    else static if (ieeePrecision!F == IEEEPrecision.doubleExtended)
+    {
+        const maskedMantissa = bits.mantissa & bits.mantissaMask;
+        if ((bits.exp == bits.expMask && maskedMantissa != 0)
+         || ((bits.exp != 0) && (bits.mantissa < (1L << 63))))
+        {
+            return true;
+        }
+        return false;
+    }
+}
+
+///
+pure nothrow @safe @nogc unittest
+{
+    assert(isNaN(float.init));
+    assert(isNaN(double.init));
+    assert(isNaN(real.init));
+}
+
+/**
+ * Determines whether $(D_PARAM x) is a positive or negative infinity.
+ *
+ * Params:
+ *  F = Type of the floating point number.
+ *  x = Floating point number.
+ *
+ * Returns: $(D_KEYWORD true) if $(D_PARAM x) is infinity, $(D_KEYWORD false)
+ *          otherwise.
+ *
+ * See_Also: $(D_PSYMBOL isFinite).
+ */
+bool isInfinity(F)(F x)
+if (isFloatingPoint!F)
+{
+    FloatBits!F bits;
+    bits.floating = abs(x);
+    static if (ieeePrecision!F == IEEEPrecision.single
+            || ieeePrecision!F == IEEEPrecision.double_)
+    {
+        return bits.integral == bits.expMask;
+    }
+    else static if (ieeePrecision!F == IEEEPrecision.doubleExtended)
+    {
+        return (bits.exp == bits.expMask)
+            && ((bits.mantissa & bits.mantissaMask) == 0);
+    }
+}
+
+///
+pure nothrow @safe @nogc unittest
+{
+    assert(isInfinity(float.infinity));
+    assert(isInfinity(-float.infinity));
+    assert(isInfinity(double.infinity));
+    assert(isInfinity(-double.infinity));
+    assert(isInfinity(real.infinity));
+    assert(isInfinity(-real.infinity));
+}
+
+/**
+ * Determines whether $(D_PARAM x) is a denormilized number or not.
+
+ * Denormalized number is a number between `0` and `1` that cannot be
+ * represented as
+ *
+ * <pre>
+ * m*2<sup>e</sup>
+ * </pre>
+ *
+ * where $(I m) is the mantissa and $(I e) is an exponent that fits into the
+ * exponent field of the type $(D_PARAM F).
+ *
+ * `0` is neither normalized nor denormalized.
+ *
+ * Params:
+ *  F = Type of the floating point number.
+ *  x = Floating point number.
+ *
+ * Returns: $(D_KEYWORD true) if $(D_PARAM x) is a denormilized number,
+ *          $(D_KEYWORD false) otherwise.
+ *
+ * See_Also: $(D_PSYMBOL isNormal).
+ */
+bool isSubnormal(F)(F x)
+if (isFloatingPoint!F)
+{
+    FloatBits!F bits;
+    bits.floating = abs(x);
+    static if (ieeePrecision!F == IEEEPrecision.single)
+    {
+        return bits.integral < (1 << 23) && bits.integral > 0;
+    }
+    else static if (ieeePrecision!F == IEEEPrecision.double_)
+    {
+        return bits.integral < (1L << 52) && bits.integral > 0;
+    }
+    else static if (ieeePrecision!F == IEEEPrecision.doubleExtended)
+    {
+        return bits.exp == 0 && bits.mantissa != 0;
+    }
+}
+
+///
+pure nothrow @safe @nogc unittest
+{
+    assert(!isSubnormal(0.0f));
+    assert(!isSubnormal(float.nan));
+    assert(!isSubnormal(float.infinity));
+    assert(!isSubnormal(0.3f));
+    assert(isSubnormal(5.87747e-38f / 10));
+
+    assert(!isSubnormal(0.0));
+    assert(!isSubnormal(double.nan));
+    assert(!isSubnormal(double.infinity));
+    assert(!isSubnormal(1.4));
+    assert(isSubnormal(1.11254e-307 / 10));
+
+    assert(!isSubnormal(0.0L));
+    assert(!isSubnormal(real.nan));
+    assert(!isSubnormal(real.infinity));
+}
+
+/**
+ * Determines whether $(D_PARAM x) is a normilized number or not.
+
+ * Normalized number is a number that can be represented as
+ *
+ * <pre>
+ * m*2<sup>e</sup>
+ * </pre>
+ *
+ * where $(I m) is the mantissa and $(I e) is an exponent that fits into the
+ * exponent field of the type $(D_PARAM F).
+ *
+ * `0` is neither normalized nor denormalized.
+ *
+ * Params:
+ *  F = Type of the floating point number.
+ *  x = Floating point number.
+ *
+ * Returns: $(D_KEYWORD true) if $(D_PARAM x) is a normilized number,
+ *          $(D_KEYWORD false) otherwise.
+ *
+ * See_Also: $(D_PSYMBOL isSubnormal).
+ */
+bool isNormal(F)(F x)
+if (isFloatingPoint!F)
+{
+    static if (ieeePrecision!F == IEEEPrecision.single
+            || ieeePrecision!F == IEEEPrecision.double_)
+    {
+        FloatBits!F bits;
+        bits.floating = x;
+        bits.integral &= bits.expMask;
+        return bits.integral != 0 && bits.integral != bits.expMask;
+    }
+    else static if (ieeePrecision!F == IEEEPrecision.doubleExtended)
+    {
+        return classify(x) == FloatingPointClass.normal;
+    }
+}
+
+///
+pure nothrow @safe @nogc unittest
+{
+    assert(!isNormal(0.0f));
+    assert(!isNormal(float.nan));
+    assert(!isNormal(float.infinity));
+    assert(isNormal(0.3f));
+    assert(!isNormal(5.87747e-38f / 10));
+
+    assert(!isNormal(0.0));
+    assert(!isNormal(double.nan));
+    assert(!isNormal(double.infinity));
+    assert(isNormal(1.4));
+    assert(!isNormal(1.11254e-307 / 10));
+
+    assert(!isNormal(0.0L));
+    assert(!isNormal(real.nan));
+    assert(!isNormal(real.infinity));
+}
+
+/**
+ * Determines whether the sign bit of $(D_PARAM x) is set or not.
+ *
+ * If the sign bit, $(D_PARAM x) is a negative number, otherwise positive.
+ *
+ * Params:
+ *  F = Type of the floating point number.
+ *  x = Floating point number.
+ *
+ * Returns: $(D_KEYWORD true) if the sign bit of $(D_PARAM x) is set,
+ *          $(D_KEYWORD false) otherwise.
+ */
+bool signBit(F)(F x)
+if (isFloatingPoint!F)
+{
+    FloatBits!F bits;
+    bits.floating = x;
+    static if (ieeePrecision!F == IEEEPrecision.single)
+    {
+        return (bits.integral & (1 << 31)) != 0;
+    }
+    else static if (ieeePrecision!F == IEEEPrecision.double_)
+    {
+        return (bits.integral & (1L << 63)) != 0;
+    }
+    else static if (ieeePrecision!F == IEEEPrecision.doubleExtended)
+    {
+        return (bits.exp & (1 << 15)) != 0;
+    }
+}
+
+///
+pure nothrow @safe @nogc unittest
+{
+    assert(signBit(-1.0f));
+    assert(!signBit(1.0f));
+
+    assert(signBit(-1.0));
+    assert(!signBit(1.0));
+
+    assert(signBit(-1.0L));
+    assert(!signBit(1.0L));
+}
+
 /**
  * Computes $(D_PARAM x) to the power $(D_PARAM y) modulo $(D_PARAM z).
  *