@@ -264,7 +264,7 @@ class AnyIntView {
264264 return digits[size () - 1 ];
265265 }
266266 double top_double () const {
267- return size () > 1 ? (double )digits[size () - 1 ] + (double )digits[size () - 2 ] * ( 1.0 / Tr::Base)
267+ return size () > 1 ? (double )digits[size () - 1 ] + (double )digits[size () - 2 ] * Tr::InvBase
268268 : (double )digits[size () - 1 ];
269269 }
270270 bool is_odd_any () const {
@@ -314,8 +314,15 @@ class BigIntG {
314314 digits[0 ] = x;
315315 }
316316 BigIntG (Normalize, word_t x) : n(1 ) {
317- digits[0 ] = x;
318- normalize_bool ();
317+ if (x >= -Tr::Half && x < Tr::Half) {
318+ digits[0 ] = x;
319+ } else if (len <= 1 ) {
320+ digits[0 ] = x;
321+ normalize_bool ();
322+ } else {
323+ digits[0 ] = ((x + Tr::Half) & (Tr::Base - 1 )) - Tr::Half;
324+ digits[n++] = (x >> Tr::word_shift) + (digits[0 ] < 0 );
325+ }
319326 }
320327 BigIntG (const BigIntG& x) : n(x.n) {
321328 std::memcpy (digits, x.digits , n * sizeof (word_t ));
@@ -757,7 +764,7 @@ bool AnyIntView<Tr>::add_pow2_any(int exponent, int factor) {
757764 while (size () <= k) {
758765 digits[inc_size ()] = 0 ;
759766 }
760- digits[k] += (factor << dm.rem );
767+ digits[k] += (( word_t ) factor << dm.rem );
761768 return true ;
762769}
763770
@@ -1087,12 +1094,16 @@ int AnyIntView<Tr>::cmp_any(const AnyIntView<Tr>& yp) const {
10871094
10881095template <class Tr >
10891096int AnyIntView<Tr>::cmp_any(word_t y) const {
1090- if (size () > 1 ) {
1091- return top_word () < 0 ? -1 : 1 ;
1092- } else if (size () == 1 ) {
1097+ if (size () == 1 ) {
10931098 return digits[0 ] < y ? -1 : (digits[0 ] > y ? 1 : 0 );
1094- } else {
1099+ } else if (! size ()) {
10951100 return 0x80000000 ;
1101+ } else if (size () == 2 && (y >= Tr::Half || y < -Tr::Half)) {
1102+ word_t x0 = digits[0 ] & (Tr::Base - 1 ), y0 = y & (Tr::Base - 1 );
1103+ word_t x1 = digits[1 ] + (digits[0 ] >> Tr::word_shift), y1 = (y >> Tr::word_shift);
1104+ return x1 < y1 ? -1 : (x1 > y1 ? 1 : (x0 < y0 ? -1 : (x0 > y0 ? 1 : 0 )));
1105+ } else {
1106+ return top_word () < 0 ? -1 : 1 ;
10961107 }
10971108}
10981109
@@ -1312,17 +1323,14 @@ bool AnyIntView<Tr>::mod_div_any(const AnyIntView<Tr>& yp, AnyIntView<Tr>& quot,
13121323 if (k > quot.max_size ()) {
13131324 return invalidate_bool ();
13141325 }
1315- quot.set_size (max (k,1 ));
1316- for (int qi=0 ; qi< max (k,1 ); qi++) {
1317- quot.digits [qi]=0 ;
1318- }
1326+ quot.set_size (std::max (k, 1 ));
1327+ quot.digits [0 ] = 0 ;
13191328 } else {
13201329 if (k >= quot.max_size ()) {
13211330 return invalidate_bool ();
13221331 }
13231332 quot.set_size (k + 1 );
1324- double x_top = top_double ();
1325- word_t q = std::llrint (x_top * y_inv * Tr::InvBase);
1333+ word_t q = std::llrint (top_double () * y_inv * Tr::InvBase);
13261334 quot.digits [k] = q;
13271335 int i = yp.size () - 1 ;
13281336 word_t hi = 0 ;
@@ -1337,24 +1345,26 @@ bool AnyIntView<Tr>::mod_div_any(const AnyIntView<Tr>& yp, AnyIntView<Tr>& quot,
13371345 quot.digits [0 ] = 0 ;
13381346 }
13391347 while (--k >= 0 ) {
1340- double x_top = top_double ();
1341- word_t q = std::llrint (x_top * y_inv);
1348+ word_t q = std::llrint (top_double () * y_inv);
13421349 quot.digits [k] = q;
13431350 for (int i = yp.size () - 1 ; i >= 0 ; --i) {
13441351 Tr::sub_mul (&digits[k + i + 1 ], &digits[k + i], q, yp.digits [i]);
13451352 }
13461353 dec_size ();
13471354 digits[size () - 1 ] += (digits[size ()] << word_shift);
13481355 }
1349- if (size () >= yp.size ()) {
1350- assert (size () = = yp.size ());
1351- double x_top = top_double ( );
1352- double t = x_top * y_inv * Tr::InvBase;
1356+ if (size () >= yp.size () - 1 ) {
1357+ assert (size () < = yp.size ());
1358+ bool grow = ( size () < yp. size () );
1359+ double t = top_double () * y_inv * (grow ? Tr::InvBase * Tr::InvBase : Tr::InvBase) ;
13531360 if (round_mode >= 0 ) {
13541361 t += (round_mode ? 1 : 0.5 );
13551362 }
13561363 word_t q = std::llrint (std::floor (t));
13571364 if (q) {
1365+ if (grow) {
1366+ digits[inc_size ()] = 0 ;
1367+ }
13581368 for (int i = 0 ; i < size (); i++) {
13591369 digits[i] -= q * yp.digits [i];
13601370 }
@@ -1411,6 +1421,7 @@ bool AnyIntView<Tr>::mod_div_any(const AnyIntView<Tr>& yp, AnyIntView<Tr>& quot,
14111421 return normalize_bool_any ();
14121422}
14131423
1424+ // works for almost-normalized numbers (digits -Base+1 .. Base-1, top non-zero), result also almost-normalized
14141425template <class Tr >
14151426bool AnyIntView<Tr>::mod_pow2_any(int exponent) {
14161427 if (!is_valid ()) {
@@ -1462,25 +1473,21 @@ bool AnyIntView<Tr>::mod_pow2_any(int exponent) {
14621473 if (exponent >= max_size () * word_shift) {
14631474 return invalidate_bool ();
14641475 }
1465- if (q - word_shift >= 0 ) {
1476+ if (q - word_shift >= 0 ) { // original top digit was a non-zero multiple of Base, impossible(?)
14661477 digits[size () - 1 ] = 0 ;
14671478 digits[inc_size ()] = ((word_t )1 << (q - word_shift));
1468- }
1469- if (q - word_shift == -1 && size () < max_size () - 1 ) {
1479+ } else if (q - word_shift == -1 && size () < max_size ()) {
14701480 digits[size () - 1 ] = -Tr::Half;
14711481 digits[inc_size ()] = 1 ;
14721482 } else {
14731483 digits[size () - 1 ] = pow;
14741484 }
14751485 return true ;
1476- } else if (v >= Tr::Half) {
1477- if (size () == max_size () - 1 ) {
1478- return invalidate_bool ();
1479- } else {
1480- digits[size () - 1 ] = v | -Tr::Half;
1481- digits[inc_size ()] = ((word_t )1 << (q - word_shift));
1482- return true ;
1483- }
1486+ } else if (v >= Tr::Half && size () < max_size ()) {
1487+ word_t w = (((v >> (word_shift - 1 )) + 1 ) >> 1 );
1488+ digits[size () - 1 ] = v - (w << word_shift);
1489+ digits[inc_size ()] = w;
1490+ return true ;
14841491 } else {
14851492 digits[size () - 1 ] = v;
14861493 return true ;
0 commit comments