|
2
|
1 /* LibTomMath, multiple-precision integer library -- Tom St Denis |
|
|
2 * |
|
|
3 * LibTomMath is a library that provides multiple-precision |
|
|
4 * integer arithmetic as well as number theoretic functionality. |
|
|
5 * |
|
|
6 * The library was designed directly after the MPI library by |
|
|
7 * Michael Fromberger but has been written from scratch with |
|
|
8 * additional optimizations in place. |
|
|
9 * |
|
|
10 * The library is free for all purposes without any express |
|
|
11 * guarantee it works. |
|
|
12 * |
|
|
13 * Tom St Denis, [email protected], http://math.libtomcrypt.org |
|
|
14 */ |
|
|
15 #include <tommath.h> |
|
|
16 |
|
|
17 /* computes the jacobi c = (a | n) (or Legendre if n is prime) |
|
|
18 * HAC pp. 73 Algorithm 2.149 |
|
|
19 */ |
|
|
20 int mp_jacobi (mp_int * a, mp_int * p, int *c) |
|
|
21 { |
|
|
22 mp_int a1, p1; |
|
|
23 int k, s, r, res; |
|
|
24 mp_digit residue; |
|
|
25 |
|
|
26 /* if p <= 0 return MP_VAL */ |
|
|
27 if (mp_cmp_d(p, 0) != MP_GT) { |
|
|
28 return MP_VAL; |
|
|
29 } |
|
|
30 |
|
|
31 /* step 1. if a == 0, return 0 */ |
|
|
32 if (mp_iszero (a) == 1) { |
|
|
33 *c = 0; |
|
|
34 return MP_OKAY; |
|
|
35 } |
|
|
36 |
|
|
37 /* step 2. if a == 1, return 1 */ |
|
|
38 if (mp_cmp_d (a, 1) == MP_EQ) { |
|
|
39 *c = 1; |
|
|
40 return MP_OKAY; |
|
|
41 } |
|
|
42 |
|
|
43 /* default */ |
|
|
44 s = 0; |
|
|
45 |
|
|
46 /* step 3. write a = a1 * 2**k */ |
|
|
47 if ((res = mp_init_copy (&a1, a)) != MP_OKAY) { |
|
|
48 return res; |
|
|
49 } |
|
|
50 |
|
|
51 if ((res = mp_init (&p1)) != MP_OKAY) { |
|
|
52 goto __A1; |
|
|
53 } |
|
|
54 |
|
|
55 /* divide out larger power of two */ |
|
|
56 k = mp_cnt_lsb(&a1); |
|
|
57 if ((res = mp_div_2d(&a1, k, &a1, NULL)) != MP_OKAY) { |
|
|
58 goto __P1; |
|
|
59 } |
|
|
60 |
|
|
61 /* step 4. if e is even set s=1 */ |
|
|
62 if ((k & 1) == 0) { |
|
|
63 s = 1; |
|
|
64 } else { |
|
|
65 /* else set s=1 if p = 1/7 (mod 8) or s=-1 if p = 3/5 (mod 8) */ |
|
|
66 residue = p->dp[0] & 7; |
|
|
67 |
|
|
68 if (residue == 1 || residue == 7) { |
|
|
69 s = 1; |
|
|
70 } else if (residue == 3 || residue == 5) { |
|
|
71 s = -1; |
|
|
72 } |
|
|
73 } |
|
|
74 |
|
|
75 /* step 5. if p == 3 (mod 4) *and* a1 == 3 (mod 4) then s = -s */ |
|
|
76 if ( ((p->dp[0] & 3) == 3) && ((a1.dp[0] & 3) == 3)) { |
|
|
77 s = -s; |
|
|
78 } |
|
|
79 |
|
|
80 /* if a1 == 1 we're done */ |
|
|
81 if (mp_cmp_d (&a1, 1) == MP_EQ) { |
|
|
82 *c = s; |
|
|
83 } else { |
|
|
84 /* n1 = n mod a1 */ |
|
|
85 if ((res = mp_mod (p, &a1, &p1)) != MP_OKAY) { |
|
|
86 goto __P1; |
|
|
87 } |
|
|
88 if ((res = mp_jacobi (&p1, &a1, &r)) != MP_OKAY) { |
|
|
89 goto __P1; |
|
|
90 } |
|
|
91 *c = s * r; |
|
|
92 } |
|
|
93 |
|
|
94 /* done */ |
|
|
95 res = MP_OKAY; |
|
|
96 __P1:mp_clear (&p1); |
|
|
97 __A1:mp_clear (&a1); |
|
|
98 return res; |
|
|
99 } |
