1 /* mpz_si_kronecker -- long+mpz Kronecker/Jacobi symbol.
3 Copyright 1999, 2000, 2001, 2002 Free Software Foundation, Inc.
5 This file is part of the GNU MP Library.
7 The GNU MP Library is free software; you can redistribute it and/or modify
8 it under the terms of the GNU Lesser General Public License as published by
9 the Free Software Foundation; either version 3 of the License, or (at your
10 option) any later version.
12 The GNU MP Library is distributed in the hope that it will be useful, but
13 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
14 or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
15 License for more details.
17 You should have received a copy of the GNU Lesser General Public License
18 along with the GNU MP Library. If not, see http://www.gnu.org/licenses/. */
26 mpz_si_kronecker (long a
, mpz_srcptr b
)
32 mp_limb_t a_limb
, b_rem
;
36 #if GMP_NUMB_BITS < BITS_PER_ULONG
37 if (a
> GMP_NUMB_MAX
|| a
< -GMP_NUMB_MAX
)
41 ALLOC(az
) = numberof (alimbs
);
44 return mpz_kronecker (az
, b
);
50 return JACOBI_S0 (a
); /* (a/0) */
52 /* account for the effect of the sign of b, then ignore it */
53 result_bit1
= JACOBI_BSGN_SS_BIT1 (a
, b_size
);
57 b_abs_size
= ABS (b_size
);
63 result_bit1
^= JACOBI_ASGN_SU_BIT1 (a
, b_low
);
64 a_limb
= (unsigned long) ABS(a
);
66 if ((a_limb
& 1) == 0)
68 /* (0/b)=1 for b=+/-1, 0 otherwise */
70 return (b_abs_size
== 1 && b_low
== 1);
73 count_trailing_zeros (twos
, a_limb
);
75 /* (a*2^n/b) = (a/b) * twos(n,a) */
76 result_bit1
^= JACOBI_TWOS_U_BIT1 (twos
, b_low
);
81 /* (even/even)=0, and (0/b)=0 for b!=+/-1 */
87 Establish shifted b_low with valid bit1 for ASGN and RECIP below.
88 Zero limbs stripped are accounted for, but zero bits on b_low are
89 not because they remain in {b_ptr,b_abs_size} for the
90 JACOBI_MOD_OR_MODEXACT_1_ODD. */
92 JACOBI_STRIP_LOW_ZEROS (result_bit1
, a
, b_ptr
, b_abs_size
, b_low
);
95 if (UNLIKELY (b_low
== GMP_NUMB_HIGHBIT
))
97 /* need b_ptr[1] to get bit1 in b_low */
100 /* (a/0x80000000) = (a/2)^(BPML-1) */
101 if ((GMP_NUMB_BITS
% 2) == 0)
102 result_bit1
^= JACOBI_TWO_U_BIT1 (a
);
103 return JACOBI_BIT1_TO_PN (result_bit1
);
107 b_low
= b_ptr
[1] << 1;
111 count_trailing_zeros (twos
, b_low
);
116 result_bit1
^= JACOBI_ASGN_SU_BIT1 (a
, b_low
);
117 a_limb
= (unsigned long) ABS(a
);
121 return JACOBI_BIT1_TO_PN (result_bit1
); /* (1/b)=1 */
123 /* (a/b*2^n) = (b*2^n mod a / a) * recip(a,b) */
124 JACOBI_MOD_OR_MODEXACT_1_ODD (result_bit1
, b_rem
, b_ptr
, b_abs_size
, a_limb
);
125 result_bit1
^= JACOBI_RECIP_UU_BIT1 (a_limb
, b_low
);
126 return mpn_jacobi_base (b_rem
, a_limb
, result_bit1
);