https://github.com/henrycg/earand
Tip revision: 362625f3847ddd2e2d5a90cf5f93b4b3ea76e059 authored by Henry Corrigan-Gibbs on 21 May 2018, 19:54:59 UTC
Merge branch 'master' of github.com:henrycg/earand
Merge branch 'master' of github.com:henrycg/earand
Tip revision: 362625f
prime.c
#include <openssl/bn.h>
#include "bn_prime.h"
#include "prime.h"
#include "util.h"
static const int bound_bits = 48;
static const int max_iterations = 4096*8;
static const double epsilon = 0.33;
BIGNUM* GeneratePrime(int n_bits,
const BIGNUM* min, const BIGNUM* max, BN_CTX* ctx)
{
BIGNUM* ret = BN_new();
CHECK_CALL(ret);
if(n_bits < bound_bits) {
CHECK_CALL(BN_generate_prime(ret, n_bits, false, NULL, NULL, NULL, NULL));
return ret;
}
if(min != NULL) {
n_bits = BN_num_bits(min);
}
const double b_min = (epsilon * (double)n_bits);
BIGNUM* factor = GeneratePrime((int)b_min, NULL, NULL, ctx);
BIGNUM* N0 = BN_new();
BIGNUM* a = BN_new();
BIGNUM* ia = BN_new();
CHECK_CALL(N0);
CHECK_CALL(a);
CHECK_CALL(ia);
const int sn = n_bits * max_iterations;
BIGNUM* t = BN_new();
BIGNUM* lower = BN_new();
BIGNUM* range = BN_new();
// Pick random t in 2^{n-2}/F, ..., ((2^{n-1}/F) - sn)
CHECK_CALL(BN_one(lower));
CHECK_CALL(BN_lshift(lower, lower, n_bits-2));
CHECK_CALL(BN_div(lower, NULL, lower, factor, ctx));
CHECK_CALL(BN_copy(range, lower));
CHECK_CALL(BN_sub_word(range, sn));
printf("$%d\n", n_bits);
puts("---");
BN_print_fp(stderr, range);
puts("<>");
BN_print_fp(stderr, lower);
puts("---");
CHECK_CALL(BN_rand_range(t, range));
CHECK_CALL(BN_add(t, t, lower));
BN_clear_free(lower);
BN_clear_free(range);
//a = 2F
CHECK_CALL(BN_lshift1(a, factor));
//N0 = ta+1
CHECK_CALL(BN_mul(N0, t, a, ctx));
CHECK_CALL(BN_add_word(N0, 1));
BN_clear_free(t);
bool success = false;
//N = N0 + ia;
//i = 0, ..., s
BIGNUM* i = BN_new();
CHECK_CALL(i);
int ntest = 0;
if(min == NULL) {
BN_zero(i);
} else {
// i = ceil((min - ta - 1)/a)
// N0 = ta+1
CHECK_CALL(BN_copy(i, min));
CHECK_CALL(BN_sub(i, i, N0));
CHECK_CALL(BN_sub_word(i, 2));
CHECK_CALL(BN_div(i, NULL, i, a, ctx));
}
// iterate until N0 is too big
while(!max || BN_cmp(max, N0) < 1) {
CHECK_CALL(BN_copy(ia, a));
CHECK_CALL(BN_mul(ia, ia, i, ctx));
CHECK_CALL(BN_add(ret, N0, ia));
bool is_comp = false;
for(int j=0; j<NUMPRIMES; j++) {
if(BN_mod_word(ret, primes[j]) == 0) {
is_comp = true;
break;
}
}
if(is_comp) continue;
printf("%d\n",ntest++);
// Check is prime and break if so
if(BN_is_prime_fasttest(ret, BN_prime_checks, NULL, ctx, NULL, 1)) {
success = true;
break;
}
BN_add_word(i, 1);
}
BN_clear_free(i);
BN_clear_free(N0);
BN_clear_free(a);
BN_clear_free(ia);
CHECK_CALL(success);
BN_clear_free(factor);
return ret;
}
int main(void) {
BN_CTX* ctx = BN_CTX_new();
GeneratePrime(512, NULL, NULL, ctx);
BN_CTX_free(ctx);
return 0;
}