From a0a57b074f79dbe5e01745da9cfbd08a52fa8cfd Mon Sep 17 00:00:00 2001
From: Iain Buclaw <ibuclaw@gcc.gnu.org>
Date: Fri, 16 Nov 2018 21:17:33 +0000
Subject: [PATCH] libphobos: Add IEEE quadruple support to
 core.internal.convert

Backport from upstream druntime 2.083 for AArch64.

Reviewed-on: https://github.com/dlang/druntime/pull/2257

From-SVN: r266222
---
 libphobos/libdruntime/core/internal/convert.d | 135 +++++++++++++++---
 1 file changed, 114 insertions(+), 21 deletions(-)

diff --git a/libphobos/libdruntime/core/internal/convert.d b/libphobos/libdruntime/core/internal/convert.d
index a612707e2672..c4745b6f5a7f 100644
--- a/libphobos/libdruntime/core/internal/convert.d
+++ b/libphobos/libdruntime/core/internal/convert.d
@@ -10,13 +10,36 @@
 module core.internal.convert;
 import core.internal.traits : Unqual;
 
+/+
+A @nogc function can allocate memory during CTFE.
++/
+@nogc nothrow pure @trusted
+private ubyte[] ctfe_alloc()(size_t n)
+{
+    if (!__ctfe)
+    {
+        assert(0, "CTFE only");
+    }
+    else
+    {
+        static ubyte[] alloc(size_t x) nothrow pure
+        {
+            if (__ctfe) // Needed to prevent _d_newarray from appearing in compiled prorgam.
+                return new ubyte[x];
+            else
+                assert(0);
+        }
+        return (cast(ubyte[] function(size_t) @nogc nothrow pure) &alloc)(n);
+    }
+}
+
 @trusted pure nothrow
-const(ubyte)[] toUbyte(T)(ref T val) if (is(Unqual!T == float) || is(Unqual!T == double) || is(Unqual!T == real) ||
+const(ubyte)[] toUbyte(T)(const ref T val) if (is(Unqual!T == float) || is(Unqual!T == double) || is(Unqual!T == real) ||
                                         is(Unqual!T == ifloat) || is(Unqual!T == idouble) || is(Unqual!T == ireal))
 {
     static const(ubyte)[] reverse_(const(ubyte)[] arr)
     {
-        ubyte[] buff = new ubyte[arr.length];
+        ubyte[] buff = ctfe_alloc(arr.length);
         foreach (k, v; arr)
         {
             buff[$-k-1] = v;
@@ -31,17 +54,35 @@ const(ubyte)[] toUbyte(T)(ref T val) if (is(Unqual!T == float) || is(Unqual!T ==
         uint exp = parsed.exponent;
         uint sign = parsed.sign;
 
-        ubyte[T.sizeof] buff;
+        ubyte[] buff = ctfe_alloc(T.sizeof);
         size_t off_bytes = 0;
         size_t off_bits  = 0;
+        // Quadruples won't fit in one ulong, so check for that.
+        enum mantissaMax = FloatTraits!T.MANTISSA < ulong.sizeof*8 ?
+                           FloatTraits!T.MANTISSA : ulong.sizeof*8;
 
-        for (; off_bytes < FloatTraits!T.MANTISSA/8; ++off_bytes)
+        for (; off_bytes < mantissaMax/8; ++off_bytes)
         {
             buff[off_bytes] = cast(ubyte)mantissa;
             mantissa >>= 8;
         }
-        off_bits = FloatTraits!T.MANTISSA%8;
-        buff[off_bytes] = cast(ubyte)mantissa;
+
+        static if (floatFormat!T == FloatFormat.Quadruple)
+        {
+            ulong mantissa2 = parsed.mantissa2;
+            off_bytes--; // go back one, since mantissa only stored data in 56
+                         // bits, ie 7 bytes
+            for(; off_bytes < FloatTraits!T.MANTISSA/8; ++off_bytes)
+            {
+                buff[off_bytes] = cast(ubyte)mantissa2;
+                mantissa2 >>= 8;
+            }
+        }
+        else
+        {
+            off_bits = FloatTraits!T.MANTISSA%8;
+            buff[off_bytes] = cast(ubyte)mantissa;
+        }
 
         for (size_t i=0; i<FloatTraits!T.EXPONENT/8; ++i)
         {
@@ -60,7 +101,7 @@ const(ubyte)[] toUbyte(T)(ref T val) if (is(Unqual!T == float) || is(Unqual!T ==
 
         version (LittleEndian)
         {
-            return buff.dup;
+            return buff;
         }
         else
         {
@@ -83,8 +124,8 @@ private Float parse(bool is_denormalized = false, T)(T x) if (is(Unqual!T == ifl
 private Float parse(bool is_denormalized = false, T:real)(T x_) if (floatFormat!T != FloatFormat.Real80)
 {
     Unqual!T x = x_;
-    assert(floatFormat!T != FloatFormat.DoubleDouble && floatFormat!T != FloatFormat.Quadruple,
-           "doubledouble and quadruple float formats are not supported in CTFE");
+    static assert(floatFormat!T != FloatFormat.DoubleDouble,
+           "doubledouble float format not supported in CTFE");
     if (x is cast(T)0.0) return FloatTraits!T.ZERO;
     if (x is cast(T)-0.0) return FloatTraits!T.NZERO;
     if (x is T.nan) return FloatTraits!T.NAN;
@@ -103,7 +144,7 @@ private Float parse(bool is_denormalized = false, T:real)(T x_) if (floatFormat!
         if (is_denormalized)
             return Float(0, 0, sign);
         else
-            return Float(denormalizedMantissa(x), 0, sign);
+            return denormalizedMantissa(x, sign);
     }
 
     x2 /= binPow2(e);
@@ -111,9 +152,30 @@ private Float parse(bool is_denormalized = false, T:real)(T x_) if (floatFormat!
     static if (!is_denormalized)
         x2 -= 1.0;
 
-    x2 *=  2UL<<(FloatTraits!T.MANTISSA);
-    ulong mant = shiftrRound(cast(ulong)x2);
-    return Float(mant, exp, sign);
+    static if (floatFormat!T == FloatFormat.Quadruple)
+    {
+        // Store the 112-bit mantissa in two ulongs, specifically the lower 56
+        // bits of each, with the most significant bits in mantissa2. There's
+        // an edge case exposed by the labeled test below, where only a subnormal
+        // with the highest bit set being the 57th bit will "overflow" to the
+        // 57th bit in mantissa2 with the following logic, but that special case
+        // is handled by an additional check in denormalizedMantissa for
+        // Quadruples below.
+
+        x2 *=  2UL<<(FloatTraits!T.MANTISSA - (ulong.sizeof - 1)*8 - 1);
+        ulong mant2 = cast(ulong) x2;
+        x2 -= mant2;
+
+        x2 *=  2UL<<((ulong.sizeof - 1)*8 - 1);
+        ulong mant = cast(ulong) x2;
+        return Float(mant, exp, sign, mant2);
+    }
+    else
+    {
+        x2 *=  2UL<<(FloatTraits!T.MANTISSA);
+        ulong mant = shiftrRound(cast(ulong)x2);
+        return Float(mant, exp, sign);
+    }
 }
 
 @safe pure nothrow
@@ -151,7 +213,7 @@ private Float parse(bool _ = false, T:real)(T x_) if (floatFormat!T == FloatForm
     uint exp = cast(uint)(e + EXPONENT_MED);
     if (!exp)
     {
-        return Float(denormalizedMantissa(x), 0, sign);
+        return denormalizedMantissa(x, sign);
     }
     int pow = (FloatTraits!T.MANTISSA-1-e);
     x *=  binPow2((pow / EXPONENT_MED)*EXPONENT_MED); //To avoid overflow in 2.0L ^^ pow
@@ -165,6 +227,7 @@ private struct Float
     ulong mantissa;
     uint exponent;
     uint sign;
+    ulong mantissa2;
 }
 
 private template FloatTraits(T) if (floatFormat!T == FloatFormat.Float)
@@ -215,14 +278,14 @@ private template FloatTraits(T) if (floatFormat!T == FloatFormat.DoubleDouble) /
     enum NINF     = Float(0, 0x7ff, 1);
 }
 
-private template FloatTraits(T) if (floatFormat!T == FloatFormat.Quadruple) //Unsupported in CTFE
+private template FloatTraits(T) if (floatFormat!T == FloatFormat.Quadruple)
 {
     enum EXPONENT = 15;
     enum MANTISSA = 112;
     enum ZERO     = Float(0, 0, 0);
     enum NZERO    = Float(0, 0, 1);
-    enum NAN      = Float(-1, 0x7fff, 0);
-    enum NNAN     = Float(-1, 0x7fff, 1);
+    enum NAN      = Float(0, 0x7fff, 0, 0x80000000000000UL);
+    enum NNAN     = Float(0, 0x7fff, 1, 0x80000000000000UL);
     enum INF      = Float(0, 0x7fff, 0);
     enum NINF     = Float(0, 0x7fff, 1);
 }
@@ -291,21 +354,49 @@ private uint binLog2(T)(T x)
 }
 
 @safe pure nothrow
-private ulong denormalizedMantissa(T)(T x) if (floatFormat!T == FloatFormat.Real80)
+private Float denormalizedMantissa(T)(T x, uint sign) if (floatFormat!T == FloatFormat.Real80)
 {
     x *= 2.0L^^FloatTraits!T.MANTISSA;
     auto fl = parse(x);
     uint pow = FloatTraits!T.MANTISSA - fl.exponent + 1;
-    return fl.mantissa >> pow;
+    return Float(fl.mantissa >> pow, 0, sign);
 }
 
 @safe pure nothrow
-private ulong denormalizedMantissa(T)(T x) if (floatFormat!T != FloatFormat.Real80)
+private Float denormalizedMantissa(T)(T x, uint sign)
+    if (floatFormat!T == FloatFormat.Float || floatFormat!T == FloatFormat.Double)
 {
     x *= 2.0L^^FloatTraits!T.MANTISSA;
     auto fl = parse!true(x);
     ulong mant = fl.mantissa >> (FloatTraits!T.MANTISSA - fl.exponent);
-    return shiftrRound(mant);
+    return Float(shiftrRound(mant), 0, sign);
+}
+
+@safe pure nothrow
+private Float denormalizedMantissa(T)(T x, uint sign) if (floatFormat!T == FloatFormat.Quadruple)
+{
+    x *= 2.0L^^FloatTraits!T.MANTISSA;
+    auto fl = parse!true(x);
+    uint offset = FloatTraits!T.MANTISSA - fl.exponent + 1;
+    enum mantissaSize = (ulong.sizeof - 1) * 8;
+
+    if (offset < mantissaSize)
+    {   // Create a new mantissa ulong with the trailing mantissa2 bits that
+        // need to be shifted into mantissa, by shifting the needed bits left,
+        // zeroing out the first byte, and then ORing it with mantissa shifted
+        // right by offset.
+
+        ulong shiftedMantissa = ((fl.mantissa2 << (mantissaSize - offset)) &
+                                 0x00FFFFFFFFFFFFFFUL) | fl.mantissa >> offset;
+        return Float(shiftedMantissa, 0, sign, fl.mantissa2 >> offset);
+    }
+    else if (offset > mantissaSize)
+        return Float(fl.mantissa2 >> offset - mantissaSize , 0, sign, 0);
+    else
+        // Handle special case mentioned in parse() above by zeroing out the
+        // 57'th bit of mantissa2, "shifting" it into mantissa, and setting the
+        // first bit of mantissa2.
+        return Float(fl.mantissa2 & 0x00FFFFFFFFFFFFFFUL , 0, sign, 1);
 }
 
 version (unittest)
@@ -403,6 +494,8 @@ version (unittest)
 
         testNumberConvert!("real.min_normal/2");
         testNumberConvert!("real.min_normal/2UL^^63");
+        // check subnormal storage edge case for Quadruple
+        testNumberConvert!("real.min_normal/2UL^^56");
         //testNumberConvert!("real.min_normal/19"); // XGDC: ct[0] == 0, rt[0] == 27
         //testNumberConvert!("real.min_normal/17"); // XGDC: ct[0= == 128, rt[0] == 136