• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1 // Copyright 2021 the V8 project authors. All rights reserved.
2 // Use of this source code is governed by a BSD-style license that can be
3 // found in the LICENSE file.
4 
5 // "Schoolbook" division. This is loosely based on Go's implementation
6 // found at https://golang.org/src/math/big/nat.go, licensed as follows:
7 //
8 // Copyright 2009 The Go Authors. All rights reserved.
9 // Use of this source code is governed by a BSD-style
10 // license that can be found in the LICENSE file [1].
11 //
12 // [1] https://golang.org/LICENSE
13 
14 #include <limits>
15 
16 #include "src/bigint/bigint-internal.h"
17 #include "src/bigint/digit-arithmetic.h"
18 #include "src/bigint/div-helpers.h"
19 #include "src/bigint/util.h"
20 #include "src/bigint/vector-arithmetic.h"
21 
22 namespace v8 {
23 namespace bigint {
24 
25 // Computes Q(uotient) and remainder for A/b, such that
26 // Q = (A - remainder) / b, with 0 <= remainder < b.
27 // If Q.len == 0, only the remainder will be returned.
28 // Q may be the same as A for an in-place division.
DivideSingle(RWDigits Q,digit_t * remainder,Digits A,digit_t b)29 void ProcessorImpl::DivideSingle(RWDigits Q, digit_t* remainder, Digits A,
30                                  digit_t b) {
31   DCHECK(b != 0);
32   DCHECK(A.len() > 0);
33   *remainder = 0;
34   int length = A.len();
35   if (Q.len() != 0) {
36     if (A[length - 1] >= b) {
37       DCHECK(Q.len() >= A.len());
38       for (int i = length - 1; i >= 0; i--) {
39         Q[i] = digit_div(*remainder, A[i], b, remainder);
40       }
41       for (int i = length; i < Q.len(); i++) Q[i] = 0;
42     } else {
43       DCHECK(Q.len() >= A.len() - 1);
44       *remainder = A[length - 1];
45       for (int i = length - 2; i >= 0; i--) {
46         Q[i] = digit_div(*remainder, A[i], b, remainder);
47       }
48       for (int i = length - 1; i < Q.len(); i++) Q[i] = 0;
49     }
50   } else {
51     for (int i = length - 1; i >= 0; i--) {
52       digit_div(*remainder, A[i], b, remainder);
53     }
54   }
55 }
56 
57 // Z += X. Returns the "carry" (0 or 1) after adding all of X's digits.
InplaceAdd(RWDigits Z,Digits X)58 inline digit_t InplaceAdd(RWDigits Z, Digits X) {
59   return AddAndReturnCarry(Z, Z, X);
60 }
61 
62 // Z -= X. Returns the "borrow" (0 or 1) after subtracting all of X's digits.
InplaceSub(RWDigits Z,Digits X)63 inline digit_t InplaceSub(RWDigits Z, Digits X) {
64   return SubtractAndReturnBorrow(Z, Z, X);
65 }
66 
67 // Returns whether (factor1 * factor2) > (high << kDigitBits) + low.
ProductGreaterThan(digit_t factor1,digit_t factor2,digit_t high,digit_t low)68 bool ProductGreaterThan(digit_t factor1, digit_t factor2, digit_t high,
69                         digit_t low) {
70   digit_t result_high;
71   digit_t result_low = digit_mul(factor1, factor2, &result_high);
72   return result_high > high || (result_high == high && result_low > low);
73 }
74 
75 #if DEBUG
QLengthOK(Digits Q,Digits A,Digits B)76 bool QLengthOK(Digits Q, Digits A, Digits B) {
77   // If A's top B.len digits are greater than or equal to B, then the division
78   // result will be greater than A.len - B.len, otherwise it will be that
79   // difference. Intuitively: 100/10 has 2 digits, 100/11 has 1.
80   if (GreaterThanOrEqual(Digits(A, A.len() - B.len(), B.len()), B)) {
81     return Q.len() >= A.len() - B.len() + 1;
82   }
83   return Q.len() >= A.len() - B.len();
84 }
85 #endif
86 
87 // Computes Q(uotient) and R(emainder) for A/B, such that
88 // Q = (A - R) / B, with 0 <= R < B.
89 // Both Q and R are optional: callers that are only interested in one of them
90 // can pass the other with len == 0.
91 // If Q is present, its length must be at least A.len - B.len + 1.
92 // If R is present, its length must be at least B.len.
93 // See Knuth, Volume 2, section 4.3.1, Algorithm D.
DivideSchoolbook(RWDigits Q,RWDigits R,Digits A,Digits B)94 void ProcessorImpl::DivideSchoolbook(RWDigits Q, RWDigits R, Digits A,
95                                      Digits B) {
96   DCHECK(B.len() >= 2);        // Use DivideSingle otherwise.
97   DCHECK(A.len() >= B.len());  // No-op otherwise.
98   DCHECK(Q.len() == 0 || QLengthOK(Q, A, B));
99   DCHECK(R.len() == 0 || R.len() >= B.len());
100   // The unusual variable names inside this function are consistent with
101   // Knuth's book, as well as with Go's implementation of this algorithm.
102   // Maintaining this consistency is probably more useful than trying to
103   // come up with more descriptive names for them.
104   const int n = B.len();
105   const int m = A.len() - n;
106 
107   // In each iteration, {qhatv} holds {divisor} * {current quotient digit}.
108   // "v" is the book's name for {divisor}, "qhat" the current quotient digit.
109   ScratchDigits qhatv(n + 1);
110 
111   // D1.
112   // Left-shift inputs so that the divisor's MSB is set. This is necessary
113   // to prevent the digit-wise divisions (see digit_div call below) from
114   // overflowing (they take a two digits wide input, and return a one digit
115   // result).
116   ShiftedDigits b_normalized(B);
117   B = b_normalized;
118   // U holds the (continuously updated) remaining part of the dividend, which
119   // eventually becomes the remainder.
120   ScratchDigits U(A.len() + 1);
121   LeftShift(U, A, b_normalized.shift());
122 
123   // D2.
124   // Iterate over the dividend's digits (like the "grad school" algorithm).
125   // {vn1} is the divisor's most significant digit.
126   digit_t vn1 = B[n - 1];
127   for (int j = m; j >= 0; j--) {
128     // D3.
129     // Estimate the current iteration's quotient digit (see Knuth for details).
130     // {qhat} is the current quotient digit.
131     digit_t qhat = std::numeric_limits<digit_t>::max();
132     // {ujn} is the dividend's most significant remaining digit.
133     digit_t ujn = U[j + n];
134     if (ujn != vn1) {
135       // {rhat} is the current iteration's remainder.
136       digit_t rhat = 0;
137       // Estimate the current quotient digit by dividing the most significant
138       // digits of dividend and divisor. The result will not be too small,
139       // but could be a bit too large.
140       qhat = digit_div(ujn, U[j + n - 1], vn1, &rhat);
141 
142       // Decrement the quotient estimate as needed by looking at the next
143       // digit, i.e. by testing whether
144       // qhat * v_{n-2} > (rhat << kDigitBits) + u_{j+n-2}.
145       digit_t vn2 = B[n - 2];
146       digit_t ujn2 = U[j + n - 2];
147       while (ProductGreaterThan(qhat, vn2, rhat, ujn2)) {
148         qhat--;
149         digit_t prev_rhat = rhat;
150         rhat += vn1;
151         // v[n-1] >= 0, so this tests for overflow.
152         if (rhat < prev_rhat) break;
153       }
154     }
155 
156     // D4.
157     // Multiply the divisor with the current quotient digit, and subtract
158     // it from the dividend. If there was "borrow", then the quotient digit
159     // was one too high, so we must correct it and undo one subtraction of
160     // the (shifted) divisor.
161     if (qhat == 0) {
162       qhatv.Clear();
163     } else {
164       MultiplySingle(qhatv, B, qhat);
165     }
166     digit_t c = InplaceSub(U + j, qhatv);
167     if (c != 0) {
168       c = InplaceAdd(U + j, B);
169       U[j + n] = U[j + n] + c;
170       qhat--;
171     }
172 
173     if (Q.len() != 0) {
174       if (j >= Q.len()) {
175         DCHECK(qhat == 0);
176       } else {
177         Q[j] = qhat;
178       }
179     }
180   }
181   if (R.len() != 0) {
182     RightShift(R, U, b_normalized.shift());
183   }
184   // If Q has extra storage, clear it.
185   for (int i = m + 1; i < Q.len(); i++) Q[i] = 0;
186 }
187 
188 }  // namespace bigint
189 }  // namespace v8
190