1 #include "os.h"
2 #include <mp.h>
3 #include "dat.h"
4
5 static mpdigit _mptwodata[1] = { 2 };
6 static mpint _mptwo =
7 {
8 1,
9 1,
10 1,
11 _mptwodata,
12 MPstatic
13 };
14 mpint *mptwo = &_mptwo;
15
16 static mpdigit _mponedata[1] = { 1 };
17 static mpint _mpone =
18 {
19 1,
20 1,
21 1,
22 _mponedata,
23 MPstatic
24 };
25 mpint *mpone = &_mpone;
26
27 static mpdigit _mpzerodata[1] = { 0 };
28 static mpint _mpzero =
29 {
30 1,
31 1,
32 0,
33 _mpzerodata,
34 MPstatic
35 };
36 mpint *mpzero = &_mpzero;
37
38 static int mpmindigits = 33;
39
40 // set minimum digit allocation
41 void
mpsetminbits(int n)42 mpsetminbits(int n)
43 {
44 if(n == 0)
45 n = 1;
46 mpmindigits = DIGITS(n);
47 }
48
49 // allocate an n bit 0'd number
50 mpint*
mpnew(int n)51 mpnew(int n)
52 {
53 mpint *b;
54
55 b = mallocz(sizeof(mpint), 1);
56 if(b == nil)
57 sysfatal("mpnew: %r");
58 n = DIGITS(n);
59 if(n < mpmindigits)
60 n = mpmindigits;
61 n = n;
62 b->p = (mpdigit*)mallocz(n*Dbytes, 1);
63 if(b->p == nil)
64 sysfatal("mpnew: %r");
65 b->size = n;
66 b->sign = 1;
67
68 return b;
69 }
70
71 // guarantee at least n significant bits
72 void
mpbits(mpint * b,int m)73 mpbits(mpint *b, int m)
74 {
75 int n;
76
77 n = DIGITS(m);
78 if(b->size >= n){
79 if(b->top >= n)
80 return;
81 memset(&b->p[b->top], 0, Dbytes*(n - b->top));
82 b->top = n;
83 return;
84 }
85 b->p = (mpdigit*)realloc(b->p, n*Dbytes);
86 if(b->p == nil)
87 sysfatal("mpbits: %r");
88 memset(&b->p[b->top], 0, Dbytes*(n - b->top));
89 b->size = n;
90 b->top = n;
91 }
92
93 void
mpfree(mpint * b)94 mpfree(mpint *b)
95 {
96 if(b == nil)
97 return;
98 if(b->flags & MPstatic)
99 sysfatal("freeing mp constant");
100 memset(b->p, 0, b->top*Dbytes); // information hiding
101 free(b->p);
102 free(b);
103 }
104
105 void
mpnorm(mpint * b)106 mpnorm(mpint *b)
107 {
108 int i;
109
110 for(i = b->top-1; i >= 0; i--)
111 if(b->p[i] != 0)
112 break;
113 b->top = i+1;
114 if(b->top == 0)
115 b->sign = 1;
116 }
117
118 mpint*
mpcopy(mpint * old)119 mpcopy(mpint *old)
120 {
121 mpint *new;
122
123 new = mpnew(Dbits*old->size);
124 new->top = old->top;
125 new->sign = old->sign;
126 memmove(new->p, old->p, Dbytes*old->top);
127 return new;
128 }
129
130 void
mpassign(mpint * old,mpint * new)131 mpassign(mpint *old, mpint *new)
132 {
133 mpbits(new, Dbits*old->top);
134 new->sign = old->sign;
135 new->top = old->top;
136 memmove(new->p, old->p, Dbytes*old->top);
137 }
138
139 // number of significant bits in mantissa
140 int
mpsignif(mpint * n)141 mpsignif(mpint *n)
142 {
143 int i, j;
144 mpdigit d;
145
146 if(n->top == 0)
147 return 0;
148 for(i = n->top-1; i >= 0; i--){
149 d = n->p[i];
150 for(j = Dbits-1; j >= 0; j--){
151 if(d & (((mpdigit)1)<<j))
152 return i*Dbits + j + 1;
153 }
154 }
155 return 0;
156 }
157
158 // k, where n = 2**k * q for odd q
159 int
mplowbits0(mpint * n)160 mplowbits0(mpint *n)
161 {
162 int k, bit, digit;
163 mpdigit d;
164
165 if(n->top==0)
166 return 0;
167 k = 0;
168 bit = 0;
169 digit = 0;
170 d = n->p[0];
171 for(;;){
172 if(d & (1<<bit))
173 break;
174 k++;
175 bit++;
176 if(bit==Dbits){
177 if(++digit >= n->top)
178 return 0;
179 d = n->p[digit];
180 bit = 0;
181 }
182 }
183 return k;
184 }
185
186