Bug Summary

File:tools/polly/lib/External/isl/imath/imrat.c
Warning:line 739, column 10
Assigned value is garbage or undefined

Annotated Source Code

1/*
2 Name: imrat.c
3 Purpose: Arbitrary precision rational arithmetic routines.
4 Author: M. J. Fromberger <http://spinning-yarns.org/michael/>
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 "imrat.h"
28#include <stdlib.h>
29#include <string.h>
30#include <ctype.h>
31#include <assert.h>
32
33#define TEMP(K)(temp + (K)) (temp + (K))
34#define SETUP(E, C)do{if((res = (E)) != MP_OK) goto CLEANUP; ++(C);}while(0) \
35do{if((res = (E)) != MP_OK) goto CLEANUP; ++(C);}while(0)
36
37/* Argument checking:
38 Use CHECK() where a return value is required; NRCHECK() elsewhere */
39#define CHECK(TEST)((TEST) ? (void) (0) : __assert_fail ("TEST", "/tmp/buildd/llvm-toolchain-snapshot-4.0~svn290870/tools/polly/lib/External/isl/imath/imrat.c"
, 39, __PRETTY_FUNCTION__))
assert(TEST)((TEST) ? (void) (0) : __assert_fail ("TEST", "/tmp/buildd/llvm-toolchain-snapshot-4.0~svn290870/tools/polly/lib/External/isl/imath/imrat.c"
, 39, __PRETTY_FUNCTION__))
40#define NRCHECK(TEST)((TEST) ? (void) (0) : __assert_fail ("TEST", "/tmp/buildd/llvm-toolchain-snapshot-4.0~svn290870/tools/polly/lib/External/isl/imath/imrat.c"
, 40, __PRETTY_FUNCTION__))
assert(TEST)((TEST) ? (void) (0) : __assert_fail ("TEST", "/tmp/buildd/llvm-toolchain-snapshot-4.0~svn290870/tools/polly/lib/External/isl/imath/imrat.c"
, 40, __PRETTY_FUNCTION__))
41
42/* Reduce the given rational, in place, to lowest terms and canonical form.
43 Zero is represented as 0/1, one as 1/1. Signs are adjusted so that the sign
44 of the numerator is definitive. */
45static mp_result s_rat_reduce(mp_rat r);
46
47/* Common code for addition and subtraction operations on rationals. */
48static mp_result s_rat_combine(mp_rat a, mp_rat b, mp_rat c,
49 mp_result (*comb_f)(mp_int, mp_int, mp_int));
50
51mp_result mp_rat_init(mp_rat r)
52{
53 mp_result res;
54
55 if ((res = mp_int_init(MP_NUMER_P(r)(&((r)->num)))) != MP_OK)
56 return res;
57 if ((res = mp_int_init(MP_DENOM_P(r)(&((r)->den)))) != MP_OK) {
58 mp_int_clear(MP_NUMER_P(r)(&((r)->num)));
59 return res;
60 }
61
62 return mp_int_set_value(MP_DENOM_P(r)(&((r)->den)), 1);
63}
64
65mp_rat mp_rat_alloc(void)
66{
67 mp_rat out = malloc(sizeof(*out));
68
69 if (out != NULL((void*)0)) {
70 if (mp_rat_init(out) != MP_OK) {
71 free(out);
72 return NULL((void*)0);
73 }
74 }
75
76 return out;
77}
78
79mp_result mp_rat_reduce(mp_rat r) {
80 return s_rat_reduce(r);
81}
82
83mp_result mp_rat_init_size(mp_rat r, mp_size n_prec, mp_size d_prec)
84{
85 mp_result res;
86
87 if ((res = mp_int_init_size(MP_NUMER_P(r)(&((r)->num)), n_prec)) != MP_OK)
88 return res;
89 if ((res = mp_int_init_size(MP_DENOM_P(r)(&((r)->den)), d_prec)) != MP_OK) {
90 mp_int_clear(MP_NUMER_P(r)(&((r)->num)));
91 return res;
92 }
93
94 return mp_int_set_value(MP_DENOM_P(r)(&((r)->den)), 1);
95}
96
97mp_result mp_rat_init_copy(mp_rat r, mp_rat old)
98{
99 mp_result res;
100
101 if ((res = mp_int_init_copy(MP_NUMER_P(r)(&((r)->num)), MP_NUMER_P(old)(&((old)->num)))) != MP_OK)
102 return res;
103 if ((res = mp_int_init_copy(MP_DENOM_P(r)(&((r)->den)), MP_DENOM_P(old)(&((old)->den)))) != MP_OK)
104 mp_int_clear(MP_NUMER_P(r)(&((r)->num)));
105
106 return res;
107}
108
109mp_result mp_rat_set_value(mp_rat r, mp_small numer, mp_small denom)
110{
111 mp_result res;
112
113 if (denom == 0)
114 return MP_UNDEF;
115
116 if ((res = mp_int_set_value(MP_NUMER_P(r)(&((r)->num)), numer)) != MP_OK)
117 return res;
118 if ((res = mp_int_set_value(MP_DENOM_P(r)(&((r)->den)), denom)) != MP_OK)
119 return res;
120
121 return s_rat_reduce(r);
122}
123
124mp_result mp_rat_set_uvalue(mp_rat r, mp_usmall numer, mp_usmall denom)
125{
126 mp_result res;
127
128 if (denom == 0)
129 return MP_UNDEF;
130
131 if ((res = mp_int_set_uvalue(MP_NUMER_P(r)(&((r)->num)), numer)) != MP_OK)
132 return res;
133 if ((res = mp_int_set_uvalue(MP_DENOM_P(r)(&((r)->den)), denom)) != MP_OK)
134 return res;
135
136 return s_rat_reduce(r);
137}
138
139void mp_rat_clear(mp_rat r)
140{
141 mp_int_clear(MP_NUMER_P(r)(&((r)->num)));
142 mp_int_clear(MP_DENOM_P(r)(&((r)->den)));
143
144}
145
146void mp_rat_free(mp_rat r)
147{
148 NRCHECK(r != NULL)((r != ((void*)0)) ? (void) (0) : __assert_fail ("r != ((void*)0)"
, "/tmp/buildd/llvm-toolchain-snapshot-4.0~svn290870/tools/polly/lib/External/isl/imath/imrat.c"
, 148, __PRETTY_FUNCTION__))
;
149
150 if (r->num.digits != NULL((void*)0))
151 mp_rat_clear(r);
152
153 free(r);
154}
155
156mp_result mp_rat_numer(mp_rat r, mp_int z)
157{
158 return mp_int_copy(MP_NUMER_P(r)(&((r)->num)), z);
159}
160
161mp_int mp_rat_numer_ref(mp_rat r)
162{
163 return MP_NUMER_P(r)(&((r)->num));
164}
165
166
167mp_result mp_rat_denom(mp_rat r, mp_int z)
168{
169 return mp_int_copy(MP_DENOM_P(r)(&((r)->den)), z);
170}
171
172mp_int mp_rat_denom_ref(mp_rat r)
173{
174 return MP_DENOM_P(r)(&((r)->den));
175}
176
177mp_sign mp_rat_sign(mp_rat r)
178{
179 return MP_SIGN(MP_NUMER_P(r))(((&((r)->num)))->sign);
180}
181
182mp_result mp_rat_copy(mp_rat a, mp_rat c)
183{
184 mp_result res;
185
186 if ((res = mp_int_copy(MP_NUMER_P(a)(&((a)->num)), MP_NUMER_P(c)(&((c)->num)))) != MP_OK)
187 return res;
188
189 res = mp_int_copy(MP_DENOM_P(a)(&((a)->den)), MP_DENOM_P(c)(&((c)->den)));
190 return res;
191}
192
193void mp_rat_zero(mp_rat r)
194{
195 mp_int_zero(MP_NUMER_P(r)(&((r)->num)));
196 mp_int_set_value(MP_DENOM_P(r)(&((r)->den)), 1);
197
198}
199
200mp_result mp_rat_abs(mp_rat a, mp_rat c)
201{
202 mp_result res;
203
204 if ((res = mp_int_abs(MP_NUMER_P(a)(&((a)->num)), MP_NUMER_P(c)(&((c)->num)))) != MP_OK)
205 return res;
206
207 res = mp_int_abs(MP_DENOM_P(a)(&((a)->den)), MP_DENOM_P(c)(&((c)->den)));
208 return res;
209}
210
211mp_result mp_rat_neg(mp_rat a, mp_rat c)
212{
213 mp_result res;
214
215 if ((res = mp_int_neg(MP_NUMER_P(a)(&((a)->num)), MP_NUMER_P(c)(&((c)->num)))) != MP_OK)
216 return res;
217
218 res = mp_int_copy(MP_DENOM_P(a)(&((a)->den)), MP_DENOM_P(c)(&((c)->den)));
219 return res;
220}
221
222mp_result mp_rat_recip(mp_rat a, mp_rat c)
223{
224 mp_result res;
225
226 if (mp_rat_compare_zero(a) == 0)
227 return MP_UNDEF;
228
229 if ((res = mp_rat_copy(a, c)) != MP_OK)
230 return res;
231
232 mp_int_swap(MP_NUMER_P(c)(&((c)->num)), MP_DENOM_P(c)(&((c)->den)));
233
234 /* Restore the signs of the swapped elements */
235 {
236 mp_sign tmp = MP_SIGN(MP_NUMER_P(c))(((&((c)->num)))->sign);
237
238 MP_SIGN(MP_NUMER_P(c))(((&((c)->num)))->sign) = MP_SIGN(MP_DENOM_P(c))(((&((c)->den)))->sign);
239 MP_SIGN(MP_DENOM_P(c))(((&((c)->den)))->sign) = tmp;
240 }
241
242 return MP_OK;
243}
244
245mp_result mp_rat_add(mp_rat a, mp_rat b, mp_rat c)
246{
247 return s_rat_combine(a, b, c, mp_int_add);
248
249}
250
251mp_result mp_rat_sub(mp_rat a, mp_rat b, mp_rat c)
252{
253 return s_rat_combine(a, b, c, mp_int_sub);
254
255}
256
257mp_result mp_rat_mul(mp_rat a, mp_rat b, mp_rat c)
258{
259 mp_result res;
260
261 if ((res = mp_int_mul(MP_NUMER_P(a)(&((a)->num)), MP_NUMER_P(b)(&((b)->num)), MP_NUMER_P(c)(&((c)->num)))) != MP_OK)
262 return res;
263
264 if (mp_int_compare_zero(MP_NUMER_P(c)(&((c)->num))) != 0) {
265 if ((res = mp_int_mul(MP_DENOM_P(a)(&((a)->den)), MP_DENOM_P(b)(&((b)->den)), MP_DENOM_P(c)(&((c)->den)))) != MP_OK)
266 return res;
267 }
268
269 return s_rat_reduce(c);
270}
271
272mp_result mp_rat_div(mp_rat a, mp_rat b, mp_rat c)
273{
274 mp_result res = MP_OK;
275
276 if (mp_rat_compare_zero(b) == 0)
277 return MP_UNDEF;
278
279 if (c == a || c == b) {
280 mpz_t tmp;
281
282 if ((res = mp_int_init(&tmp)) != MP_OK) return res;
283 if ((res = mp_int_mul(MP_NUMER_P(a)(&((a)->num)), MP_DENOM_P(b)(&((b)->den)), &tmp)) != MP_OK)
284 goto CLEANUP;
285 if ((res = mp_int_mul(MP_DENOM_P(a)(&((a)->den)), MP_NUMER_P(b)(&((b)->num)), MP_DENOM_P(c)(&((c)->den)))) != MP_OK)
286 goto CLEANUP;
287 res = mp_int_copy(&tmp, MP_NUMER_P(c)(&((c)->num)));
288
289 CLEANUP:
290 mp_int_clear(&tmp);
291 }
292 else {
293 if ((res = mp_int_mul(MP_NUMER_P(a)(&((a)->num)), MP_DENOM_P(b)(&((b)->den)), MP_NUMER_P(c)(&((c)->num)))) != MP_OK)
294 return res;
295 if ((res = mp_int_mul(MP_DENOM_P(a)(&((a)->den)), MP_NUMER_P(b)(&((b)->num)), MP_DENOM_P(c)(&((c)->den)))) != MP_OK)
296 return res;
297 }
298
299 if (res != MP_OK)
300 return res;
301 else
302 return s_rat_reduce(c);
303}
304
305mp_result mp_rat_add_int(mp_rat a, mp_int b, mp_rat c)
306{
307 mpz_t tmp;
308 mp_result res;
309
310 if ((res = mp_int_init_copy(&tmp, b)) != MP_OK)
311 return res;
312
313 if ((res = mp_int_mul(&tmp, MP_DENOM_P(a)(&((a)->den)), &tmp)) != MP_OK)
314 goto CLEANUP;
315
316 if ((res = mp_rat_copy(a, c)) != MP_OK)
317 goto CLEANUP;
318
319 if ((res = mp_int_add(MP_NUMER_P(c)(&((c)->num)), &tmp, MP_NUMER_P(c)(&((c)->num)))) != MP_OK)
320 goto CLEANUP;
321
322 res = s_rat_reduce(c);
323
324 CLEANUP:
325 mp_int_clear(&tmp);
326 return res;
327}
328
329mp_result mp_rat_sub_int(mp_rat a, mp_int b, mp_rat c)
330{
331 mpz_t tmp;
332 mp_result res;
333
334 if ((res = mp_int_init_copy(&tmp, b)) != MP_OK)
335 return res;
336
337 if ((res = mp_int_mul(&tmp, MP_DENOM_P(a)(&((a)->den)), &tmp)) != MP_OK)
338 goto CLEANUP;
339
340 if ((res = mp_rat_copy(a, c)) != MP_OK)
341 goto CLEANUP;
342
343 if ((res = mp_int_sub(MP_NUMER_P(c)(&((c)->num)), &tmp, MP_NUMER_P(c)(&((c)->num)))) != MP_OK)
344 goto CLEANUP;
345
346 res = s_rat_reduce(c);
347
348 CLEANUP:
349 mp_int_clear(&tmp);
350 return res;
351}
352
353mp_result mp_rat_mul_int(mp_rat a, mp_int b, mp_rat c)
354{
355 mp_result res;
356
357 if ((res = mp_rat_copy(a, c)) != MP_OK)
358 return res;
359
360 if ((res = mp_int_mul(MP_NUMER_P(c)(&((c)->num)), b, MP_NUMER_P(c)(&((c)->num)))) != MP_OK)
361 return res;
362
363 return s_rat_reduce(c);
364}
365
366mp_result mp_rat_div_int(mp_rat a, mp_int b, mp_rat c)
367{
368 mp_result res;
369
370 if (mp_int_compare_zero(b) == 0)
371 return MP_UNDEF;
372
373 if ((res = mp_rat_copy(a, c)) != MP_OK)
374 return res;
375
376 if ((res = mp_int_mul(MP_DENOM_P(c)(&((c)->den)), b, MP_DENOM_P(c)(&((c)->den)))) != MP_OK)
377 return res;
378
379 return s_rat_reduce(c);
380}
381
382mp_result mp_rat_expt(mp_rat a, mp_small b, mp_rat c)
383{
384 mp_result res;
385
386 /* Special cases for easy powers. */
387 if (b == 0)
388 return mp_rat_set_value(c, 1, 1);
389 else if(b == 1)
390 return mp_rat_copy(a, c);
391
392 /* Since rationals are always stored in lowest terms, it is not necessary to
393 reduce again when raising to an integer power. */
394 if ((res = mp_int_expt(MP_NUMER_P(a)(&((a)->num)), b, MP_NUMER_P(c)(&((c)->num)))) != MP_OK)
395 return res;
396
397 return mp_int_expt(MP_DENOM_P(a)(&((a)->den)), b, MP_DENOM_P(c)(&((c)->den)));
398}
399
400int mp_rat_compare(mp_rat a, mp_rat b)
401{
402 /* Quick check for opposite signs. Works because the sign of the numerator
403 is always definitive. */
404 if (MP_SIGN(MP_NUMER_P(a))(((&((a)->num)))->sign) != MP_SIGN(MP_NUMER_P(b))(((&((b)->num)))->sign)) {
405 if (MP_SIGN(MP_NUMER_P(a))(((&((a)->num)))->sign) == MP_ZPOS)
406 return 1;
407 else
408 return -1;
409 }
410 else {
411 /* Compare absolute magnitudes; if both are positive, the answer stands,
412 otherwise it needs to be reflected about zero. */
413 int cmp = mp_rat_compare_unsigned(a, b);
414
415 if (MP_SIGN(MP_NUMER_P(a))(((&((a)->num)))->sign) == MP_ZPOS)
416 return cmp;
417 else
418 return -cmp;
419 }
420}
421
422int mp_rat_compare_unsigned(mp_rat a, mp_rat b)
423{
424 /* If the denominators are equal, we can quickly compare numerators without
425 multiplying. Otherwise, we actually have to do some work. */
426 if (mp_int_compare_unsigned(MP_DENOM_P(a)(&((a)->den)), MP_DENOM_P(b)(&((b)->den))) == 0)
427 return mp_int_compare_unsigned(MP_NUMER_P(a)(&((a)->num)), MP_NUMER_P(b)(&((b)->num)));
428
429 else {
430 mpz_t temp[2];
431 mp_result res;
432 int cmp = INT_MAX2147483647, last = 0;
433
434 /* t0 = num(a) * den(b), t1 = num(b) * den(a) */
435 SETUP(mp_int_init_copy(TEMP(last), MP_NUMER_P(a)), last)do{if((res = (mp_int_init_copy((temp + (last)), (&((a)->
num))))) != MP_OK) goto CLEANUP; ++(last);}while(0)
;
436 SETUP(mp_int_init_copy(TEMP(last), MP_NUMER_P(b)), last)do{if((res = (mp_int_init_copy((temp + (last)), (&((b)->
num))))) != MP_OK) goto CLEANUP; ++(last);}while(0)
;
437
438 if ((res = mp_int_mul(TEMP(0)(temp + (0)), MP_DENOM_P(b)(&((b)->den)), TEMP(0)(temp + (0)))) != MP_OK ||
439 (res = mp_int_mul(TEMP(1)(temp + (1)), MP_DENOM_P(a)(&((a)->den)), TEMP(1)(temp + (1)))) != MP_OK)
440 goto CLEANUP;
441
442 cmp = mp_int_compare_unsigned(TEMP(0)(temp + (0)), TEMP(1)(temp + (1)));
443
444 CLEANUP:
445 while (--last >= 0)
446 mp_int_clear(TEMP(last)(temp + (last)));
447
448 return cmp;
449 }
450}
451
452int mp_rat_compare_zero(mp_rat r)
453{
454 return mp_int_compare_zero(MP_NUMER_P(r)(&((r)->num)));
455}
456
457int mp_rat_compare_value(mp_rat r, mp_small n, mp_small d)
458{
459 mpq_t tmp;
460 mp_result res;
461 int out = INT_MAX2147483647;
462
463 if ((res = mp_rat_init(&tmp)) != MP_OK)
464 return out;
465 if ((res = mp_rat_set_value(&tmp, n, d)) != MP_OK)
466 goto CLEANUP;
467
468 out = mp_rat_compare(r, &tmp);
469
470 CLEANUP:
471 mp_rat_clear(&tmp);
472 return out;
473}
474
475int mp_rat_is_integer(mp_rat r)
476{
477 return (mp_int_compare_value(MP_DENOM_P(r)(&((r)->den)), 1) == 0);
478}
479
480mp_result mp_rat_to_ints(mp_rat r, mp_small *num, mp_small *den)
481{
482 mp_result res;
483
484 if ((res = mp_int_to_int(MP_NUMER_P(r)(&((r)->num)), num)) != MP_OK)
485 return res;
486
487 res = mp_int_to_int(MP_DENOM_P(r)(&((r)->den)), den);
488 return res;
489}
490
491mp_result mp_rat_to_string(mp_rat r, mp_size radix, char *str, int limit)
492{
493 char *start;
494 int len;
495 mp_result res;
496
497 /* Write the numerator. The sign of the rational number is written by the
498 underlying integer implementation. */
499 if ((res = mp_int_to_string(MP_NUMER_P(r)(&((r)->num)), radix, str, limit)) != MP_OK)
500 return res;
501
502 /* If the value is zero, don't bother writing any denominator */
503 if (mp_int_compare_zero(MP_NUMER_P(r)(&((r)->num))) == 0)
504 return MP_OK;
505
506 /* Locate the end of the numerator, and make sure we are not going to exceed
507 the limit by writing a slash. */
508 len = strlen(str);
509 start = str + len;
510 limit -= len;
511 if(limit == 0)
512 return MP_TRUNC;
513
514 *start++ = '/';
515 limit -= 1;
516
517 res = mp_int_to_string(MP_DENOM_P(r)(&((r)->den)), radix, start, limit);
518 return res;
519}
520
521mp_result mp_rat_to_decimal(mp_rat r, mp_size radix, mp_size prec,
522 mp_round_mode round, char *str, int limit)
523{
524 mpz_t temp[3];
525 mp_result res;
526 char *start = str;
527 int len, lead_0, left = limit, last = 0;
528
529 SETUP(mp_int_init_copy(TEMP(last), MP_NUMER_P(r)), last)do{if((res = (mp_int_init_copy((temp + (last)), (&((r)->
num))))) != MP_OK) goto CLEANUP; ++(last);}while(0)
;
530 SETUP(mp_int_init(TEMP(last)), last)do{if((res = (mp_int_init((temp + (last))))) != MP_OK) goto CLEANUP
; ++(last);}while(0)
;
531 SETUP(mp_int_init(TEMP(last)), last)do{if((res = (mp_int_init((temp + (last))))) != MP_OK) goto CLEANUP
; ++(last);}while(0)
;
532
533 /* Get the unsigned integer part by dividing denominator into the absolute
534 value of the numerator. */
535 mp_int_abs(TEMP(0)(temp + (0)), TEMP(0)(temp + (0)));
536 if ((res = mp_int_div(TEMP(0)(temp + (0)), MP_DENOM_P(r)(&((r)->den)), TEMP(0)(temp + (0)), TEMP(1)(temp + (1)))) != MP_OK)
537 goto CLEANUP;
538
539 /* Now: T0 = integer portion, unsigned;
540 T1 = remainder, from which fractional part is computed. */
541
542 /* Count up leading zeroes after the radix point. */
543 for (lead_0 = 0; lead_0 < prec && mp_int_compare(TEMP(1)(temp + (1)), MP_DENOM_P(r)(&((r)->den))) < 0;
544 ++lead_0) {
545 if ((res = mp_int_mul_value(TEMP(1)(temp + (1)), radix, TEMP(1)(temp + (1)))) != MP_OK)
546 goto CLEANUP;
547 }
548
549 /* Multiply remainder by a power of the radix sufficient to get the right
550 number of significant figures. */
551 if (prec > lead_0) {
552 if ((res = mp_int_expt_value(radix, prec - lead_0, TEMP(2)(temp + (2)))) != MP_OK)
553 goto CLEANUP;
554 if ((res = mp_int_mul(TEMP(1)(temp + (1)), TEMP(2)(temp + (2)), TEMP(1)(temp + (1)))) != MP_OK)
555 goto CLEANUP;
556 }
557 if ((res = mp_int_div(TEMP(1)(temp + (1)), MP_DENOM_P(r)(&((r)->den)), TEMP(1)(temp + (1)), TEMP(2)(temp + (2)))) != MP_OK)
558 goto CLEANUP;
559
560 /* Now: T1 = significant digits of fractional part;
561 T2 = leftovers, to use for rounding.
562
563 At this point, what we do depends on the rounding mode. The default is
564 MP_ROUND_DOWN, for which everything is as it should be already.
565 */
566 switch (round) {
567 int cmp;
568
569 case MP_ROUND_UP:
570 if (mp_int_compare_zero(TEMP(2)(temp + (2))) != 0) {
571 if (prec == 0)
572 res = mp_int_add_value(TEMP(0)(temp + (0)), 1, TEMP(0)(temp + (0)));
573 else
574 res = mp_int_add_value(TEMP(1)(temp + (1)), 1, TEMP(1)(temp + (1)));
575 }
576 break;
577
578 case MP_ROUND_HALF_UP:
579 case MP_ROUND_HALF_DOWN:
580 if ((res = mp_int_mul_pow2(TEMP(2)(temp + (2)), 1, TEMP(2)(temp + (2)))) != MP_OK)
581 goto CLEANUP;
582
583 cmp = mp_int_compare(TEMP(2)(temp + (2)), MP_DENOM_P(r)(&((r)->den)));
584
585 if (round == MP_ROUND_HALF_UP)
586 cmp += 1;
587
588 if (cmp > 0) {
589 if (prec == 0)
590 res = mp_int_add_value(TEMP(0)(temp + (0)), 1, TEMP(0)(temp + (0)));
591 else
592 res = mp_int_add_value(TEMP(1)(temp + (1)), 1, TEMP(1)(temp + (1)));
593 }
594 break;
595
596 case MP_ROUND_DOWN:
597 break; /* No action required */
598
599 default:
600 return MP_BADARG; /* Invalid rounding specifier */
601 }
602
603 /* The sign of the output should be the sign of the numerator, but if all the
604 displayed digits will be zero due to the precision, a negative shouldn't
605 be shown. */
606 if (MP_SIGN(MP_NUMER_P(r))(((&((r)->num)))->sign) == MP_NEG &&
607 (mp_int_compare_zero(TEMP(0)(temp + (0))) != 0 ||
608 mp_int_compare_zero(TEMP(1)(temp + (1))) != 0)) {
609 *start++ = '-';
610 left -= 1;
611 }
612
613 if ((res = mp_int_to_string(TEMP(0)(temp + (0)), radix, start, left)) != MP_OK)
614 goto CLEANUP;
615
616 len = strlen(start);
617 start += len;
618 left -= len;
619
620 if (prec == 0)
621 goto CLEANUP;
622
623 *start++ = '.';
624 left -= 1;
625
626 if (left < prec + 1) {
627 res = MP_TRUNC;
628 goto CLEANUP;
629 }
630
631 memset(start, '0', lead_0 - 1);
632 left -= lead_0;
633 start += lead_0 - 1;
634
635 res = mp_int_to_string(TEMP(1)(temp + (1)), radix, start, left);
636
637 CLEANUP:
638 while (--last >= 0)
639 mp_int_clear(TEMP(last)(temp + (last)));
640
641 return res;
642}
643
644mp_result mp_rat_string_len(mp_rat r, mp_size radix)
645{
646 mp_result n_len, d_len = 0;
647
648 n_len = mp_int_string_len(MP_NUMER_P(r)(&((r)->num)), radix);
649
650 if (mp_int_compare_zero(MP_NUMER_P(r)(&((r)->num))) != 0)
651 d_len = mp_int_string_len(MP_DENOM_P(r)(&((r)->den)), radix);
652
653 /* Though simplistic, this formula is correct. Space for the sign flag is
654 included in n_len, and the space for the NUL that is counted in n_len
655 counts for the separator here. The space for the NUL counted in d_len
656 counts for the final terminator here. */
657
658 return n_len + d_len;
659
660}
661
662mp_result mp_rat_decimal_len(mp_rat r, mp_size radix, mp_size prec)
663{
664 int z_len, f_len;
665
666 z_len = mp_int_string_len(MP_NUMER_P(r)(&((r)->num)), radix);
667
668 if (prec == 0)
669 f_len = 1; /* terminator only */
670 else
671 f_len = 1 + prec + 1; /* decimal point, digits, terminator */
672
673 return z_len + f_len;
674}
675
676mp_result mp_rat_read_string(mp_rat r, mp_size radix, const char *str)
677{
678 return mp_rat_read_cstring(r, radix, str, NULL((void*)0));
679}
680
681mp_result mp_rat_read_cstring(mp_rat r, mp_size radix, const char *str,
682 char **end)
683{
684 mp_result res;
685 char *endp;
686
687 if ((res = mp_int_read_cstring(MP_NUMER_P(r)(&((r)->num)), radix, str, &endp)) != MP_OK &&
5
Taking true branch
688 (res != MP_TRUNC))
689 return res;
690
691 /* Skip whitespace between numerator and (possible) separator */
692 while (isspace((unsigned char) *endp)((*__ctype_b_loc ())[(int) (((unsigned char) *endp))] & (
unsigned short int) _ISspace)
)
693 ++endp;
694
695 /* If there is no separator, we will stop reading at this point. */
696 if (*endp != '/') {
697 mp_int_set_value(MP_DENOM_P(r)(&((r)->den)), 1);
698 if (end != NULL((void*)0))
699 *end = endp;
700 return res;
701 }
702
703 ++endp; /* skip separator */
704 if ((res = mp_int_read_cstring(MP_DENOM_P(r)(&((r)->den)), radix, endp, end)) != MP_OK)
705 return res;
706
707 /* Make sure the value is well-defined */
708 if (mp_int_compare_zero(MP_DENOM_P(r)(&((r)->den))) == 0)
709 return MP_UNDEF;
710
711 /* Reduce to lowest terms */
712 return s_rat_reduce(r);
713}
714
715/* Read a string and figure out what format it's in. The radix may be supplied
716 as zero to use "default" behaviour.
717
718 This function will accept either a/b notation or decimal notation.
719 */
720mp_result mp_rat_read_ustring(mp_rat r, mp_size radix, const char *str,
721 char **end)
722{
723 char *endp;
1
'endp' declared without an initial value
724 mp_result res;
725
726 if (radix == 0)
2
Assuming 'radix' is not equal to 0
3
Taking false branch
727 radix = 10; /* default to decimal input */
728
729 if ((res = mp_rat_read_cstring(r, radix, str, &endp)) != MP_OK) {
4
Calling 'mp_rat_read_cstring'
6
Returning from 'mp_rat_read_cstring'
7
Taking false branch
730 if (res == MP_TRUNC) {
731 if (*endp == '.')
732 res = mp_rat_read_cdecimal(r, radix, str, &endp);
733 }
734 else
735 return res;
736 }
737
738 if (end != NULL((void*)0))
8
Assuming 'end' is not equal to NULL
9
Taking true branch
739 *end = endp;
10
Assigned value is garbage or undefined
740
741 return res;
742}
743
744mp_result mp_rat_read_decimal(mp_rat r, mp_size radix, const char *str)
745{
746 return mp_rat_read_cdecimal(r, radix, str, NULL((void*)0));
747}
748
749mp_result mp_rat_read_cdecimal(mp_rat r, mp_size radix, const char *str,
750 char **end)
751{
752 mp_result res;
753 mp_sign osign;
754 char *endp;
755
756 while (isspace((unsigned char) *str)((*__ctype_b_loc ())[(int) (((unsigned char) *str))] & (unsigned
short int) _ISspace)
)
757 ++str;
758
759 switch (*str) {
760 case '-':
761 osign = MP_NEG;
762 break;
763 default:
764 osign = MP_ZPOS;
765 }
766
767 if ((res = mp_int_read_cstring(MP_NUMER_P(r)(&((r)->num)), radix, str, &endp)) != MP_OK &&
768 (res != MP_TRUNC))
769 return res;
770
771 /* This needs to be here. */
772 (void) mp_int_set_value(MP_DENOM_P(r)(&((r)->den)), 1);
773
774 if (*endp != '.') {
775 if (end != NULL((void*)0))
776 *end = endp;
777 return res;
778 }
779
780 /* If the character following the decimal point is whitespace or a sign flag,
781 we will consider this a truncated value. This special case is because
782 mp_int_read_string() will consider whitespace or sign flags to be valid
783 starting characters for a value, and we do not want them following the
784 decimal point.
785
786 Once we have done this check, it is safe to read in the value of the
787 fractional piece as a regular old integer.
788 */
789 ++endp;
790 if (*endp == '\0') {
791 if (end != NULL((void*)0))
792 *end = endp;
793 return MP_OK;
794 }
795 else if(isspace((unsigned char) *endp)((*__ctype_b_loc ())[(int) (((unsigned char) *endp))] & (
unsigned short int) _ISspace)
|| *endp == '-' || *endp == '+') {
796 return MP_TRUNC;
797 }
798 else {
799 mpz_t frac;
800 mp_result save_res;
801 char *save = endp;
802 int num_lz = 0;
803
804 /* Make a temporary to hold the part after the decimal point. */
805 if ((res = mp_int_init(&frac)) != MP_OK)
806 return res;
807
808 if ((res = mp_int_read_cstring(&frac, radix, endp, &endp)) != MP_OK &&
809 (res != MP_TRUNC))
810 goto CLEANUP;
811
812 /* Save this response for later. */
813 save_res = res;
814
815 if (mp_int_compare_zero(&frac) == 0)
816 goto FINISHED;
817
818 /* Discard trailing zeroes (somewhat inefficiently) */
819 while (mp_int_divisible_value(&frac, radix))
820 if ((res = mp_int_div_value(&frac, radix, &frac, NULL((void*)0))) != MP_OK)
821 goto CLEANUP;
822
823 /* Count leading zeros after the decimal point */
824 while (save[num_lz] == '0')
825 ++num_lz;
826
827 /* Find the least power of the radix that is at least as large as the
828 significant value of the fractional part, ignoring leading zeroes. */
829 (void) mp_int_set_value(MP_DENOM_P(r)(&((r)->den)), radix);
830
831 while (mp_int_compare(MP_DENOM_P(r)(&((r)->den)), &frac) < 0) {
832 if ((res = mp_int_mul_value(MP_DENOM_P(r)(&((r)->den)), radix, MP_DENOM_P(r)(&((r)->den)))) != MP_OK)
833 goto CLEANUP;
834 }
835
836 /* Also shift by enough to account for leading zeroes */
837 while (num_lz > 0) {
838 if ((res = mp_int_mul_value(MP_DENOM_P(r)(&((r)->den)), radix, MP_DENOM_P(r)(&((r)->den)))) != MP_OK)
839 goto CLEANUP;
840
841 --num_lz;
842 }
843
844 /* Having found this power, shift the numerator leftward that many, digits,
845 and add the nonzero significant digits of the fractional part to get the
846 result. */
847 if ((res = mp_int_mul(MP_NUMER_P(r)(&((r)->num)), MP_DENOM_P(r)(&((r)->den)), MP_NUMER_P(r)(&((r)->num)))) != MP_OK)
848 goto CLEANUP;
849
850 { /* This addition needs to be unsigned. */
851 MP_SIGN(MP_NUMER_P(r))(((&((r)->num)))->sign) = MP_ZPOS;
852 if ((res = mp_int_add(MP_NUMER_P(r)(&((r)->num)), &frac, MP_NUMER_P(r)(&((r)->num)))) != MP_OK)
853 goto CLEANUP;
854
855 MP_SIGN(MP_NUMER_P(r))(((&((r)->num)))->sign) = osign;
856 }
857 if ((res = s_rat_reduce(r)) != MP_OK)
858 goto CLEANUP;
859
860 /* At this point, what we return depends on whether reading the fractional
861 part was truncated or not. That information is saved from when we
862 called mp_int_read_string() above. */
863 FINISHED:
864 res = save_res;
865 if (end != NULL((void*)0))
866 *end = endp;
867
868 CLEANUP:
869 mp_int_clear(&frac);
870
871 return res;
872 }
873}
874
875/* Private functions for internal use. Make unchecked assumptions about format
876 and validity of inputs. */
877
878static mp_result s_rat_reduce(mp_rat r)
879{
880 mpz_t gcd;
881 mp_result res = MP_OK;
882
883 if (mp_int_compare_zero(MP_NUMER_P(r)(&((r)->num))) == 0) {
884 mp_int_set_value(MP_DENOM_P(r)(&((r)->den)), 1);
885 return MP_OK;
886 }
887
888 /* If the greatest common divisor of the numerator and denominator is greater
889 than 1, divide it out. */
890 if ((res = mp_int_init(&gcd)) != MP_OK)
891 return res;
892
893 if ((res = mp_int_gcd(MP_NUMER_P(r)(&((r)->num)), MP_DENOM_P(r)(&((r)->den)), &gcd)) != MP_OK)
894 goto CLEANUP;
895
896 if (mp_int_compare_value(&gcd, 1) != 0) {
897 if ((res = mp_int_div(MP_NUMER_P(r)(&((r)->num)), &gcd, MP_NUMER_P(r)(&((r)->num)), NULL((void*)0))) != MP_OK)
898 goto CLEANUP;
899 if ((res = mp_int_div(MP_DENOM_P(r)(&((r)->den)), &gcd, MP_DENOM_P(r)(&((r)->den)), NULL((void*)0))) != MP_OK)
900 goto CLEANUP;
901 }
902
903 /* Fix up the signs of numerator and denominator */
904 if (MP_SIGN(MP_NUMER_P(r))(((&((r)->num)))->sign) == MP_SIGN(MP_DENOM_P(r))(((&((r)->den)))->sign))
905 MP_SIGN(MP_NUMER_P(r))(((&((r)->num)))->sign) = MP_SIGN(MP_DENOM_P(r))(((&((r)->den)))->sign) = MP_ZPOS;
906 else {
907 MP_SIGN(MP_NUMER_P(r))(((&((r)->num)))->sign) = MP_NEG;
908 MP_SIGN(MP_DENOM_P(r))(((&((r)->den)))->sign) = MP_ZPOS;
909 }
910
911 CLEANUP:
912 mp_int_clear(&gcd);
913
914 return res;
915}
916
917static mp_result s_rat_combine(mp_rat a, mp_rat b, mp_rat c,
918 mp_result (*comb_f)(mp_int, mp_int, mp_int))
919{
920 mp_result res;
921
922 /* Shortcut when denominators are already common */
923 if (mp_int_compare(MP_DENOM_P(a)(&((a)->den)), MP_DENOM_P(b)(&((b)->den))) == 0) {
924 if ((res = (comb_f)(MP_NUMER_P(a)(&((a)->num)), MP_NUMER_P(b)(&((b)->num)), MP_NUMER_P(c)(&((c)->num)))) != MP_OK)
925 return res;
926 if ((res = mp_int_copy(MP_DENOM_P(a)(&((a)->den)), MP_DENOM_P(c)(&((c)->den)))) != MP_OK)
927 return res;
928
929 return s_rat_reduce(c);
930 }
931 else {
932 mpz_t temp[2];
933 int last = 0;
934
935 SETUP(mp_int_init_copy(TEMP(last), MP_NUMER_P(a)), last)do{if((res = (mp_int_init_copy((temp + (last)), (&((a)->
num))))) != MP_OK) goto CLEANUP; ++(last);}while(0)
;
936 SETUP(mp_int_init_copy(TEMP(last), MP_NUMER_P(b)), last)do{if((res = (mp_int_init_copy((temp + (last)), (&((b)->
num))))) != MP_OK) goto CLEANUP; ++(last);}while(0)
;
937
938 if ((res = mp_int_mul(TEMP(0)(temp + (0)), MP_DENOM_P(b)(&((b)->den)), TEMP(0)(temp + (0)))) != MP_OK)
939 goto CLEANUP;
940 if ((res = mp_int_mul(TEMP(1)(temp + (1)), MP_DENOM_P(a)(&((a)->den)), TEMP(1)(temp + (1)))) != MP_OK)
941 goto CLEANUP;
942 if ((res = (comb_f)(TEMP(0)(temp + (0)), TEMP(1)(temp + (1)), MP_NUMER_P(c)(&((c)->num)))) != MP_OK)
943 goto CLEANUP;
944
945 res = mp_int_mul(MP_DENOM_P(a)(&((a)->den)), MP_DENOM_P(b)(&((b)->den)), MP_DENOM_P(c)(&((c)->den)));
946
947 CLEANUP:
948 while (--last >= 0)
949 mp_int_clear(TEMP(last)(temp + (last)));
950
951 if (res == MP_OK)
952 return s_rat_reduce(c);
953 else
954 return res;
955 }
956}
957
958/* Here there be dragons */