📄 t-jac.c
字号:
/* Exercise mpz_*_kronecker_*() and mpz_jacobi() functions. *//*Copyright 1999, 2000, 2001, 2002, 2003 Free Software Foundation, Inc.This file is part of the GNU MP Library.The GNU MP Library is free software; you can redistribute it and/or modifyit under the terms of the GNU Lesser General Public License as published bythe Free Software Foundation; either version 2.1 of the License, or (at youroption) any later version.The GNU MP Library is distributed in the hope that it will be useful, butWITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITYor FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General PublicLicense for more details.You should have received a copy of the GNU Lesser General Public Licensealong with the GNU MP Library; see the file COPYING.LIB. If not, write tothe Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,MA 02111-1307, USA.*//* With no arguments the various Kronecker/Jacobi symbol routines are checked against some test data and a lot of derived data. To check the test data against PARI-GP, run t-jac -p | gp -q It takes a while because the output from "t-jac -p" is big. Enhancements: More big test cases than those given by check_squares_zi would be good. */#include <stdio.h>#include <stdlib.h>#include <string.h>#include "gmp.h"#include "gmp-impl.h"#include "tests.h"#ifdef _LONG_LONG_LIMB#define LL(l,ll) ll#else#define LL(l,ll) l#endifint option_pari = 0;unsigned longmpz_mod4 (mpz_srcptr z){ mpz_t m; unsigned long ret; mpz_init (m); mpz_fdiv_r_2exp (m, z, 2); ret = mpz_get_ui (m); mpz_clear (m); return ret;}intmpz_fits_ulimb_p (mpz_srcptr z){ return (SIZ(z) == 1 || SIZ(z) == 0);}mp_limb_tmpz_get_ulimb (mpz_srcptr z){ if (SIZ(z) == 0) return 0; else return PTR(z)[0];}voidtry_base (mp_limb_t a, mp_limb_t b, int answer){ int got; if ((b & 1) == 0 || b == 1 || a > b) return; got = mpn_jacobi_base (a, b, 0); if (got != answer) { printf (LL("mpn_jacobi_base (%lu, %lu) is %d should be %d\n", "mpn_jacobi_base (%llu, %llu) is %d should be %d\n"), a, b, got, answer); abort (); }}voidtry_zi_ui (mpz_srcptr a, unsigned long b, int answer){ int got; got = mpz_kronecker_ui (a, b); if (got != answer) { printf ("mpz_kronecker_ui ("); mpz_out_str (stdout, 10, a); printf (", %lu) is %d should be %d\n", b, got, answer); abort (); }}voidtry_zi_si (mpz_srcptr a, long b, int answer){ int got; got = mpz_kronecker_si (a, b); if (got != answer) { printf ("mpz_kronecker_si ("); mpz_out_str (stdout, 10, a); printf (", %ld) is %d should be %d\n", b, got, answer); abort (); }}voidtry_ui_zi (unsigned long a, mpz_srcptr b, int answer){ int got; got = mpz_ui_kronecker (a, b); if (got != answer) { printf ("mpz_ui_kronecker (%lu, ", a); mpz_out_str (stdout, 10, b); printf (") is %d should be %d\n", got, answer); abort (); }}voidtry_si_zi (long a, mpz_srcptr b, int answer){ int got; got = mpz_si_kronecker (a, b); if (got != answer) { printf ("mpz_si_kronecker (%ld, ", a); mpz_out_str (stdout, 10, b); printf (") is %d should be %d\n", got, answer); abort (); }}/* Don't bother checking mpz_jacobi, since it only differs for b even, and we don't have an actual expected answer for it. tests/devel/try.c does some checks though. */voidtry_zi_zi (mpz_srcptr a, mpz_srcptr b, int answer){ int got; got = mpz_kronecker (a, b); if (got != answer) { printf ("mpz_kronecker ("); mpz_out_str (stdout, 10, a); printf (", "); mpz_out_str (stdout, 10, b); printf (") is %d should be %d\n", got, answer); abort (); }}voidtry_pari (mpz_srcptr a, mpz_srcptr b, int answer){ printf ("try("); mpz_out_str (stdout, 10, a); printf (","); mpz_out_str (stdout, 10, b); printf (",%d)\n", answer);}voidtry_each (mpz_srcptr a, mpz_srcptr b, int answer){ if (option_pari) { try_pari (a, b, answer); return; } if (mpz_fits_ulimb_p (a) && mpz_fits_ulimb_p (b)) try_base (mpz_get_ulimb (a), mpz_get_ulimb (b), answer); if (mpz_fits_ulong_p (b)) try_zi_ui (a, mpz_get_ui (b), answer); if (mpz_fits_slong_p (b)) try_zi_si (a, mpz_get_si (b), answer); if (mpz_fits_ulong_p (a)) try_ui_zi (mpz_get_ui (a), b, answer); if (mpz_fits_sint_p (a)) try_si_zi (mpz_get_si (a), b, answer); try_zi_zi (a, b, answer);} /* Try (a/b) and (a/-b). */voidtry_pn (mpz_srcptr a, mpz_srcptr b_orig, int answer){ mpz_t b; mpz_init_set (b, b_orig); try_each (a, b, answer); mpz_neg (b, b); if (mpz_sgn (a) < 0) answer = -answer; try_each (a, b, answer); mpz_clear (b);}/* Try (a+k*p/b) for various k, using the fact (a/b) is periodic in a with period p. For b>0, p=b if b!=2mod4 or p=4*b if b==2mod4. */voidtry_periodic_num (mpz_srcptr a_orig, mpz_srcptr b, int answer){ mpz_t a, a_period; int i; if (mpz_sgn (b) <= 0) return; mpz_init_set (a, a_orig); mpz_init_set (a_period, b); if (mpz_mod4 (b) == 2) mpz_mul_ui (a_period, a_period, 4); /* don't bother with these tests if they're only going to produce even/even */ if (mpz_even_p (a) && mpz_even_p (b) && mpz_even_p (a_period)) goto done; for (i = 0; i < 10; i++) { mpz_add (a, a, a_period); try_pn (a, b, answer); } mpz_set (a, a_orig); for (i = 0; i < 10; i++) { mpz_sub (a, a, a_period); try_pn (a, b, answer); } done: mpz_clear (a); mpz_clear (a_period);}/* Try (a/b+k*p) for various k, using the fact (a/b) is periodic in b of period p. period p a==0,1mod4 a a==2mod4 4*a a==3mod4 and b odd 4*a a==3mod4 and b even 8*a In Henri Cohen's book the period is given as 4*a for all a==2,3mod4, but a counterexample would seem to be (3/2)=-1 which with (3/14)=+1 doesn't have period 4*a (but rather 8*a with (3/26)=-1). Maybe the plain 4*a is to be read as applying to a plain Jacobi symbol with b odd, rather than the Kronecker extension to b even. */voidtry_periodic_den (mpz_srcptr a, mpz_srcptr b_orig, int answer){ mpz_t b, b_period; int i; if (mpz_sgn (a) == 0 || mpz_sgn (b_orig) == 0) return; mpz_init_set (b, b_orig); mpz_init_set (b_period, a); if (mpz_mod4 (a) == 3 && mpz_even_p (b)) mpz_mul_ui (b_period, b_period, 8); else if (mpz_mod4 (a) >= 2) mpz_mul_ui (b_period, b_period, 4); /* don't bother with these tests if they're only going to produce even/even */ if (mpz_even_p (a) && mpz_even_p (b) && mpz_even_p (b_period)) goto done; for (i = 0; i < 10; i++) { mpz_add (b, b, b_period); try_pn (a, b, answer); } mpz_set (b, b_orig); for (i = 0; i < 10; i++) { mpz_sub (b, b, b_period); try_pn (a, b, answer); } done: mpz_clear (b); mpz_clear (b_period);}/* Try (a/b*2^k) for various k. */voidtry_2den (mpz_srcptr a, mpz_srcptr b_orig, int answer){ mpz_t b; int i; int answer_two; /* don't bother when b==0 */ if (mpz_sgn (b_orig) == 0) return; mpz_init_set (b, b_orig); /* answer_two = (a/2) */ switch (mpz_fdiv_ui (a, 8)) { case 1: case 7: answer_two = 1; break; case 3: case 5: answer_two = -1; break;
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -