1 /* Test mpz_perfect_square_p. 2 3 Copyright 2000, 2001, 2002 Free Software Foundation, Inc. 4 5 This file is part of the GNU MP Library test suite. 6 7 The GNU MP Library test suite is free software; you can redistribute it 8 and/or modify it under the terms of the GNU General Public License as 9 published by the Free Software Foundation; either version 3 of the License, 10 or (at your option) any later version. 11 12 The GNU MP Library test suite is distributed in the hope that it will be 13 useful, but WITHOUT ANY WARRANTY; without even the implied warranty of 14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General 15 Public License for more details. 16 17 You should have received a copy of the GNU General Public License along with 18 the GNU MP Library test suite. If not, see http://www.gnu.org/licenses/. */ 19 20 #include <stdio.h> 21 #include <stdlib.h> 22 23 #include "gmp.h" 24 #include "gmp-impl.h" 25 #include "tests.h" 26 27 #include "mpn/perfsqr.h" 28 29 30 /* check_modulo() exercises mpz_perfect_square_p on squares which cover each 31 possible quadratic residue to each divisor used within 32 mpn_perfect_square_p, ensuring those residues aren't incorrectly claimed 33 to be non-residues. 34 35 Each divisor is taken separately. It's arranged that n is congruent to 0 36 modulo the other divisors, 0 of course being a quadratic residue to any 37 modulus. 38 39 The values "(j*others)^2" cover all quadratic residues mod divisor[i], 40 but in no particular order. j is run from 1<=j<=divisor[i] so that zero 41 is excluded. A literal n==0 doesn't reach the residue tests. */ 42 43 void 44 check_modulo (void) 45 { 46 static const unsigned long divisor[] = PERFSQR_DIVISORS; 47 unsigned long i, j; 48 49 mpz_t alldiv, others, n; 50 51 mpz_init (alldiv); 52 mpz_init (others); 53 mpz_init (n); 54 55 /* product of all divisors */ 56 mpz_set_ui (alldiv, 1L); 57 for (i = 0; i < numberof (divisor); i++) 58 mpz_mul_ui (alldiv, alldiv, divisor[i]); 59 60 for (i = 0; i < numberof (divisor); i++) 61 { 62 /* product of all divisors except i */ 63 mpz_set_ui (others, 1L); 64 for (j = 0; j < numberof (divisor); j++) 65 if (i != j) 66 mpz_mul_ui (others, others, divisor[j]); 67 68 for (j = 1; j <= divisor[i]; j++) 69 { 70 /* square */ 71 mpz_mul_ui (n, others, j); 72 mpz_mul (n, n, n); 73 if (! mpz_perfect_square_p (n)) 74 { 75 printf ("mpz_perfect_square_p got 0, want 1\n"); 76 mpz_trace (" n", n); 77 abort (); 78 } 79 } 80 } 81 82 mpz_clear (alldiv); 83 mpz_clear (others); 84 mpz_clear (n); 85 } 86 87 88 /* Exercise mpz_perfect_square_p compared to what mpz_sqrt says. */ 89 void 90 check_sqrt (int reps) 91 { 92 mpz_t x2, x2t, x; 93 mp_size_t x2n; 94 int res; 95 int i; 96 /* int cnt = 0; */ 97 gmp_randstate_ptr rands = RANDS; 98 mpz_t bs; 99 100 mpz_init (bs); 101 102 mpz_init (x2); 103 mpz_init (x); 104 mpz_init (x2t); 105 106 for (i = 0; i < reps; i++) 107 { 108 mpz_urandomb (bs, rands, 9); 109 x2n = mpz_get_ui (bs); 110 mpz_rrandomb (x2, rands, x2n); 111 /* mpz_out_str (stdout, -16, x2); puts (""); */ 112 113 res = mpz_perfect_square_p (x2); 114 mpz_sqrt (x, x2); 115 mpz_mul (x2t, x, x); 116 117 if (res != (mpz_cmp (x2, x2t) == 0)) 118 { 119 printf ("mpz_perfect_square_p and mpz_sqrt differ\n"); 120 mpz_trace (" x ", x); 121 mpz_trace (" x2 ", x2); 122 mpz_trace (" x2t", x2t); 123 printf (" mpz_perfect_square_p %d\n", res); 124 printf (" mpz_sqrt %d\n", mpz_cmp (x2, x2t) == 0); 125 abort (); 126 } 127 128 /* cnt += res != 0; */ 129 } 130 /* printf ("%d/%d perfect squares\n", cnt, reps); */ 131 132 mpz_clear (bs); 133 mpz_clear (x2); 134 mpz_clear (x); 135 mpz_clear (x2t); 136 } 137 138 139 int 140 main (int argc, char **argv) 141 { 142 int reps = 200000; 143 144 tests_start (); 145 mp_trace_base = -16; 146 147 if (argc == 2) 148 reps = atoi (argv[1]); 149 150 check_modulo (); 151 check_sqrt (reps); 152 153 tests_end (); 154 exit (0); 155 } 156