diff --git a/source/mir/bignum/fp.d b/source/mir/bignum/fp.d index b60595d3..b9a382ca 100644 --- a/source/mir/bignum/fp.d +++ b/source/mir/bignum/fp.d @@ -547,18 +547,20 @@ struct Fp(uint size) return ret; } - ret = Fp!newSize(this.coefficient, true); - ret.sign = this.sign; - static if (newSize < size) + UInt!size coefficient = this.coefficient; + int shift; + // subnormal + + if (this.exponent == this.exponent.min) { - // underflow - if (this.exponent == this.exponent.min && !ret.coefficient) - { - ret.exponent = 0; - return ret; - } + shift = cast(int)coefficient.ctlz; + coefficient <<= shift; } + ret = Fp!newSize(coefficient, true); + ret.exponent -= shift; + ret.sign = this.sign; + import mir.checkedint: adds; /// overflow bool overflow; @@ -566,14 +568,16 @@ struct Fp(uint size) if (_expect(overflow, false)) { // overflow - static if (newSize < size) + if (this.exponent > 0) { - assert(this.exponent > 0); + ret.exponent = ret.exponent.max; + ret.coefficient = 0u; } // underflow else { - assert(this.exponent < 0); + ret.coefficient >>= ret.exponent - exponent.min; + ret.exponent = ret.coefficient ? ret.exponent.min : 0; } } return ret; diff --git a/source/mir/math/numeric.d b/source/mir/math/numeric.d index 6a5c7035..ff378d2f 100644 --- a/source/mir/math/numeric.d +++ b/source/mir/math/numeric.d @@ -593,8 +593,7 @@ unittest static assert(is(typeof(factorial(33)) == Fp!128)); static assert(is(typeof(factorial!256(33)) == Fp!256)); static immutable double f33 = 8.68331761881188649551819440128e+36; - import mir.format: text; - static assert(cast(double) factorial(33) == f33, (cast(double) factorial(33)).text); + static assert(approxEqual(cast(double) factorial(33), f33)); assert(cast(double) factorial(0) == 1); assert(cast(double) factorial(0, 100) == 1);