1 /* Test mpz_perfect_square_p. 2 3 Copyright 2000-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 https://www.gnu.org/licenses/. */ 19 20 #include <stdio.h> 21 #include <stdlib.h> 22 23 #include "gmp-impl.h" 24 #include "tests.h" 25 26 #include "mpn/perfsqr.h" 27 28 29 /* check_modulo() exercises mpz_perfect_square_p on squares which cover each 30 possible quadratic residue to each divisor used within 31 mpn_perfect_square_p, ensuring those residues aren't incorrectly claimed 32 to be non-residues. 33 34 Each divisor is taken separately. It's arranged that n is congruent to 0 35 modulo the other divisors, 0 of course being a quadratic residue to any 36 modulus. 37 38 The values "(j*others)^2" cover all quadratic residues mod divisor[i], 39 but in no particular order. j is run from 1<=j<=divisor[i] so that zero 40 is excluded. A literal n==0 doesn't reach the residue tests. */ 41 42 void 43 check_modulo (void) 44 { 45 static const unsigned long divisor[] = PERFSQR_DIVISORS; 46 unsigned long i, j; 47 48 mpz_t alldiv, others, n; 49 50 mpz_init (alldiv); 51 mpz_init (others); 52 mpz_init (n); 53 54 /* product of all divisors */ 55 mpz_set_ui (alldiv, 1L); 56 for (i = 0; i < numberof (divisor); i++) 57 mpz_mul_ui (alldiv, alldiv, divisor[i]); 58 59 for (i = 0; i < numberof (divisor); i++) 60 { 61 /* product of all divisors except i */ 62 mpz_set_ui (others, 1L); 63 for (j = 0; j < numberof (divisor); j++) 64 if (i != j) 65 mpz_mul_ui (others, others, divisor[j]); 66 67 for (j = 1; j <= divisor[i]; j++) 68 { 69 /* square */ 70 mpz_mul_ui (n, others, j); 71 mpz_mul (n, n, n); 72 if (! mpz_perfect_square_p (n)) 73 { 74 printf ("mpz_perfect_square_p got 0, want 1\n"); 75 mpz_trace (" n", n); 76 abort (); 77 } 78 } 79 } 80 81 mpz_clear (alldiv); 82 mpz_clear (others); 83 mpz_clear (n); 84 } 85 86 87 /* Exercise mpz_perfect_square_p compared to what mpz_sqrt says. */ 88 void 89 check_sqrt (int reps) 90 { 91 mpz_t x2, x2t, x; 92 mp_size_t x2n; 93 int res; 94 int i; 95 /* int cnt = 0; */ 96 gmp_randstate_ptr rands = RANDS; 97 mpz_t bs; 98 99 mpz_init (bs); 100 101 mpz_init (x2); 102 mpz_init (x); 103 mpz_init (x2t); 104 105 for (i = 0; i < reps; i++) 106 { 107 mpz_urandomb (bs, rands, 9); 108 x2n = mpz_get_ui (bs); 109 mpz_rrandomb (x2, rands, x2n); 110 /* mpz_out_str (stdout, -16, x2); puts (""); */ 111 112 res = mpz_perfect_square_p (x2); 113 mpz_sqrt (x, x2); 114 mpz_mul (x2t, x, x); 115 116 if (res != (mpz_cmp (x2, x2t) == 0)) 117 { 118 printf ("mpz_perfect_square_p and mpz_sqrt differ\n"); 119 mpz_trace (" x ", x); 120 mpz_trace (" x2 ", x2); 121 mpz_trace (" x2t", x2t); 122 printf (" mpz_perfect_square_p %d\n", res); 123 printf (" mpz_sqrt %d\n", mpz_cmp (x2, x2t) == 0); 124 abort (); 125 } 126 127 /* cnt += res != 0; */ 128 } 129 /* printf ("%d/%d perfect squares\n", cnt, reps); */ 130 131 mpz_clear (bs); 132 mpz_clear (x2); 133 mpz_clear (x); 134 mpz_clear (x2t); 135 } 136 137 138 int 139 main (int argc, char **argv) 140 { 141 int reps = 200000; 142 143 tests_start (); 144 mp_trace_base = -16; 145 146 if (argc == 2) 147 reps = atoi (argv[1]); 148 149 check_modulo (); 150 check_sqrt (reps); 151 152 tests_end (); 153 exit (0); 154 } 155