Reduce array allocations within BigIntegerCalculator

- An allocation for BigInteger.DivRem with small divisor is unnecessary,
  which gets fixed. BigIngeger handles these cases better now.
- The square and multiply code submitted with dotnet#1436 used stackalloc,
  which led to stack overflows for huge numbers. With dotnet#1618 these has
  been changed to ordinary array allocations, which put some pressure
  on the managed heap. Thus, a hybrid solution should be better.

src/System.Runtime.Numerics/src/System/Numerics/BigInteger.cs

@@ -945,10 +945,10 @@ public static BigInteger DivRem(BigInteger dividend, BigInteger divisor, out Big
 
             if (trivialDivisor)
             {
-                uint[] rest;
+                uint rest;
                 uint[] bits = BigIntegerCalculator.Divide(dividend._bits, NumericsHelpers.Abs(divisor._sign), out rest);
 
-                remainder = new BigInteger(rest, dividend._sign < 0);
+                remainder = dividend._sign < 0 ? -1 * (long)rest : rest;
                 return new BigInteger(bits, (dividend._sign < 0) ^ (divisor._sign < 0));
             }
 
@@ -1065,12 +1065,12 @@ public static BigInteger ModPow(BigInteger value, BigInteger exponent, BigIntege
 
             if (trivialModulus)
             {
-                long bits = trivialValue && trivialExponent ? BigIntegerCalculator.Pow(NumericsHelpers.Abs(value._sign), NumericsHelpers.Abs(exponent._sign), NumericsHelpers.Abs(modulus._sign)) :
+                uint bits = trivialValue && trivialExponent ? BigIntegerCalculator.Pow(NumericsHelpers.Abs(value._sign), NumericsHelpers.Abs(exponent._sign), NumericsHelpers.Abs(modulus._sign)) :
                             trivialValue ? BigIntegerCalculator.Pow(NumericsHelpers.Abs(value._sign), exponent._bits, NumericsHelpers.Abs(modulus._sign)) :
                             trivialExponent ? BigIntegerCalculator.Pow(value._bits, NumericsHelpers.Abs(exponent._sign), NumericsHelpers.Abs(modulus._sign)) :
                             BigIntegerCalculator.Pow(value._bits, exponent._bits, NumericsHelpers.Abs(modulus._sign));
 
-                return value._sign < 0 && !exponent.IsEven ? -1 * bits : bits;
+                return value._sign < 0 && !exponent.IsEven ? -1 * (long)bits : bits;
             }
             else
             {
@@ -1718,8 +1718,8 @@ private static BigInteger Subtract(uint[] leftBits, int leftSign, uint[] rightBi
 
             if (trivialDivisor)
             {
-                long bits = BigIntegerCalculator.Remainder(dividend._bits, NumericsHelpers.Abs(divisor._sign));
-                return dividend._sign < 0 ? -1 * bits : bits;
+                uint bits = BigIntegerCalculator.Remainder(dividend._bits, NumericsHelpers.Abs(divisor._sign));
+                return dividend._sign < 0 ? -1 * (long)bits : bits;
             }
 
             if (dividend._bits.Length < divisor._bits.Length)

src/System.Runtime.Numerics/src/System/Numerics/BigIntegerCalculator.DivRem.cs

@@ -6,7 +6,7 @@ namespace System.Numerics
     internal static partial class BigIntegerCalculator
     {
         public static uint[] Divide(uint[] left, uint right,
-                                    out uint[] remainder)
+                                    out uint remainder)
         {
             Debug.Assert(left != null);
             Debug.Assert(left.Length >= 1);
@@ -24,7 +24,7 @@ public static uint[] Divide(uint[] left, uint right,
                 quotient[i] = (uint)(value / right);
                 carry = value % right;
             }
-            remainder = new uint[] { (uint)carry };
+            remainder = (uint)carry;
 
             return quotient;
         }

src/System.Runtime.Numerics/src/System/Numerics/BigIntegerCalculator.SquMul.cs

@@ -26,6 +26,7 @@ public unsafe static uint[] Square(uint[] value)
 
         // Mutable for unit testing...
         private static int SquareThreshold = 32;
+        private static int AllocationThreshold = 4096;
 
         [SecuritySafeCritical]
         private unsafe static void Square(uint* value, int valueLength,
@@ -111,9 +112,12 @@ private unsafe static void Square(uint* value, int valueLength,
 
                 int foldLength = valueHighLength + 1;
                 int coreLength = foldLength + foldLength;
-                fixed (uint* fold = new uint[foldLength],
-                             core = new uint[coreLength])
+
+                if (coreLength < AllocationThreshold)
                 {
+                    uint* fold = stackalloc uint[foldLength];
+                    uint* core = stackalloc uint[coreLength];
+
                     // ... compute z_a = a_1 + a_0 (call it fold...)
                     Add(valueHigh, valueHighLength,
                         valueLow, valueLowLength,
@@ -129,6 +133,27 @@ private unsafe static void Square(uint* value, int valueLength,
                     // ... and finally merge the result! :-)
                     AddSelf(bits + n, bitsLength - n, core, coreLength);
                 }
+                else
+                {
+                    fixed (uint* fold = new uint[foldLength],
+                                 core = new uint[coreLength])
+                    {
+                        // ... compute z_a = a_1 + a_0 (call it fold...)
+                        Add(valueHigh, valueHighLength,
+                            valueLow, valueLowLength,
+                            fold, foldLength);
+
+                        // ... compute z_1 = z_a * z_a - z_0 - z_2
+                        Square(fold, foldLength,
+                               core, coreLength);
+                        SubtractCore(bitsHigh, bitsHighLength,
+                                     bitsLow, bitsLowLength,
+                                     core, coreLength);
+
+                        // ... and finally merge the result! :-)
+                        AddSelf(bits + n, bitsLength - n, core, coreLength);
+                    }
+                }
             }
         }
 
@@ -270,10 +295,13 @@ private unsafe static void Multiply(uint* left, int leftLength,
                 int leftFoldLength = leftHighLength + 1;
                 int rightFoldLength = rightHighLength + 1;
                 int coreLength = leftFoldLength + rightFoldLength;
-                fixed (uint* leftFold = new uint[leftFoldLength],
-                             rightFold = new uint[rightFoldLength],
-                             core = new uint[coreLength])
+
+                if (coreLength < AllocationThreshold)
                 {
+                    uint* leftFold = stackalloc uint[leftFoldLength];
+                    uint* rightFold = stackalloc uint[rightFoldLength];
+                    uint* core = stackalloc uint[coreLength];
+
                     // ... compute z_a = a_1 + a_0 (call it fold...)
                     Add(leftHigh, leftHighLength,
                         leftLow, leftLowLength,
@@ -295,6 +323,34 @@ private unsafe static void Multiply(uint* left, int leftLength,
                     // ... and finally merge the result! :-)
                     AddSelf(bits + n, bitsLength - n, core, coreLength);
                 }
+                else
+                {
+                    fixed (uint* leftFold = new uint[leftFoldLength],
+                                 rightFold = new uint[rightFoldLength],
+                                 core = new uint[coreLength])
+                    {
+                        // ... compute z_a = a_1 + a_0 (call it fold...)
+                        Add(leftHigh, leftHighLength,
+                            leftLow, leftLowLength,
+                            leftFold, leftFoldLength);
+
+                        // ... compute z_b = b_1 + b_0 (call it fold...)
+                        Add(rightHigh, rightHighLength,
+                            rightLow, rightLowLength,
+                            rightFold, rightFoldLength);
+
+                        // ... compute z_1 = z_a * z_b - z_0 - z_2
+                        Multiply(leftFold, leftFoldLength,
+                                 rightFold, rightFoldLength,
+                                 core, coreLength);
+                        SubtractCore(bitsHigh, bitsHighLength,
+                                     bitsLow, bitsLowLength,
+                                     core, coreLength);
+
+                        // ... and finally merge the result! :-)
+                        AddSelf(bits + n, bitsLength - n, core, coreLength);
+                    }
+                }
             }
         }
 

src/System.Runtime.Numerics/tests/BigInteger/multiply.cs

@@ -38,8 +38,16 @@ public static void RunMultiply_TwoLargeBigIntegers()
         public static void RunMultiply_TwoLargeBigIntegers_Threshold()
         {
             // Again, with lower threshold
-            BigIntTools.Utils.RunWithFakeThreshold("SquareThreshold", 8, RunMultiply_TwoLargeBigIntegers);
-            BigIntTools.Utils.RunWithFakeThreshold("MultiplyThreshold", 8, RunMultiply_TwoLargeBigIntegers);
+            BigIntTools.Utils.RunWithFakeThreshold("SquareThreshold", 8, () =>
+                BigIntTools.Utils.RunWithFakeThreshold("MultiplyThreshold", 8, RunMultiply_TwoLargeBigIntegers)
+            );
+
+            // Again, with lower threshold
+            BigIntTools.Utils.RunWithFakeThreshold("SquareThreshold", 8, () =>
+                BigIntTools.Utils.RunWithFakeThreshold("MultiplyThreshold", 8, () =>
+                    BigIntTools.Utils.RunWithFakeThreshold("AllocationThreshold", 8, RunMultiply_TwoLargeBigIntegers)
+                )
+            );
         }
 
         [Fact]