sebastienros.jint/Jint/Native/Number/Dtoa/Bignum.cs

#nullable disable

using System.Diagnostics;
using Jint.Runtime;

namespace Jint.Native.Number.Dtoa;

internal sealed class Bignum
{
    // 3584 = 128 * 28. We can represent 2^3584 > 10^1000 accurately.
    // This bignum can encode much bigger numbers, since it contains an
    // exponent.
    private const int kMaxSignificantBits = 3584;

    private const int kChunkSize = sizeof(uint) * 8;
    private const int kDoubleChunkSize = sizeof(ulong) * 8;

    // With bigit size of 28 we loose some bits, but a double still fits easily
    // into two chunks, and more importantly we can use the Comba multiplication.
    private const int kBigitSize = 28;
    private const uint kBigitMask = (1 << kBigitSize) - 1;

    // Every instance allocates kBigitLength chunks on the stack. Bignums cannot
    // grow. There are no checks if the stack-allocated space is sufficient.
    private const int kBigitCapacity = kMaxSignificantBits / kBigitSize;

    private readonly uint[] bigits_ = new uint[kBigitCapacity];

    // The Bignum's value equals value(bigits_) * 2^(exponent_ * kBigitSize).
    private int exponent_;
    private int used_digits_;

    private int BigitLength()
    {
        return used_digits_ + exponent_;
    }

    // Precondition: this/other < 16bit.
    public uint DivideModuloIntBignum(Bignum other)
    {
        Debug.Assert(IsClamped());
        Debug.Assert(other.IsClamped());
        Debug.Assert(other.used_digits_ > 0);

        // Easy case: if we have less digits than the divisor than the result is 0.
        // Note: this handles the case where this == 0, too.
        if (BigitLength() < other.BigitLength()) return 0;

        Align(other);

        uint result = 0;

        // Start by removing multiples of 'other' until both numbers have the same
        // number of digits.
        while (BigitLength() > other.BigitLength())
        {
            // This naive approach is extremely inefficient if the this divided other
            // might be big. This function is implemented for doubleToString where
            // the result should be small (less than 10).
            Debug.Assert(other.bigits_[other.used_digits_ - 1] >= (1 << kBigitSize) / 16);
            // Remove the multiples of the first digit.
            // Example this = 23 and other equals 9. -> Remove 2 multiples.
            result += bigits_[used_digits_ - 1];
            SubtractTimes(other, bigits_[used_digits_ - 1]);
        }

        Debug.Assert(BigitLength() == other.BigitLength());

        // Both bignums are at the same length now.
        // Since other has more than 0 digits we know that the access to
        // bigits_[used_digits_ - 1] is safe.
        var this_bigit = bigits_[used_digits_ - 1];
        var other_bigit = other.bigits_[other.used_digits_ - 1];

        if (other.used_digits_ == 1)
        {
            // Shortcut for easy (and common) case.
            uint quotient = this_bigit / other_bigit;
            bigits_[used_digits_ - 1] = this_bigit - other_bigit * quotient;
            result += quotient;
            Clamp();
            return result;
        }

        uint division_estimate = this_bigit / (other_bigit + 1);
        result += division_estimate;
        SubtractTimes(other, division_estimate);

        if (other_bigit * (division_estimate + 1) > this_bigit) return result;

        while (LessEqual(other, this))
        {
            SubtractBignum(other);
            result++;
        }

        return result;
    }

    void Align(Bignum other)
    {
        if (exponent_ > other.exponent_)
        {
            // If "X" represents a "hidden" digit (by the exponent) then we are in the
            // following case (a == this, b == other):
            // a:  aaaaaaXXXX   or a:   aaaaaXXX
            // b:     bbbbbbX      b: bbbbbbbbXX
            // We replace some of the hidden digits (X) of a with 0 digits.
            // a:  aaaaaa000X   or a:   aaaaa0XX
            int zero_digits = exponent_ - other.exponent_;
            ValidateCapacity(used_digits_ + zero_digits);
            for (int i = used_digits_ - 1; i >= 0; --i)
            {
                bigits_[i + zero_digits] = bigits_[i];
            }

            for (int i = 0; i < zero_digits; ++i)
            {
                bigits_[i] = 0;
            }

            used_digits_ += zero_digits;
            exponent_ -= zero_digits;
            Debug.Assert(used_digits_ >= 0);
            Debug.Assert(exponent_ >= 0);
        }
    }

    private static void ValidateCapacity(int size)
    {
        if (size > kBigitCapacity)
        {
            ExceptionHelper.ThrowInvalidOperationException();
        }
    }

    private void Clamp()
    {
        while (used_digits_ > 0 && bigits_[used_digits_ - 1] == 0) used_digits_--;
        if (used_digits_ == 0) exponent_ = 0;
    }


    private bool IsClamped()
    {
        return used_digits_ == 0 || bigits_[used_digits_ - 1] != 0;
    }


    private void Zero()
    {
        for (var i = 0; i < used_digits_; ++i) bigits_[i] = 0;
        used_digits_ = 0;
        exponent_ = 0;
    }

    // Guaranteed to lie in one Bigit.
    internal void AssignUInt16(uint value)
    {
        Debug.Assert(kBigitSize <= 8 * sizeof(uint));
        Zero();
        if (value == 0) return;

        ValidateCapacity(1);
        bigits_[0] = value;
        used_digits_ = 1;
    }

    internal void AssignUInt64(ulong value)
    {
        const int kUInt64Size = 64;

        Zero();
        if (value == 0) return;

        int needed_bigits = kUInt64Size / kBigitSize + 1;
        ValidateCapacity(needed_bigits);
        for (int i = 0; i < needed_bigits; ++i)
        {
            bigits_[i] = (uint) (value & kBigitMask);
            value = value >> kBigitSize;
        }

        used_digits_ = needed_bigits;
        Clamp();
    }


    internal void AssignBignum(Bignum other)
    {
        exponent_ = other.exponent_;
        for (int i = 0; i < other.used_digits_; ++i)
        {
            bigits_[i] = other.bigits_[i];
        }

        // Clear the excess digits (if there were any).
        for (int i = other.used_digits_; i < used_digits_; ++i)
        {
            bigits_[i] = 0;
        }

        used_digits_ = other.used_digits_;
    }


    void SubtractTimes(Bignum other, uint factor)
    {
#if DEBUG
        var a = new Bignum();
        var b = new Bignum();
        a.AssignBignum(this);
        b.AssignBignum(other);
        b.MultiplyByUInt32(factor);
        a.SubtractBignum(b);
#endif
        Debug.Assert(exponent_ <= other.exponent_);
        if (factor < 3)
        {
            for (int i = 0; i < factor; ++i)
            {
                SubtractBignum(other);
            }

            return;
        }

        uint borrow = 0;
        int exponent_diff = other.exponent_ - exponent_;
        for (int i = 0; i < other.used_digits_; ++i)
        {
            ulong product = factor * other.bigits_[i];
            ulong remove = borrow + product;
            uint difference = bigits_[i + exponent_diff] - (uint) (remove & kBigitMask);
            bigits_[i + exponent_diff] = difference & kBigitMask;
            borrow = (uint) ((difference >> (kChunkSize - 1)) + (remove >> kBigitSize));
        }

        for (int i = other.used_digits_ + exponent_diff; i < used_digits_; ++i)
        {
            if (borrow == 0) return;
            uint difference = bigits_[i] - borrow;
            bigits_[i] = difference & kBigitMask;
            borrow = difference >> (kChunkSize - 1);
        }

        Clamp();

#if DEBUG
        Debug.Assert(Equal(a, this));
#endif
    }


    void SubtractBignum(Bignum other)
    {
        Debug.Assert(IsClamped());
        Debug.Assert(other.IsClamped());
        // We require this to be bigger than other.
        Debug.Assert(LessEqual(other, this));

        Align(other);

        int offset = other.exponent_ - exponent_;
        uint borrow = 0;
        int i;
        for (i = 0; i < other.used_digits_; ++i)
        {
            Debug.Assert((borrow == 0) || (borrow == 1));
            uint difference = bigits_[i + offset] - other.bigits_[i] - borrow;
            bigits_[i + offset] = difference & kBigitMask;
            borrow = difference >> (kChunkSize - 1);
        }

        while (borrow != 0)
        {
            uint difference = bigits_[i + offset] - borrow;
            bigits_[i + offset] = difference & kBigitMask;
            borrow = difference >> (kChunkSize - 1);
            ++i;
        }

        Clamp();
    }

    internal static bool Equal(Bignum a, Bignum b)
    {
        return Compare(a, b) == 0;
    }

    internal static bool LessEqual(Bignum a, Bignum b)
    {
        return Compare(a, b) <= 0;
    }

    internal static bool Less(Bignum a, Bignum b)
    {
        return Compare(a, b) < 0;
    }

    // Returns a + b == c
    static bool PlusEqual(Bignum a, Bignum b, Bignum c)
    {
        return PlusCompare(a, b, c) == 0;
    }

    // Returns a + b <= c
    static bool PlusLessEqual(Bignum a, Bignum b, Bignum c)
    {
        return PlusCompare(a, b, c) <= 0;
    }

    // Returns a + b < c
    static bool PlusLess(Bignum a, Bignum b, Bignum c)
    {
        return PlusCompare(a, b, c) < 0;
    }

    uint BigitAt(int index)
    {
        if (index >= BigitLength()) return 0;
        if (index < exponent_) return 0;
        return bigits_[index - exponent_];
    }


    static int Compare(Bignum a, Bignum b)
    {
        Debug.Assert(a.IsClamped());
        Debug.Assert(b.IsClamped());
        int bigit_length_a = a.BigitLength();
        int bigit_length_b = b.BigitLength();
        if (bigit_length_a < bigit_length_b) return -1;
        if (bigit_length_a > bigit_length_b) return +1;
        for (int i = bigit_length_a - 1; i >= System.Math.Min(a.exponent_, b.exponent_); --i)
        {
            uint bigit_a = a.BigitAt(i);
            uint bigit_b = b.BigitAt(i);
            if (bigit_a < bigit_b) return -1;
            if (bigit_a > bigit_b) return +1;
            // Otherwise they are equal up to this digit. Try the next digit.
        }

        return 0;
    }


    internal static int PlusCompare(Bignum a, Bignum b, Bignum c)
    {
        Debug.Assert(a.IsClamped());
        Debug.Assert(b.IsClamped());
        Debug.Assert(c.IsClamped());
        if (a.BigitLength() < b.BigitLength())
        {
            return PlusCompare(b, a, c);
        }

        if (a.BigitLength() + 1 < c.BigitLength()) return -1;
        if (a.BigitLength() > c.BigitLength()) return +1;
        // The exponent encodes 0-bigits. So if there are more 0-digits in 'a' than
        // 'b' has digits, then the bigit-length of 'a'+'b' must be equal to the one
        // of 'a'.
        if (a.exponent_ >= b.BigitLength() && a.BigitLength() < c.BigitLength())
        {
            return -1;
        }

        uint borrow = 0;
        // Starting at min_exponent all digits are == 0. So no need to compare them.
        int min_exponent = System.Math.Min(System.Math.Min(a.exponent_, b.exponent_), c.exponent_);
        for (int i = c.BigitLength() - 1; i >= min_exponent; --i)
        {
            uint chunk_a = a.BigitAt(i);
            uint chunk_b = b.BigitAt(i);
            uint chunk_c = c.BigitAt(i);
            uint sum = chunk_a + chunk_b;
            if (sum > chunk_c + borrow)
            {
                return +1;
            }
            else
            {
                borrow = chunk_c + borrow - sum;
                if (borrow > 1) return -1;
                borrow <<= kBigitSize;
            }
        }

        if (borrow == 0) return 0;
        return -1;
    }

    internal void Times10()
    {
        MultiplyByUInt32(10);
    }

    internal void MultiplyByUInt32(uint factor)
    {
        if (factor == 1) return;
        if (factor == 0)
        {
            Zero();
            return;
        }

        if (used_digits_ == 0) return;

        // The product of a bigit with the factor is of size kBigitSize + 32.
        // Assert that this number + 1 (for the carry) fits into double chunk.
        Debug.Assert(kDoubleChunkSize >= kBigitSize + 32 + 1);
        ulong carry = 0;
        for (int i = 0; i < used_digits_; ++i)
        {
            ulong product = ((ulong) factor) * bigits_[i] + carry;
            bigits_[i] = (uint) (product & kBigitMask);
            carry = (product >> kBigitSize);
        }

        while (carry != 0)
        {
            ValidateCapacity(used_digits_ + 1);
            bigits_[used_digits_] = (uint) (carry & kBigitMask);
            used_digits_++;
            carry >>= kBigitSize;
        }
    }

    internal void MultiplyByUInt64(ulong factor)
    {
        if (factor == 1) return;
        if (factor == 0)
        {
            Zero();
            return;
        }
        Debug.Assert(kBigitSize < 32);
        ulong carry = 0;
        ulong low = factor & 0xFFFFFFFF;
        ulong high = factor >> 32;
        for (int i = 0; i < used_digits_; ++i)
        {
            ulong product_low = low * bigits_[i];
            ulong product_high = high * bigits_[i];
            ulong tmp = (carry & kBigitMask) + product_low;
            bigits_[i] = (uint) (tmp & kBigitMask);
            carry = (carry >> kBigitSize) + (tmp >> kBigitSize) +
                    (product_high << (32 - kBigitSize));
        }
        while (carry != 0)
        {
            ValidateCapacity(used_digits_ + 1);
            bigits_[used_digits_] = (uint) (carry & kBigitMask);
            used_digits_++;
            carry >>= kBigitSize;
        }
    }

    internal void ShiftLeft(int shift_amount)
    {
        if (used_digits_ == 0) return;
        exponent_ += shift_amount / kBigitSize;
        int local_shift = shift_amount % kBigitSize;
        ValidateCapacity(used_digits_ + 1);
        BigitsShiftLeft(local_shift);
    }

    void BigitsShiftLeft(int shift_amount)
    {
        Debug.Assert(shift_amount < kBigitSize);
        Debug.Assert(shift_amount >= 0);
        uint carry = 0;
        for (int i = 0; i < used_digits_; ++i)
        {
            uint new_carry = bigits_[i] >> (kBigitSize - shift_amount);
            bigits_[i] = ((bigits_[i] << shift_amount) + carry) & kBigitMask;
            carry = new_carry;
        }

        if (carry != 0)
        {
            bigits_[used_digits_] = carry;
            used_digits_++;
        }
    }


    internal void AssignPowerUInt16(uint baseValue, int power_exponent)
    {
        Debug.Assert(baseValue != 0);
        Debug.Assert(power_exponent >= 0);
        if (power_exponent == 0)
        {
            AssignUInt16(1);
            return;
        }

        Zero();
        int shifts = 0;
        // We expect baseValue to be in range 2-32, and most often to be 10.
        // It does not make much sense to implement different algorithms for counting
        // the bits.
        while ((baseValue & 1) == 0)
        {
            baseValue >>= 1;
            shifts++;
        }

        int bit_size = 0;
        uint tmp_base = baseValue;
        while (tmp_base != 0)
        {
            tmp_base >>= 1;
            bit_size++;
        }

        int final_size = bit_size * power_exponent;
        // 1 extra bigit for the shifting, and one for rounded final_size.
        ValidateCapacity(final_size / kBigitSize + 2);

        // Left to Right exponentiation.
        int mask = 1;
        while (power_exponent >= mask) mask <<= 1;

        // The mask is now pointing to the bit above the most significant 1-bit of
        // power_exponent.
        // Get rid of first 1-bit;
        mask >>= 2;
        ulong this_value = baseValue;

        bool delayed_multipliciation = false;
        const ulong max_32bits = 0xFFFFFFFF;
        while (mask != 0 && this_value <= max_32bits)
        {
            this_value = this_value * this_value;
            // Verify that there is enough space in this_value to perform the
            // multiplication.  The first bit_size bits must be 0.
            if ((power_exponent & mask) != 0)
            {
                ulong base_bits_mask = ~((((ulong) 1) << (64 - bit_size)) - 1);
                bool high_bits_zero = (this_value & base_bits_mask) == 0;
                if (high_bits_zero)
                {
                    this_value *= baseValue;
                }
                else
                {
                    delayed_multipliciation = true;
                }
            }

            mask >>= 1;
        }

        AssignUInt64(this_value);
        if (delayed_multipliciation)
        {
            MultiplyByUInt32(baseValue);
        }

        // Now do the same thing as a bignum.
        while (mask != 0)
        {
            Square();
            if ((power_exponent & mask) != 0)
            {
                MultiplyByUInt32(baseValue);
            }

            mask >>= 1;
        }

        // And finally add the saved shifts.
        ShiftLeft(shifts * power_exponent);
    }

    void Square()
    {
        Debug.Assert(IsClamped());
        int product_length = 2 * used_digits_;
        ValidateCapacity(product_length);

        // Comba multiplication: compute each column separately.
        // Example: r = a2a1a0 * b2b1b0.
        //    r =  1    * a0b0 +
        //        10    * (a1b0 + a0b1) +
        //        100   * (a2b0 + a1b1 + a0b2) +
        //        1000  * (a2b1 + a1b2) +
        //        10000 * a2b2
        //
        // In the worst case we have to accumulate nb-digits products of digit*digit.
        //
        // Assert that the additional number of bits in a DoubleChunk are enough to
        // sum up used_digits of Bigit*Bigit.
        if ((1 << (2 * (kChunkSize - kBigitSize))) <= used_digits_)
        {
            ExceptionHelper.ThrowNotImplementedException();
        }

        ulong accumulator = 0;
        // First shift the digits so we don't overwrite them.
        int copy_offset = used_digits_;
        for (int i = 0; i < used_digits_; ++i)
        {
            bigits_[copy_offset + i] = bigits_[i];
        }

        // We have two loops to avoid some 'if's in the loop.
        for (int i = 0; i < used_digits_; ++i)
        {
            // Process temporary digit i with power i.
            // The sum of the two indices must be equal to i.
            int bigit_index1 = i;
            int bigit_index2 = 0;
            // Sum all of the sub-products.
            while (bigit_index1 >= 0)
            {
                uint chunk1 = bigits_[copy_offset + bigit_index1];
                uint chunk2 = bigits_[copy_offset + bigit_index2];
                accumulator += (ulong) chunk1 * chunk2;
                bigit_index1--;
                bigit_index2++;
            }

            bigits_[i] = (uint) accumulator & kBigitMask;
            accumulator >>= kBigitSize;
        }

        for (int i = used_digits_; i < product_length; ++i)
        {
            int bigit_index1 = used_digits_ - 1;
            int bigit_index2 = i - bigit_index1;
            // Invariant: sum of both indices is again equal to i.
            // Inner loop runs 0 times on last iteration, emptying accumulator.
            while (bigit_index2 < used_digits_)
            {
                uint chunk1 = bigits_[copy_offset + bigit_index1];
                uint chunk2 = bigits_[copy_offset + bigit_index2];
                accumulator += (ulong) chunk1 * chunk2;
                bigit_index1--;
                bigit_index2++;
            }

            // The overwritten bigits_[i] will never be read in further loop iterations,
            // because bigit_index1 and bigit_index2 are always greater
            // than i - used_digits_.
            bigits_[i] = (uint) accumulator & kBigitMask;
            accumulator >>= kBigitSize;
        }

        // Since the result was guaranteed to lie inside the number the
        // accumulator must be 0 now.
        Debug.Assert(accumulator == 0);

        // Don't forget to update the used_digits and the exponent.
        used_digits_ = product_length;
        exponent_ *= 2;
        Clamp();
    }
}