File: | build/source/polly/lib/External/isl/imath/imath.c |
Warning: | line 1145, column 3 Potential leak of memory pointed to by field 'digits' |
Press '?' to see keyboard shortcuts
Keyboard shortcuts:
1 | /* | |||
2 | Name: imath.c | |||
3 | Purpose: Arbitrary precision integer arithmetic routines. | |||
4 | Author: M. J. Fromberger | |||
5 | ||||
6 | Copyright (C) 2002-2007 Michael J. Fromberger, All Rights Reserved. | |||
7 | ||||
8 | Permission is hereby granted, free of charge, to any person obtaining a copy | |||
9 | of this software and associated documentation files (the "Software"), to deal | |||
10 | in the Software without restriction, including without limitation the rights | |||
11 | to use, copy, modify, merge, publish, distribute, sublicense, and/or sell | |||
12 | copies of the Software, and to permit persons to whom the Software is | |||
13 | furnished to do so, subject to the following conditions: | |||
14 | ||||
15 | The above copyright notice and this permission notice shall be included in | |||
16 | all copies or substantial portions of the Software. | |||
17 | ||||
18 | THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR | |||
19 | IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, | |||
20 | FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE | |||
21 | AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER | |||
22 | LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, | |||
23 | OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE | |||
24 | SOFTWARE. | |||
25 | */ | |||
26 | ||||
27 | #include "imath.h" | |||
28 | ||||
29 | #include <assert.h> | |||
30 | #include <ctype.h> | |||
31 | #include <stdlib.h> | |||
32 | #include <string.h> | |||
33 | ||||
34 | const mp_result MP_OK = 0; /* no error, all is well */ | |||
35 | const mp_result MP_FALSE = 0; /* boolean false */ | |||
36 | const mp_result MP_TRUE = -1; /* boolean true */ | |||
37 | const mp_result MP_MEMORY = -2; /* out of memory */ | |||
38 | const mp_result MP_RANGE = -3; /* argument out of range */ | |||
39 | const mp_result MP_UNDEF = -4; /* result undefined */ | |||
40 | const mp_result MP_TRUNC = -5; /* output truncated */ | |||
41 | const mp_result MP_BADARG = -6; /* invalid null argument */ | |||
42 | const mp_result MP_MINERR = -6; | |||
43 | ||||
44 | const mp_sign MP_NEG = 1; /* value is strictly negative */ | |||
45 | const mp_sign MP_ZPOS = 0; /* value is non-negative */ | |||
46 | ||||
47 | static const char *s_unknown_err = "unknown result code"; | |||
48 | static const char *s_error_msg[] = {"error code 0", "boolean true", | |||
49 | "out of memory", "argument out of range", | |||
50 | "result undefined", "output truncated", | |||
51 | "invalid argument", NULL((void*)0)}; | |||
52 | ||||
53 | /* The ith entry of this table gives the value of log_i(2). | |||
54 | ||||
55 | An integer value n requires ceil(log_i(n)) digits to be represented | |||
56 | in base i. Since it is easy to compute lg(n), by counting bits, we | |||
57 | can compute log_i(n) = lg(n) * log_i(2). | |||
58 | ||||
59 | The use of this table eliminates a dependency upon linkage against | |||
60 | the standard math libraries. | |||
61 | ||||
62 | If MP_MAX_RADIX is increased, this table should be expanded too. | |||
63 | */ | |||
64 | static const double s_log2[] = { | |||
65 | 0.000000000, 0.000000000, 1.000000000, 0.630929754, /* (D)(D) 2 3 */ | |||
66 | 0.500000000, 0.430676558, 0.386852807, 0.356207187, /* 4 5 6 7 */ | |||
67 | 0.333333333, 0.315464877, 0.301029996, 0.289064826, /* 8 9 10 11 */ | |||
68 | 0.278942946, 0.270238154, 0.262649535, 0.255958025, /* 12 13 14 15 */ | |||
69 | 0.250000000, 0.244650542, 0.239812467, 0.235408913, /* 16 17 18 19 */ | |||
70 | 0.231378213, 0.227670249, 0.224243824, 0.221064729, /* 20 21 22 23 */ | |||
71 | 0.218104292, 0.215338279, 0.212746054, 0.210309918, /* 24 25 26 27 */ | |||
72 | 0.208014598, 0.205846832, 0.203795047, 0.201849087, /* 28 29 30 31 */ | |||
73 | 0.200000000, 0.198239863, 0.196561632, 0.194959022, /* 32 33 34 35 */ | |||
74 | 0.193426404, /* 36 */ | |||
75 | }; | |||
76 | ||||
77 | /* Return the number of digits needed to represent a static value */ | |||
78 | #define MP_VALUE_DIGITS(V)((sizeof(V) + (sizeof(mp_digit) - 1)) / sizeof(mp_digit)) \ | |||
79 | ((sizeof(V) + (sizeof(mp_digit) - 1)) / sizeof(mp_digit)) | |||
80 | ||||
81 | /* Round precision P to nearest word boundary */ | |||
82 | static inline mp_size s_round_prec(mp_size P) { return 2 * ((P + 1) / 2); } | |||
83 | ||||
84 | /* Set array P of S digits to zero */ | |||
85 | static inline void ZERO(mp_digit *P, mp_size S) { | |||
86 | mp_size i__ = S * sizeof(mp_digit); | |||
87 | mp_digit *p__ = P; | |||
88 | memset(p__, 0, i__); | |||
89 | } | |||
90 | ||||
91 | /* Copy S digits from array P to array Q */ | |||
92 | static inline void COPY(mp_digit *P, mp_digit *Q, mp_size S) { | |||
93 | mp_size i__ = S * sizeof(mp_digit); | |||
94 | mp_digit *p__ = P; | |||
95 | mp_digit *q__ = Q; | |||
96 | memcpy(q__, p__, i__); | |||
97 | } | |||
98 | ||||
99 | /* Reverse N elements of unsigned char in A. */ | |||
100 | static inline void REV(unsigned char *A, int N) { | |||
101 | unsigned char *u_ = A; | |||
102 | unsigned char *v_ = u_ + N - 1; | |||
103 | while (u_ < v_) { | |||
104 | unsigned char xch = *u_; | |||
105 | *u_++ = *v_; | |||
106 | *v_-- = xch; | |||
107 | } | |||
108 | } | |||
109 | ||||
110 | /* Strip leading zeroes from z_ in-place. */ | |||
111 | static inline void CLAMP(mp_int z_) { | |||
112 | mp_size uz_ = MP_USED(z_); | |||
113 | mp_digit *dz_ = MP_DIGITS(z_) + uz_ - 1; | |||
114 | while (uz_ > 1 && (*dz_-- == 0)) --uz_; | |||
115 | z_->used = uz_; | |||
116 | } | |||
117 | ||||
118 | /* Select min/max. */ | |||
119 | static inline int MIN(int A, int B) { return (B < A ? B : A); } | |||
120 | static inline mp_size MAX(mp_size A, mp_size B) { return (B > A ? B : A); } | |||
121 | ||||
122 | /* Exchange lvalues A and B of type T, e.g. | |||
123 | SWAP(int, x, y) where x and y are variables of type int. */ | |||
124 | #define SWAP(T, A, B)do { T t_ = (A); A = (B); B = t_; } while (0) \ | |||
125 | do { \ | |||
126 | T t_ = (A); \ | |||
127 | A = (B); \ | |||
128 | B = t_; \ | |||
129 | } while (0) | |||
130 | ||||
131 | /* Declare a block of N temporary mpz_t values. | |||
132 | These values are initialized to zero. | |||
133 | You must add CLEANUP_TEMP() at the end of the function. | |||
134 | Use TEMP(i) to access a pointer to the ith value. | |||
135 | */ | |||
136 | #define DECLARE_TEMP(N)struct { mpz_t value[(N)]; int len; mp_result err; } temp_ = { .len = (N), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0) \ | |||
137 | struct { \ | |||
138 | mpz_t value[(N)]; \ | |||
139 | int len; \ | |||
140 | mp_result err; \ | |||
141 | } temp_ = { \ | |||
142 | .len = (N), \ | |||
143 | .err = MP_OK, \ | |||
144 | }; \ | |||
145 | do { \ | |||
146 | for (int i = 0; i < temp_.len; i++) { \ | |||
147 | mp_int_init(TEMP(i)(temp_.value + (i))); \ | |||
148 | } \ | |||
149 | } while (0) | |||
150 | ||||
151 | /* Clear all allocated temp values. */ | |||
152 | #define CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0) \ | |||
153 | CLEANUP: \ | |||
154 | do { \ | |||
155 | for (int i = 0; i < temp_.len; i++) { \ | |||
156 | mp_int_clear(TEMP(i)(temp_.value + (i))); \ | |||
157 | } \ | |||
158 | if (temp_.err != MP_OK) { \ | |||
159 | return temp_.err; \ | |||
160 | } \ | |||
161 | } while (0) | |||
162 | ||||
163 | /* A pointer to the kth temp value. */ | |||
164 | #define TEMP(K)(temp_.value + (K)) (temp_.value + (K)) | |||
165 | ||||
166 | /* Evaluate E, an expression of type mp_result expected to return MP_OK. If | |||
167 | the value is not MP_OK, the error is cached and control resumes at the | |||
168 | cleanup handler, which returns it. | |||
169 | */ | |||
170 | #define REQUIRE(E)do { temp_.err = (E); if (temp_.err != MP_OK) goto CLEANUP; } while (0) \ | |||
171 | do { \ | |||
172 | temp_.err = (E); \ | |||
173 | if (temp_.err != MP_OK) goto CLEANUP; \ | |||
174 | } while (0) | |||
175 | ||||
176 | /* Compare value to zero. */ | |||
177 | static inline int CMPZ(mp_int Z) { | |||
178 | if (Z->used == 1 && Z->digits[0] == 0) return 0; | |||
179 | return (Z->sign == MP_NEG) ? -1 : 1; | |||
180 | } | |||
181 | ||||
182 | static inline mp_word UPPER_HALF(mp_word W) { return (W >> MP_DIGIT_BIT(sizeof(mp_digit) * 8)); } | |||
183 | static inline mp_digit LOWER_HALF(mp_word W) { return (mp_digit)(W); } | |||
184 | ||||
185 | /* Report whether the highest-order bit of W is 1. */ | |||
186 | static inline bool_Bool HIGH_BIT_SET(mp_word W) { | |||
187 | return (W >> (MP_WORD_BIT(sizeof(mp_word) * 8) - 1)) != 0; | |||
188 | } | |||
189 | ||||
190 | /* Report whether adding W + V will carry out. */ | |||
191 | static inline bool_Bool ADD_WILL_OVERFLOW(mp_word W, mp_word V) { | |||
192 | return ((MP_WORD_MAX((18446744073709551615UL)) - V) < W); | |||
193 | } | |||
194 | ||||
195 | /* Default number of digits allocated to a new mp_int */ | |||
196 | static mp_size default_precision = 8; | |||
197 | ||||
198 | void mp_int_default_precision(mp_size size) { | |||
199 | assert(size > 0)((void) sizeof ((size > 0) ? 1 : 0), __extension__ ({ if ( size > 0) ; else __assert_fail ("size > 0", "polly/lib/External/isl/imath/imath.c" , 199, __extension__ __PRETTY_FUNCTION__); })); | |||
200 | default_precision = size; | |||
201 | } | |||
202 | ||||
203 | /* Minimum number of digits to invoke recursive multiply */ | |||
204 | static mp_size multiply_threshold = 32; | |||
205 | ||||
206 | void mp_int_multiply_threshold(mp_size thresh) { | |||
207 | assert(thresh >= sizeof(mp_word))((void) sizeof ((thresh >= sizeof(mp_word)) ? 1 : 0), __extension__ ({ if (thresh >= sizeof(mp_word)) ; else __assert_fail ("thresh >= sizeof(mp_word)" , "polly/lib/External/isl/imath/imath.c", 207, __extension__ __PRETTY_FUNCTION__ ); })); | |||
208 | multiply_threshold = thresh; | |||
209 | } | |||
210 | ||||
211 | /* Allocate a buffer of (at least) num digits, or return | |||
212 | NULL if that couldn't be done. */ | |||
213 | static mp_digit *s_alloc(mp_size num); | |||
214 | ||||
215 | /* Release a buffer of digits allocated by s_alloc(). */ | |||
216 | static void s_free(void *ptr); | |||
217 | ||||
218 | /* Insure that z has at least min digits allocated, resizing if | |||
219 | necessary. Returns true if successful, false if out of memory. */ | |||
220 | static bool_Bool s_pad(mp_int z, mp_size min); | |||
221 | ||||
222 | /* Ensure Z has at least N digits allocated. */ | |||
223 | static inline mp_result GROW(mp_int Z, mp_size N) { | |||
224 | return s_pad(Z, N) ? MP_OK : MP_MEMORY; | |||
225 | } | |||
226 | ||||
227 | /* Fill in a "fake" mp_int on the stack with a given value */ | |||
228 | static void s_fake(mp_int z, mp_small value, mp_digit vbuf[]); | |||
229 | static void s_ufake(mp_int z, mp_usmall value, mp_digit vbuf[]); | |||
230 | ||||
231 | /* Compare two runs of digits of given length, returns <0, 0, >0 */ | |||
232 | static int s_cdig(mp_digit *da, mp_digit *db, mp_size len); | |||
233 | ||||
234 | /* Pack the unsigned digits of v into array t */ | |||
235 | static int s_uvpack(mp_usmall v, mp_digit t[]); | |||
236 | ||||
237 | /* Compare magnitudes of a and b, returns <0, 0, >0 */ | |||
238 | static int s_ucmp(mp_int a, mp_int b); | |||
239 | ||||
240 | /* Compare magnitudes of a and v, returns <0, 0, >0 */ | |||
241 | static int s_vcmp(mp_int a, mp_small v); | |||
242 | static int s_uvcmp(mp_int a, mp_usmall uv); | |||
243 | ||||
244 | /* Unsigned magnitude addition; assumes dc is big enough. | |||
245 | Carry out is returned (no memory allocated). */ | |||
246 | static mp_digit s_uadd(mp_digit *da, mp_digit *db, mp_digit *dc, mp_size size_a, | |||
247 | mp_size size_b); | |||
248 | ||||
249 | /* Unsigned magnitude subtraction. Assumes dc is big enough. */ | |||
250 | static void s_usub(mp_digit *da, mp_digit *db, mp_digit *dc, mp_size size_a, | |||
251 | mp_size size_b); | |||
252 | ||||
253 | /* Unsigned recursive multiplication. Assumes dc is big enough. */ | |||
254 | static int s_kmul(mp_digit *da, mp_digit *db, mp_digit *dc, mp_size size_a, | |||
255 | mp_size size_b); | |||
256 | ||||
257 | /* Unsigned magnitude multiplication. Assumes dc is big enough. */ | |||
258 | static void s_umul(mp_digit *da, mp_digit *db, mp_digit *dc, mp_size size_a, | |||
259 | mp_size size_b); | |||
260 | ||||
261 | /* Unsigned recursive squaring. Assumes dc is big enough. */ | |||
262 | static int s_ksqr(mp_digit *da, mp_digit *dc, mp_size size_a); | |||
263 | ||||
264 | /* Unsigned magnitude squaring. Assumes dc is big enough. */ | |||
265 | static void s_usqr(mp_digit *da, mp_digit *dc, mp_size size_a); | |||
266 | ||||
267 | /* Single digit addition. Assumes a is big enough. */ | |||
268 | static void s_dadd(mp_int a, mp_digit b); | |||
269 | ||||
270 | /* Single digit multiplication. Assumes a is big enough. */ | |||
271 | static void s_dmul(mp_int a, mp_digit b); | |||
272 | ||||
273 | /* Single digit multiplication on buffers; assumes dc is big enough. */ | |||
274 | static void s_dbmul(mp_digit *da, mp_digit b, mp_digit *dc, mp_size size_a); | |||
275 | ||||
276 | /* Single digit division. Replaces a with the quotient, | |||
277 | returns the remainder. */ | |||
278 | static mp_digit s_ddiv(mp_int a, mp_digit b); | |||
279 | ||||
280 | /* Quick division by a power of 2, replaces z (no allocation) */ | |||
281 | static void s_qdiv(mp_int z, mp_size p2); | |||
282 | ||||
283 | /* Quick remainder by a power of 2, replaces z (no allocation) */ | |||
284 | static void s_qmod(mp_int z, mp_size p2); | |||
285 | ||||
286 | /* Quick multiplication by a power of 2, replaces z. | |||
287 | Allocates if necessary; returns false in case this fails. */ | |||
288 | static int s_qmul(mp_int z, mp_size p2); | |||
289 | ||||
290 | /* Quick subtraction from a power of 2, replaces z. | |||
291 | Allocates if necessary; returns false in case this fails. */ | |||
292 | static int s_qsub(mp_int z, mp_size p2); | |||
293 | ||||
294 | /* Return maximum k such that 2^k divides z. */ | |||
295 | static int s_dp2k(mp_int z); | |||
296 | ||||
297 | /* Return k >= 0 such that z = 2^k, or -1 if there is no such k. */ | |||
298 | static int s_isp2(mp_int z); | |||
299 | ||||
300 | /* Set z to 2^k. May allocate; returns false in case this fails. */ | |||
301 | static int s_2expt(mp_int z, mp_small k); | |||
302 | ||||
303 | /* Normalize a and b for division, returns normalization constant */ | |||
304 | static int s_norm(mp_int a, mp_int b); | |||
305 | ||||
306 | /* Compute constant mu for Barrett reduction, given modulus m, result | |||
307 | replaces z, m is untouched. */ | |||
308 | static mp_result s_brmu(mp_int z, mp_int m); | |||
309 | ||||
310 | /* Reduce a modulo m, using Barrett's algorithm. */ | |||
311 | static int s_reduce(mp_int x, mp_int m, mp_int mu, mp_int q1, mp_int q2); | |||
312 | ||||
313 | /* Modular exponentiation, using Barrett reduction */ | |||
314 | static mp_result s_embar(mp_int a, mp_int b, mp_int m, mp_int mu, mp_int c); | |||
315 | ||||
316 | /* Unsigned magnitude division. Assumes |a| > |b|. Allocates temporaries; | |||
317 | overwrites a with quotient, b with remainder. */ | |||
318 | static mp_result s_udiv_knuth(mp_int a, mp_int b); | |||
319 | ||||
320 | /* Compute the number of digits in radix r required to represent the given | |||
321 | value. Does not account for sign flags, terminators, etc. */ | |||
322 | static int s_outlen(mp_int z, mp_size r); | |||
323 | ||||
324 | /* Guess how many digits of precision will be needed to represent a radix r | |||
325 | value of the specified number of digits. Returns a value guaranteed to be | |||
326 | no smaller than the actual number required. */ | |||
327 | static mp_size s_inlen(int len, mp_size r); | |||
328 | ||||
329 | /* Convert a character to a digit value in radix r, or | |||
330 | -1 if out of range */ | |||
331 | static int s_ch2val(char c, int r); | |||
332 | ||||
333 | /* Convert a digit value to a character */ | |||
334 | static char s_val2ch(int v, int caps); | |||
335 | ||||
336 | /* Take 2's complement of a buffer in place */ | |||
337 | static void s_2comp(unsigned char *buf, int len); | |||
338 | ||||
339 | /* Convert a value to binary, ignoring sign. On input, *limpos is the bound on | |||
340 | how many bytes should be written to buf; on output, *limpos is set to the | |||
341 | number of bytes actually written. */ | |||
342 | static mp_result s_tobin(mp_int z, unsigned char *buf, int *limpos, int pad); | |||
343 | ||||
344 | /* Multiply X by Y into Z, ignoring signs. Requires that Z have enough storage | |||
345 | preallocated to hold the result. */ | |||
346 | static inline void UMUL(mp_int X, mp_int Y, mp_int Z) { | |||
347 | mp_size ua_ = MP_USED(X); | |||
348 | mp_size ub_ = MP_USED(Y); | |||
349 | mp_size o_ = ua_ + ub_; | |||
350 | ZERO(MP_DIGITS(Z), o_); | |||
351 | (void)s_kmul(MP_DIGITS(X), MP_DIGITS(Y), MP_DIGITS(Z), ua_, ub_); | |||
352 | Z->used = o_; | |||
353 | CLAMP(Z); | |||
354 | } | |||
355 | ||||
356 | /* Square X into Z. Requires that Z have enough storage to hold the result. */ | |||
357 | static inline void USQR(mp_int X, mp_int Z) { | |||
358 | mp_size ua_ = MP_USED(X); | |||
359 | mp_size o_ = ua_ + ua_; | |||
360 | ZERO(MP_DIGITS(Z), o_); | |||
361 | (void)s_ksqr(MP_DIGITS(X), MP_DIGITS(Z), ua_); | |||
362 | Z->used = o_; | |||
363 | CLAMP(Z); | |||
364 | } | |||
365 | ||||
366 | mp_result mp_int_init(mp_int z) { | |||
367 | if (z == NULL((void*)0)) return MP_BADARG; | |||
368 | ||||
369 | z->single = 0; | |||
370 | z->digits = &(z->single); | |||
371 | z->alloc = 1; | |||
372 | z->used = 1; | |||
373 | z->sign = MP_ZPOS; | |||
374 | ||||
375 | return MP_OK; | |||
376 | } | |||
377 | ||||
378 | mp_int mp_int_alloc(void) { | |||
379 | mp_int out = malloc(sizeof(mpz_t)); | |||
380 | ||||
381 | if (out != NULL((void*)0)) mp_int_init(out); | |||
382 | ||||
383 | return out; | |||
384 | } | |||
385 | ||||
386 | mp_result mp_int_init_size(mp_int z, mp_size prec) { | |||
387 | assert(z != NULL)((void) sizeof ((z != ((void*)0)) ? 1 : 0), __extension__ ({ if (z != ((void*)0)) ; else __assert_fail ("z != NULL", "polly/lib/External/isl/imath/imath.c" , 387, __extension__ __PRETTY_FUNCTION__); })); | |||
388 | ||||
389 | if (prec == 0) { | |||
390 | prec = default_precision; | |||
391 | } else if (prec == 1) { | |||
392 | return mp_int_init(z); | |||
393 | } else { | |||
394 | prec = s_round_prec(prec); | |||
395 | } | |||
396 | ||||
397 | z->digits = s_alloc(prec); | |||
398 | if (MP_DIGITS(z) == NULL((void*)0)) return MP_MEMORY; | |||
399 | ||||
400 | z->digits[0] = 0; | |||
401 | z->used = 1; | |||
402 | z->alloc = prec; | |||
403 | z->sign = MP_ZPOS; | |||
404 | ||||
405 | return MP_OK; | |||
406 | } | |||
407 | ||||
408 | mp_result mp_int_init_copy(mp_int z, mp_int old) { | |||
409 | assert(z != NULL && old != NULL)((void) sizeof ((z != ((void*)0) && old != ((void*)0) ) ? 1 : 0), __extension__ ({ if (z != ((void*)0) && old != ((void*)0)) ; else __assert_fail ("z != NULL && old != NULL" , "polly/lib/External/isl/imath/imath.c", 409, __extension__ __PRETTY_FUNCTION__ ); })); | |||
410 | ||||
411 | mp_size uold = MP_USED(old); | |||
412 | if (uold == 1) { | |||
413 | mp_int_init(z); | |||
414 | } else { | |||
415 | mp_size target = MAX(uold, default_precision); | |||
416 | mp_result res = mp_int_init_size(z, target); | |||
417 | if (res != MP_OK) return res; | |||
418 | } | |||
419 | ||||
420 | z->used = uold; | |||
421 | z->sign = old->sign; | |||
422 | COPY(MP_DIGITS(old), MP_DIGITS(z), uold); | |||
423 | ||||
424 | return MP_OK; | |||
425 | } | |||
426 | ||||
427 | mp_result mp_int_init_value(mp_int z, mp_small value) { | |||
428 | mpz_t vtmp; | |||
429 | mp_digit vbuf[MP_VALUE_DIGITS(value)((sizeof(value) + (sizeof(mp_digit) - 1)) / sizeof(mp_digit))]; | |||
430 | ||||
431 | s_fake(&vtmp, value, vbuf); | |||
432 | return mp_int_init_copy(z, &vtmp); | |||
433 | } | |||
434 | ||||
435 | mp_result mp_int_init_uvalue(mp_int z, mp_usmall uvalue) { | |||
436 | mpz_t vtmp; | |||
437 | mp_digit vbuf[MP_VALUE_DIGITS(uvalue)((sizeof(uvalue) + (sizeof(mp_digit) - 1)) / sizeof(mp_digit) )]; | |||
438 | ||||
439 | s_ufake(&vtmp, uvalue, vbuf); | |||
440 | return mp_int_init_copy(z, &vtmp); | |||
441 | } | |||
442 | ||||
443 | mp_result mp_int_set_value(mp_int z, mp_small value) { | |||
444 | mpz_t vtmp; | |||
445 | mp_digit vbuf[MP_VALUE_DIGITS(value)((sizeof(value) + (sizeof(mp_digit) - 1)) / sizeof(mp_digit))]; | |||
446 | ||||
447 | s_fake(&vtmp, value, vbuf); | |||
448 | return mp_int_copy(&vtmp, z); | |||
449 | } | |||
450 | ||||
451 | mp_result mp_int_set_uvalue(mp_int z, mp_usmall uvalue) { | |||
452 | mpz_t vtmp; | |||
453 | mp_digit vbuf[MP_VALUE_DIGITS(uvalue)((sizeof(uvalue) + (sizeof(mp_digit) - 1)) / sizeof(mp_digit) )]; | |||
454 | ||||
455 | s_ufake(&vtmp, uvalue, vbuf); | |||
456 | return mp_int_copy(&vtmp, z); | |||
457 | } | |||
458 | ||||
459 | void mp_int_clear(mp_int z) { | |||
460 | if (z == NULL((void*)0)) return; | |||
461 | ||||
462 | if (MP_DIGITS(z) != NULL((void*)0)) { | |||
463 | if (MP_DIGITS(z) != &(z->single)) s_free(MP_DIGITS(z)); | |||
464 | ||||
465 | z->digits = NULL((void*)0); | |||
466 | } | |||
467 | } | |||
468 | ||||
469 | void mp_int_free(mp_int z) { | |||
470 | assert(z != NULL)((void) sizeof ((z != ((void*)0)) ? 1 : 0), __extension__ ({ if (z != ((void*)0)) ; else __assert_fail ("z != NULL", "polly/lib/External/isl/imath/imath.c" , 470, __extension__ __PRETTY_FUNCTION__); })); | |||
471 | ||||
472 | mp_int_clear(z); | |||
473 | free(z); /* note: NOT s_free() */ | |||
474 | } | |||
475 | ||||
476 | mp_result mp_int_copy(mp_int a, mp_int c) { | |||
477 | assert(a != NULL && c != NULL)((void) sizeof ((a != ((void*)0) && c != ((void*)0)) ? 1 : 0), __extension__ ({ if (a != ((void*)0) && c != ((void*)0)) ; else __assert_fail ("a != NULL && c != NULL" , "polly/lib/External/isl/imath/imath.c", 477, __extension__ __PRETTY_FUNCTION__ ); })); | |||
478 | ||||
479 | if (a != c) { | |||
480 | mp_size ua = MP_USED(a); | |||
481 | mp_digit *da, *dc; | |||
482 | ||||
483 | if (!s_pad(c, ua)) return MP_MEMORY; | |||
484 | ||||
485 | da = MP_DIGITS(a); | |||
486 | dc = MP_DIGITS(c); | |||
487 | COPY(da, dc, ua); | |||
488 | ||||
489 | c->used = ua; | |||
490 | c->sign = a->sign; | |||
491 | } | |||
492 | ||||
493 | return MP_OK; | |||
494 | } | |||
495 | ||||
496 | void mp_int_swap(mp_int a, mp_int c) { | |||
497 | if (a != c) { | |||
498 | mpz_t tmp = *a; | |||
499 | ||||
500 | *a = *c; | |||
501 | *c = tmp; | |||
502 | ||||
503 | if (MP_DIGITS(a) == &(c->single)) a->digits = &(a->single); | |||
504 | if (MP_DIGITS(c) == &(a->single)) c->digits = &(c->single); | |||
505 | } | |||
506 | } | |||
507 | ||||
508 | void mp_int_zero(mp_int z) { | |||
509 | assert(z != NULL)((void) sizeof ((z != ((void*)0)) ? 1 : 0), __extension__ ({ if (z != ((void*)0)) ; else __assert_fail ("z != NULL", "polly/lib/External/isl/imath/imath.c" , 509, __extension__ __PRETTY_FUNCTION__); })); | |||
510 | ||||
511 | z->digits[0] = 0; | |||
512 | z->used = 1; | |||
513 | z->sign = MP_ZPOS; | |||
514 | } | |||
515 | ||||
516 | mp_result mp_int_abs(mp_int a, mp_int c) { | |||
517 | assert(a != NULL && c != NULL)((void) sizeof ((a != ((void*)0) && c != ((void*)0)) ? 1 : 0), __extension__ ({ if (a != ((void*)0) && c != ((void*)0)) ; else __assert_fail ("a != NULL && c != NULL" , "polly/lib/External/isl/imath/imath.c", 517, __extension__ __PRETTY_FUNCTION__ ); })); | |||
518 | ||||
519 | mp_result res; | |||
520 | if ((res = mp_int_copy(a, c)) != MP_OK) return res; | |||
521 | ||||
522 | c->sign = MP_ZPOS; | |||
523 | return MP_OK; | |||
524 | } | |||
525 | ||||
526 | mp_result mp_int_neg(mp_int a, mp_int c) { | |||
527 | assert(a != NULL && c != NULL)((void) sizeof ((a != ((void*)0) && c != ((void*)0)) ? 1 : 0), __extension__ ({ if (a != ((void*)0) && c != ((void*)0)) ; else __assert_fail ("a != NULL && c != NULL" , "polly/lib/External/isl/imath/imath.c", 527, __extension__ __PRETTY_FUNCTION__ ); })); | |||
528 | ||||
529 | mp_result res; | |||
530 | if ((res = mp_int_copy(a, c)) != MP_OK) return res; | |||
531 | ||||
532 | if (CMPZ(c) != 0) c->sign = 1 - MP_SIGN(a); | |||
533 | ||||
534 | return MP_OK; | |||
535 | } | |||
536 | ||||
537 | mp_result mp_int_add(mp_int a, mp_int b, mp_int c) { | |||
538 | assert(a != NULL && b != NULL && c != NULL)((void) sizeof ((a != ((void*)0) && b != ((void*)0) && c != ((void*)0)) ? 1 : 0), __extension__ ({ if (a != ((void* )0) && b != ((void*)0) && c != ((void*)0)) ; else __assert_fail ("a != NULL && b != NULL && c != NULL" , "polly/lib/External/isl/imath/imath.c", 538, __extension__ __PRETTY_FUNCTION__ ); })); | |||
539 | ||||
540 | mp_size ua = MP_USED(a); | |||
541 | mp_size ub = MP_USED(b); | |||
542 | mp_size max = MAX(ua, ub); | |||
543 | ||||
544 | if (MP_SIGN(a) == MP_SIGN(b)) { | |||
545 | /* Same sign -- add magnitudes, preserve sign of addends */ | |||
546 | if (!s_pad(c, max)) return MP_MEMORY; | |||
547 | ||||
548 | mp_digit carry = s_uadd(MP_DIGITS(a), MP_DIGITS(b), MP_DIGITS(c), ua, ub); | |||
549 | mp_size uc = max; | |||
550 | ||||
551 | if (carry) { | |||
552 | if (!s_pad(c, max + 1)) return MP_MEMORY; | |||
553 | ||||
554 | c->digits[max] = carry; | |||
555 | ++uc; | |||
556 | } | |||
557 | ||||
558 | c->used = uc; | |||
559 | c->sign = a->sign; | |||
560 | ||||
561 | } else { | |||
562 | /* Different signs -- subtract magnitudes, preserve sign of greater */ | |||
563 | int cmp = s_ucmp(a, b); /* magnitude comparison, sign ignored */ | |||
564 | ||||
565 | /* Set x to max(a, b), y to min(a, b) to simplify later code. | |||
566 | A special case yields zero for equal magnitudes. | |||
567 | */ | |||
568 | mp_int x, y; | |||
569 | if (cmp == 0) { | |||
570 | mp_int_zero(c); | |||
571 | return MP_OK; | |||
572 | } else if (cmp < 0) { | |||
573 | x = b; | |||
574 | y = a; | |||
575 | } else { | |||
576 | x = a; | |||
577 | y = b; | |||
578 | } | |||
579 | ||||
580 | if (!s_pad(c, MP_USED(x))) return MP_MEMORY; | |||
581 | ||||
582 | /* Subtract smaller from larger */ | |||
583 | s_usub(MP_DIGITS(x), MP_DIGITS(y), MP_DIGITS(c), MP_USED(x), MP_USED(y)); | |||
584 | c->used = x->used; | |||
585 | CLAMP(c); | |||
586 | ||||
587 | /* Give result the sign of the larger */ | |||
588 | c->sign = x->sign; | |||
589 | } | |||
590 | ||||
591 | return MP_OK; | |||
592 | } | |||
593 | ||||
594 | mp_result mp_int_add_value(mp_int a, mp_small value, mp_int c) { | |||
595 | mpz_t vtmp; | |||
596 | mp_digit vbuf[MP_VALUE_DIGITS(value)((sizeof(value) + (sizeof(mp_digit) - 1)) / sizeof(mp_digit))]; | |||
597 | ||||
598 | s_fake(&vtmp, value, vbuf); | |||
599 | ||||
600 | return mp_int_add(a, &vtmp, c); | |||
601 | } | |||
602 | ||||
603 | mp_result mp_int_sub(mp_int a, mp_int b, mp_int c) { | |||
604 | assert(a != NULL && b != NULL && c != NULL)((void) sizeof ((a != ((void*)0) && b != ((void*)0) && c != ((void*)0)) ? 1 : 0), __extension__ ({ if (a != ((void* )0) && b != ((void*)0) && c != ((void*)0)) ; else __assert_fail ("a != NULL && b != NULL && c != NULL" , "polly/lib/External/isl/imath/imath.c", 604, __extension__ __PRETTY_FUNCTION__ ); })); | |||
605 | ||||
606 | mp_size ua = MP_USED(a); | |||
607 | mp_size ub = MP_USED(b); | |||
608 | mp_size max = MAX(ua, ub); | |||
609 | ||||
610 | if (MP_SIGN(a) != MP_SIGN(b)) { | |||
611 | /* Different signs -- add magnitudes and keep sign of a */ | |||
612 | if (!s_pad(c, max)) return MP_MEMORY; | |||
613 | ||||
614 | mp_digit carry = s_uadd(MP_DIGITS(a), MP_DIGITS(b), MP_DIGITS(c), ua, ub); | |||
615 | mp_size uc = max; | |||
616 | ||||
617 | if (carry) { | |||
618 | if (!s_pad(c, max + 1)) return MP_MEMORY; | |||
619 | ||||
620 | c->digits[max] = carry; | |||
621 | ++uc; | |||
622 | } | |||
623 | ||||
624 | c->used = uc; | |||
625 | c->sign = a->sign; | |||
626 | ||||
627 | } else { | |||
628 | /* Same signs -- subtract magnitudes */ | |||
629 | if (!s_pad(c, max)) return MP_MEMORY; | |||
630 | mp_int x, y; | |||
631 | mp_sign osign; | |||
632 | ||||
633 | int cmp = s_ucmp(a, b); | |||
634 | if (cmp >= 0) { | |||
635 | x = a; | |||
636 | y = b; | |||
637 | osign = MP_ZPOS; | |||
638 | } else { | |||
639 | x = b; | |||
640 | y = a; | |||
641 | osign = MP_NEG; | |||
642 | } | |||
643 | ||||
644 | if (MP_SIGN(a) == MP_NEG && cmp != 0) osign = 1 - osign; | |||
645 | ||||
646 | s_usub(MP_DIGITS(x), MP_DIGITS(y), MP_DIGITS(c), MP_USED(x), MP_USED(y)); | |||
647 | c->used = x->used; | |||
648 | CLAMP(c); | |||
649 | ||||
650 | c->sign = osign; | |||
651 | } | |||
652 | ||||
653 | return MP_OK; | |||
654 | } | |||
655 | ||||
656 | mp_result mp_int_sub_value(mp_int a, mp_small value, mp_int c) { | |||
657 | mpz_t vtmp; | |||
658 | mp_digit vbuf[MP_VALUE_DIGITS(value)((sizeof(value) + (sizeof(mp_digit) - 1)) / sizeof(mp_digit))]; | |||
659 | ||||
660 | s_fake(&vtmp, value, vbuf); | |||
661 | ||||
662 | return mp_int_sub(a, &vtmp, c); | |||
663 | } | |||
664 | ||||
665 | mp_result mp_int_mul(mp_int a, mp_int b, mp_int c) { | |||
666 | assert(a != NULL && b != NULL && c != NULL)((void) sizeof ((a != ((void*)0) && b != ((void*)0) && c != ((void*)0)) ? 1 : 0), __extension__ ({ if (a != ((void* )0) && b != ((void*)0) && c != ((void*)0)) ; else __assert_fail ("a != NULL && b != NULL && c != NULL" , "polly/lib/External/isl/imath/imath.c", 666, __extension__ __PRETTY_FUNCTION__ ); })); | |||
667 | ||||
668 | /* If either input is zero, we can shortcut multiplication */ | |||
669 | if (mp_int_compare_zero(a) == 0 || mp_int_compare_zero(b) == 0) { | |||
670 | mp_int_zero(c); | |||
671 | return MP_OK; | |||
672 | } | |||
673 | ||||
674 | /* Output is positive if inputs have same sign, otherwise negative */ | |||
675 | mp_sign osign = (MP_SIGN(a) == MP_SIGN(b)) ? MP_ZPOS : MP_NEG; | |||
676 | ||||
677 | /* If the output is not identical to any of the inputs, we'll write the | |||
678 | results directly; otherwise, allocate a temporary space. */ | |||
679 | mp_size ua = MP_USED(a); | |||
680 | mp_size ub = MP_USED(b); | |||
681 | mp_size osize = MAX(ua, ub); | |||
682 | osize = 4 * ((osize + 1) / 2); | |||
683 | ||||
684 | mp_digit *out; | |||
685 | mp_size p = 0; | |||
686 | if (c == a || c == b) { | |||
687 | p = MAX(s_round_prec(osize), default_precision); | |||
688 | ||||
689 | if ((out = s_alloc(p)) == NULL((void*)0)) return MP_MEMORY; | |||
690 | } else { | |||
691 | if (!s_pad(c, osize)) return MP_MEMORY; | |||
692 | ||||
693 | out = MP_DIGITS(c); | |||
694 | } | |||
695 | ZERO(out, osize); | |||
696 | ||||
697 | if (!s_kmul(MP_DIGITS(a), MP_DIGITS(b), out, ua, ub)) return MP_MEMORY; | |||
698 | ||||
699 | /* If we allocated a new buffer, get rid of whatever memory c was already | |||
700 | using, and fix up its fields to reflect that. | |||
701 | */ | |||
702 | if (out != MP_DIGITS(c)) { | |||
703 | if ((void *)MP_DIGITS(c) != (void *)c) s_free(MP_DIGITS(c)); | |||
704 | c->digits = out; | |||
705 | c->alloc = p; | |||
706 | } | |||
707 | ||||
708 | c->used = osize; /* might not be true, but we'll fix it ... */ | |||
709 | CLAMP(c); /* ... right here */ | |||
710 | c->sign = osign; | |||
711 | ||||
712 | return MP_OK; | |||
713 | } | |||
714 | ||||
715 | mp_result mp_int_mul_value(mp_int a, mp_small value, mp_int c) { | |||
716 | mpz_t vtmp; | |||
717 | mp_digit vbuf[MP_VALUE_DIGITS(value)((sizeof(value) + (sizeof(mp_digit) - 1)) / sizeof(mp_digit))]; | |||
718 | ||||
719 | s_fake(&vtmp, value, vbuf); | |||
720 | ||||
721 | return mp_int_mul(a, &vtmp, c); | |||
722 | } | |||
723 | ||||
724 | mp_result mp_int_mul_pow2(mp_int a, mp_small p2, mp_int c) { | |||
725 | assert(a != NULL && c != NULL && p2 >= 0)((void) sizeof ((a != ((void*)0) && c != ((void*)0) && p2 >= 0) ? 1 : 0), __extension__ ({ if (a != ((void*)0) && c != ((void*)0) && p2 >= 0) ; else __assert_fail ( "a != NULL && c != NULL && p2 >= 0", "polly/lib/External/isl/imath/imath.c" , 725, __extension__ __PRETTY_FUNCTION__); })); | |||
726 | ||||
727 | mp_result res = mp_int_copy(a, c); | |||
728 | if (res != MP_OK) return res; | |||
729 | ||||
730 | if (s_qmul(c, (mp_size)p2)) { | |||
731 | return MP_OK; | |||
732 | } else { | |||
733 | return MP_MEMORY; | |||
734 | } | |||
735 | } | |||
736 | ||||
737 | mp_result mp_int_sqr(mp_int a, mp_int c) { | |||
738 | assert(a != NULL && c != NULL)((void) sizeof ((a != ((void*)0) && c != ((void*)0)) ? 1 : 0), __extension__ ({ if (a != ((void*)0) && c != ((void*)0)) ; else __assert_fail ("a != NULL && c != NULL" , "polly/lib/External/isl/imath/imath.c", 738, __extension__ __PRETTY_FUNCTION__ ); })); | |||
739 | ||||
740 | /* Get a temporary buffer big enough to hold the result */ | |||
741 | mp_size osize = (mp_size)4 * ((MP_USED(a) + 1) / 2); | |||
742 | mp_size p = 0; | |||
743 | mp_digit *out; | |||
744 | if (a == c) { | |||
745 | p = s_round_prec(osize); | |||
746 | p = MAX(p, default_precision); | |||
747 | ||||
748 | if ((out = s_alloc(p)) == NULL((void*)0)) return MP_MEMORY; | |||
749 | } else { | |||
750 | if (!s_pad(c, osize)) return MP_MEMORY; | |||
751 | ||||
752 | out = MP_DIGITS(c); | |||
753 | } | |||
754 | ZERO(out, osize); | |||
755 | ||||
756 | s_ksqr(MP_DIGITS(a), out, MP_USED(a)); | |||
757 | ||||
758 | /* Get rid of whatever memory c was already using, and fix up its fields to | |||
759 | reflect the new digit array it's using | |||
760 | */ | |||
761 | if (out != MP_DIGITS(c)) { | |||
762 | if ((void *)MP_DIGITS(c) != (void *)c) s_free(MP_DIGITS(c)); | |||
763 | c->digits = out; | |||
764 | c->alloc = p; | |||
765 | } | |||
766 | ||||
767 | c->used = osize; /* might not be true, but we'll fix it ... */ | |||
768 | CLAMP(c); /* ... right here */ | |||
769 | c->sign = MP_ZPOS; | |||
770 | ||||
771 | return MP_OK; | |||
772 | } | |||
773 | ||||
774 | mp_result mp_int_div(mp_int a, mp_int b, mp_int q, mp_int r) { | |||
775 | assert(a != NULL && b != NULL && q != r)((void) sizeof ((a != ((void*)0) && b != ((void*)0) && q != r) ? 1 : 0), __extension__ ({ if (a != ((void*)0) && b != ((void*)0) && q != r) ; else __assert_fail ("a != NULL && b != NULL && q != r" , "polly/lib/External/isl/imath/imath.c", 775, __extension__ __PRETTY_FUNCTION__ ); })); | |||
776 | ||||
777 | int cmp; | |||
778 | mp_result res = MP_OK; | |||
779 | mp_int qout, rout; | |||
780 | mp_sign sa = MP_SIGN(a); | |||
781 | mp_sign sb = MP_SIGN(b); | |||
782 | if (CMPZ(b) == 0) { | |||
783 | return MP_UNDEF; | |||
784 | } else if ((cmp = s_ucmp(a, b)) < 0) { | |||
785 | /* If |a| < |b|, no division is required: | |||
786 | q = 0, r = a | |||
787 | */ | |||
788 | if (r && (res = mp_int_copy(a, r)) != MP_OK) return res; | |||
789 | ||||
790 | if (q) mp_int_zero(q); | |||
791 | ||||
792 | return MP_OK; | |||
793 | } else if (cmp == 0) { | |||
794 | /* If |a| = |b|, no division is required: | |||
795 | q = 1 or -1, r = 0 | |||
796 | */ | |||
797 | if (r) mp_int_zero(r); | |||
798 | ||||
799 | if (q) { | |||
800 | mp_int_zero(q); | |||
801 | q->digits[0] = 1; | |||
802 | ||||
803 | if (sa != sb) q->sign = MP_NEG; | |||
804 | } | |||
805 | ||||
806 | return MP_OK; | |||
807 | } | |||
808 | ||||
809 | /* When |a| > |b|, real division is required. We need someplace to store | |||
810 | quotient and remainder, but q and r are allowed to be NULL or to overlap | |||
811 | with the inputs. | |||
812 | */ | |||
813 | DECLARE_TEMP(2)struct { mpz_t value[(2)]; int len; mp_result err; } temp_ = { .len = (2), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0); | |||
814 | int lg; | |||
815 | if ((lg = s_isp2(b)) < 0) { | |||
816 | if (q && b != q) { | |||
817 | REQUIRE(mp_int_copy(a, q))do { temp_.err = (mp_int_copy(a, q)); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
818 | qout = q; | |||
819 | } else { | |||
820 | REQUIRE(mp_int_copy(a, TEMP(0)))do { temp_.err = (mp_int_copy(a, (temp_.value + (0)))); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
821 | qout = TEMP(0)(temp_.value + (0)); | |||
822 | } | |||
823 | ||||
824 | if (r && a != r) { | |||
825 | REQUIRE(mp_int_copy(b, r))do { temp_.err = (mp_int_copy(b, r)); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
826 | rout = r; | |||
827 | } else { | |||
828 | REQUIRE(mp_int_copy(b, TEMP(1)))do { temp_.err = (mp_int_copy(b, (temp_.value + (1)))); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
829 | rout = TEMP(1)(temp_.value + (1)); | |||
830 | } | |||
831 | ||||
832 | REQUIRE(s_udiv_knuth(qout, rout))do { temp_.err = (s_udiv_knuth(qout, rout)); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
833 | } else { | |||
834 | if (q) REQUIRE(mp_int_copy(a, q))do { temp_.err = (mp_int_copy(a, q)); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
835 | if (r) REQUIRE(mp_int_copy(a, r))do { temp_.err = (mp_int_copy(a, r)); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
836 | ||||
837 | if (q) s_qdiv(q, (mp_size)lg); | |||
838 | qout = q; | |||
839 | if (r) s_qmod(r, (mp_size)lg); | |||
840 | rout = r; | |||
841 | } | |||
842 | ||||
843 | /* Recompute signs for output */ | |||
844 | if (rout) { | |||
845 | rout->sign = sa; | |||
846 | if (CMPZ(rout) == 0) rout->sign = MP_ZPOS; | |||
847 | } | |||
848 | if (qout) { | |||
849 | qout->sign = (sa == sb) ? MP_ZPOS : MP_NEG; | |||
850 | if (CMPZ(qout) == 0) qout->sign = MP_ZPOS; | |||
851 | } | |||
852 | ||||
853 | if (q) REQUIRE(mp_int_copy(qout, q))do { temp_.err = (mp_int_copy(qout, q)); if (temp_.err != MP_OK ) goto CLEANUP; } while (0); | |||
854 | if (r) REQUIRE(mp_int_copy(rout, r))do { temp_.err = (mp_int_copy(rout, r)); if (temp_.err != MP_OK ) goto CLEANUP; } while (0); | |||
855 | CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0); | |||
856 | return res; | |||
857 | } | |||
858 | ||||
859 | mp_result mp_int_mod(mp_int a, mp_int m, mp_int c) { | |||
860 | DECLARE_TEMP(1)struct { mpz_t value[(1)]; int len; mp_result err; } temp_ = { .len = (1), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0); | |||
861 | mp_int out = (m
| |||
862 | REQUIRE(mp_int_div(a, m, NULL, out))do { temp_.err = (mp_int_div(a, m, ((void*)0), out)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
863 | if (CMPZ(out) < 0) { | |||
864 | REQUIRE(mp_int_add(out, m, c))do { temp_.err = (mp_int_add(out, m, c)); if (temp_.err != MP_OK ) goto CLEANUP; } while (0); | |||
865 | } else { | |||
866 | REQUIRE(mp_int_copy(out, c))do { temp_.err = (mp_int_copy(out, c)); if (temp_.err != MP_OK ) goto CLEANUP; } while (0); | |||
867 | } | |||
868 | CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0); | |||
869 | return MP_OK; | |||
870 | } | |||
871 | ||||
872 | mp_result mp_int_div_value(mp_int a, mp_small value, mp_int q, mp_small *r) { | |||
873 | mpz_t vtmp; | |||
874 | mp_digit vbuf[MP_VALUE_DIGITS(value)((sizeof(value) + (sizeof(mp_digit) - 1)) / sizeof(mp_digit))]; | |||
875 | s_fake(&vtmp, value, vbuf); | |||
876 | ||||
877 | DECLARE_TEMP(1)struct { mpz_t value[(1)]; int len; mp_result err; } temp_ = { .len = (1), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0); | |||
878 | REQUIRE(mp_int_div(a, &vtmp, q, TEMP(0)))do { temp_.err = (mp_int_div(a, &vtmp, q, (temp_.value + ( 0)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
879 | ||||
880 | if (r) (void)mp_int_to_int(TEMP(0)(temp_.value + (0)), r); /* can't fail */ | |||
881 | ||||
882 | CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0); | |||
883 | return MP_OK; | |||
884 | } | |||
885 | ||||
886 | mp_result mp_int_div_pow2(mp_int a, mp_small p2, mp_int q, mp_int r) { | |||
887 | assert(a != NULL && p2 >= 0 && q != r)((void) sizeof ((a != ((void*)0) && p2 >= 0 && q != r) ? 1 : 0), __extension__ ({ if (a != ((void*)0) && p2 >= 0 && q != r) ; else __assert_fail ("a != NULL && p2 >= 0 && q != r" , "polly/lib/External/isl/imath/imath.c", 887, __extension__ __PRETTY_FUNCTION__ ); })); | |||
888 | ||||
889 | mp_result res = MP_OK; | |||
890 | if (q != NULL((void*)0) && (res = mp_int_copy(a, q)) == MP_OK) { | |||
891 | s_qdiv(q, (mp_size)p2); | |||
892 | } | |||
893 | ||||
894 | if (res == MP_OK && r != NULL((void*)0) && (res = mp_int_copy(a, r)) == MP_OK) { | |||
895 | s_qmod(r, (mp_size)p2); | |||
896 | } | |||
897 | ||||
898 | return res; | |||
899 | } | |||
900 | ||||
901 | mp_result mp_int_expt(mp_int a, mp_small b, mp_int c) { | |||
902 | assert(c != NULL)((void) sizeof ((c != ((void*)0)) ? 1 : 0), __extension__ ({ if (c != ((void*)0)) ; else __assert_fail ("c != NULL", "polly/lib/External/isl/imath/imath.c" , 902, __extension__ __PRETTY_FUNCTION__); })); | |||
903 | if (b < 0) return MP_RANGE; | |||
904 | ||||
905 | DECLARE_TEMP(1)struct { mpz_t value[(1)]; int len; mp_result err; } temp_ = { .len = (1), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0); | |||
906 | REQUIRE(mp_int_copy(a, TEMP(0)))do { temp_.err = (mp_int_copy(a, (temp_.value + (0)))); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
907 | ||||
908 | (void)mp_int_set_value(c, 1); | |||
909 | unsigned int v = labs(b); | |||
910 | while (v != 0) { | |||
911 | if (v & 1) { | |||
912 | REQUIRE(mp_int_mul(c, TEMP(0), c))do { temp_.err = (mp_int_mul(c, (temp_.value + (0)), c)); if ( temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
913 | } | |||
914 | ||||
915 | v >>= 1; | |||
916 | if (v == 0) break; | |||
917 | ||||
918 | REQUIRE(mp_int_sqr(TEMP(0), TEMP(0)))do { temp_.err = (mp_int_sqr((temp_.value + (0)), (temp_.value + (0)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
919 | } | |||
920 | ||||
921 | CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0); | |||
922 | return MP_OK; | |||
923 | } | |||
924 | ||||
925 | mp_result mp_int_expt_value(mp_small a, mp_small b, mp_int c) { | |||
926 | assert(c != NULL)((void) sizeof ((c != ((void*)0)) ? 1 : 0), __extension__ ({ if (c != ((void*)0)) ; else __assert_fail ("c != NULL", "polly/lib/External/isl/imath/imath.c" , 926, __extension__ __PRETTY_FUNCTION__); })); | |||
927 | if (b < 0) return MP_RANGE; | |||
928 | ||||
929 | DECLARE_TEMP(1)struct { mpz_t value[(1)]; int len; mp_result err; } temp_ = { .len = (1), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0); | |||
930 | REQUIRE(mp_int_set_value(TEMP(0), a))do { temp_.err = (mp_int_set_value((temp_.value + (0)), a)); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
931 | ||||
932 | (void)mp_int_set_value(c, 1); | |||
933 | unsigned int v = labs(b); | |||
934 | while (v != 0) { | |||
935 | if (v & 1) { | |||
936 | REQUIRE(mp_int_mul(c, TEMP(0), c))do { temp_.err = (mp_int_mul(c, (temp_.value + (0)), c)); if ( temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
937 | } | |||
938 | ||||
939 | v >>= 1; | |||
940 | if (v == 0) break; | |||
941 | ||||
942 | REQUIRE(mp_int_sqr(TEMP(0), TEMP(0)))do { temp_.err = (mp_int_sqr((temp_.value + (0)), (temp_.value + (0)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
943 | } | |||
944 | ||||
945 | CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0); | |||
946 | return MP_OK; | |||
947 | } | |||
948 | ||||
949 | mp_result mp_int_expt_full(mp_int a, mp_int b, mp_int c) { | |||
950 | assert(a != NULL && b != NULL && c != NULL)((void) sizeof ((a != ((void*)0) && b != ((void*)0) && c != ((void*)0)) ? 1 : 0), __extension__ ({ if (a != ((void* )0) && b != ((void*)0) && c != ((void*)0)) ; else __assert_fail ("a != NULL && b != NULL && c != NULL" , "polly/lib/External/isl/imath/imath.c", 950, __extension__ __PRETTY_FUNCTION__ ); })); | |||
951 | if (MP_SIGN(b) == MP_NEG) return MP_RANGE; | |||
952 | ||||
953 | DECLARE_TEMP(1)struct { mpz_t value[(1)]; int len; mp_result err; } temp_ = { .len = (1), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0); | |||
954 | REQUIRE(mp_int_copy(a, TEMP(0)))do { temp_.err = (mp_int_copy(a, (temp_.value + (0)))); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
955 | ||||
956 | (void)mp_int_set_value(c, 1); | |||
957 | for (unsigned ix = 0; ix < MP_USED(b); ++ix) { | |||
958 | mp_digit d = b->digits[ix]; | |||
959 | ||||
960 | for (unsigned jx = 0; jx < MP_DIGIT_BIT(sizeof(mp_digit) * 8); ++jx) { | |||
961 | if (d & 1) { | |||
962 | REQUIRE(mp_int_mul(c, TEMP(0), c))do { temp_.err = (mp_int_mul(c, (temp_.value + (0)), c)); if ( temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
963 | } | |||
964 | ||||
965 | d >>= 1; | |||
966 | if (d == 0 && ix + 1 == MP_USED(b)) break; | |||
967 | REQUIRE(mp_int_sqr(TEMP(0), TEMP(0)))do { temp_.err = (mp_int_sqr((temp_.value + (0)), (temp_.value + (0)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
968 | } | |||
969 | } | |||
970 | ||||
971 | CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0); | |||
972 | return MP_OK; | |||
973 | } | |||
974 | ||||
975 | int mp_int_compare(mp_int a, mp_int b) { | |||
976 | assert(a != NULL && b != NULL)((void) sizeof ((a != ((void*)0) && b != ((void*)0)) ? 1 : 0), __extension__ ({ if (a != ((void*)0) && b != ((void*)0)) ; else __assert_fail ("a != NULL && b != NULL" , "polly/lib/External/isl/imath/imath.c", 976, __extension__ __PRETTY_FUNCTION__ ); })); | |||
977 | ||||
978 | mp_sign sa = MP_SIGN(a); | |||
979 | if (sa == MP_SIGN(b)) { | |||
980 | int cmp = s_ucmp(a, b); | |||
981 | ||||
982 | /* If they're both zero or positive, the normal comparison applies; if both | |||
983 | negative, the sense is reversed. */ | |||
984 | if (sa == MP_ZPOS) { | |||
985 | return cmp; | |||
986 | } else { | |||
987 | return -cmp; | |||
988 | } | |||
989 | } else if (sa == MP_ZPOS) { | |||
990 | return 1; | |||
991 | } else { | |||
992 | return -1; | |||
993 | } | |||
994 | } | |||
995 | ||||
996 | int mp_int_compare_unsigned(mp_int a, mp_int b) { | |||
997 | assert(a != NULL && b != NULL)((void) sizeof ((a != ((void*)0) && b != ((void*)0)) ? 1 : 0), __extension__ ({ if (a != ((void*)0) && b != ((void*)0)) ; else __assert_fail ("a != NULL && b != NULL" , "polly/lib/External/isl/imath/imath.c", 997, __extension__ __PRETTY_FUNCTION__ ); })); | |||
998 | ||||
999 | return s_ucmp(a, b); | |||
1000 | } | |||
1001 | ||||
1002 | int mp_int_compare_zero(mp_int z) { | |||
1003 | assert(z != NULL)((void) sizeof ((z != ((void*)0)) ? 1 : 0), __extension__ ({ if (z != ((void*)0)) ; else __assert_fail ("z != NULL", "polly/lib/External/isl/imath/imath.c" , 1003, __extension__ __PRETTY_FUNCTION__); })); | |||
1004 | ||||
1005 | if (MP_USED(z) == 1 && z->digits[0] == 0) { | |||
1006 | return 0; | |||
1007 | } else if (MP_SIGN(z) == MP_ZPOS) { | |||
1008 | return 1; | |||
1009 | } else { | |||
1010 | return -1; | |||
1011 | } | |||
1012 | } | |||
1013 | ||||
1014 | int mp_int_compare_value(mp_int z, mp_small value) { | |||
1015 | assert(z != NULL)((void) sizeof ((z != ((void*)0)) ? 1 : 0), __extension__ ({ if (z != ((void*)0)) ; else __assert_fail ("z != NULL", "polly/lib/External/isl/imath/imath.c" , 1015, __extension__ __PRETTY_FUNCTION__); })); | |||
1016 | ||||
1017 | mp_sign vsign = (value < 0) ? MP_NEG : MP_ZPOS; | |||
1018 | if (vsign == MP_SIGN(z)) { | |||
1019 | int cmp = s_vcmp(z, value); | |||
1020 | ||||
1021 | return (vsign == MP_ZPOS) ? cmp : -cmp; | |||
1022 | } else { | |||
1023 | return (value < 0) ? 1 : -1; | |||
1024 | } | |||
1025 | } | |||
1026 | ||||
1027 | int mp_int_compare_uvalue(mp_int z, mp_usmall uv) { | |||
1028 | assert(z != NULL)((void) sizeof ((z != ((void*)0)) ? 1 : 0), __extension__ ({ if (z != ((void*)0)) ; else __assert_fail ("z != NULL", "polly/lib/External/isl/imath/imath.c" , 1028, __extension__ __PRETTY_FUNCTION__); })); | |||
1029 | ||||
1030 | if (MP_SIGN(z) == MP_NEG) { | |||
1031 | return -1; | |||
1032 | } else { | |||
1033 | return s_uvcmp(z, uv); | |||
1034 | } | |||
1035 | } | |||
1036 | ||||
1037 | mp_result mp_int_exptmod(mp_int a, mp_int b, mp_int m, mp_int c) { | |||
1038 | assert(a != NULL && b != NULL && c != NULL && m != NULL)((void) sizeof ((a != ((void*)0) && b != ((void*)0) && c != ((void*)0) && m != ((void*)0)) ? 1 : 0), __extension__ ({ if (a != ((void*)0) && b != ((void*)0) && c != ((void*)0) && m != ((void*)0)) ; else __assert_fail ("a != NULL && b != NULL && c != NULL && m != NULL" , "polly/lib/External/isl/imath/imath.c", 1038, __extension__ __PRETTY_FUNCTION__); })); | |||
1039 | ||||
1040 | /* Zero moduli and negative exponents are not considered. */ | |||
1041 | if (CMPZ(m) == 0) return MP_UNDEF; | |||
1042 | if (CMPZ(b) < 0) return MP_RANGE; | |||
1043 | ||||
1044 | mp_size um = MP_USED(m); | |||
1045 | DECLARE_TEMP(3)struct { mpz_t value[(3)]; int len; mp_result err; } temp_ = { .len = (3), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0); | |||
1046 | REQUIRE(GROW(TEMP(0), 2 * um))do { temp_.err = (GROW((temp_.value + (0)), 2 * um)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1047 | REQUIRE(GROW(TEMP(1), 2 * um))do { temp_.err = (GROW((temp_.value + (1)), 2 * um)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1048 | ||||
1049 | mp_int s; | |||
1050 | if (c == b || c == m) { | |||
1051 | REQUIRE(GROW(TEMP(2), 2 * um))do { temp_.err = (GROW((temp_.value + (2)), 2 * um)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1052 | s = TEMP(2)(temp_.value + (2)); | |||
1053 | } else { | |||
1054 | s = c; | |||
1055 | } | |||
1056 | ||||
1057 | REQUIRE(mp_int_mod(a, m, TEMP(0)))do { temp_.err = (mp_int_mod(a, m, (temp_.value + (0)))); if ( temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1058 | REQUIRE(s_brmu(TEMP(1), m))do { temp_.err = (s_brmu((temp_.value + (1)), m)); if (temp_. err != MP_OK) goto CLEANUP; } while (0); | |||
1059 | REQUIRE(s_embar(TEMP(0), b, m, TEMP(1), s))do { temp_.err = (s_embar((temp_.value + (0)), b, m, (temp_.value + (1)), s)); if (temp_.err != MP_OK) goto CLEANUP; } while ( 0); | |||
1060 | REQUIRE(mp_int_copy(s, c))do { temp_.err = (mp_int_copy(s, c)); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1061 | ||||
1062 | CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0); | |||
1063 | return MP_OK; | |||
1064 | } | |||
1065 | ||||
1066 | mp_result mp_int_exptmod_evalue(mp_int a, mp_small value, mp_int m, mp_int c) { | |||
1067 | mpz_t vtmp; | |||
1068 | mp_digit vbuf[MP_VALUE_DIGITS(value)((sizeof(value) + (sizeof(mp_digit) - 1)) / sizeof(mp_digit))]; | |||
1069 | ||||
1070 | s_fake(&vtmp, value, vbuf); | |||
1071 | ||||
1072 | return mp_int_exptmod(a, &vtmp, m, c); | |||
1073 | } | |||
1074 | ||||
1075 | mp_result mp_int_exptmod_bvalue(mp_small value, mp_int b, mp_int m, mp_int c) { | |||
1076 | mpz_t vtmp; | |||
1077 | mp_digit vbuf[MP_VALUE_DIGITS(value)((sizeof(value) + (sizeof(mp_digit) - 1)) / sizeof(mp_digit))]; | |||
1078 | ||||
1079 | s_fake(&vtmp, value, vbuf); | |||
1080 | ||||
1081 | return mp_int_exptmod(&vtmp, b, m, c); | |||
1082 | } | |||
1083 | ||||
1084 | mp_result mp_int_exptmod_known(mp_int a, mp_int b, mp_int m, mp_int mu, | |||
1085 | mp_int c) { | |||
1086 | assert(a && b && m && c)((void) sizeof ((a && b && m && c) ? 1 : 0), __extension__ ({ if (a && b && m && c) ; else __assert_fail ("a && b && m && c" , "polly/lib/External/isl/imath/imath.c", 1086, __extension__ __PRETTY_FUNCTION__); })); | |||
1087 | ||||
1088 | /* Zero moduli and negative exponents are not considered. */ | |||
1089 | if (CMPZ(m) == 0) return MP_UNDEF; | |||
1090 | if (CMPZ(b) < 0) return MP_RANGE; | |||
1091 | ||||
1092 | DECLARE_TEMP(2)struct { mpz_t value[(2)]; int len; mp_result err; } temp_ = { .len = (2), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0); | |||
1093 | mp_size um = MP_USED(m); | |||
1094 | REQUIRE(GROW(TEMP(0), 2 * um))do { temp_.err = (GROW((temp_.value + (0)), 2 * um)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1095 | ||||
1096 | mp_int s; | |||
1097 | if (c == b || c == m) { | |||
1098 | REQUIRE(GROW(TEMP(1), 2 * um))do { temp_.err = (GROW((temp_.value + (1)), 2 * um)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1099 | s = TEMP(1)(temp_.value + (1)); | |||
1100 | } else { | |||
1101 | s = c; | |||
1102 | } | |||
1103 | ||||
1104 | REQUIRE(mp_int_mod(a, m, TEMP(0)))do { temp_.err = (mp_int_mod(a, m, (temp_.value + (0)))); if ( temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1105 | REQUIRE(s_embar(TEMP(0), b, m, mu, s))do { temp_.err = (s_embar((temp_.value + (0)), b, m, mu, s)); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1106 | REQUIRE(mp_int_copy(s, c))do { temp_.err = (mp_int_copy(s, c)); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1107 | ||||
1108 | CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0); | |||
1109 | return MP_OK; | |||
1110 | } | |||
1111 | ||||
1112 | mp_result mp_int_redux_const(mp_int m, mp_int c) { | |||
1113 | assert(m != NULL && c != NULL && m != c)((void) sizeof ((m != ((void*)0) && c != ((void*)0) && m != c) ? 1 : 0), __extension__ ({ if (m != ((void*)0) && c != ((void*)0) && m != c) ; else __assert_fail ("m != NULL && c != NULL && m != c" , "polly/lib/External/isl/imath/imath.c", 1113, __extension__ __PRETTY_FUNCTION__); })); | |||
1114 | ||||
1115 | return s_brmu(c, m); | |||
1116 | } | |||
1117 | ||||
1118 | mp_result mp_int_invmod(mp_int a, mp_int m, mp_int c) { | |||
1119 | assert(a != NULL && m != NULL && c != NULL)((void) sizeof ((a != ((void*)0) && m != ((void*)0) && c != ((void*)0)) ? 1 : 0), __extension__ ({ if (a != ((void* )0) && m != ((void*)0) && c != ((void*)0)) ; else __assert_fail ("a != NULL && m != NULL && c != NULL" , "polly/lib/External/isl/imath/imath.c", 1119, __extension__ __PRETTY_FUNCTION__); })); | |||
| ||||
1120 | ||||
1121 | if (CMPZ(a) == 0 || CMPZ(m) <= 0) return MP_RANGE; | |||
1122 | ||||
1123 | DECLARE_TEMP(2)struct { mpz_t value[(2)]; int len; mp_result err; } temp_ = { .len = (2), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0); | |||
1124 | ||||
1125 | REQUIRE(mp_int_egcd(a, m, TEMP(0), TEMP(1), NULL))do { temp_.err = (mp_int_egcd(a, m, (temp_.value + (0)), (temp_ .value + (1)), ((void*)0))); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1126 | ||||
1127 | if (mp_int_compare_value(TEMP(0)(temp_.value + (0)), 1) != 0) { | |||
1128 | REQUIRE(MP_UNDEF)do { temp_.err = (MP_UNDEF); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1129 | } | |||
1130 | ||||
1131 | /* It is first necessary to constrain the value to the proper range */ | |||
1132 | REQUIRE(mp_int_mod(TEMP(1), m, TEMP(1)))do { temp_.err = (mp_int_mod((temp_.value + (1)), m, (temp_.value + (1)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1133 | ||||
1134 | /* Now, if 'a' was originally negative, the value we have is actually the | |||
1135 | magnitude of the negative representative; to get the positive value we | |||
1136 | have to subtract from the modulus. Otherwise, the value is okay as it | |||
1137 | stands. | |||
1138 | */ | |||
1139 | if (MP_SIGN(a) == MP_NEG) { | |||
1140 | REQUIRE(mp_int_sub(m, TEMP(1), c))do { temp_.err = (mp_int_sub(m, (temp_.value + (1)), c)); if ( temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1141 | } else { | |||
1142 | REQUIRE(mp_int_copy(TEMP(1), c))do { temp_.err = (mp_int_copy((temp_.value + (1)), c)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1143 | } | |||
1144 | ||||
1145 | CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0); | |||
| ||||
1146 | return MP_OK; | |||
1147 | } | |||
1148 | ||||
1149 | /* Binary GCD algorithm due to Josef Stein, 1961 */ | |||
1150 | mp_result mp_int_gcd(mp_int a, mp_int b, mp_int c) { | |||
1151 | assert(a != NULL && b != NULL && c != NULL)((void) sizeof ((a != ((void*)0) && b != ((void*)0) && c != ((void*)0)) ? 1 : 0), __extension__ ({ if (a != ((void* )0) && b != ((void*)0) && c != ((void*)0)) ; else __assert_fail ("a != NULL && b != NULL && c != NULL" , "polly/lib/External/isl/imath/imath.c", 1151, __extension__ __PRETTY_FUNCTION__); })); | |||
1152 | ||||
1153 | int ca = CMPZ(a); | |||
1154 | int cb = CMPZ(b); | |||
1155 | if (ca == 0 && cb == 0) { | |||
1156 | return MP_UNDEF; | |||
1157 | } else if (ca == 0) { | |||
1158 | return mp_int_abs(b, c); | |||
1159 | } else if (cb == 0) { | |||
1160 | return mp_int_abs(a, c); | |||
1161 | } | |||
1162 | ||||
1163 | DECLARE_TEMP(3)struct { mpz_t value[(3)]; int len; mp_result err; } temp_ = { .len = (3), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0); | |||
1164 | REQUIRE(mp_int_copy(a, TEMP(0)))do { temp_.err = (mp_int_copy(a, (temp_.value + (0)))); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1165 | REQUIRE(mp_int_copy(b, TEMP(1)))do { temp_.err = (mp_int_copy(b, (temp_.value + (1)))); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1166 | ||||
1167 | TEMP(0)(temp_.value + (0))->sign = MP_ZPOS; | |||
1168 | TEMP(1)(temp_.value + (1))->sign = MP_ZPOS; | |||
1169 | ||||
1170 | int k = 0; | |||
1171 | { /* Divide out common factors of 2 from u and v */ | |||
1172 | int div2_u = s_dp2k(TEMP(0)(temp_.value + (0))); | |||
1173 | int div2_v = s_dp2k(TEMP(1)(temp_.value + (1))); | |||
1174 | ||||
1175 | k = MIN(div2_u, div2_v); | |||
1176 | s_qdiv(TEMP(0)(temp_.value + (0)), (mp_size)k); | |||
1177 | s_qdiv(TEMP(1)(temp_.value + (1)), (mp_size)k); | |||
1178 | } | |||
1179 | ||||
1180 | if (mp_int_is_odd(TEMP(0)(temp_.value + (0)))) { | |||
1181 | REQUIRE(mp_int_neg(TEMP(1), TEMP(2)))do { temp_.err = (mp_int_neg((temp_.value + (1)), (temp_.value + (2)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1182 | } else { | |||
1183 | REQUIRE(mp_int_copy(TEMP(0), TEMP(2)))do { temp_.err = (mp_int_copy((temp_.value + (0)), (temp_.value + (2)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1184 | } | |||
1185 | ||||
1186 | for (;;) { | |||
1187 | s_qdiv(TEMP(2)(temp_.value + (2)), s_dp2k(TEMP(2)(temp_.value + (2)))); | |||
1188 | ||||
1189 | if (CMPZ(TEMP(2)(temp_.value + (2))) > 0) { | |||
1190 | REQUIRE(mp_int_copy(TEMP(2), TEMP(0)))do { temp_.err = (mp_int_copy((temp_.value + (2)), (temp_.value + (0)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1191 | } else { | |||
1192 | REQUIRE(mp_int_neg(TEMP(2), TEMP(1)))do { temp_.err = (mp_int_neg((temp_.value + (2)), (temp_.value + (1)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1193 | } | |||
1194 | ||||
1195 | REQUIRE(mp_int_sub(TEMP(0), TEMP(1), TEMP(2)))do { temp_.err = (mp_int_sub((temp_.value + (0)), (temp_.value + (1)), (temp_.value + (2)))); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1196 | ||||
1197 | if (CMPZ(TEMP(2)(temp_.value + (2))) == 0) break; | |||
1198 | } | |||
1199 | ||||
1200 | REQUIRE(mp_int_abs(TEMP(0), c))do { temp_.err = (mp_int_abs((temp_.value + (0)), c)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1201 | if (!s_qmul(c, (mp_size)k)) REQUIRE(MP_MEMORY)do { temp_.err = (MP_MEMORY); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1202 | ||||
1203 | CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0); | |||
1204 | return MP_OK; | |||
1205 | } | |||
1206 | ||||
1207 | /* This is the binary GCD algorithm again, but this time we keep track of the | |||
1208 | elementary matrix operations as we go, so we can get values x and y | |||
1209 | satisfying c = ax + by. | |||
1210 | */ | |||
1211 | mp_result mp_int_egcd(mp_int a, mp_int b, mp_int c, mp_int x, mp_int y) { | |||
1212 | assert(a != NULL && b != NULL && c != NULL && (x != NULL || y != NULL))((void) sizeof ((a != ((void*)0) && b != ((void*)0) && c != ((void*)0) && (x != ((void*)0) || y != ((void*) 0))) ? 1 : 0), __extension__ ({ if (a != ((void*)0) && b != ((void*)0) && c != ((void*)0) && (x != ( (void*)0) || y != ((void*)0))) ; else __assert_fail ("a != NULL && b != NULL && c != NULL && (x != NULL || y != NULL)" , "polly/lib/External/isl/imath/imath.c", 1212, __extension__ __PRETTY_FUNCTION__); })); | |||
1213 | ||||
1214 | mp_result res = MP_OK; | |||
1215 | int ca = CMPZ(a); | |||
1216 | int cb = CMPZ(b); | |||
1217 | if (ca == 0 && cb == 0) { | |||
1218 | return MP_UNDEF; | |||
1219 | } else if (ca == 0) { | |||
1220 | if ((res = mp_int_abs(b, c)) != MP_OK) return res; | |||
1221 | mp_int_zero(x); | |||
1222 | (void)mp_int_set_value(y, 1); | |||
1223 | return MP_OK; | |||
1224 | } else if (cb == 0) { | |||
1225 | if ((res = mp_int_abs(a, c)) != MP_OK) return res; | |||
1226 | (void)mp_int_set_value(x, 1); | |||
1227 | mp_int_zero(y); | |||
1228 | return MP_OK; | |||
1229 | } | |||
1230 | ||||
1231 | /* Initialize temporaries: | |||
1232 | A:0, B:1, C:2, D:3, u:4, v:5, ou:6, ov:7 */ | |||
1233 | DECLARE_TEMP(8)struct { mpz_t value[(8)]; int len; mp_result err; } temp_ = { .len = (8), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0); | |||
1234 | REQUIRE(mp_int_set_value(TEMP(0), 1))do { temp_.err = (mp_int_set_value((temp_.value + (0)), 1)); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1235 | REQUIRE(mp_int_set_value(TEMP(3), 1))do { temp_.err = (mp_int_set_value((temp_.value + (3)), 1)); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1236 | REQUIRE(mp_int_copy(a, TEMP(4)))do { temp_.err = (mp_int_copy(a, (temp_.value + (4)))); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1237 | REQUIRE(mp_int_copy(b, TEMP(5)))do { temp_.err = (mp_int_copy(b, (temp_.value + (5)))); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1238 | ||||
1239 | /* We will work with absolute values here */ | |||
1240 | TEMP(4)(temp_.value + (4))->sign = MP_ZPOS; | |||
1241 | TEMP(5)(temp_.value + (5))->sign = MP_ZPOS; | |||
1242 | ||||
1243 | int k = 0; | |||
1244 | { /* Divide out common factors of 2 from u and v */ | |||
1245 | int div2_u = s_dp2k(TEMP(4)(temp_.value + (4))), div2_v = s_dp2k(TEMP(5)(temp_.value + (5))); | |||
1246 | ||||
1247 | k = MIN(div2_u, div2_v); | |||
1248 | s_qdiv(TEMP(4)(temp_.value + (4)), k); | |||
1249 | s_qdiv(TEMP(5)(temp_.value + (5)), k); | |||
1250 | } | |||
1251 | ||||
1252 | REQUIRE(mp_int_copy(TEMP(4), TEMP(6)))do { temp_.err = (mp_int_copy((temp_.value + (4)), (temp_.value + (6)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1253 | REQUIRE(mp_int_copy(TEMP(5), TEMP(7)))do { temp_.err = (mp_int_copy((temp_.value + (5)), (temp_.value + (7)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1254 | ||||
1255 | for (;;) { | |||
1256 | while (mp_int_is_even(TEMP(4)(temp_.value + (4)))) { | |||
1257 | s_qdiv(TEMP(4)(temp_.value + (4)), 1); | |||
1258 | ||||
1259 | if (mp_int_is_odd(TEMP(0)(temp_.value + (0))) || mp_int_is_odd(TEMP(1)(temp_.value + (1)))) { | |||
1260 | REQUIRE(mp_int_add(TEMP(0), TEMP(7), TEMP(0)))do { temp_.err = (mp_int_add((temp_.value + (0)), (temp_.value + (7)), (temp_.value + (0)))); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1261 | REQUIRE(mp_int_sub(TEMP(1), TEMP(6), TEMP(1)))do { temp_.err = (mp_int_sub((temp_.value + (1)), (temp_.value + (6)), (temp_.value + (1)))); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1262 | } | |||
1263 | ||||
1264 | s_qdiv(TEMP(0)(temp_.value + (0)), 1); | |||
1265 | s_qdiv(TEMP(1)(temp_.value + (1)), 1); | |||
1266 | } | |||
1267 | ||||
1268 | while (mp_int_is_even(TEMP(5)(temp_.value + (5)))) { | |||
1269 | s_qdiv(TEMP(5)(temp_.value + (5)), 1); | |||
1270 | ||||
1271 | if (mp_int_is_odd(TEMP(2)(temp_.value + (2))) || mp_int_is_odd(TEMP(3)(temp_.value + (3)))) { | |||
1272 | REQUIRE(mp_int_add(TEMP(2), TEMP(7), TEMP(2)))do { temp_.err = (mp_int_add((temp_.value + (2)), (temp_.value + (7)), (temp_.value + (2)))); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1273 | REQUIRE(mp_int_sub(TEMP(3), TEMP(6), TEMP(3)))do { temp_.err = (mp_int_sub((temp_.value + (3)), (temp_.value + (6)), (temp_.value + (3)))); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1274 | } | |||
1275 | ||||
1276 | s_qdiv(TEMP(2)(temp_.value + (2)), 1); | |||
1277 | s_qdiv(TEMP(3)(temp_.value + (3)), 1); | |||
1278 | } | |||
1279 | ||||
1280 | if (mp_int_compare(TEMP(4)(temp_.value + (4)), TEMP(5)(temp_.value + (5))) >= 0) { | |||
1281 | REQUIRE(mp_int_sub(TEMP(4), TEMP(5), TEMP(4)))do { temp_.err = (mp_int_sub((temp_.value + (4)), (temp_.value + (5)), (temp_.value + (4)))); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1282 | REQUIRE(mp_int_sub(TEMP(0), TEMP(2), TEMP(0)))do { temp_.err = (mp_int_sub((temp_.value + (0)), (temp_.value + (2)), (temp_.value + (0)))); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1283 | REQUIRE(mp_int_sub(TEMP(1), TEMP(3), TEMP(1)))do { temp_.err = (mp_int_sub((temp_.value + (1)), (temp_.value + (3)), (temp_.value + (1)))); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1284 | } else { | |||
1285 | REQUIRE(mp_int_sub(TEMP(5), TEMP(4), TEMP(5)))do { temp_.err = (mp_int_sub((temp_.value + (5)), (temp_.value + (4)), (temp_.value + (5)))); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1286 | REQUIRE(mp_int_sub(TEMP(2), TEMP(0), TEMP(2)))do { temp_.err = (mp_int_sub((temp_.value + (2)), (temp_.value + (0)), (temp_.value + (2)))); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1287 | REQUIRE(mp_int_sub(TEMP(3), TEMP(1), TEMP(3)))do { temp_.err = (mp_int_sub((temp_.value + (3)), (temp_.value + (1)), (temp_.value + (3)))); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1288 | } | |||
1289 | ||||
1290 | if (CMPZ(TEMP(4)(temp_.value + (4))) == 0) { | |||
1291 | if (x) REQUIRE(mp_int_copy(TEMP(2), x))do { temp_.err = (mp_int_copy((temp_.value + (2)), x)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1292 | if (y) REQUIRE(mp_int_copy(TEMP(3), y))do { temp_.err = (mp_int_copy((temp_.value + (3)), y)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1293 | if (c) { | |||
1294 | if (!s_qmul(TEMP(5)(temp_.value + (5)), k)) { | |||
1295 | REQUIRE(MP_MEMORY)do { temp_.err = (MP_MEMORY); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1296 | } | |||
1297 | REQUIRE(mp_int_copy(TEMP(5), c))do { temp_.err = (mp_int_copy((temp_.value + (5)), c)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1298 | } | |||
1299 | ||||
1300 | break; | |||
1301 | } | |||
1302 | } | |||
1303 | ||||
1304 | CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0); | |||
1305 | return MP_OK; | |||
1306 | } | |||
1307 | ||||
1308 | mp_result mp_int_lcm(mp_int a, mp_int b, mp_int c) { | |||
1309 | assert(a != NULL && b != NULL && c != NULL)((void) sizeof ((a != ((void*)0) && b != ((void*)0) && c != ((void*)0)) ? 1 : 0), __extension__ ({ if (a != ((void* )0) && b != ((void*)0) && c != ((void*)0)) ; else __assert_fail ("a != NULL && b != NULL && c != NULL" , "polly/lib/External/isl/imath/imath.c", 1309, __extension__ __PRETTY_FUNCTION__); })); | |||
1310 | ||||
1311 | /* Since a * b = gcd(a, b) * lcm(a, b), we can compute | |||
1312 | lcm(a, b) = (a / gcd(a, b)) * b. | |||
1313 | ||||
1314 | This formulation insures everything works even if the input | |||
1315 | variables share space. | |||
1316 | */ | |||
1317 | DECLARE_TEMP(1)struct { mpz_t value[(1)]; int len; mp_result err; } temp_ = { .len = (1), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0); | |||
1318 | REQUIRE(mp_int_gcd(a, b, TEMP(0)))do { temp_.err = (mp_int_gcd(a, b, (temp_.value + (0)))); if ( temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1319 | REQUIRE(mp_int_div(a, TEMP(0), TEMP(0), NULL))do { temp_.err = (mp_int_div(a, (temp_.value + (0)), (temp_.value + (0)), ((void*)0))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1320 | REQUIRE(mp_int_mul(TEMP(0), b, TEMP(0)))do { temp_.err = (mp_int_mul((temp_.value + (0)), b, (temp_.value + (0)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1321 | REQUIRE(mp_int_copy(TEMP(0), c))do { temp_.err = (mp_int_copy((temp_.value + (0)), c)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1322 | ||||
1323 | CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0); | |||
1324 | return MP_OK; | |||
1325 | } | |||
1326 | ||||
1327 | bool_Bool mp_int_divisible_value(mp_int a, mp_small v) { | |||
1328 | mp_small rem = 0; | |||
1329 | ||||
1330 | if (mp_int_div_value(a, v, NULL((void*)0), &rem) != MP_OK) { | |||
1331 | return false0; | |||
1332 | } | |||
1333 | return rem == 0; | |||
1334 | } | |||
1335 | ||||
1336 | int mp_int_is_pow2(mp_int z) { | |||
1337 | assert(z != NULL)((void) sizeof ((z != ((void*)0)) ? 1 : 0), __extension__ ({ if (z != ((void*)0)) ; else __assert_fail ("z != NULL", "polly/lib/External/isl/imath/imath.c" , 1337, __extension__ __PRETTY_FUNCTION__); })); | |||
1338 | ||||
1339 | return s_isp2(z); | |||
1340 | } | |||
1341 | ||||
1342 | /* Implementation of Newton's root finding method, based loosely on a patch | |||
1343 | contributed by Hal Finkel <half@halssoftware.com> | |||
1344 | modified by M. J. Fromberger. | |||
1345 | */ | |||
1346 | mp_result mp_int_root(mp_int a, mp_small b, mp_int c) { | |||
1347 | assert(a != NULL && c != NULL && b > 0)((void) sizeof ((a != ((void*)0) && c != ((void*)0) && b > 0) ? 1 : 0), __extension__ ({ if (a != ((void*)0) && c != ((void*)0) && b > 0) ; else __assert_fail ("a != NULL && c != NULL && b > 0" , "polly/lib/External/isl/imath/imath.c", 1347, __extension__ __PRETTY_FUNCTION__); })); | |||
1348 | ||||
1349 | if (b == 1) { | |||
1350 | return mp_int_copy(a, c); | |||
1351 | } | |||
1352 | bool_Bool flips = false0; | |||
1353 | if (MP_SIGN(a) == MP_NEG) { | |||
1354 | if (b % 2 == 0) { | |||
1355 | return MP_UNDEF; /* root does not exist for negative a with even b */ | |||
1356 | } else { | |||
1357 | flips = true1; | |||
1358 | } | |||
1359 | } | |||
1360 | ||||
1361 | DECLARE_TEMP(5)struct { mpz_t value[(5)]; int len; mp_result err; } temp_ = { .len = (5), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0); | |||
1362 | REQUIRE(mp_int_copy(a, TEMP(0)))do { temp_.err = (mp_int_copy(a, (temp_.value + (0)))); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1363 | REQUIRE(mp_int_copy(a, TEMP(1)))do { temp_.err = (mp_int_copy(a, (temp_.value + (1)))); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1364 | TEMP(0)(temp_.value + (0))->sign = MP_ZPOS; | |||
1365 | TEMP(1)(temp_.value + (1))->sign = MP_ZPOS; | |||
1366 | ||||
1367 | for (;;) { | |||
1368 | REQUIRE(mp_int_expt(TEMP(1), b, TEMP(2)))do { temp_.err = (mp_int_expt((temp_.value + (1)), b, (temp_. value + (2)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1369 | ||||
1370 | if (mp_int_compare_unsigned(TEMP(2)(temp_.value + (2)), TEMP(0)(temp_.value + (0))) <= 0) break; | |||
1371 | ||||
1372 | REQUIRE(mp_int_sub(TEMP(2), TEMP(0), TEMP(2)))do { temp_.err = (mp_int_sub((temp_.value + (2)), (temp_.value + (0)), (temp_.value + (2)))); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1373 | REQUIRE(mp_int_expt(TEMP(1), b - 1, TEMP(3)))do { temp_.err = (mp_int_expt((temp_.value + (1)), b - 1, (temp_ .value + (3)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1374 | REQUIRE(mp_int_mul_value(TEMP(3), b, TEMP(3)))do { temp_.err = (mp_int_mul_value((temp_.value + (3)), b, (temp_ .value + (3)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1375 | REQUIRE(mp_int_div(TEMP(2), TEMP(3), TEMP(4), NULL))do { temp_.err = (mp_int_div((temp_.value + (2)), (temp_.value + (3)), (temp_.value + (4)), ((void*)0))); if (temp_.err != MP_OK ) goto CLEANUP; } while (0); | |||
1376 | REQUIRE(mp_int_sub(TEMP(1), TEMP(4), TEMP(4)))do { temp_.err = (mp_int_sub((temp_.value + (1)), (temp_.value + (4)), (temp_.value + (4)))); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
1377 | ||||
1378 | if (mp_int_compare_unsigned(TEMP(1)(temp_.value + (1)), TEMP(4)(temp_.value + (4))) == 0) { | |||
1379 | REQUIRE(mp_int_sub_value(TEMP(4), 1, TEMP(4)))do { temp_.err = (mp_int_sub_value((temp_.value + (4)), 1, (temp_ .value + (4)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1380 | } | |||
1381 | REQUIRE(mp_int_copy(TEMP(4), TEMP(1)))do { temp_.err = (mp_int_copy((temp_.value + (4)), (temp_.value + (1)))); if (temp_.err != MP_OK) goto CLEANUP; } while (0); | |||
1382 | } | |||
1383 | ||||
1384 | REQUIRE(mp_int_copy(TEMP(1), c))do { temp_.err = (mp_int_copy((temp_.value + (1)), c)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
1385 | ||||
1386 | /* If the original value of a was negative, flip the output sign. */ | |||
1387 | if (flips) (void)mp_int_neg(c, c); /* cannot fail */ | |||
1388 | ||||
1389 | CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0); | |||
1390 | return MP_OK; | |||
1391 | } | |||
1392 | ||||
1393 | mp_result mp_int_to_int(mp_int z, mp_small *out) { | |||
1394 | assert(z != NULL)((void) sizeof ((z != ((void*)0)) ? 1 : 0), __extension__ ({ if (z != ((void*)0)) ; else __assert_fail ("z != NULL", "polly/lib/External/isl/imath/imath.c" , 1394, __extension__ __PRETTY_FUNCTION__); })); | |||
1395 | ||||
1396 | /* Make sure the value is representable as a small integer */ | |||
1397 | mp_sign sz = MP_SIGN(z); | |||
1398 | if ((sz == MP_ZPOS && mp_int_compare_value(z, MP_SMALL_MAX9223372036854775807L) > 0) || | |||
1399 | mp_int_compare_value(z, MP_SMALL_MIN(-9223372036854775807L -1L)) < 0) { | |||
1400 | return MP_RANGE; | |||
1401 | } | |||
1402 | ||||
1403 | mp_usmall uz = MP_USED(z); | |||
1404 | mp_digit *dz = MP_DIGITS(z) + uz - 1; | |||
1405 | mp_small uv = 0; | |||
1406 | while (uz > 0) { | |||
1407 | uv <<= MP_DIGIT_BIT(sizeof(mp_digit) * 8) / 2; | |||
1408 | uv = (uv << (MP_DIGIT_BIT(sizeof(mp_digit) * 8) / 2)) | *dz--; | |||
1409 | --uz; | |||
1410 | } | |||
1411 | ||||
1412 | if (out) *out = (mp_small)((sz == MP_NEG) ? -uv : uv); | |||
1413 | ||||
1414 | return MP_OK; | |||
1415 | } | |||
1416 | ||||
1417 | mp_result mp_int_to_uint(mp_int z, mp_usmall *out) { | |||
1418 | assert(z != NULL)((void) sizeof ((z != ((void*)0)) ? 1 : 0), __extension__ ({ if (z != ((void*)0)) ; else __assert_fail ("z != NULL", "polly/lib/External/isl/imath/imath.c" , 1418, __extension__ __PRETTY_FUNCTION__); })); | |||
1419 | ||||
1420 | /* Make sure the value is representable as an unsigned small integer */ | |||
1421 | mp_size sz = MP_SIGN(z); | |||
1422 | if (sz == MP_NEG || mp_int_compare_uvalue(z, MP_USMALL_MAX(9223372036854775807L *2UL+1UL)) > 0) { | |||
1423 | return MP_RANGE; | |||
1424 | } | |||
1425 | ||||
1426 | mp_size uz = MP_USED(z); | |||
1427 | mp_digit *dz = MP_DIGITS(z) + uz - 1; | |||
1428 | mp_usmall uv = 0; | |||
1429 | ||||
1430 | while (uz > 0) { | |||
1431 | uv <<= MP_DIGIT_BIT(sizeof(mp_digit) * 8) / 2; | |||
1432 | uv = (uv << (MP_DIGIT_BIT(sizeof(mp_digit) * 8) / 2)) | *dz--; | |||
1433 | --uz; | |||
1434 | } | |||
1435 | ||||
1436 | if (out) *out = uv; | |||
1437 | ||||
1438 | return MP_OK; | |||
1439 | } | |||
1440 | ||||
1441 | mp_result mp_int_to_string(mp_int z, mp_size radix, char *str, int limit) { | |||
1442 | assert(z != NULL && str != NULL && limit >= 2)((void) sizeof ((z != ((void*)0) && str != ((void*)0) && limit >= 2) ? 1 : 0), __extension__ ({ if (z != ((void*)0) && str != ((void*)0) && limit >= 2) ; else __assert_fail ("z != NULL && str != NULL && limit >= 2" , "polly/lib/External/isl/imath/imath.c", 1442, __extension__ __PRETTY_FUNCTION__); })); | |||
1443 | assert(radix >= MP_MIN_RADIX && radix <= MP_MAX_RADIX)((void) sizeof ((radix >= 2 && radix <= 36) ? 1 : 0), __extension__ ({ if (radix >= 2 && radix <= 36) ; else __assert_fail ("radix >= MP_MIN_RADIX && radix <= MP_MAX_RADIX" , "polly/lib/External/isl/imath/imath.c", 1443, __extension__ __PRETTY_FUNCTION__); })); | |||
1444 | ||||
1445 | int cmp = 0; | |||
1446 | if (CMPZ(z) == 0) { | |||
1447 | *str++ = s_val2ch(0, 1); | |||
1448 | } else { | |||
1449 | mp_result res; | |||
1450 | mpz_t tmp; | |||
1451 | char *h, *t; | |||
1452 | ||||
1453 | if ((res = mp_int_init_copy(&tmp, z)) != MP_OK) return res; | |||
1454 | ||||
1455 | if (MP_SIGN(z) == MP_NEG) { | |||
1456 | *str++ = '-'; | |||
1457 | --limit; | |||
1458 | } | |||
1459 | h = str; | |||
1460 | ||||
1461 | /* Generate digits in reverse order until finished or limit reached */ | |||
1462 | for (/* */; limit > 0; --limit) { | |||
1463 | mp_digit d; | |||
1464 | ||||
1465 | if ((cmp = CMPZ(&tmp)) == 0) break; | |||
1466 | ||||
1467 | d = s_ddiv(&tmp, (mp_digit)radix); | |||
1468 | *str++ = s_val2ch(d, 1); | |||
1469 | } | |||
1470 | t = str - 1; | |||
1471 | ||||
1472 | /* Put digits back in correct output order */ | |||
1473 | while (h < t) { | |||
1474 | char tc = *h; | |||
1475 | *h++ = *t; | |||
1476 | *t-- = tc; | |||
1477 | } | |||
1478 | ||||
1479 | mp_int_clear(&tmp); | |||
1480 | } | |||
1481 | ||||
1482 | *str = '\0'; | |||
1483 | if (cmp == 0) { | |||
1484 | return MP_OK; | |||
1485 | } else { | |||
1486 | return MP_TRUNC; | |||
1487 | } | |||
1488 | } | |||
1489 | ||||
1490 | mp_result mp_int_string_len(mp_int z, mp_size radix) { | |||
1491 | assert(z != NULL)((void) sizeof ((z != ((void*)0)) ? 1 : 0), __extension__ ({ if (z != ((void*)0)) ; else __assert_fail ("z != NULL", "polly/lib/External/isl/imath/imath.c" , 1491, __extension__ __PRETTY_FUNCTION__); })); | |||
1492 | assert(radix >= MP_MIN_RADIX && radix <= MP_MAX_RADIX)((void) sizeof ((radix >= 2 && radix <= 36) ? 1 : 0), __extension__ ({ if (radix >= 2 && radix <= 36) ; else __assert_fail ("radix >= MP_MIN_RADIX && radix <= MP_MAX_RADIX" , "polly/lib/External/isl/imath/imath.c", 1492, __extension__ __PRETTY_FUNCTION__); })); | |||
1493 | ||||
1494 | int len = s_outlen(z, radix) + 1; /* for terminator */ | |||
1495 | ||||
1496 | /* Allow for sign marker on negatives */ | |||
1497 | if (MP_SIGN(z) == MP_NEG) len += 1; | |||
1498 | ||||
1499 | return len; | |||
1500 | } | |||
1501 | ||||
1502 | /* Read zero-terminated string into z */ | |||
1503 | mp_result mp_int_read_string(mp_int z, mp_size radix, const char *str) { | |||
1504 | return mp_int_read_cstring(z, radix, str, NULL((void*)0)); | |||
1505 | } | |||
1506 | ||||
1507 | mp_result mp_int_read_cstring(mp_int z, mp_size radix, const char *str, | |||
1508 | char **end) { | |||
1509 | assert(z != NULL && str != NULL)((void) sizeof ((z != ((void*)0) && str != ((void*)0) ) ? 1 : 0), __extension__ ({ if (z != ((void*)0) && str != ((void*)0)) ; else __assert_fail ("z != NULL && str != NULL" , "polly/lib/External/isl/imath/imath.c", 1509, __extension__ __PRETTY_FUNCTION__); })); | |||
1510 | assert(radix >= MP_MIN_RADIX && radix <= MP_MAX_RADIX)((void) sizeof ((radix >= 2 && radix <= 36) ? 1 : 0), __extension__ ({ if (radix >= 2 && radix <= 36) ; else __assert_fail ("radix >= MP_MIN_RADIX && radix <= MP_MAX_RADIX" , "polly/lib/External/isl/imath/imath.c", 1510, __extension__ __PRETTY_FUNCTION__); })); | |||
1511 | ||||
1512 | /* Skip leading whitespace */ | |||
1513 | while (isspace((unsigned char)*str)((*__ctype_b_loc ())[(int) (((unsigned char)*str))] & (unsigned short int) _ISspace)) ++str; | |||
1514 | ||||
1515 | /* Handle leading sign tag (+/-, positive default) */ | |||
1516 | switch (*str) { | |||
1517 | case '-': | |||
1518 | z->sign = MP_NEG; | |||
1519 | ++str; | |||
1520 | break; | |||
1521 | case '+': | |||
1522 | ++str; /* fallthrough */ | |||
1523 | default: | |||
1524 | z->sign = MP_ZPOS; | |||
1525 | break; | |||
1526 | } | |||
1527 | ||||
1528 | /* Skip leading zeroes */ | |||
1529 | int ch; | |||
1530 | while ((ch = s_ch2val(*str, radix)) == 0) ++str; | |||
1531 | ||||
1532 | /* Make sure there is enough space for the value */ | |||
1533 | if (!s_pad(z, s_inlen(strlen(str), radix))) return MP_MEMORY; | |||
1534 | ||||
1535 | z->used = 1; | |||
1536 | z->digits[0] = 0; | |||
1537 | ||||
1538 | while (*str != '\0' && ((ch = s_ch2val(*str, radix)) >= 0)) { | |||
1539 | s_dmul(z, (mp_digit)radix); | |||
1540 | s_dadd(z, (mp_digit)ch); | |||
1541 | ++str; | |||
1542 | } | |||
1543 | ||||
1544 | CLAMP(z); | |||
1545 | ||||
1546 | /* Override sign for zero, even if negative specified. */ | |||
1547 | if (CMPZ(z) == 0) z->sign = MP_ZPOS; | |||
1548 | ||||
1549 | if (end != NULL((void*)0)) *end = (char *)str; | |||
1550 | ||||
1551 | /* Return a truncation error if the string has unprocessed characters | |||
1552 | remaining, so the caller can tell if the whole string was done */ | |||
1553 | if (*str != '\0') { | |||
1554 | return MP_TRUNC; | |||
1555 | } else { | |||
1556 | return MP_OK; | |||
1557 | } | |||
1558 | } | |||
1559 | ||||
1560 | mp_result mp_int_count_bits(mp_int z) { | |||
1561 | assert(z != NULL)((void) sizeof ((z != ((void*)0)) ? 1 : 0), __extension__ ({ if (z != ((void*)0)) ; else __assert_fail ("z != NULL", "polly/lib/External/isl/imath/imath.c" , 1561, __extension__ __PRETTY_FUNCTION__); })); | |||
1562 | ||||
1563 | mp_size uz = MP_USED(z); | |||
1564 | if (uz == 1 && z->digits[0] == 0) return 1; | |||
1565 | ||||
1566 | --uz; | |||
1567 | mp_size nbits = uz * MP_DIGIT_BIT(sizeof(mp_digit) * 8); | |||
1568 | mp_digit d = z->digits[uz]; | |||
1569 | ||||
1570 | while (d != 0) { | |||
1571 | d >>= 1; | |||
1572 | ++nbits; | |||
1573 | } | |||
1574 | ||||
1575 | return nbits; | |||
1576 | } | |||
1577 | ||||
1578 | mp_result mp_int_to_binary(mp_int z, unsigned char *buf, int limit) { | |||
1579 | static const int PAD_FOR_2C = 1; | |||
1580 | ||||
1581 | assert(z != NULL && buf != NULL)((void) sizeof ((z != ((void*)0) && buf != ((void*)0) ) ? 1 : 0), __extension__ ({ if (z != ((void*)0) && buf != ((void*)0)) ; else __assert_fail ("z != NULL && buf != NULL" , "polly/lib/External/isl/imath/imath.c", 1581, __extension__ __PRETTY_FUNCTION__); })); | |||
1582 | ||||
1583 | int limpos = limit; | |||
1584 | mp_result res = s_tobin(z, buf, &limpos, PAD_FOR_2C); | |||
1585 | ||||
1586 | if (MP_SIGN(z) == MP_NEG) s_2comp(buf, limpos); | |||
1587 | ||||
1588 | return res; | |||
1589 | } | |||
1590 | ||||
1591 | mp_result mp_int_read_binary(mp_int z, unsigned char *buf, int len) { | |||
1592 | assert(z != NULL && buf != NULL && len > 0)((void) sizeof ((z != ((void*)0) && buf != ((void*)0) && len > 0) ? 1 : 0), __extension__ ({ if (z != ( (void*)0) && buf != ((void*)0) && len > 0) ; else __assert_fail ("z != NULL && buf != NULL && len > 0" , "polly/lib/External/isl/imath/imath.c", 1592, __extension__ __PRETTY_FUNCTION__); })); | |||
1593 | ||||
1594 | /* Figure out how many digits are needed to represent this value */ | |||
1595 | mp_size need = ((len * CHAR_BIT8) + (MP_DIGIT_BIT(sizeof(mp_digit) * 8) - 1)) / MP_DIGIT_BIT(sizeof(mp_digit) * 8); | |||
1596 | if (!s_pad(z, need)) return MP_MEMORY; | |||
1597 | ||||
1598 | mp_int_zero(z); | |||
1599 | ||||
1600 | /* If the high-order bit is set, take the 2's complement before reading the | |||
1601 | value (it will be restored afterward) */ | |||
1602 | if (buf[0] >> (CHAR_BIT8 - 1)) { | |||
1603 | z->sign = MP_NEG; | |||
1604 | s_2comp(buf, len); | |||
1605 | } | |||
1606 | ||||
1607 | mp_digit *dz = MP_DIGITS(z); | |||
1608 | unsigned char *tmp = buf; | |||
1609 | for (int i = len; i > 0; --i, ++tmp) { | |||
1610 | s_qmul(z, (mp_size)CHAR_BIT8); | |||
1611 | *dz |= *tmp; | |||
1612 | } | |||
1613 | ||||
1614 | /* Restore 2's complement if we took it before */ | |||
1615 | if (MP_SIGN(z) == MP_NEG) s_2comp(buf, len); | |||
1616 | ||||
1617 | return MP_OK; | |||
1618 | } | |||
1619 | ||||
1620 | mp_result mp_int_binary_len(mp_int z) { | |||
1621 | mp_result res = mp_int_count_bits(z); | |||
1622 | if (res <= 0) return res; | |||
1623 | ||||
1624 | int bytes = mp_int_unsigned_len(z); | |||
1625 | ||||
1626 | /* If the highest-order bit falls exactly on a byte boundary, we need to pad | |||
1627 | with an extra byte so that the sign will be read correctly when reading it | |||
1628 | back in. */ | |||
1629 | if (bytes * CHAR_BIT8 == res) ++bytes; | |||
1630 | ||||
1631 | return bytes; | |||
1632 | } | |||
1633 | ||||
1634 | mp_result mp_int_to_unsigned(mp_int z, unsigned char *buf, int limit) { | |||
1635 | static const int NO_PADDING = 0; | |||
1636 | ||||
1637 | assert(z != NULL && buf != NULL)((void) sizeof ((z != ((void*)0) && buf != ((void*)0) ) ? 1 : 0), __extension__ ({ if (z != ((void*)0) && buf != ((void*)0)) ; else __assert_fail ("z != NULL && buf != NULL" , "polly/lib/External/isl/imath/imath.c", 1637, __extension__ __PRETTY_FUNCTION__); })); | |||
1638 | ||||
1639 | return s_tobin(z, buf, &limit, NO_PADDING); | |||
1640 | } | |||
1641 | ||||
1642 | mp_result mp_int_read_unsigned(mp_int z, unsigned char *buf, int len) { | |||
1643 | assert(z != NULL && buf != NULL && len > 0)((void) sizeof ((z != ((void*)0) && buf != ((void*)0) && len > 0) ? 1 : 0), __extension__ ({ if (z != ( (void*)0) && buf != ((void*)0) && len > 0) ; else __assert_fail ("z != NULL && buf != NULL && len > 0" , "polly/lib/External/isl/imath/imath.c", 1643, __extension__ __PRETTY_FUNCTION__); })); | |||
1644 | ||||
1645 | /* Figure out how many digits are needed to represent this value */ | |||
1646 | mp_size need = ((len * CHAR_BIT8) + (MP_DIGIT_BIT(sizeof(mp_digit) * 8) - 1)) / MP_DIGIT_BIT(sizeof(mp_digit) * 8); | |||
1647 | if (!s_pad(z, need)) return MP_MEMORY; | |||
1648 | ||||
1649 | mp_int_zero(z); | |||
1650 | ||||
1651 | unsigned char *tmp = buf; | |||
1652 | for (int i = len; i > 0; --i, ++tmp) { | |||
1653 | (void)s_qmul(z, CHAR_BIT8); | |||
1654 | *MP_DIGITS(z) |= *tmp; | |||
1655 | } | |||
1656 | ||||
1657 | return MP_OK; | |||
1658 | } | |||
1659 | ||||
1660 | mp_result mp_int_unsigned_len(mp_int z) { | |||
1661 | mp_result res = mp_int_count_bits(z); | |||
1662 | if (res <= 0) return res; | |||
1663 | ||||
1664 | int bytes = (res + (CHAR_BIT8 - 1)) / CHAR_BIT8; | |||
1665 | return bytes; | |||
1666 | } | |||
1667 | ||||
1668 | const char *mp_error_string(mp_result res) { | |||
1669 | if (res > 0) return s_unknown_err; | |||
1670 | ||||
1671 | res = -res; | |||
1672 | int ix; | |||
1673 | for (ix = 0; ix < res && s_error_msg[ix] != NULL((void*)0); ++ix) | |||
1674 | ; | |||
1675 | ||||
1676 | if (s_error_msg[ix] != NULL((void*)0)) { | |||
1677 | return s_error_msg[ix]; | |||
1678 | } else { | |||
1679 | return s_unknown_err; | |||
1680 | } | |||
1681 | } | |||
1682 | ||||
1683 | /*------------------------------------------------------------------------*/ | |||
1684 | /* Private functions for internal use. These make assumptions. */ | |||
1685 | ||||
1686 | #if DEBUG | |||
1687 | static const mp_digit fill = (mp_digit)0xdeadbeefabad1dea; | |||
1688 | #endif | |||
1689 | ||||
1690 | static mp_digit *s_alloc(mp_size num) { | |||
1691 | mp_digit *out = malloc(num * sizeof(mp_digit)); | |||
1692 | assert(out != NULL)((void) sizeof ((out != ((void*)0)) ? 1 : 0), __extension__ ( { if (out != ((void*)0)) ; else __assert_fail ("out != NULL", "polly/lib/External/isl/imath/imath.c", 1692, __extension__ __PRETTY_FUNCTION__ ); })); | |||
1693 | ||||
1694 | #if DEBUG | |||
1695 | for (mp_size ix = 0; ix < num; ++ix) out[ix] = fill; | |||
1696 | #endif | |||
1697 | return out; | |||
1698 | } | |||
1699 | ||||
1700 | static mp_digit *s_realloc(mp_digit *old, mp_size osize, mp_size nsize) { | |||
1701 | #if DEBUG | |||
1702 | mp_digit *new = s_alloc(nsize); | |||
1703 | assert(new != NULL)((void) sizeof ((new != ((void*)0)) ? 1 : 0), __extension__ ( { if (new != ((void*)0)) ; else __assert_fail ("new != NULL", "polly/lib/External/isl/imath/imath.c", 1703, __extension__ __PRETTY_FUNCTION__ ); })); | |||
1704 | ||||
1705 | for (mp_size ix = 0; ix < nsize; ++ix) new[ix] = fill; | |||
1706 | memcpy(new, old, osize * sizeof(mp_digit)); | |||
1707 | #else | |||
1708 | mp_digit *new = realloc(old, nsize * sizeof(mp_digit)); | |||
1709 | assert(new != NULL)((void) sizeof ((new != ((void*)0)) ? 1 : 0), __extension__ ( { if (new != ((void*)0)) ; else __assert_fail ("new != NULL", "polly/lib/External/isl/imath/imath.c", 1709, __extension__ __PRETTY_FUNCTION__ ); })); | |||
1710 | #endif | |||
1711 | ||||
1712 | return new; | |||
1713 | } | |||
1714 | ||||
1715 | static void s_free(void *ptr) { free(ptr); } | |||
1716 | ||||
1717 | static bool_Bool s_pad(mp_int z, mp_size min) { | |||
1718 | if (MP_ALLOC(z) < min) { | |||
1719 | mp_size nsize = s_round_prec(min); | |||
1720 | mp_digit *tmp; | |||
1721 | ||||
1722 | if (z->digits == &(z->single)) { | |||
1723 | if ((tmp = s_alloc(nsize)) == NULL((void*)0)) return false0; | |||
1724 | tmp[0] = z->single; | |||
1725 | } else if ((tmp = s_realloc(MP_DIGITS(z), MP_ALLOC(z), nsize)) == NULL((void*)0)) { | |||
1726 | return false0; | |||
1727 | } | |||
1728 | ||||
1729 | z->digits = tmp; | |||
1730 | z->alloc = nsize; | |||
1731 | } | |||
1732 | ||||
1733 | return true1; | |||
1734 | } | |||
1735 | ||||
1736 | /* Note: This will not work correctly when value == MP_SMALL_MIN */ | |||
1737 | static void s_fake(mp_int z, mp_small value, mp_digit vbuf[]) { | |||
1738 | mp_usmall uv = (mp_usmall)(value < 0) ? -value : value; | |||
1739 | s_ufake(z, uv, vbuf); | |||
1740 | if (value < 0) z->sign = MP_NEG; | |||
1741 | } | |||
1742 | ||||
1743 | static void s_ufake(mp_int z, mp_usmall value, mp_digit vbuf[]) { | |||
1744 | mp_size ndig = (mp_size)s_uvpack(value, vbuf); | |||
1745 | ||||
1746 | z->used = ndig; | |||
1747 | z->alloc = MP_VALUE_DIGITS(value)((sizeof(value) + (sizeof(mp_digit) - 1)) / sizeof(mp_digit)); | |||
1748 | z->sign = MP_ZPOS; | |||
1749 | z->digits = vbuf; | |||
1750 | } | |||
1751 | ||||
1752 | static int s_cdig(mp_digit *da, mp_digit *db, mp_size len) { | |||
1753 | mp_digit *dat = da + len - 1, *dbt = db + len - 1; | |||
1754 | ||||
1755 | for (/* */; len != 0; --len, --dat, --dbt) { | |||
1756 | if (*dat > *dbt) { | |||
1757 | return 1; | |||
1758 | } else if (*dat < *dbt) { | |||
1759 | return -1; | |||
1760 | } | |||
1761 | } | |||
1762 | ||||
1763 | return 0; | |||
1764 | } | |||
1765 | ||||
1766 | static int s_uvpack(mp_usmall uv, mp_digit t[]) { | |||
1767 | int ndig = 0; | |||
1768 | ||||
1769 | if (uv == 0) | |||
1770 | t[ndig++] = 0; | |||
1771 | else { | |||
1772 | while (uv != 0) { | |||
1773 | t[ndig++] = (mp_digit)uv; | |||
1774 | uv >>= MP_DIGIT_BIT(sizeof(mp_digit) * 8) / 2; | |||
1775 | uv >>= MP_DIGIT_BIT(sizeof(mp_digit) * 8) / 2; | |||
1776 | } | |||
1777 | } | |||
1778 | ||||
1779 | return ndig; | |||
1780 | } | |||
1781 | ||||
1782 | static int s_ucmp(mp_int a, mp_int b) { | |||
1783 | mp_size ua = MP_USED(a), ub = MP_USED(b); | |||
1784 | ||||
1785 | if (ua > ub) { | |||
1786 | return 1; | |||
1787 | } else if (ub > ua) { | |||
1788 | return -1; | |||
1789 | } else { | |||
1790 | return s_cdig(MP_DIGITS(a), MP_DIGITS(b), ua); | |||
1791 | } | |||
1792 | } | |||
1793 | ||||
1794 | static int s_vcmp(mp_int a, mp_small v) { | |||
1795 | mp_usmall uv = (v < 0) ? -(mp_usmall)v : (mp_usmall)v; | |||
1796 | return s_uvcmp(a, uv); | |||
1797 | } | |||
1798 | ||||
1799 | static int s_uvcmp(mp_int a, mp_usmall uv) { | |||
1800 | mpz_t vtmp; | |||
1801 | mp_digit vdig[MP_VALUE_DIGITS(uv)((sizeof(uv) + (sizeof(mp_digit) - 1)) / sizeof(mp_digit))]; | |||
1802 | ||||
1803 | s_ufake(&vtmp, uv, vdig); | |||
1804 | return s_ucmp(a, &vtmp); | |||
1805 | } | |||
1806 | ||||
1807 | static mp_digit s_uadd(mp_digit *da, mp_digit *db, mp_digit *dc, mp_size size_a, | |||
1808 | mp_size size_b) { | |||
1809 | mp_size pos; | |||
1810 | mp_word w = 0; | |||
1811 | ||||
1812 | /* Insure that da is the longer of the two to simplify later code */ | |||
1813 | if (size_b > size_a) { | |||
1814 | SWAP(mp_digit *, da, db)do { mp_digit * t_ = (da); da = (db); db = t_; } while (0); | |||
1815 | SWAP(mp_size, size_a, size_b)do { mp_size t_ = (size_a); size_a = (size_b); size_b = t_; } while (0); | |||
1816 | } | |||
1817 | ||||
1818 | /* Add corresponding digits until the shorter number runs out */ | |||
1819 | for (pos = 0; pos < size_b; ++pos, ++da, ++db, ++dc) { | |||
1820 | w = w + (mp_word)*da + (mp_word)*db; | |||
1821 | *dc = LOWER_HALF(w); | |||
1822 | w = UPPER_HALF(w); | |||
1823 | } | |||
1824 | ||||
1825 | /* Propagate carries as far as necessary */ | |||
1826 | for (/* */; pos < size_a; ++pos, ++da, ++dc) { | |||
1827 | w = w + *da; | |||
1828 | ||||
1829 | *dc = LOWER_HALF(w); | |||
1830 | w = UPPER_HALF(w); | |||
1831 | } | |||
1832 | ||||
1833 | /* Return carry out */ | |||
1834 | return (mp_digit)w; | |||
1835 | } | |||
1836 | ||||
1837 | static void s_usub(mp_digit *da, mp_digit *db, mp_digit *dc, mp_size size_a, | |||
1838 | mp_size size_b) { | |||
1839 | mp_size pos; | |||
1840 | mp_word w = 0; | |||
1841 | ||||
1842 | /* We assume that |a| >= |b| so this should definitely hold */ | |||
1843 | assert(size_a >= size_b)((void) sizeof ((size_a >= size_b) ? 1 : 0), __extension__ ({ if (size_a >= size_b) ; else __assert_fail ("size_a >= size_b" , "polly/lib/External/isl/imath/imath.c", 1843, __extension__ __PRETTY_FUNCTION__); })); | |||
1844 | ||||
1845 | /* Subtract corresponding digits and propagate borrow */ | |||
1846 | for (pos = 0; pos < size_b; ++pos, ++da, ++db, ++dc) { | |||
1847 | w = ((mp_word)MP_DIGIT_MAX((4294967295U) * 1UL) + 1 + /* MP_RADIX */ | |||
1848 | (mp_word)*da) - | |||
1849 | w - (mp_word)*db; | |||
1850 | ||||
1851 | *dc = LOWER_HALF(w); | |||
1852 | w = (UPPER_HALF(w) == 0); | |||
1853 | } | |||
1854 | ||||
1855 | /* Finish the subtraction for remaining upper digits of da */ | |||
1856 | for (/* */; pos < size_a; ++pos, ++da, ++dc) { | |||
1857 | w = ((mp_word)MP_DIGIT_MAX((4294967295U) * 1UL) + 1 + /* MP_RADIX */ | |||
1858 | (mp_word)*da) - | |||
1859 | w; | |||
1860 | ||||
1861 | *dc = LOWER_HALF(w); | |||
1862 | w = (UPPER_HALF(w) == 0); | |||
1863 | } | |||
1864 | ||||
1865 | /* If there is a borrow out at the end, it violates the precondition */ | |||
1866 | assert(w == 0)((void) sizeof ((w == 0) ? 1 : 0), __extension__ ({ if (w == 0 ) ; else __assert_fail ("w == 0", "polly/lib/External/isl/imath/imath.c" , 1866, __extension__ __PRETTY_FUNCTION__); })); | |||
1867 | } | |||
1868 | ||||
1869 | static int s_kmul(mp_digit *da, mp_digit *db, mp_digit *dc, mp_size size_a, | |||
1870 | mp_size size_b) { | |||
1871 | mp_size bot_size; | |||
1872 | ||||
1873 | /* Make sure b is the smaller of the two input values */ | |||
1874 | if (size_b > size_a) { | |||
1875 | SWAP(mp_digit *, da, db)do { mp_digit * t_ = (da); da = (db); db = t_; } while (0); | |||
1876 | SWAP(mp_size, size_a, size_b)do { mp_size t_ = (size_a); size_a = (size_b); size_b = t_; } while (0); | |||
1877 | } | |||
1878 | ||||
1879 | /* Insure that the bottom is the larger half in an odd-length split; the code | |||
1880 | below relies on this being true. | |||
1881 | */ | |||
1882 | bot_size = (size_a + 1) / 2; | |||
1883 | ||||
1884 | /* If the values are big enough to bother with recursion, use the Karatsuba | |||
1885 | algorithm to compute the product; otherwise use the normal multiplication | |||
1886 | algorithm | |||
1887 | */ | |||
1888 | if (multiply_threshold && size_a >= multiply_threshold && size_b > bot_size) { | |||
1889 | mp_digit *t1, *t2, *t3, carry; | |||
1890 | ||||
1891 | mp_digit *a_top = da + bot_size; | |||
1892 | mp_digit *b_top = db + bot_size; | |||
1893 | ||||
1894 | mp_size at_size = size_a - bot_size; | |||
1895 | mp_size bt_size = size_b - bot_size; | |||
1896 | mp_size buf_size = 2 * bot_size; | |||
1897 | ||||
1898 | /* Do a single allocation for all three temporary buffers needed; each | |||
1899 | buffer must be big enough to hold the product of two bottom halves, and | |||
1900 | one buffer needs space for the completed product; twice the space is | |||
1901 | plenty. | |||
1902 | */ | |||
1903 | if ((t1 = s_alloc(4 * buf_size)) == NULL((void*)0)) return 0; | |||
1904 | t2 = t1 + buf_size; | |||
1905 | t3 = t2 + buf_size; | |||
1906 | ZERO(t1, 4 * buf_size); | |||
1907 | ||||
1908 | /* t1 and t2 are initially used as temporaries to compute the inner product | |||
1909 | (a1 + a0)(b1 + b0) = a1b1 + a1b0 + a0b1 + a0b0 | |||
1910 | */ | |||
1911 | carry = s_uadd(da, a_top, t1, bot_size, at_size); /* t1 = a1 + a0 */ | |||
1912 | t1[bot_size] = carry; | |||
1913 | ||||
1914 | carry = s_uadd(db, b_top, t2, bot_size, bt_size); /* t2 = b1 + b0 */ | |||
1915 | t2[bot_size] = carry; | |||
1916 | ||||
1917 | (void)s_kmul(t1, t2, t3, bot_size + 1, bot_size + 1); /* t3 = t1 * t2 */ | |||
1918 | ||||
1919 | /* Now we'll get t1 = a0b0 and t2 = a1b1, and subtract them out so that | |||
1920 | we're left with only the pieces we want: t3 = a1b0 + a0b1 | |||
1921 | */ | |||
1922 | ZERO(t1, buf_size); | |||
1923 | ZERO(t2, buf_size); | |||
1924 | (void)s_kmul(da, db, t1, bot_size, bot_size); /* t1 = a0 * b0 */ | |||
1925 | (void)s_kmul(a_top, b_top, t2, at_size, bt_size); /* t2 = a1 * b1 */ | |||
1926 | ||||
1927 | /* Subtract out t1 and t2 to get the inner product */ | |||
1928 | s_usub(t3, t1, t3, buf_size + 2, buf_size); | |||
1929 | s_usub(t3, t2, t3, buf_size + 2, buf_size); | |||
1930 | ||||
1931 | /* Assemble the output value */ | |||
1932 | COPY(t1, dc, buf_size); | |||
1933 | carry = s_uadd(t3, dc + bot_size, dc + bot_size, buf_size + 1, buf_size); | |||
1934 | assert(carry == 0)((void) sizeof ((carry == 0) ? 1 : 0), __extension__ ({ if (carry == 0) ; else __assert_fail ("carry == 0", "polly/lib/External/isl/imath/imath.c" , 1934, __extension__ __PRETTY_FUNCTION__); })); | |||
1935 | ||||
1936 | carry = | |||
1937 | s_uadd(t2, dc + 2 * bot_size, dc + 2 * bot_size, buf_size, buf_size); | |||
1938 | assert(carry == 0)((void) sizeof ((carry == 0) ? 1 : 0), __extension__ ({ if (carry == 0) ; else __assert_fail ("carry == 0", "polly/lib/External/isl/imath/imath.c" , 1938, __extension__ __PRETTY_FUNCTION__); })); | |||
1939 | ||||
1940 | s_free(t1); /* note t2 and t3 are just internal pointers to t1 */ | |||
1941 | } else { | |||
1942 | s_umul(da, db, dc, size_a, size_b); | |||
1943 | } | |||
1944 | ||||
1945 | return 1; | |||
1946 | } | |||
1947 | ||||
1948 | static void s_umul(mp_digit *da, mp_digit *db, mp_digit *dc, mp_size size_a, | |||
1949 | mp_size size_b) { | |||
1950 | mp_size a, b; | |||
1951 | mp_word w; | |||
1952 | ||||
1953 | for (a = 0; a < size_a; ++a, ++dc, ++da) { | |||
1954 | mp_digit *dct = dc; | |||
1955 | mp_digit *dbt = db; | |||
1956 | ||||
1957 | if (*da == 0) continue; | |||
1958 | ||||
1959 | w = 0; | |||
1960 | for (b = 0; b < size_b; ++b, ++dbt, ++dct) { | |||
1961 | w = (mp_word)*da * (mp_word)*dbt + w + (mp_word)*dct; | |||
1962 | ||||
1963 | *dct = LOWER_HALF(w); | |||
1964 | w = UPPER_HALF(w); | |||
1965 | } | |||
1966 | ||||
1967 | *dct = (mp_digit)w; | |||
1968 | } | |||
1969 | } | |||
1970 | ||||
1971 | static int s_ksqr(mp_digit *da, mp_digit *dc, mp_size size_a) { | |||
1972 | if (multiply_threshold && size_a > multiply_threshold) { | |||
1973 | mp_size bot_size = (size_a + 1) / 2; | |||
1974 | mp_digit *a_top = da + bot_size; | |||
1975 | mp_digit *t1, *t2, *t3, carry; | |||
1976 | mp_size at_size = size_a - bot_size; | |||
1977 | mp_size buf_size = 2 * bot_size; | |||
1978 | ||||
1979 | if ((t1 = s_alloc(4 * buf_size)) == NULL((void*)0)) return 0; | |||
1980 | t2 = t1 + buf_size; | |||
1981 | t3 = t2 + buf_size; | |||
1982 | ZERO(t1, 4 * buf_size); | |||
1983 | ||||
1984 | (void)s_ksqr(da, t1, bot_size); /* t1 = a0 ^ 2 */ | |||
1985 | (void)s_ksqr(a_top, t2, at_size); /* t2 = a1 ^ 2 */ | |||
1986 | ||||
1987 | (void)s_kmul(da, a_top, t3, bot_size, at_size); /* t3 = a0 * a1 */ | |||
1988 | ||||
1989 | /* Quick multiply t3 by 2, shifting left (can't overflow) */ | |||
1990 | { | |||
1991 | int i, top = bot_size + at_size; | |||
1992 | mp_word w, save = 0; | |||
1993 | ||||
1994 | for (i = 0; i < top; ++i) { | |||
1995 | w = t3[i]; | |||
1996 | w = (w << 1) | save; | |||
1997 | t3[i] = LOWER_HALF(w); | |||
1998 | save = UPPER_HALF(w); | |||
1999 | } | |||
2000 | t3[i] = LOWER_HALF(save); | |||
2001 | } | |||
2002 | ||||
2003 | /* Assemble the output value */ | |||
2004 | COPY(t1, dc, 2 * bot_size); | |||
2005 | carry = s_uadd(t3, dc + bot_size, dc + bot_size, buf_size + 1, buf_size); | |||
2006 | assert(carry == 0)((void) sizeof ((carry == 0) ? 1 : 0), __extension__ ({ if (carry == 0) ; else __assert_fail ("carry == 0", "polly/lib/External/isl/imath/imath.c" , 2006, __extension__ __PRETTY_FUNCTION__); })); | |||
2007 | ||||
2008 | carry = | |||
2009 | s_uadd(t2, dc + 2 * bot_size, dc + 2 * bot_size, buf_size, buf_size); | |||
2010 | assert(carry == 0)((void) sizeof ((carry == 0) ? 1 : 0), __extension__ ({ if (carry == 0) ; else __assert_fail ("carry == 0", "polly/lib/External/isl/imath/imath.c" , 2010, __extension__ __PRETTY_FUNCTION__); })); | |||
2011 | ||||
2012 | s_free(t1); /* note that t2 and t2 are internal pointers only */ | |||
2013 | ||||
2014 | } else { | |||
2015 | s_usqr(da, dc, size_a); | |||
2016 | } | |||
2017 | ||||
2018 | return 1; | |||
2019 | } | |||
2020 | ||||
2021 | static void s_usqr(mp_digit *da, mp_digit *dc, mp_size size_a) { | |||
2022 | mp_size i, j; | |||
2023 | mp_word w; | |||
2024 | ||||
2025 | for (i = 0; i < size_a; ++i, dc += 2, ++da) { | |||
2026 | mp_digit *dct = dc, *dat = da; | |||
2027 | ||||
2028 | if (*da == 0) continue; | |||
2029 | ||||
2030 | /* Take care of the first digit, no rollover */ | |||
2031 | w = (mp_word)*dat * (mp_word)*dat + (mp_word)*dct; | |||
2032 | *dct = LOWER_HALF(w); | |||
2033 | w = UPPER_HALF(w); | |||
2034 | ++dat; | |||
2035 | ++dct; | |||
2036 | ||||
2037 | for (j = i + 1; j < size_a; ++j, ++dat, ++dct) { | |||
2038 | mp_word t = (mp_word)*da * (mp_word)*dat; | |||
2039 | mp_word u = w + (mp_word)*dct, ov = 0; | |||
2040 | ||||
2041 | /* Check if doubling t will overflow a word */ | |||
2042 | if (HIGH_BIT_SET(t)) ov = 1; | |||
2043 | ||||
2044 | w = t + t; | |||
2045 | ||||
2046 | /* Check if adding u to w will overflow a word */ | |||
2047 | if (ADD_WILL_OVERFLOW(w, u)) ov = 1; | |||
2048 | ||||
2049 | w += u; | |||
2050 | ||||
2051 | *dct = LOWER_HALF(w); | |||
2052 | w = UPPER_HALF(w); | |||
2053 | if (ov) { | |||
2054 | w += MP_DIGIT_MAX((4294967295U) * 1UL); /* MP_RADIX */ | |||
2055 | ++w; | |||
2056 | } | |||
2057 | } | |||
2058 | ||||
2059 | w = w + *dct; | |||
2060 | *dct = (mp_digit)w; | |||
2061 | while ((w = UPPER_HALF(w)) != 0) { | |||
2062 | ++dct; | |||
2063 | w = w + *dct; | |||
2064 | *dct = LOWER_HALF(w); | |||
2065 | } | |||
2066 | ||||
2067 | assert(w == 0)((void) sizeof ((w == 0) ? 1 : 0), __extension__ ({ if (w == 0 ) ; else __assert_fail ("w == 0", "polly/lib/External/isl/imath/imath.c" , 2067, __extension__ __PRETTY_FUNCTION__); })); | |||
2068 | } | |||
2069 | } | |||
2070 | ||||
2071 | static void s_dadd(mp_int a, mp_digit b) { | |||
2072 | mp_word w = 0; | |||
2073 | mp_digit *da = MP_DIGITS(a); | |||
2074 | mp_size ua = MP_USED(a); | |||
2075 | ||||
2076 | w = (mp_word)*da + b; | |||
2077 | *da++ = LOWER_HALF(w); | |||
2078 | w = UPPER_HALF(w); | |||
2079 | ||||
2080 | for (ua -= 1; ua > 0; --ua, ++da) { | |||
2081 | w = (mp_word)*da + w; | |||
2082 | ||||
2083 | *da = LOWER_HALF(w); | |||
2084 | w = UPPER_HALF(w); | |||
2085 | } | |||
2086 | ||||
2087 | if (w) { | |||
2088 | *da = (mp_digit)w; | |||
2089 | a->used += 1; | |||
2090 | } | |||
2091 | } | |||
2092 | ||||
2093 | static void s_dmul(mp_int a, mp_digit b) { | |||
2094 | mp_word w = 0; | |||
2095 | mp_digit *da = MP_DIGITS(a); | |||
2096 | mp_size ua = MP_USED(a); | |||
2097 | ||||
2098 | while (ua > 0) { | |||
2099 | w = (mp_word)*da * b + w; | |||
2100 | *da++ = LOWER_HALF(w); | |||
2101 | w = UPPER_HALF(w); | |||
2102 | --ua; | |||
2103 | } | |||
2104 | ||||
2105 | if (w) { | |||
2106 | *da = (mp_digit)w; | |||
2107 | a->used += 1; | |||
2108 | } | |||
2109 | } | |||
2110 | ||||
2111 | static void s_dbmul(mp_digit *da, mp_digit b, mp_digit *dc, mp_size size_a) { | |||
2112 | mp_word w = 0; | |||
2113 | ||||
2114 | while (size_a > 0) { | |||
2115 | w = (mp_word)*da++ * (mp_word)b + w; | |||
2116 | ||||
2117 | *dc++ = LOWER_HALF(w); | |||
2118 | w = UPPER_HALF(w); | |||
2119 | --size_a; | |||
2120 | } | |||
2121 | ||||
2122 | if (w) *dc = LOWER_HALF(w); | |||
2123 | } | |||
2124 | ||||
2125 | static mp_digit s_ddiv(mp_int a, mp_digit b) { | |||
2126 | mp_word w = 0, qdigit; | |||
2127 | mp_size ua = MP_USED(a); | |||
2128 | mp_digit *da = MP_DIGITS(a) + ua - 1; | |||
2129 | ||||
2130 | for (/* */; ua > 0; --ua, --da) { | |||
2131 | w = (w << MP_DIGIT_BIT(sizeof(mp_digit) * 8)) | *da; | |||
2132 | ||||
2133 | if (w >= b) { | |||
2134 | qdigit = w / b; | |||
2135 | w = w % b; | |||
2136 | } else { | |||
2137 | qdigit = 0; | |||
2138 | } | |||
2139 | ||||
2140 | *da = (mp_digit)qdigit; | |||
2141 | } | |||
2142 | ||||
2143 | CLAMP(a); | |||
2144 | return (mp_digit)w; | |||
2145 | } | |||
2146 | ||||
2147 | static void s_qdiv(mp_int z, mp_size p2) { | |||
2148 | mp_size ndig = p2 / MP_DIGIT_BIT(sizeof(mp_digit) * 8), nbits = p2 % MP_DIGIT_BIT(sizeof(mp_digit) * 8); | |||
2149 | mp_size uz = MP_USED(z); | |||
2150 | ||||
2151 | if (ndig) { | |||
2152 | mp_size mark; | |||
2153 | mp_digit *to, *from; | |||
2154 | ||||
2155 | if (ndig >= uz) { | |||
2156 | mp_int_zero(z); | |||
2157 | return; | |||
2158 | } | |||
2159 | ||||
2160 | to = MP_DIGITS(z); | |||
2161 | from = to + ndig; | |||
2162 | ||||
2163 | for (mark = ndig; mark < uz; ++mark) { | |||
2164 | *to++ = *from++; | |||
2165 | } | |||
2166 | ||||
2167 | z->used = uz - ndig; | |||
2168 | } | |||
2169 | ||||
2170 | if (nbits) { | |||
2171 | mp_digit d = 0, *dz, save; | |||
2172 | mp_size up = MP_DIGIT_BIT(sizeof(mp_digit) * 8) - nbits; | |||
2173 | ||||
2174 | uz = MP_USED(z); | |||
2175 | dz = MP_DIGITS(z) + uz - 1; | |||
2176 | ||||
2177 | for (/* */; uz > 0; --uz, --dz) { | |||
2178 | save = *dz; | |||
2179 | ||||
2180 | *dz = (*dz >> nbits) | (d << up); | |||
2181 | d = save; | |||
2182 | } | |||
2183 | ||||
2184 | CLAMP(z); | |||
2185 | } | |||
2186 | ||||
2187 | if (MP_USED(z) == 1 && z->digits[0] == 0) z->sign = MP_ZPOS; | |||
2188 | } | |||
2189 | ||||
2190 | static void s_qmod(mp_int z, mp_size p2) { | |||
2191 | mp_size start = p2 / MP_DIGIT_BIT(sizeof(mp_digit) * 8) + 1, rest = p2 % MP_DIGIT_BIT(sizeof(mp_digit) * 8); | |||
2192 | mp_size uz = MP_USED(z); | |||
2193 | mp_digit mask = (1u << rest) - 1; | |||
2194 | ||||
2195 | if (start <= uz) { | |||
2196 | z->used = start; | |||
2197 | z->digits[start - 1] &= mask; | |||
2198 | CLAMP(z); | |||
2199 | } | |||
2200 | } | |||
2201 | ||||
2202 | static int s_qmul(mp_int z, mp_size p2) { | |||
2203 | mp_size uz, need, rest, extra, i; | |||
2204 | mp_digit *from, *to, d; | |||
2205 | ||||
2206 | if (p2 == 0) return 1; | |||
2207 | ||||
2208 | uz = MP_USED(z); | |||
2209 | need = p2 / MP_DIGIT_BIT(sizeof(mp_digit) * 8); | |||
2210 | rest = p2 % MP_DIGIT_BIT(sizeof(mp_digit) * 8); | |||
2211 | ||||
2212 | /* Figure out if we need an extra digit at the top end; this occurs if the | |||
2213 | topmost `rest' bits of the high-order digit of z are not zero, meaning | |||
2214 | they will be shifted off the end if not preserved */ | |||
2215 | extra = 0; | |||
2216 | if (rest != 0) { | |||
2217 | mp_digit *dz = MP_DIGITS(z) + uz - 1; | |||
2218 | ||||
2219 | if ((*dz >> (MP_DIGIT_BIT(sizeof(mp_digit) * 8) - rest)) != 0) extra = 1; | |||
2220 | } | |||
2221 | ||||
2222 | if (!s_pad(z, uz + need + extra)) return 0; | |||
2223 | ||||
2224 | /* If we need to shift by whole digits, do that in one pass, then | |||
2225 | to back and shift by partial digits. | |||
2226 | */ | |||
2227 | if (need > 0) { | |||
2228 | from = MP_DIGITS(z) + uz - 1; | |||
2229 | to = from + need; | |||
2230 | ||||
2231 | for (i = 0; i < uz; ++i) *to-- = *from--; | |||
2232 | ||||
2233 | ZERO(MP_DIGITS(z), need); | |||
2234 | uz += need; | |||
2235 | } | |||
2236 | ||||
2237 | if (rest) { | |||
2238 | d = 0; | |||
2239 | for (i = need, from = MP_DIGITS(z) + need; i < uz; ++i, ++from) { | |||
2240 | mp_digit save = *from; | |||
2241 | ||||
2242 | *from = (*from << rest) | (d >> (MP_DIGIT_BIT(sizeof(mp_digit) * 8) - rest)); | |||
2243 | d = save; | |||
2244 | } | |||
2245 | ||||
2246 | d >>= (MP_DIGIT_BIT(sizeof(mp_digit) * 8) - rest); | |||
2247 | if (d != 0) { | |||
2248 | *from = d; | |||
2249 | uz += extra; | |||
2250 | } | |||
2251 | } | |||
2252 | ||||
2253 | z->used = uz; | |||
2254 | CLAMP(z); | |||
2255 | ||||
2256 | return 1; | |||
2257 | } | |||
2258 | ||||
2259 | /* Compute z = 2^p2 - |z|; requires that 2^p2 >= |z| | |||
2260 | The sign of the result is always zero/positive. | |||
2261 | */ | |||
2262 | static int s_qsub(mp_int z, mp_size p2) { | |||
2263 | mp_digit hi = (1u << (p2 % MP_DIGIT_BIT(sizeof(mp_digit) * 8))), *zp; | |||
2264 | mp_size tdig = (p2 / MP_DIGIT_BIT(sizeof(mp_digit) * 8)), pos; | |||
2265 | mp_word w = 0; | |||
2266 | ||||
2267 | if (!s_pad(z, tdig + 1)) return 0; | |||
2268 | ||||
2269 | for (pos = 0, zp = MP_DIGITS(z); pos < tdig; ++pos, ++zp) { | |||
2270 | w = ((mp_word)MP_DIGIT_MAX((4294967295U) * 1UL) + 1) - w - (mp_word)*zp; | |||
2271 | ||||
2272 | *zp = LOWER_HALF(w); | |||
2273 | w = UPPER_HALF(w) ? 0 : 1; | |||
2274 | } | |||
2275 | ||||
2276 | w = ((mp_word)MP_DIGIT_MAX((4294967295U) * 1UL) + 1 + hi) - w - (mp_word)*zp; | |||
2277 | *zp = LOWER_HALF(w); | |||
2278 | ||||
2279 | assert(UPPER_HALF(w) != 0)((void) sizeof ((UPPER_HALF(w) != 0) ? 1 : 0), __extension__ ( { if (UPPER_HALF(w) != 0) ; else __assert_fail ("UPPER_HALF(w) != 0" , "polly/lib/External/isl/imath/imath.c", 2279, __extension__ __PRETTY_FUNCTION__); })); /* no borrow out should be possible */ | |||
2280 | ||||
2281 | z->sign = MP_ZPOS; | |||
2282 | CLAMP(z); | |||
2283 | ||||
2284 | return 1; | |||
2285 | } | |||
2286 | ||||
2287 | static int s_dp2k(mp_int z) { | |||
2288 | int k = 0; | |||
2289 | mp_digit *dp = MP_DIGITS(z), d; | |||
2290 | ||||
2291 | if (MP_USED(z) == 1 && *dp == 0) return 1; | |||
2292 | ||||
2293 | while (*dp == 0) { | |||
2294 | k += MP_DIGIT_BIT(sizeof(mp_digit) * 8); | |||
2295 | ++dp; | |||
2296 | } | |||
2297 | ||||
2298 | d = *dp; | |||
2299 | while ((d & 1) == 0) { | |||
2300 | d >>= 1; | |||
2301 | ++k; | |||
2302 | } | |||
2303 | ||||
2304 | return k; | |||
2305 | } | |||
2306 | ||||
2307 | static int s_isp2(mp_int z) { | |||
2308 | mp_size uz = MP_USED(z), k = 0; | |||
2309 | mp_digit *dz = MP_DIGITS(z), d; | |||
2310 | ||||
2311 | while (uz > 1) { | |||
2312 | if (*dz++ != 0) return -1; | |||
2313 | k += MP_DIGIT_BIT(sizeof(mp_digit) * 8); | |||
2314 | --uz; | |||
2315 | } | |||
2316 | ||||
2317 | d = *dz; | |||
2318 | while (d > 1) { | |||
2319 | if (d & 1) return -1; | |||
2320 | ++k; | |||
2321 | d >>= 1; | |||
2322 | } | |||
2323 | ||||
2324 | return (int)k; | |||
2325 | } | |||
2326 | ||||
2327 | static int s_2expt(mp_int z, mp_small k) { | |||
2328 | mp_size ndig, rest; | |||
2329 | mp_digit *dz; | |||
2330 | ||||
2331 | ndig = (k + MP_DIGIT_BIT(sizeof(mp_digit) * 8)) / MP_DIGIT_BIT(sizeof(mp_digit) * 8); | |||
2332 | rest = k % MP_DIGIT_BIT(sizeof(mp_digit) * 8); | |||
2333 | ||||
2334 | if (!s_pad(z, ndig)) return 0; | |||
2335 | ||||
2336 | dz = MP_DIGITS(z); | |||
2337 | ZERO(dz, ndig); | |||
2338 | *(dz + ndig - 1) = (1u << rest); | |||
2339 | z->used = ndig; | |||
2340 | ||||
2341 | return 1; | |||
2342 | } | |||
2343 | ||||
2344 | static int s_norm(mp_int a, mp_int b) { | |||
2345 | mp_digit d = b->digits[MP_USED(b) - 1]; | |||
2346 | int k = 0; | |||
2347 | ||||
2348 | while (d < (1u << (mp_digit)(MP_DIGIT_BIT(sizeof(mp_digit) * 8) - 1))) { /* d < (MP_RADIX / 2) */ | |||
2349 | d <<= 1; | |||
2350 | ++k; | |||
2351 | } | |||
2352 | ||||
2353 | /* These multiplications can't fail */ | |||
2354 | if (k != 0) { | |||
2355 | (void)s_qmul(a, (mp_size)k); | |||
2356 | (void)s_qmul(b, (mp_size)k); | |||
2357 | } | |||
2358 | ||||
2359 | return k; | |||
2360 | } | |||
2361 | ||||
2362 | static mp_result s_brmu(mp_int z, mp_int m) { | |||
2363 | mp_size um = MP_USED(m) * 2; | |||
2364 | ||||
2365 | if (!s_pad(z, um)) return MP_MEMORY; | |||
2366 | ||||
2367 | s_2expt(z, MP_DIGIT_BIT(sizeof(mp_digit) * 8) * um); | |||
2368 | return mp_int_div(z, m, z, NULL((void*)0)); | |||
2369 | } | |||
2370 | ||||
2371 | static int s_reduce(mp_int x, mp_int m, mp_int mu, mp_int q1, mp_int q2) { | |||
2372 | mp_size um = MP_USED(m), umb_p1, umb_m1; | |||
2373 | ||||
2374 | umb_p1 = (um + 1) * MP_DIGIT_BIT(sizeof(mp_digit) * 8); | |||
2375 | umb_m1 = (um - 1) * MP_DIGIT_BIT(sizeof(mp_digit) * 8); | |||
2376 | ||||
2377 | if (mp_int_copy(x, q1) != MP_OK) return 0; | |||
2378 | ||||
2379 | /* Compute q2 = floor((floor(x / b^(k-1)) * mu) / b^(k+1)) */ | |||
2380 | s_qdiv(q1, umb_m1); | |||
2381 | UMUL(q1, mu, q2); | |||
2382 | s_qdiv(q2, umb_p1); | |||
2383 | ||||
2384 | /* Set x = x mod b^(k+1) */ | |||
2385 | s_qmod(x, umb_p1); | |||
2386 | ||||
2387 | /* Now, q is a guess for the quotient a / m. | |||
2388 | Compute x - q * m mod b^(k+1), replacing x. This may be off | |||
2389 | by a factor of 2m, but no more than that. | |||
2390 | */ | |||
2391 | UMUL(q2, m, q1); | |||
2392 | s_qmod(q1, umb_p1); | |||
2393 | (void)mp_int_sub(x, q1, x); /* can't fail */ | |||
2394 | ||||
2395 | /* The result may be < 0; if it is, add b^(k+1) to pin it in the proper | |||
2396 | range. */ | |||
2397 | if ((CMPZ(x) < 0) && !s_qsub(x, umb_p1)) return 0; | |||
2398 | ||||
2399 | /* If x > m, we need to back it off until it is in range. This will be | |||
2400 | required at most twice. */ | |||
2401 | if (mp_int_compare(x, m) >= 0) { | |||
2402 | (void)mp_int_sub(x, m, x); | |||
2403 | if (mp_int_compare(x, m) >= 0) { | |||
2404 | (void)mp_int_sub(x, m, x); | |||
2405 | } | |||
2406 | } | |||
2407 | ||||
2408 | /* At this point, x has been properly reduced. */ | |||
2409 | return 1; | |||
2410 | } | |||
2411 | ||||
2412 | /* Perform modular exponentiation using Barrett's method, where mu is the | |||
2413 | reduction constant for m. Assumes a < m, b > 0. */ | |||
2414 | static mp_result s_embar(mp_int a, mp_int b, mp_int m, mp_int mu, mp_int c) { | |||
2415 | mp_digit umu = MP_USED(mu); | |||
2416 | mp_digit *db = MP_DIGITS(b); | |||
2417 | mp_digit *dbt = db + MP_USED(b) - 1; | |||
2418 | ||||
2419 | DECLARE_TEMP(3)struct { mpz_t value[(3)]; int len; mp_result err; } temp_ = { .len = (3), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0); | |||
2420 | REQUIRE(GROW(TEMP(0), 4 * umu))do { temp_.err = (GROW((temp_.value + (0)), 4 * umu)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
2421 | REQUIRE(GROW(TEMP(1), 4 * umu))do { temp_.err = (GROW((temp_.value + (1)), 4 * umu)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
2422 | REQUIRE(GROW(TEMP(2), 4 * umu))do { temp_.err = (GROW((temp_.value + (2)), 4 * umu)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
2423 | ZERO(TEMP(0)(temp_.value + (0))->digits, TEMP(0)(temp_.value + (0))->alloc); | |||
2424 | ZERO(TEMP(1)(temp_.value + (1))->digits, TEMP(1)(temp_.value + (1))->alloc); | |||
2425 | ZERO(TEMP(2)(temp_.value + (2))->digits, TEMP(2)(temp_.value + (2))->alloc); | |||
2426 | ||||
2427 | (void)mp_int_set_value(c, 1); | |||
2428 | ||||
2429 | /* Take care of low-order digits */ | |||
2430 | while (db < dbt) { | |||
2431 | mp_digit d = *db; | |||
2432 | ||||
2433 | for (int i = MP_DIGIT_BIT(sizeof(mp_digit) * 8); i > 0; --i, d >>= 1) { | |||
2434 | if (d & 1) { | |||
2435 | /* The use of a second temporary avoids allocation */ | |||
2436 | UMUL(c, a, TEMP(0)(temp_.value + (0))); | |||
2437 | if (!s_reduce(TEMP(0)(temp_.value + (0)), m, mu, TEMP(1)(temp_.value + (1)), TEMP(2)(temp_.value + (2)))) { | |||
2438 | REQUIRE(MP_MEMORY)do { temp_.err = (MP_MEMORY); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
2439 | } | |||
2440 | mp_int_copy(TEMP(0)(temp_.value + (0)), c); | |||
2441 | } | |||
2442 | ||||
2443 | USQR(a, TEMP(0)(temp_.value + (0))); | |||
2444 | assert(MP_SIGN(TEMP(0)) == MP_ZPOS)((void) sizeof ((MP_SIGN((temp_.value + (0))) == MP_ZPOS) ? 1 : 0), __extension__ ({ if (MP_SIGN((temp_.value + (0))) == MP_ZPOS ) ; else __assert_fail ("MP_SIGN(TEMP(0)) == MP_ZPOS", "polly/lib/External/isl/imath/imath.c" , 2444, __extension__ __PRETTY_FUNCTION__); })); | |||
2445 | if (!s_reduce(TEMP(0)(temp_.value + (0)), m, mu, TEMP(1)(temp_.value + (1)), TEMP(2)(temp_.value + (2)))) { | |||
2446 | REQUIRE(MP_MEMORY)do { temp_.err = (MP_MEMORY); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
2447 | } | |||
2448 | assert(MP_SIGN(TEMP(0)) == MP_ZPOS)((void) sizeof ((MP_SIGN((temp_.value + (0))) == MP_ZPOS) ? 1 : 0), __extension__ ({ if (MP_SIGN((temp_.value + (0))) == MP_ZPOS ) ; else __assert_fail ("MP_SIGN(TEMP(0)) == MP_ZPOS", "polly/lib/External/isl/imath/imath.c" , 2448, __extension__ __PRETTY_FUNCTION__); })); | |||
2449 | mp_int_copy(TEMP(0)(temp_.value + (0)), a); | |||
2450 | } | |||
2451 | ||||
2452 | ++db; | |||
2453 | } | |||
2454 | ||||
2455 | /* Take care of highest-order digit */ | |||
2456 | mp_digit d = *dbt; | |||
2457 | for (;;) { | |||
2458 | if (d & 1) { | |||
2459 | UMUL(c, a, TEMP(0)(temp_.value + (0))); | |||
2460 | if (!s_reduce(TEMP(0)(temp_.value + (0)), m, mu, TEMP(1)(temp_.value + (1)), TEMP(2)(temp_.value + (2)))) { | |||
2461 | REQUIRE(MP_MEMORY)do { temp_.err = (MP_MEMORY); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
2462 | } | |||
2463 | mp_int_copy(TEMP(0)(temp_.value + (0)), c); | |||
2464 | } | |||
2465 | ||||
2466 | d >>= 1; | |||
2467 | if (!d) break; | |||
2468 | ||||
2469 | USQR(a, TEMP(0)(temp_.value + (0))); | |||
2470 | if (!s_reduce(TEMP(0)(temp_.value + (0)), m, mu, TEMP(1)(temp_.value + (1)), TEMP(2)(temp_.value + (2)))) { | |||
2471 | REQUIRE(MP_MEMORY)do { temp_.err = (MP_MEMORY); if (temp_.err != MP_OK) goto CLEANUP ; } while (0); | |||
2472 | } | |||
2473 | (void)mp_int_copy(TEMP(0)(temp_.value + (0)), a); | |||
2474 | } | |||
2475 | ||||
2476 | CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0); | |||
2477 | return MP_OK; | |||
2478 | } | |||
2479 | ||||
2480 | /* Division of nonnegative integers | |||
2481 | ||||
2482 | This function implements division algorithm for unsigned multi-precision | |||
2483 | integers. The algorithm is based on Algorithm D from Knuth's "The Art of | |||
2484 | Computer Programming", 3rd ed. 1998, pg 272-273. | |||
2485 | ||||
2486 | We diverge from Knuth's algorithm in that we do not perform the subtraction | |||
2487 | from the remainder until we have determined that we have the correct | |||
2488 | quotient digit. This makes our algorithm less efficient that Knuth because | |||
2489 | we might have to perform multiple multiplication and comparison steps before | |||
2490 | the subtraction. The advantage is that it is easy to implement and ensure | |||
2491 | correctness without worrying about underflow from the subtraction. | |||
2492 | ||||
2493 | inputs: u a n+m digit integer in base b (b is 2^MP_DIGIT_BIT) | |||
2494 | v a n digit integer in base b (b is 2^MP_DIGIT_BIT) | |||
2495 | n >= 1 | |||
2496 | m >= 0 | |||
2497 | outputs: u / v stored in u | |||
2498 | u % v stored in v | |||
2499 | */ | |||
2500 | static mp_result s_udiv_knuth(mp_int u, mp_int v) { | |||
2501 | /* Force signs to positive */ | |||
2502 | u->sign = MP_ZPOS; | |||
2503 | v->sign = MP_ZPOS; | |||
2504 | ||||
2505 | /* Use simple division algorithm when v is only one digit long */ | |||
2506 | if (MP_USED(v) == 1) { | |||
2507 | mp_digit d, rem; | |||
2508 | d = v->digits[0]; | |||
2509 | rem = s_ddiv(u, d); | |||
2510 | mp_int_set_value(v, rem); | |||
2511 | return MP_OK; | |||
2512 | } | |||
2513 | ||||
2514 | /* Algorithm D | |||
2515 | ||||
2516 | The n and m variables are defined as used by Knuth. | |||
2517 | u is an n digit number with digits u_{n-1}..u_0. | |||
2518 | v is an n+m digit number with digits from v_{m+n-1}..v_0. | |||
2519 | We require that n > 1 and m >= 0 | |||
2520 | */ | |||
2521 | mp_size n = MP_USED(v); | |||
2522 | mp_size m = MP_USED(u) - n; | |||
2523 | assert(n > 1)((void) sizeof ((n > 1) ? 1 : 0), __extension__ ({ if (n > 1) ; else __assert_fail ("n > 1", "polly/lib/External/isl/imath/imath.c" , 2523, __extension__ __PRETTY_FUNCTION__); })); | |||
2524 | /* assert(m >= 0) follows because m is unsigned. */ | |||
2525 | ||||
2526 | /* D1: Normalize. | |||
2527 | The normalization step provides the necessary condition for Theorem B, | |||
2528 | which states that the quotient estimate for q_j, call it qhat | |||
2529 | ||||
2530 | qhat = u_{j+n}u_{j+n-1} / v_{n-1} | |||
2531 | ||||
2532 | is bounded by | |||
2533 | ||||
2534 | qhat - 2 <= q_j <= qhat. | |||
2535 | ||||
2536 | That is, qhat is always greater than the actual quotient digit q, | |||
2537 | and it is never more than two larger than the actual quotient digit. | |||
2538 | */ | |||
2539 | int k = s_norm(u, v); | |||
2540 | ||||
2541 | /* Extend size of u by one if needed. | |||
2542 | ||||
2543 | The algorithm begins with a value of u that has one more digit of input. | |||
2544 | The normalization step sets u_{m+n}..u_0 = 2^k * u_{m+n-1}..u_0. If the | |||
2545 | multiplication did not increase the number of digits of u, we need to add | |||
2546 | a leading zero here. | |||
2547 | */ | |||
2548 | if (k == 0 || MP_USED(u) != m + n + 1) { | |||
2549 | if (!s_pad(u, m + n + 1)) return MP_MEMORY; | |||
2550 | u->digits[m + n] = 0; | |||
2551 | u->used = m + n + 1; | |||
2552 | } | |||
2553 | ||||
2554 | /* Add a leading 0 to v. | |||
2555 | ||||
2556 | The multiplication in step D4 multiplies qhat * 0v_{n-1}..v_0. We need to | |||
2557 | add the leading zero to v here to ensure that the multiplication will | |||
2558 | produce the full n+1 digit result. | |||
2559 | */ | |||
2560 | if (!s_pad(v, n + 1)) return MP_MEMORY; | |||
2561 | v->digits[n] = 0; | |||
2562 | ||||
2563 | /* Initialize temporary variables q and t. | |||
2564 | q allocates space for m+1 digits to store the quotient digits | |||
2565 | t allocates space for n+1 digits to hold the result of q_j*v | |||
2566 | */ | |||
2567 | DECLARE_TEMP(2)struct { mpz_t value[(2)]; int len; mp_result err; } temp_ = { .len = (2), .err = MP_OK, }; do { for (int i = 0; i < temp_ .len; i++) { mp_int_init((temp_.value + (i))); } } while (0); | |||
2568 | REQUIRE(GROW(TEMP(0), m + 1))do { temp_.err = (GROW((temp_.value + (0)), m + 1)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
2569 | REQUIRE(GROW(TEMP(1), n + 1))do { temp_.err = (GROW((temp_.value + (1)), n + 1)); if (temp_ .err != MP_OK) goto CLEANUP; } while (0); | |||
2570 | ||||
2571 | /* D2: Initialize j */ | |||
2572 | int j = m; | |||
2573 | mpz_t r; | |||
2574 | r.digits = MP_DIGITS(u) + j; /* The contents of r are shared with u */ | |||
2575 | r.used = n + 1; | |||
2576 | r.sign = MP_ZPOS; | |||
2577 | r.alloc = MP_ALLOC(u); | |||
2578 | ZERO(TEMP(1)(temp_.value + (1))->digits, TEMP(1)(temp_.value + (1))->alloc); | |||
2579 | ||||
2580 | /* Calculate the m+1 digits of the quotient result */ | |||
2581 | for (; j >= 0; j--) { | |||
2582 | /* D3: Calculate q' */ | |||
2583 | /* r->digits is aligned to position j of the number u */ | |||
2584 | mp_word pfx, qhat; | |||
2585 | pfx = r.digits[n]; | |||
2586 | pfx <<= MP_DIGIT_BIT(sizeof(mp_digit) * 8) / 2; | |||
2587 | pfx <<= MP_DIGIT_BIT(sizeof(mp_digit) * 8) / 2; | |||
2588 | pfx |= r.digits[n - 1]; /* pfx = u_{j+n}{j+n-1} */ | |||
2589 | ||||
2590 | qhat = pfx / v->digits[n - 1]; | |||
2591 | /* Check to see if qhat > b, and decrease qhat if so. | |||
2592 | Theorem B guarantess that qhat is at most 2 larger than the | |||
2593 | actual value, so it is possible that qhat is greater than | |||
2594 | the maximum value that will fit in a digit */ | |||
2595 | if (qhat > MP_DIGIT_MAX((4294967295U) * 1UL)) qhat = MP_DIGIT_MAX((4294967295U) * 1UL); | |||
2596 | ||||
2597 | /* D4,D5,D6: Multiply qhat * v and test for a correct value of q | |||
2598 | ||||
2599 | We proceed a bit different than the way described by Knuth. This way is | |||
2600 | simpler but less efficent. Instead of doing the multiply and subtract | |||
2601 | then checking for underflow, we first do the multiply of qhat * v and | |||
2602 | see if it is larger than the current remainder r. If it is larger, we | |||
2603 | decrease qhat by one and try again. We may need to decrease qhat one | |||
2604 | more time before we get a value that is smaller than r. | |||
2605 | ||||
2606 | This way is less efficent than Knuth because we do more multiplies, but | |||
2607 | we do not need to worry about underflow this way. | |||
2608 | */ | |||
2609 | /* t = qhat * v */ | |||
2610 | s_dbmul(MP_DIGITS(v), (mp_digit)qhat, TEMP(1)(temp_.value + (1))->digits, n + 1); | |||
2611 | TEMP(1)(temp_.value + (1))->used = n + 1; | |||
2612 | CLAMP(TEMP(1)(temp_.value + (1))); | |||
2613 | ||||
2614 | /* Clamp r for the comparison. Comparisons do not like leading zeros. */ | |||
2615 | CLAMP(&r); | |||
2616 | if (s_ucmp(TEMP(1)(temp_.value + (1)), &r) > 0) { /* would the remainder be negative? */ | |||
2617 | qhat -= 1; /* try a smaller q */ | |||
2618 | s_dbmul(MP_DIGITS(v), (mp_digit)qhat, TEMP(1)(temp_.value + (1))->digits, n + 1); | |||
2619 | TEMP(1)(temp_.value + (1))->used = n + 1; | |||
2620 | CLAMP(TEMP(1)(temp_.value + (1))); | |||
2621 | if (s_ucmp(TEMP(1)(temp_.value + (1)), &r) > 0) { /* would the remainder be negative? */ | |||
2622 | assert(qhat > 0)((void) sizeof ((qhat > 0) ? 1 : 0), __extension__ ({ if ( qhat > 0) ; else __assert_fail ("qhat > 0", "polly/lib/External/isl/imath/imath.c" , 2622, __extension__ __PRETTY_FUNCTION__); })); | |||
2623 | qhat -= 1; /* try a smaller q */ | |||
2624 | s_dbmul(MP_DIGITS(v), (mp_digit)qhat, TEMP(1)(temp_.value + (1))->digits, n + 1); | |||
2625 | TEMP(1)(temp_.value + (1))->used = n + 1; | |||
2626 | CLAMP(TEMP(1)(temp_.value + (1))); | |||
2627 | } | |||
2628 | assert(s_ucmp(TEMP(1), &r) <= 0 && "The mathematics failed us.")((void) sizeof ((s_ucmp((temp_.value + (1)), &r) <= 0 && "The mathematics failed us.") ? 1 : 0), __extension__ ({ if ( s_ucmp((temp_.value + (1)), &r) <= 0 && "The mathematics failed us." ) ; else __assert_fail ("s_ucmp(TEMP(1), &r) <= 0 && \"The mathematics failed us.\"" , "polly/lib/External/isl/imath/imath.c", 2628, __extension__ __PRETTY_FUNCTION__); })); | |||
2629 | } | |||
2630 | /* Unclamp r. The D algorithm expects r = u_{j+n}..u_j to always be n+1 | |||
2631 | digits long. */ | |||
2632 | r.used = n + 1; | |||
2633 | ||||
2634 | /* D4: Multiply and subtract | |||
2635 | ||||
2636 | Note: The multiply was completed above so we only need to subtract here. | |||
2637 | */ | |||
2638 | s_usub(r.digits, TEMP(1)(temp_.value + (1))->digits, r.digits, r.used, TEMP(1)(temp_.value + (1))->used); | |||
2639 | ||||
2640 | /* D5: Test remainder | |||
2641 | ||||
2642 | Note: Not needed because we always check that qhat is the correct value | |||
2643 | before performing the subtract. Value cast to mp_digit to prevent | |||
2644 | warning, qhat has been clamped to MP_DIGIT_MAX | |||
2645 | */ | |||
2646 | TEMP(0)(temp_.value + (0))->digits[j] = (mp_digit)qhat; | |||
2647 | ||||
2648 | /* D6: Add back | |||
2649 | Note: Not needed because we always check that qhat is the correct value | |||
2650 | before performing the subtract. | |||
2651 | */ | |||
2652 | ||||
2653 | /* D7: Loop on j */ | |||
2654 | r.digits--; | |||
2655 | ZERO(TEMP(1)(temp_.value + (1))->digits, TEMP(1)(temp_.value + (1))->alloc); | |||
2656 | } | |||
2657 | ||||
2658 | /* Get rid of leading zeros in q */ | |||
2659 | TEMP(0)(temp_.value + (0))->used = m + 1; | |||
2660 | CLAMP(TEMP(0)(temp_.value + (0))); | |||
2661 | ||||
2662 | /* Denormalize the remainder */ | |||
2663 | CLAMP(u); /* use u here because the r.digits pointer is off-by-one */ | |||
2664 | if (k != 0) s_qdiv(u, k); | |||
2665 | ||||
2666 | mp_int_copy(u, v); /* ok: 0 <= r < v */ | |||
2667 | mp_int_copy(TEMP(0)(temp_.value + (0)), u); /* ok: q <= u */ | |||
2668 | ||||
2669 | CLEANUP_TEMP()CLEANUP: do { for (int i = 0; i < temp_.len; i++) { mp_int_clear ((temp_.value + (i))); } if (temp_.err != MP_OK) { return temp_ .err; } } while (0); | |||
2670 | return MP_OK; | |||
2671 | } | |||
2672 | ||||
2673 | static int s_outlen(mp_int z, mp_size r) { | |||
2674 | assert(r >= MP_MIN_RADIX && r <= MP_MAX_RADIX)((void) sizeof ((r >= 2 && r <= 36) ? 1 : 0), __extension__ ({ if (r >= 2 && r <= 36) ; else __assert_fail ("r >= MP_MIN_RADIX && r <= MP_MAX_RADIX", "polly/lib/External/isl/imath/imath.c" , 2674, __extension__ __PRETTY_FUNCTION__); })); | |||
2675 | ||||
2676 | mp_result bits = mp_int_count_bits(z); | |||
2677 | double raw = (double)bits * s_log2[r]; | |||
2678 | ||||
2679 | return (int)(raw + 0.999999); | |||
2680 | } | |||
2681 | ||||
2682 | static mp_size s_inlen(int len, mp_size r) { | |||
2683 | double raw = (double)len / s_log2[r]; | |||
2684 | mp_size bits = (mp_size)(raw + 0.5); | |||
2685 | ||||
2686 | return (mp_size)((bits + (MP_DIGIT_BIT(sizeof(mp_digit) * 8) - 1)) / MP_DIGIT_BIT(sizeof(mp_digit) * 8)) + 1; | |||
2687 | } | |||
2688 | ||||
2689 | static int s_ch2val(char c, int r) { | |||
2690 | int out; | |||
2691 | ||||
2692 | /* | |||
2693 | * In some locales, isalpha() accepts characters outside the range A-Z, | |||
2694 | * producing out<0 or out>=36. The "out >= r" check will always catch | |||
2695 | * out>=36. Though nothing explicitly catches out<0, our caller reacts the | |||
2696 | * same way to every negative return value. | |||
2697 | */ | |||
2698 | if (isdigit((unsigned char)c)((*__ctype_b_loc ())[(int) (((unsigned char)c))] & (unsigned short int) _ISdigit)) | |||
2699 | out = c - '0'; | |||
2700 | else if (r > 10 && isalpha((unsigned char)c)((*__ctype_b_loc ())[(int) (((unsigned char)c))] & (unsigned short int) _ISalpha)) | |||
2701 | out = toupper((unsigned char)c)(__extension__ ({ int __res; if (sizeof ((unsigned char)c) > 1) { if (__builtin_constant_p ((unsigned char)c)) { int __c = ((unsigned char)c); __res = __c < -128 || __c > 255 ? __c : (*__ctype_toupper_loc ())[__c]; } else __res = toupper ((unsigned char)c); } else __res = (*__ctype_toupper_loc ())[(int) ((unsigned char)c)]; __res; })) - 'A' + 10; | |||
2702 | else | |||
2703 | return -1; | |||
2704 | ||||
2705 | return (out >= r) ? -1 : out; | |||
2706 | } | |||
2707 | ||||
2708 | static char s_val2ch(int v, int caps) { | |||
2709 | assert(v >= 0)((void) sizeof ((v >= 0) ? 1 : 0), __extension__ ({ if (v >= 0) ; else __assert_fail ("v >= 0", "polly/lib/External/isl/imath/imath.c" , 2709, __extension__ __PRETTY_FUNCTION__); })); | |||
2710 | ||||
2711 | if (v < 10) { | |||
2712 | return v + '0'; | |||
2713 | } else { | |||
2714 | char out = (v - 10) + 'a'; | |||
2715 | ||||
2716 | if (caps) { | |||
2717 | return toupper((unsigned char)out)(__extension__ ({ int __res; if (sizeof ((unsigned char)out) > 1) { if (__builtin_constant_p ((unsigned char)out)) { int __c = ((unsigned char)out); __res = __c < -128 || __c > 255 ? __c : (*__ctype_toupper_loc ())[__c]; } else __res = toupper ((unsigned char)out); } else __res = (*__ctype_toupper_loc ( ))[(int) ((unsigned char)out)]; __res; })); | |||
2718 | } else { | |||
2719 | return out; | |||
2720 | } | |||
2721 | } | |||
2722 | } | |||
2723 | ||||
2724 | static void s_2comp(unsigned char *buf, int len) { | |||
2725 | unsigned short s = 1; | |||
2726 | ||||
2727 | for (int i = len - 1; i >= 0; --i) { | |||
2728 | unsigned char c = ~buf[i]; | |||
2729 | ||||
2730 | s = c + s; | |||
2731 | c = s & UCHAR_MAX(127*2 +1); | |||
2732 | s >>= CHAR_BIT8; | |||
2733 | ||||
2734 | buf[i] = c; | |||
2735 | } | |||
2736 | ||||
2737 | /* last carry out is ignored */ | |||
2738 | } | |||
2739 | ||||
2740 | static mp_result s_tobin(mp_int z, unsigned char *buf, int *limpos, int pad) { | |||
2741 | int pos = 0, limit = *limpos; | |||
2742 | mp_size uz = MP_USED(z); | |||
2743 | mp_digit *dz = MP_DIGITS(z); | |||
2744 | ||||
2745 | while (uz > 0 && pos < limit) { | |||
2746 | mp_digit d = *dz++; | |||
2747 | int i; | |||
2748 | ||||
2749 | for (i = sizeof(mp_digit); i > 0 && pos < limit; --i) { | |||
2750 | buf[pos++] = (unsigned char)d; | |||
2751 | d >>= CHAR_BIT8; | |||
2752 | ||||
2753 | /* Don't write leading zeroes */ | |||
2754 | if (d == 0 && uz == 1) i = 0; /* exit loop without signaling truncation */ | |||
2755 | } | |||
2756 | ||||
2757 | /* Detect truncation (loop exited with pos >= limit) */ | |||
2758 | if (i > 0) break; | |||
2759 | ||||
2760 | --uz; | |||
2761 | } | |||
2762 | ||||
2763 | if (pad != 0 && (buf[pos - 1] >> (CHAR_BIT8 - 1))) { | |||
2764 | if (pos < limit) { | |||
2765 | buf[pos++] = 0; | |||
2766 | } else { | |||
2767 | uz = 1; | |||
2768 | } | |||
2769 | } | |||
2770 | ||||
2771 | /* Digits are in reverse order, fix that */ | |||
2772 | REV(buf, pos); | |||
2773 | ||||
2774 | /* Return the number of bytes actually written */ | |||
2775 | *limpos = pos; | |||
2776 | ||||
2777 | return (uz == 0) ? MP_OK : MP_TRUNC; | |||
2778 | } | |||
2779 | ||||
2780 | /* Here there be dragons */ |