1 /* TomsFastMath, a fast ISO C bignum library.
3 * This project is meant to fill in where LibTomMath
4 * falls short. That is speed ;-)
6 * This project is public domain and free for all purposes.
8 * Tom St Denis, tomstdenis@gmail.com
10 #include <tfm_private.h>
12 static int fp_invmod_slow (fp_int * a, fp_int * b, fp_int * c)
14 fp_int x, y, u, v, A, B, C, D;
17 /* b cannot be negative */
18 if (b->sign == FP_NEG || fp_iszero(b) == 1) {
23 fp_init(&x); fp_init(&y);
24 fp_init(&u); fp_init(&v);
25 fp_init(&A); fp_init(&B);
26 fp_init(&C); fp_init(&D);
29 if ((res = fp_mod(a, b, &x)) != FP_OKAY) {
34 /* 2. [modified] if x,y are both even then return an error! */
35 if (fp_iseven (&x) == 1 && fp_iseven (&y) == 1) {
39 /* 3. u=x, v=y, A=1, B=0, C=0,D=1 */
46 /* 4. while u is even do */
47 while (fp_iseven (&u) == 1) {
51 /* 4.2 if A or B is odd then */
52 if (fp_isodd (&A) == 1 || fp_isodd (&B) == 1) {
53 /* A = (A+y)/2, B = (B-x)/2 */
57 /* A = A/2, B = B/2 */
62 /* 5. while v is even do */
63 while (fp_iseven (&v) == 1) {
67 /* 5.2 if C or D is odd then */
68 if (fp_isodd (&C) == 1 || fp_isodd (&D) == 1) {
69 /* C = (C+y)/2, D = (D-x)/2 */
73 /* C = C/2, D = D/2 */
78 /* 6. if u >= v then */
79 if (fp_cmp (&u, &v) != FP_LT) {
80 /* u = u - v, A = A - C, B = B - D */
85 /* v - v - u, C = C - A, D = D - B */
91 /* if not zero goto step 4 */
92 if (fp_iszero (&u) == 0)
95 /* now a = C, b = D, gcd == g*v */
97 /* if v != 1 then there is no inverse */
98 if (fp_cmp_d (&v, 1) != FP_EQ) {
103 while (fp_cmp_d(&C, 0) == FP_LT) {
108 while (fp_cmp_mag(&C, b) != FP_LT) {
112 /* C is now the inverse */
117 /* c = 1/a (mod b) for odd b only */
118 int fp_invmod(fp_int *a, fp_int *b, fp_int *c)
120 fp_int x, y, u, v, B, D;
123 /* 2. [modified] b must be odd */
124 if (fp_iseven (b) == FP_YES) {
125 return fp_invmod_slow(a,b,c);
128 /* init all our temps */
129 fp_init(&x); fp_init(&y);
130 fp_init(&u); fp_init(&v);
131 fp_init(&B); fp_init(&D);
133 /* x == modulus, y == value to invert */
136 /* we need y = |a| */
139 /* 3. u=x, v=y, A=1, B=0, C=0,D=1 */
145 /* 4. while u is even do */
146 while (fp_iseven (&u) == FP_YES) {
150 /* 4.2 if B is odd then */
151 if (fp_isodd (&B) == FP_YES) {
158 /* 5. while v is even do */
159 while (fp_iseven (&v) == FP_YES) {
163 /* 5.2 if D is odd then */
164 if (fp_isodd (&D) == FP_YES) {
172 /* 6. if u >= v then */
173 if (fp_cmp (&u, &v) != FP_LT) {
174 /* u = u - v, B = B - D */
178 /* v - v - u, D = D - B */
183 /* if not zero goto step 4 */
184 if (fp_iszero (&u) == FP_NO) {
188 /* now a = C, b = D, gcd == g*v */
190 /* if v != 1 then there is no inverse */
191 if (fp_cmp_d (&v, 1) != FP_EQ) {
195 /* b is now the inverse */
197 while (D.sign == FP_NEG) {