botan.math.bigint.divide source code

1 /**
2 * Division
3 * 
4 * Copyright:
5 * (C) 1999-2007 Jack Lloyd
6 * (C) 2014-2015 Etienne Cimon
7 *
8 * License:
9 * Botan is released under the Simplified BSD License (see LICENSE.md)
10 */
11 module botan.math.bigint.divide;
12 
13 import botan.constants;
14 import botan.math.bigint.bigint;
15 import botan.math.mp.mp_core;
16 import botan.constants;
17 import std.algorithm : max;
18 /**
19 * BigInt Division
20 * Params:
21 *  x = an integer
22 *  y_arg = a non-zero integer
23 *  q = will be set to x / y
24 *  r = will be set to x % y
25 */
26 void divide(const(BigInt)* x, const(BigInt)* y_arg, BigInt* q_out, BigInt* r_out, bool trace = false)
27 {
28     /*
29     * Solve x = q * y + r
30     */
31     if (y_arg.isZero())
32         throw new BigInt.DivideByZero();
33     BigInt y = y_arg.clone;
34     const size_t y_words = y.sigWords();
35     
36     BigInt r = x.clone;
37     BigInt q = 0;
38     
39     r.setSign(BigInt.Positive);
40     y.setSign(BigInt.Positive);
41     
42     int compare = r.cmp(y);
43     
44     if (compare == 0)
45     {
46         q = 1;
47         r = 0;
48     }
49     else if (compare > 0)
50     {
51         size_t shifts = 0;
52         word y_top = y.wordAt(y.sigWords()-1);
53         while (y_top < MP_WORD_TOP_BIT) { y_top <<= 1; ++shifts; }
54         y <<= shifts;
55         r <<= shifts;
56         
57         const size_t n = r.sigWords() - 1, t = y_words - 1;
58         if (n < t)
59             throw new InternalError("BigInt division word sizes");
60         
61         q.growTo(n - t + 1);
62         
63         word* q_words = q.mutablePtr();
64         
65         if (n <= t)
66         {
67             while (r > y) { r -= y; ++(q); }
68             r >>= shifts;
69             goto ret;
70         }
71         
72         BigInt shifted_y = y << (MP_WORD_BITS * (n-t));
73         
74         while (r >= shifted_y) { r -= shifted_y; q_words[n-t] += 1; }
75         
76         for (size_t j = n; j != t; --j)
77         {
78             const word x_j0  = r.wordAt(j);
79             const word x_j1 = r.wordAt(j-1);
80             const word x_j2 = r.wordAt(j-2);
81             const word y_t0  = y.wordAt(t);
82             const word y_t1  = y.wordAt(t-1);
83 
84             word qjt = (x_j0 == y_t0) ? MP_WORD_MAX : bigint_divop(x_j0, x_j1, y_t0);
85             
86             while (divisionCheck(qjt, y_t0, y_t1, x_j0, x_j1, x_j2))
87             {
88                 qjt -= 1;
89             }
90 
91             shifted_y >>= BOTAN_MP_WORD_BITS;
92             // Now shifted_y == y << (BOTAN_MP_WORD_BITS * (j-t-1))
93             r -= shifted_y * qjt;
94             if (r.isNegative()) {
95                 // overcorrected
96                 qjt -= 1;
97                 r += shifted_y;
98                 
99             }
100             q_words[j-t-1] = qjt;
101             
102         }
103         r >>= shifts;
104     }
105     
106 ret:
107     signFixup(x, y_arg, &q, &r);
108     *r_out = r.move();
109     *q_out = q.move();
110 }
111 
112 private:
113 /*
114 * Handle signed operands, if necessary
115 */
116 void signFixup(const(BigInt)* x, const(BigInt)* y, BigInt* q, BigInt* r)
117 {
118     if (x.sign() == BigInt.Negative)
119     {
120         q.flipSign();
121         if (r.isNonzero()) { --*q; *r = y.abs() - *r; }
122     }
123     if (y.sign() == BigInt.Negative)
124         q.flipSign();
125 }
126 
127 bool divisionCheck(word q, word y2, word y1, word x3, word x2, word x1)
128 {
129     // Compute (y3,y2,y1) = (y2,y1) * q
130     
131     word y3 = 0;
132     y1 = word_madd2(q, y1, &y3);
133     y2 = word_madd2(q, y2, &y3);
134 
135     // Return (y3,y2,y1) >? (x3,x2,x1)
136     
137     if (y3 > x3) return true;
138     if (y3 < x3) return false;
139     
140     if (y2 > x2) return true;
141     if (y2 < x2) return false;
142     
143     if (y1 > x1) return true;
144     if (y1 < x1) return false;
145     
146     return false;
147 }