2 #ifdef BN_MP_PRIME_NEXT_PRIME_C
3 /* LibTomMath, multiple-precision integer library -- Tom St Denis
5 * LibTomMath is a library that provides multiple-precision
6 * integer arithmetic as well as number theoretic functionality.
8 * The library was designed directly after the MPI library by
9 * Michael Fromberger but has been written from scratch with
10 * additional optimizations in place.
12 * The library is free for all purposes without any express
15 * Tom St Denis, tomstdenis@iahu.ca, http://math.libtomcrypt.org
18 /* finds the next prime after the number "a" using "t" trials
21 * bbs_style = 1 means the prime must be congruent to 3 mod 4
23 int mp_prime_next_prime(mp_int *a, int t, int bbs_style)
26 mp_digit res_tab[PRIME_SIZE], step, kstep;
29 /* ensure t is valid */
30 if (t <= 0 || t > PRIME_SIZE) {
37 /* simple algo if a is less than the largest prime in the table */
38 if (mp_cmp_d(a, ltm_prime_tab[PRIME_SIZE-1]) == MP_LT) {
39 /* find which prime it is bigger than */
40 for (x = PRIME_SIZE - 2; x >= 0; x--) {
41 if (mp_cmp_d(a, ltm_prime_tab[x]) != MP_LT) {
43 /* ok we found a prime smaller or
44 * equal [so the next is larger]
46 * however, the prime must be
47 * congruent to 3 mod 4
49 if ((ltm_prime_tab[x + 1] & 3) != 3) {
50 /* scan upwards for a prime congruent to 3 mod 4 */
51 for (y = x + 1; y < PRIME_SIZE; y++) {
52 if ((ltm_prime_tab[y] & 3) == 3) {
53 mp_set(a, ltm_prime_tab[y]);
59 mp_set(a, ltm_prime_tab[x + 1]);
64 /* at this point a maybe 1 */
65 if (mp_cmp_d(a, 1) == MP_EQ) {
69 /* fall through to the sieve */
72 /* generate a prime congruent to 3 mod 4 or 1/3 mod 4? */
79 /* at this point we will use a combination of a sieve and Miller-Rabin */
82 /* if a mod 4 != 3 subtract the correct value to make it so */
83 if ((a->dp[0] & 3) != 3) {
84 if ((err = mp_sub_d(a, (a->dp[0] & 3) + 1, a)) != MP_OKAY) { return err; };
87 if (mp_iseven(a) == 1) {
89 if ((err = mp_sub_d(a, 1, a)) != MP_OKAY) {
95 /* generate the restable */
96 for (x = 1; x < PRIME_SIZE; x++) {
97 if ((err = mp_mod_d(a, ltm_prime_tab[x], res_tab + x)) != MP_OKAY) {
102 /* init temp used for Miller-Rabin Testing */
103 if ((err = mp_init(&b)) != MP_OKAY) {
108 /* skip to the next non-trivially divisible candidate */
111 /* y == 1 if any residue was zero [e.g. cannot be prime] */
114 /* increase step to next candidate */
117 /* compute the new residue without using division */
118 for (x = 1; x < PRIME_SIZE; x++) {
119 /* add the step to each residue */
122 /* subtract the modulus [instead of using division] */
123 if (res_tab[x] >= ltm_prime_tab[x]) {
124 res_tab[x] -= ltm_prime_tab[x];
127 /* set flag if zero */
128 if (res_tab[x] == 0) {
132 } while (y == 1 && step < ((((mp_digit)1)<<DIGIT_BIT) - kstep));
135 if ((err = mp_add_d(a, step, a)) != MP_OKAY) {
139 /* if didn't pass sieve and step == MAX then skip test */
140 if (y == 1 && step >= ((((mp_digit)1)<<DIGIT_BIT) - kstep)) {
145 for (x = 0; x < t; x++) {
146 mp_set(&b, ltm_prime_tab[t]);
147 if ((err = mp_prime_miller_rabin(a, &b, &res)) != MP_OKAY) {