1 /* Exercise mpz_probab_prime_p.
2
3 Copyright 2002, 2018-2019 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 #include "gmp-impl.h"
23 #include "tests.h"
24
25
26 /* Enhancements:
27
28 - Test some big primes don't come back claimed to be composite.
29 - Test some big composites don't come back claimed to be certainly prime.
30 - Test some big composites with small factors are identified as certainly
31 composite. */
32
33
34 /* return 1 if prime, 0 if composite */
35 int
isprime(long n)36 isprime (long n)
37 {
38 long i;
39
40 n = ABS(n);
41
42 if (n < 2)
43 return 0;
44 if (n < 4)
45 return 1;
46 if ((n & 1) == 0)
47 return 0;
48
49 for (i = 3; i*i <= n; i+=2)
50 if ((n % i) == 0)
51 return 0;
52
53 return 1;
54 }
55
56 void
check_one(mpz_srcptr n,int want)57 check_one (mpz_srcptr n, int want)
58 {
59 int got;
60
61 got = mpz_probab_prime_p (n, 25);
62
63 /* "definitely prime" is fine if we only wanted "probably prime" */
64 if (got == 2 && want == 1)
65 want = 2;
66
67 if (got != want)
68 {
69 printf ("mpz_probab_prime_p\n");
70 mpz_trace (" n ", n);
71 printf (" got =%d", got);
72 printf (" want=%d", want);
73 abort ();
74 }
75 }
76
77 void
check_pn(mpz_ptr n,int want)78 check_pn (mpz_ptr n, int want)
79 {
80 check_one (n, want);
81 mpz_neg (n, n);
82 check_one (n, want);
83 }
84
85 /* expect certainty for small n */
86 void
check_small(void)87 check_small (void)
88 {
89 mpz_t n;
90 long i;
91
92 mpz_init (n);
93
94 for (i = 0; i < 300; i++)
95 {
96 mpz_set_si (n, i);
97 check_pn (n, isprime (i));
98 }
99
100 mpz_clear (n);
101 }
102
103 void
check_composites(int count)104 check_composites (int count)
105 {
106 int i;
107 mpz_t a, b, n, bs;
108 unsigned long size_range, size;
109 gmp_randstate_ptr rands = RANDS;
110
111 mpz_init (a);
112 mpz_init (b);
113 mpz_init (n);
114 mpz_init (bs);
115
116 for (i = 0; i < count; i++)
117 {
118 mpz_urandomb (bs, rands, 32);
119 size_range = mpz_get_ui (bs) % 13 + 1; /* 0..8192 bit operands */
120
121 mpz_urandomb (bs, rands, size_range);
122 size = mpz_get_ui (bs);
123 mpz_rrandomb (a, rands, size);
124
125 mpz_urandomb (bs, rands, 32);
126 size_range = mpz_get_ui (bs) % 13 + 1; /* 0..8192 bit operands */
127 mpz_rrandomb (b, rands, size);
128
129 /* Exclude trivial factors */
130 if (mpz_cmp_ui (a, 1) == 0)
131 mpz_set_ui (a, 2);
132 if (mpz_cmp_ui (b, 1) == 0)
133 mpz_set_ui (b, 2);
134
135 mpz_mul (n, a, b);
136
137 check_pn (n, 0);
138 }
139 mpz_clear (a);
140 mpz_clear (b);
141 mpz_clear (n);
142 mpz_clear (bs);
143 }
144
145 static void
check_primes(void)146 check_primes (void)
147 {
148 static const char * const primes[] = {
149 "2", "53", "1234567891",
150 "2055693949", "1125899906842597", "16412292043871650369",
151 /* diffie-hellman-group1-sha1, also "Well known group 2" in RFC
152 2412, 2^1024 - 2^960 - 1 + 2^64 * { [2^894 pi] + 129093 } */
153 "0xFFFFFFFFFFFFFFFFC90FDAA22168C234C4C6628B80DC1CD1"
154 "29024E088A67CC74020BBEA63B139B22514A08798E3404DD"
155 "EF9519B3CD3A431B302B0A6DF25F14374FE1356D6D51C245"
156 "E485B576625E7EC6F44C42E9A637ED6B0BFF5CB6F406B7ED"
157 "EE386BFB5A899FA5AE9F24117C4B1FE649286651ECE65381"
158 "FFFFFFFFFFFFFFFF",
159 NULL
160 };
161
162 mpz_t n;
163 int i;
164
165 mpz_init (n);
166
167 for (i = 0; primes[i]; i++)
168 {
169 mpz_set_str_or_abort (n, primes[i], 0);
170 check_one (n, 1);
171 }
172 mpz_clear (n);
173 }
174
175 static void
check_fermat_mersenne(int count)176 check_fermat_mersenne (int count)
177 {
178 int fermat_exponents [] = {1, 2, 4, 8, 16};
179 int mersenne_exponents [] = {2, 3, 5, 7, 13, 17, 19, 31, 61, 89,
180 107, 127, 521, 607, 1279, 2203, 2281,
181 3217, 4253, 4423, 9689, 9941, 11213,
182 19937, 21701, 23209, 44497, 86243};
183 mpz_t pp;
184 int i, j, want;
185
186 mpz_init (pp);
187 count = MIN (110000, count);
188
189 for (i=1; i<count; ++i)
190 {
191 mpz_set_ui (pp, 1);
192 mpz_setbit (pp, i); /* 2^i + 1 */
193 want = 0;
194 for (j = 0; j < numberof (fermat_exponents); j++)
195 if (fermat_exponents[j] == i)
196 {
197 want = 1;
198 break;
199 }
200 check_one (pp, want);
201
202 mpz_sub_ui (pp, pp, 2); /* 2^i - 1 */
203 want = 0;
204 for (j = 0; j < numberof (mersenne_exponents); j++)
205 if (mersenne_exponents[j] == i)
206 {
207 want = 1;
208 break;
209 }
210 check_one (pp, want);
211 }
212 mpz_clear (pp);
213 }
214
215 int
main(int argc,char ** argv)216 main (int argc, char **argv)
217 {
218 int count = 1000;
219
220 TESTS_REPS (count, argv, argc);
221
222 tests_start ();
223
224 check_small ();
225 check_fermat_mersenne (count >> 3);
226 check_composites (count);
227 check_primes ();
228
229 tests_end ();
230 exit (0);
231 }
232