File tree Expand file tree Collapse file tree 2 files changed +17
-14
lines changed Expand file tree Collapse file tree 2 files changed +17
-14
lines changed Original file line number Diff line number Diff line change @@ -547,33 +547,37 @@ struct Fp(uint size)
547547 return ret;
548548 }
549549
550- ret = Fp! newSize(this .coefficient, true );
551- ret.sign = this .sign;
552- static if (newSize < size)
550+ UInt! size coefficient = this .coefficient;
551+ int shift;
552+ // subnormal
553+
554+ if (this .exponent == this .exponent.min)
553555 {
554- // underflow
555- if (this .exponent == this .exponent.min && ! ret.coefficient)
556- {
557- ret.exponent = 0 ;
558- return ret;
559- }
556+ shift = cast (int )coefficient.ctlz;
557+ coefficient <<= shift;
560558 }
561559
560+ ret = Fp! newSize(coefficient, true );
561+ ret.exponent -= shift;
562+ ret.sign = this .sign;
563+
562564 import mir.checkedint: adds;
563565 // / overflow
564566 bool overflow;
565567 ret.exponent = adds(ret.exponent, this .exponent, overflow);
566568 if (_expect(overflow, false ))
567569 {
568570 // overflow
569- static if (newSize < size )
571+ if (this .exponent > 0 )
570572 {
571- assert (this .exponent > 0 );
573+ ret.exponent = ret.exponent.max;
574+ ret.coefficient = 0u ;
572575 }
573576 // underflow
574577 else
575578 {
576- assert (this .exponent < 0 );
579+ ret.coefficient >>= ret.exponent - exponent.min;
580+ ret.exponent = ret.coefficient ? ret.exponent.min : 0 ;
577581 }
578582 }
579583 return ret;
Original file line number Diff line number Diff line change @@ -593,8 +593,7 @@ unittest
593593 static assert (is (typeof (factorial(33 )) == Fp! 128 ));
594594 static assert (is (typeof (factorial! 256 (33 )) == Fp! 256 ));
595595 static immutable double f33 = 8.68331761881188649551819440128e+36 ;
596- import mir.format: text;
597- static assert (cast (double ) factorial(33 ) == f33, (cast (double ) factorial(33 )).text);
596+ static assert (approxEqual(cast (double ) factorial(33 ), f33));
598597
599598 assert (cast (double ) factorial(0 ) == 1 );
600599 assert (cast (double ) factorial(0 , 100 ) == 1 );
You can’t perform that action at this time.
0 commit comments