1 /*
2     Copyright (C) 2012 William Hart
3 
4     This file is part of FLINT.
5 
6     FLINT is free software: you can redistribute it and/or modify it under
7     the terms of the GNU Lesser General Public License (LGPL) as published
8     by the Free Software Foundation; either version 2.1 of the License, or
9     (at your option) any later version.  See <https://www.gnu.org/licenses/>.
10 */
11 
12 #include <stdio.h>
13 #include <stdlib.h>
14 #include <gmp.h>
15 #include "flint.h"
16 #include "ulong_extras.h"
17 #include "fmpz.h"
18 
19 /* Composite strong pseudoprimes from https://oeis.org/A014233 */
20 static const char * composites[] = {
21     "2047",
22     "1373653",
23     "25326001",
24     "3215031751",
25     "2152302898747",
26     "3474749660383",
27     "341550071728321",
28     "341550071728321",
29     "3825123056546413051",
30     "3825123056546413051",
31     "3825123056546413051",
32     "318665857834031151167461",
33     "3317044064679887385961981",
34     NULL,
35 };
36 
37 int
main(void)38 main(void)
39 {
40     int i, result, count = 0;
41     FLINT_TEST_INIT(state);
42 
43     flint_printf("is_strong_probabprime....");
44     fflush(stdout);
45 
46     /* test table */
47     {
48         for (i = 0; composites[i] != NULL; i++)
49         {
50             int j;
51             fmpz_t n, a;
52             fmpz_init(n);
53             fmpz_init(a);
54 
55             fmpz_set_str(n, composites[i], 10);
56 
57             for (j = 0; j <= i; j++)
58             {
59                 fmpz_set_ui(a, n_nth_prime(j + 1));
60 
61                 if (!fmpz_is_strong_probabprime(n, a))
62                 {
63                     flint_printf("FAIL (composite expected to pass test):\n");
64                     fmpz_print(n); printf("\n");
65                     fmpz_print(a); printf("\n");
66                     abort();
67                 }
68             }
69 
70             fmpz_set_ui(a, n_nth_prime(i + 2
71                     + (i == 6) + 2*(i==8) + (i==9))  /* because of repeated entries */
72                 );
73 
74             if (fmpz_is_strong_probabprime(n, a))
75             {
76                 flint_printf("FAIL (composite expected to fail test):\n");
77                 fmpz_print(n); printf("\n");
78                 fmpz_print(a); printf("\n");
79                 abort();
80             }
81 
82             fmpz_clear(n);
83             fmpz_clear(a);
84         }
85     }
86 
87     /* test primes always pass */
88     for (i = 0; i < 100 * flint_test_multiplier(); i++)
89     {
90         fmpz_t p, b, F;
91 
92         fmpz_init(p);
93         fmpz_init(b);
94         fmpz_init(F);
95 
96         do {
97            fmpz_randbits(p, state, n_randint(state, 330) + 2);
98            fmpz_abs(p, p);
99         } while (!fmpz_is_probabprime(p) || fmpz_cmp_ui(p, 2) == 0);
100 
101         do {
102            fmpz_randbits(b, state, n_randint(state, 100) + 1);
103            fmpz_abs(b, b);
104         } while (fmpz_is_zero(b) || fmpz_is_one(b));
105 
106         result = fmpz_is_strong_probabprime(p, b);
107         if (!result)
108         {
109             flint_printf("FAIL:\n");
110             fmpz_print(p); printf("\n");
111             fmpz_print(b); printf("\n");
112             abort();
113         }
114 
115         fmpz_clear(p);
116         fmpz_clear(b);
117         fmpz_clear(F);
118     }
119 
120     /* test composites rarely pass */
121     for (i = 0; i < 100 * flint_test_multiplier(); i++)
122     {
123         fmpz_t p, a, b, F;
124 
125         fmpz_init(p);
126         fmpz_init(a);
127         fmpz_init(b);
128         fmpz_init(F);
129 
130         do {
131            fmpz_randbits(p, state, n_randint(state, 100) + 2);
132         } while (fmpz_cmp_ui(p, 2) < 0);
133         do {
134            fmpz_randbits(a, state, n_randint(state, 100) + 2);
135         } while (fmpz_cmp_ui(a, 2) < 0);
136 
137         do {
138            fmpz_randbits(b, state, n_randint(state, 100) + 1);
139            fmpz_abs(b, b);
140         } while (fmpz_cmp_ui(b, 2) < 0);
141 
142         fmpz_mul(p, p, a);
143 
144         if (fmpz_is_strong_probabprime(p, b))
145            count++;
146 
147         fmpz_clear(p);
148         fmpz_clear(a);
149         fmpz_clear(b);
150         fmpz_clear(F);
151     }
152 
153     result = (count < flint_test_multiplier());
154     if (!result)
155     {
156         flint_printf("FAIL:\n");
157         flint_printf("count = %ld\n", count);
158         abort();
159     }
160 
161     FLINT_TEST_CLEANUP(state);
162 
163     flint_printf("PASS\n");
164     return 0;
165 }
166