diff options
author | Matt A. Tobin <mattatobin@localhost.localdomain> | 2018-02-02 04:16:08 -0500 |
---|---|---|
committer | Matt A. Tobin <mattatobin@localhost.localdomain> | 2018-02-02 04:16:08 -0500 |
commit | 5f8de423f190bbb79a62f804151bc24824fa32d8 (patch) | |
tree | 10027f336435511475e392454359edea8e25895d /security/nss/lib/freebl/mpi/utils/sieve.c | |
parent | 49ee0794b5d912db1f95dce6eb52d781dc210db5 (diff) | |
download | UXP-5f8de423f190bbb79a62f804151bc24824fa32d8.tar UXP-5f8de423f190bbb79a62f804151bc24824fa32d8.tar.gz UXP-5f8de423f190bbb79a62f804151bc24824fa32d8.tar.lz UXP-5f8de423f190bbb79a62f804151bc24824fa32d8.tar.xz UXP-5f8de423f190bbb79a62f804151bc24824fa32d8.zip |
Add m-esr52 at 52.6.0
Diffstat (limited to 'security/nss/lib/freebl/mpi/utils/sieve.c')
-rw-r--r-- | security/nss/lib/freebl/mpi/utils/sieve.c | 243 |
1 files changed, 243 insertions, 0 deletions
diff --git a/security/nss/lib/freebl/mpi/utils/sieve.c b/security/nss/lib/freebl/mpi/utils/sieve.c new file mode 100644 index 000000000..57768af9e --- /dev/null +++ b/security/nss/lib/freebl/mpi/utils/sieve.c @@ -0,0 +1,243 @@ +/* + * sieve.c + * + * Finds prime numbers using the Sieve of Eratosthenes + * + * This implementation uses a bitmap to represent all odd integers in a + * given range. We iterate over this bitmap, crossing off the + * multiples of each prime we find. At the end, all the remaining set + * bits correspond to prime integers. + * + * Here, we make two passes -- once we have generated a sieve-ful of + * primes, we copy them out, reset the sieve using the highest + * generated prime from the first pass as a base. Then we cross out + * all the multiples of all the primes we found the first time through, + * and re-sieve. In this way, we get double use of the memory we + * allocated for the sieve the first time though. Since we also + * implicitly ignore multiples of 2, this amounts to 4 times the + * values. + * + * This could (and probably will) be generalized to re-use the sieve a + * few more times. + * + * This Source Code Form is subject to the terms of the Mozilla Public + * License, v. 2.0. If a copy of the MPL was not distributed with this + * file, You can obtain one at http://mozilla.org/MPL/2.0/. */ + +#include <stdio.h> +#include <stdlib.h> +#include <limits.h> + +typedef unsigned char byte; + +typedef struct { + int size; + byte *bits; + long base; + int next; + int nbits; +} sieve; + +void sieve_init(sieve *sp, long base, int nbits); +void sieve_grow(sieve *sp, int nbits); +long sieve_next(sieve *sp); +void sieve_reset(sieve *sp, long base); +void sieve_cross(sieve *sp, long val); +void sieve_clear(sieve *sp); + +#define S_ISSET(S, B) (((S)->bits[(B) / CHAR_BIT] >> ((B) % CHAR_BIT)) & 1) +#define S_SET(S, B) ((S)->bits[(B) / CHAR_BIT] |= (1 << ((B) % CHAR_BIT))) +#define S_CLR(S, B) ((S)->bits[(B) / CHAR_BIT] &= ~(1 << ((B) % CHAR_BIT))) +#define S_VAL(S, B) ((S)->base + (2 * (B))) +#define S_BIT(S, V) (((V) - ((S)->base)) / 2) + +int +main(int argc, char *argv[]) +{ + sieve s; + long pr, *p; + int c, ix, cur = 0; + + if (argc < 2) { + fprintf(stderr, "Usage: %s <width>\n", argv[0]); + return 1; + } + + c = atoi(argv[1]); + if (c < 0) + c = -c; + + fprintf(stderr, "%s: sieving to %d positions\n", argv[0], c); + + sieve_init(&s, 3, c); + + c = 0; + while ((pr = sieve_next(&s)) > 0) { + ++c; + } + + p = calloc(c, sizeof(long)); + if (!p) { + fprintf(stderr, "%s: out of memory after first half\n", argv[0]); + sieve_clear(&s); + exit(1); + } + + fprintf(stderr, "%s: half done ... \n", argv[0]); + + for (ix = 0; ix < s.nbits; ix++) { + if (S_ISSET(&s, ix)) { + p[cur] = S_VAL(&s, ix); + printf("%ld\n", p[cur]); + ++cur; + } + } + + sieve_reset(&s, p[cur - 1]); + fprintf(stderr, "%s: crossing off %d found primes ... \n", argv[0], cur); + for (ix = 0; ix < cur; ix++) { + sieve_cross(&s, p[ix]); + if (!(ix % 1000)) + fputc('.', stderr); + } + fputc('\n', stderr); + + free(p); + + fprintf(stderr, "%s: sieving again from %ld ... \n", argv[0], p[cur - 1]); + c = 0; + while ((pr = sieve_next(&s)) > 0) { + ++c; + } + + fprintf(stderr, "%s: done!\n", argv[0]); + for (ix = 0; ix < s.nbits; ix++) { + if (S_ISSET(&s, ix)) { + printf("%ld\n", S_VAL(&s, ix)); + } + } + + sieve_clear(&s); + + return 0; +} + +void +sieve_init(sieve *sp, long base, int nbits) +{ + sp->size = (nbits / CHAR_BIT); + + if (nbits % CHAR_BIT) + ++sp->size; + + sp->bits = calloc(sp->size, sizeof(byte)); + memset(sp->bits, UCHAR_MAX, sp->size); + if (!(base & 1)) + ++base; + sp->base = base; + + sp->next = 0; + sp->nbits = sp->size * CHAR_BIT; +} + +void +sieve_grow(sieve *sp, int nbits) +{ + int ns = (nbits / CHAR_BIT); + + if (nbits % CHAR_BIT) + ++ns; + + if (ns > sp->size) { + byte *tmp; + int ix; + + tmp = calloc(ns, sizeof(byte)); + if (tmp == NULL) { + fprintf(stderr, "Error: out of memory in sieve_grow\n"); + return; + } + + memcpy(tmp, sp->bits, sp->size); + for (ix = sp->size; ix < ns; ix++) { + tmp[ix] = UCHAR_MAX; + } + + free(sp->bits); + sp->bits = tmp; + sp->size = ns; + + sp->nbits = sp->size * CHAR_BIT; + } +} + +long +sieve_next(sieve *sp) +{ + long out; + int ix = 0; + long val; + + if (sp->next > sp->nbits) + return -1; + + out = S_VAL(sp, sp->next); +#ifdef DEBUG + fprintf(stderr, "Sieving %ld\n", out); +#endif + + /* Sieve out all multiples of the current prime */ + val = out; + while (ix < sp->nbits) { + val += out; + ix = S_BIT(sp, val); + if ((val & 1) && ix < sp->nbits) { /* && S_ISSET(sp, ix)) { */ + S_CLR(sp, ix); +#ifdef DEBUG + fprintf(stderr, "Crossing out %ld (bit %d)\n", val, ix); +#endif + } + } + + /* Scan ahead to the next prime */ + ++sp->next; + while (sp->next < sp->nbits && !S_ISSET(sp, sp->next)) + ++sp->next; + + return out; +} + +void +sieve_cross(sieve *sp, long val) +{ + int ix = 0; + long cur = val; + + while (cur < sp->base) + cur += val; + + ix = S_BIT(sp, cur); + while (ix < sp->nbits) { + if (cur & 1) + S_CLR(sp, ix); + cur += val; + ix = S_BIT(sp, cur); + } +} + +void +sieve_reset(sieve *sp, long base) +{ + memset(sp->bits, UCHAR_MAX, sp->size); + sp->base = base; + sp->next = 0; +} + +void +sieve_clear(sieve *sp) +{ + if (sp->bits) + free(sp->bits); + + sp->bits = NULL; +} |