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
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
|
#include <cstring>
#include <cmath>
#include <cassert>
#include "numalgo.h"
#include "primes.h"
using namespace std;
vector<int> smallprimes;
bool smallprimes_inited=false;
void fillsmallprimes(){
smallprimes_inited=true;
//TODO: reserve expected amount of space in smallprimes
smallprimes.push_back(2);
const int highbound=65000;
bool composite[highbound/2]; //entries for 3, 5, 7, 9, etc.
memset(composite,0,highbound/2*sizeof(bool));
int roothighbound=sqrt(highbound);
for(int i=3;i<=highbound;i+=2){
if(composite[i/2-1])continue;
smallprimes.push_back(i);
if(i>roothighbound)continue;
for(int j=i*i;j<=highbound;j+=2*i){
composite[j/2-1]=true;
}
}
//for(int p : smallprimes)cerr<<p<<' ';
//cerr<<endl;
}
pair<Bigint,Bigint> genprimepair(int nbits){
// for x = nbits/2:
// (2^x)^2 = 2^(2x)
// (2^x + 2^(x-2))^2 = 2^(2x) + 2^(2x-1) + 2^(2x-4)
// ergo: (2^x + lambda*2^(x-2))^2 \in [2^(2x), 2^(2x+1)), for lambda \in [0,1]
// To make sure the primes "differ in length by a few digits" [RSA78], we use x1=x-2 in the first
// prime and x2-x+2 in the second
int x1=nbits/2-2,x2=(nbits+1)/2+2;
assert(x1+x2==nbits);
return make_pair(
randprime(Bigint::one<<x1,(Bigint::one<<x1)+(Bigint::one<<(x1-2))),
randprime(Bigint::one<<x2,(Bigint::one<<x2)+(Bigint::one<<(x2-2))));
}
Bigint randprime(const Bigint &biglow,const Bigint &bighigh){
if(!smallprimes_inited)fillsmallprimes();
assert(bighigh>biglow);
static const int maxrangesize=100001;
Bigint diff(bighigh-biglow);
Bigint low,high; //inclusive
if(diff<=maxrangesize){
low=biglow;
high=bighigh;
} else {
high=low=cryptrandom_big(diff-maxrangesize);
high+=maxrangesize;
}
if(low.even())low+=1;
if(high.even())high-=1;
Bigint sizeb(high-low+1);
assert(sizeb>=0&&sizeb<=maxrangesize);
int nnums=(sizeb.lowdigits()+1)/2;
bool composite[nnums]; //low, low+2, low+4, ..., high (all odd numbers)
memset(composite,0,nnums*sizeof(bool));
int nsmallprimes=smallprimes.size();
for(int i=1;i<nsmallprimes;i++){
int pr=smallprimes[i];
int lowoffset=low.divmod(Bigint(pr)).second.lowdigits();
int startat;
if(lowoffset==0)startat=0;
else if((pr-lowoffset)%2==0)startat=(pr-lowoffset)/2;
else startat=pr-lowoffset;
for(int i=startat;i<nnums;i+=pr){ //skips ahead `2*pr` each time (so `pr` array elements)
composite[i]=true;
}
}
}
bool bailliePSW(const Bigint&);
|