30  * each stored in a 32-bit slot (top bit is zero) in little-endian
33  * some cases, the top word is allowed to have a 32th bit.
63  *   if neg = 1, then -m <= a < 0
64  *   if neg = 0, then 0 <= a < 2*m
66  * If neg = 0, then the top word of a[] may use 32 bits.
68  * Also, modulus m must be odd.
71 finish_mod(uint32_t *a, size_t len, const uint32_t *m, uint32_t neg)  in finish_mod()  argument
77 	 * First pass: compare a (assumed nonnegative) with m.  in finish_mod()
79 	 * subtracting m must yield a value less than 2^31, since we  in finish_mod()
80 	 * assumed that a < 2*m.  in finish_mod()
87 		mw = m[k];  in finish_mod()
93 	 *   if neg = 1, then we must add m (regardless of cc)  in finish_mod()
94 	 *   if neg = 0 and cc = 0, then we must subtract m  in finish_mod()
104 		mw = (m[k] ^ xm) & ym;  in finish_mod()
123  * to contain an extra 32th bit.
146 		 *   |cca| <= 2^32 - 1  in co_reduce()
148 		 *   |za| <= (2^31-1)*(2^32) + (2^32-1) = 2^63 - 1  in co_reduce()
150 		 * Thus, the new value of cca is such that |cca| <= 2^32 - 1.  in co_reduce()
168 #define M   ((uint64_t)1 << 32)  in co_reduce()  macro
171 		tta = (tta ^ M) - M;  in co_reduce()
172 		ttb = (ttb ^ M) - M;  in co_reduce()
175 #undef M  in co_reduce()
189  *   a <- (a*pa+b*pb)/(2^31) mod m
190  *   b <- (a*qa+b*qb)/(2^31) mod m
192  * m0i is equal to -1/m[0] mod 2^31.
196  * to contain an extra 32th bit.
201 	const uint32_t *m, uint32_t m0i)  in co_reduce_mod()  argument
223 			+ m[k] * (uint64_t)fa + (uint64_t)cca;  in co_reduce_mod()
225 			+ m[k] * (uint64_t)fb + (uint64_t)ccb;  in co_reduce_mod()
231 #define M   ((uint64_t)1 << 32)  in co_reduce_mod()  macro
234 		tta = (tta ^ M) - M;  in co_reduce_mod()
235 		ttb = (ttb ^ M) - M;  in co_reduce_mod()
238 #undef M  in co_reduce_mod()
245 	 *   -m <= a < 2*m  in co_reduce_mod()
246 	 *   -m <= b < 2*m  in co_reduce_mod()
248 	 * The top word of 'a' and 'b' may have a 32-th bit set.  in co_reduce_mod()
251 	finish_mod(a, len, m, (uint32_t)((uint64_t)cca >> 63));  in co_reduce_mod()
252 	finish_mod(b, len, m, (uint32_t)((uint64_t)ccb >> 63));  in co_reduce_mod()
257 br_i31_moddiv(uint32_t *x, const uint32_t *y, const uint32_t *m, uint32_t m0i,  in br_i31_moddiv()  argument
264 	 *   a * x = y * u mod m  in br_i31_moddiv()
265 	 *   b * x = y * v mod m  in br_i31_moddiv()
270 	 *   b = m  in br_i31_moddiv()
276 	 *   - If a is even, then a <- a/2 and u <- u/2 mod m.  in br_i31_moddiv()
277 	 *   - Otherwise, if b is even, then b <- b/2 and v <- v/2 mod m.  in br_i31_moddiv()
278 	 *   - Otherwise, if a > b, then a <- (a-b)/2 and u <- (u-v)/2 mod m.  in br_i31_moddiv()
279 	 *   - Otherwise, b <- (b-a)/2 and v <- (v-u)/2 mod m.  in br_i31_moddiv()
282 	 * to GCD(y,m); the modular division succeeds if that value is 1.  in br_i31_moddiv()
287 	 * if m has bit length k bits, then 2k-2 steps are sufficient.  in br_i31_moddiv()
290 	 * Though complexity is quadratic in the size of m, the bit-by-bit  in br_i31_moddiv()
311 	 * and b become equal, they must be both odd (since m is odd,  in br_i31_moddiv()
321 	len = (m[0] + 31) >> 5;  in br_i31_moddiv()
327 	memcpy(b, m + 1, len * sizeof *m);  in br_i31_moddiv()
334 	for (num = ((m[0] - (m[0] >> 5)) << 1) + 30; num >= 30; num -= 30) {  in br_i31_moddiv()
472 		co_reduce_mod(u, v, len, pa, pb, qa, qb, m + 1, m0i);  in br_i31_moddiv()