5 * Compute Jacobi symbol
7 * (c) 1999 Straylight/Edgeware
10 /*----- Licensing notice --------------------------------------------------*
12 * This file is part of Catacomb.
14 * Catacomb is free software; you can redistribute it and/or modify
15 * it under the terms of the GNU Library General Public License as
16 * published by the Free Software Foundation; either version 2 of the
17 * License, or (at your option) any later version.
19 * Catacomb is distributed in the hope that it will be useful,
20 * but WITHOUT ANY WARRANTY; without even the implied warranty of
21 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
22 * GNU Library General Public License for more details.
24 * You should have received a copy of the GNU Library General Public
25 * License along with Catacomb; if not, write to the Free
26 * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
30 /*----- Header files ------------------------------------------------------*/
34 /*----- Main code ---------------------------------------------------------*/
36 /* --- @mp_jacobi@ --- *
38 * Arguments: @mp *a@ = an integer
39 * @mp *n@ = another integer
41 * Returns: @-1@, @0@ or @1@ -- the Jacobi symbol %$J(a, n)$%.
43 * Use: Computes the Kronecker symbol %$\jacobi{a}{n}$%. If @n@ is
44 * prime, this is the Legendre symbol and is equal to 1 if and
45 * only if @a@ is a quadratic residue mod @n@. The result is
46 * zero if and only if @a@ and @n@ have a common factor greater
49 * If @n@ is composite, then this computes the Kronecker symbol
51 * %$\jacobi{a}{n}=\jacobi{a}{u}\prod_i\jacobi{a}{p_i}^{e_i}$%
53 * where %$n = u p_0^{e_0} \ldots p_{n-1}^{e_{n-1}}$% is the
54 * prime factorization of %$n$%. The missing bits are:
56 * * %$\jacobi{a}{1} = 1$%;
57 * * %$\jacobi{a}{-1} = 1$% if @a@ is negative, or 1 if
59 * * %$\jacobi{a}{0} = 0$%;
60 * * %$\jacobi{a}{2}$ is 0 if @a@ is even, 1 if @a@ is
61 * congruent to 1 or 7 (mod 8), or %$-1$% otherwise.
63 * If %$n$% is positive and odd, then this is the Jacobi
64 * symbol. (The Kronecker symbol is a consistant domain
65 * extension; the Jacobi symbol was implemented first, and the
69 int mp_jacobi(mp
*a
, mp
*n
)
74 /* --- Handle zero specially --- *
76 * I can't find any specific statement for what to do when %$n = 0$%; PARI
77 * opts to set %$\jacobi{\pm1}{0} = \pm 1$% and %$\jacobi{a}{0} = 0$% for
82 if (MP_EQ(a
, MP_ONE
)) return (+1);
83 else if (MP_EQ(a
, MP_MONE
)) return (-1);
87 /* --- Deal with powers of two --- *
89 * This implicitly takes a copy of %$n$%. Copy %$a$% at the same time to
90 * make cleanup easier.
94 n
= mp_odd(MP_NEW
, n
, &p2
);
99 } else if ((p2
& 1) && ((a
->v
[0] & 7) == 3 || (a
->v
[0] & 7) == 5))
103 /* --- Deal with negative %$n$% --- */
111 /* --- Check for unit %$n$% --- */
113 if (MP_EQ(n
, MP_ONE
))
116 /* --- Reduce %$a$% modulo %$n$% --- */
118 if (MP_NEGP(a
) || MP_CMP(a
, >=, n
))
121 /* --- Main recursive mess, flattened out into something nice --- */
127 /* --- Some simple special cases --- */
135 /* --- Main case with powers of two --- */
137 a
= mp_odd(a
, a
, &e
);
139 if ((e
& 1) && (nn
== 3 || nn
== 5))
141 if (MP_LEN(a
) == 1 && a
->v
[0] == 1)
143 if ((nn
& 3) == 3 && (a
->v
[0] & 3) == 3)
146 /* --- Reduce and swap --- */
149 { mp
*t
= n
; n
= a
; a
= t
; }
152 /* --- Wrap everything up --- */
160 /*----- Test rig ----------------------------------------------------------*/
164 #include <mLib/testrig.h>
166 static int verify(dstr
*v
)
168 mp
*a
= *(mp
**)v
[0].buf
;
169 mp
*n
= *(mp
**)v
[1].buf
;
170 int s
= *(int *)v
[2].buf
;
171 int j
= mp_jacobi(a
, n
);
175 fputs("\n*** fail", stderr
);
176 fputs("a = ", stderr
); mp_writefile(a
, stderr
, 10); fputc('\n', stderr
);
177 fputs("n = ", stderr
); mp_writefile(n
, stderr
, 10); fputc('\n', stderr
);
178 fprintf(stderr
, "s = %i\n", s
);
179 fprintf(stderr
, "j = %i\n", j
);
185 assert(mparena_count(MPARENA_GLOBAL
) == 0);
189 static test_chunk tests
[] = {
190 { "jacobi", verify
, { &type_mp
, &type_mp
, &type_int
, 0 } },
194 int main(int argc
, char *argv
[])
197 test_run(argc
, argv
, tests
, SRCDIR
"/tests/mp");
203 /*----- That's all, folks -------------------------------------------------*/