@@ -444,6 +444,95 @@ def IEEEContext(bits, /):
444444
445445##### Decimal class #######################################################
446446
447+ # Observation: For all q >= 0 and a >= 1, q < 10**a iff len(str(q)) <= a.
448+ #
449+ # The constants below are used to speed-up "q < 10 ** a" checks to avoid
450+ # computing len(str(q)) as much as possible. Those speed-ups are based on
451+ # the following claims.
452+ #
453+ # See https://github.com/python/cpython/issues/140036 for details.
454+
455+ _LOG_10_BASE_2_LO = float .fromhex ('0x1.a934f0979a371p+1' )
456+ assert pow (2 , _LOG_10_BASE_2_LO ) < 10
457+
458+ _LOG_10_BASE_2_HI = float .fromhex ('0x1.a934f0979a372p+1' )
459+ assert pow (2 , _LOG_10_BASE_2_HI ) > 10
460+
461+
462+ def _tento (n ):
463+ """Compute 10 ** n with 1 base-5 exponentiation and 1 bit-shift."""
464+ return (5 ** n ) << n
465+
466+
467+ def _is_less_than_pow10a_use_str (q , a ):
468+ """Try to efficiently check len(str(q)) <= a, or equivalently q < 10**a.
469+
470+ If the comparison cannot be obtained from q.bit_length(),
471+ then str(q) is explicitly computed and may raise ValueError.
472+
473+ Return (len(str(q)) <= a, None) or (len(str(q)) <= a, str(q)).
474+ """
475+ if q .bit_length () < a * _LOG_10_BASE_2_LO :
476+ # Claim: If 0 < z <= log2(10) and q.bit_length() < a*z, then q < 10**a.
477+ # Proof: By contradiction, q >= 10**a. By definition,
478+ # log2(q) >= a*log2(10) >= a*z > q.bit_length().
479+ # In particular, q > 2**q.bit_length(), which is impossible.
480+
481+ # assert q < 10 ** a
482+ return True , None
483+ elif q .bit_length () >= 1 + a * _LOG_10_BASE_2_HI :
484+ # Claim: If z > log2(10) and q.bit_length() >= 1 + a*z, then q > 10**a.
485+ # Proof: Since q >= 2**(q.bit_length()-1), we have
486+ # q >= 2**(q.bit_length()-1) >= 2**(a*z) > 2**(a*log2(10)) = 10**a.
487+
488+ # assert q > 10 ** a
489+ return False , None
490+ # Handle cases that fail due to floating point precision loss
491+ # when computing _LOG_10_BASE_2_LO and _LOG_10_BASE_2_HI, or
492+ # that cannot be distinguished with (q.bit_length(), a) only.
493+ #
494+ # For instance, (q1, a) = (95, 2) and (q2, a) = (105, 2) produce
495+ # different results but q1.bit_length() == q2.bit_length() == 7.
496+ str_q = str (q ) # can raise a ValueError
497+ is_valid = len (str_q ) <= a
498+ return is_valid , str_q
499+
500+
501+ def _is_less_than_pow10a (q , a , * , exact = True , ulp_order = 20 ):
502+ """Check that len(str(q)) <= a without computing str(q).
503+
504+ When *exact* is false, computing len(str(q)) is replaced by f(q):
505+
506+ f(q) = floor(log10(q) + ulp(log10(q)) * ulp_order + 1.0)
507+
508+ Most of the time, f(q) = len(str(q)) but in some cases, it may
509+ happen that f(q) > len(str(q)).
510+
511+ When *exact* is true, computing len(str(q)) requires one bigint
512+ exponentiation that only depends on q.
513+ """
514+ if q < 10 :
515+ return a >= 1
516+
517+ z = _math .log10 (q )
518+ t = _math .ulp (z ) * ulp_order
519+
520+ if exact :
521+ intlo = int (z - t )
522+ inthi = int (z + t )
523+ diff = inthi - intlo
524+ assert diff in (0 , 1 )
525+ if diff == 1 :
526+ lo = _tento (inthi ) # may be slow
527+ if q < lo :
528+ inthi -= 1
529+ assert q >= (lo // 10 )
530+ ndigits = inthi + 1
531+ else :
532+ ndigits = int (z + t + 1.0 )
533+ return ndigits <= a
534+
535+
447536# Do not subclass Decimal from numbers.Real and do not register it as such
448537# (because Decimals are not interoperable with floats). See the notes in
449538# numbers.py for more detail.
@@ -1357,8 +1446,13 @@ def _divide(self, other, context):
13571446 else :
13581447 op2 .int *= 10 ** (op2 .exp - op1 .exp )
13591448 q , r = divmod (op1 .int , op2 .int )
1360- if q < 10 ** context .prec :
1361- return (_dec_from_triple (sign , str (q ), 0 ),
1449+ is_valid , str_q = _is_less_than_pow10a_use_str (q , context .prec )
1450+ if is_valid :
1451+ if str_q is None :
1452+ str_q = str (q )
1453+ # assert q < 10 ** context.prec
1454+ # assert len(str(q)) <= context.prec
1455+ return (_dec_from_triple (sign , str_q , 0 ),
13621456 _dec_from_triple (self ._sign , str (r ), ideal_exp ))
13631457
13641458 # Here the quotient is too large to be representable
@@ -1513,7 +1607,9 @@ def remainder_near(self, other, context=None):
15131607 r -= op2 .int
15141608 q += 1
15151609
1516- if q >= 10 ** context .prec :
1610+ if not _is_less_than_pow10a (q , context .prec ):
1611+ # assert q >= 10 ** context.prec
1612+ # assert len(str(q)) > context.prec
15171613 return context ._raise_error (DivisionImpossible )
15181614
15191615 # result has same sign as self unless r is negative
0 commit comments