-
Notifications
You must be signed in to change notification settings - Fork 0
/
sieve-of-pritchard.c
57 lines (48 loc) · 1.53 KB
/
sieve-of-pritchard.c
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
#include <stdio.h>
#include <gmp.h>
#include <stdlib.h>
void sieve_of_pritchard(mpz_t limit) {
mpz_t n, m, a, b, c, i, j, p, q, x, y, z;
mpz_inits(n, m, a, b, c, i, j, p, q, x, y, z, NULL);
mpz_set_ui(n, 1);
mpz_set_ui(m, 1);
mpz_set_ui(a, 1);
mpz_set_ui(b, 1);
mpz_set_ui(c, 1);
mpz_t *primes = malloc(mpz_get_ui(limit) * sizeof(mpz_t));
int *is_prime = malloc(mpz_get_ui(limit) * sizeof(int));
for (mpz_set_ui(i, 2); mpz_cmp(i, limit) <= 0; mpz_add_ui(i, i, 1)) {
mpz_init(primes[mpz_get_ui(i)]);
is_prime[mpz_get_ui(i)] = 1;
}
for (mpz_set_ui(i, 2); mpz_cmp(i, limit) <= 0; mpz_add_ui(i, i, 1)) {
if (is_prime[mpz_get_ui(i)]) {
mpz_set(primes[mpz_get_ui(m)], i);
mpz_add_ui(m, m, 1);
for (mpz_mul(j, i, i); mpz_cmp(j, limit) <= 0; mpz_add(j, j, i)) {
is_prime[mpz_get_ui(j)] = 0;
}
}
}
for (mpz_set_ui(i, 0); mpz_cmp(i, m) < 0; mpz_add_ui(i, i, 1)) {
gmp_printf("%Zd\n", primes[mpz_get_ui(i)]);
}
for (mpz_set_ui(i, 2); mpz_cmp(i, limit) <= 0; mpz_add_ui(i, i, 1)) {
mpz_clear(primes[mpz_get_ui(i)]);
}
free(primes);
free(is_prime);
mpz_clears(n, m, a, b, c, i, j, p, q, x, y, z, NULL);
}
int main(int argc, char *argv[]) {
if (argc != 2) {
printf("Usage: %s <limit>\n", argv[0]);
return 1;
}
mpz_t limit;
mpz_init(limit);
mpz_set_str(limit, argv[1], 10);
sieve_of_pritchard(limit);
mpz_clear(limit);
return 0;
}