1 /* Helper functions for checked numbers.
2    Copyright (C) 2001-2010 Roberto Bagnara <bagnara@cs.unipr.it>
3    Copyright (C) 2010-2016 BUGSENG srl (http://bugseng.com)
4 
5 This file is part of the Parma Polyhedra Library (PPL).
6 
7 The PPL is free software; you can redistribute it and/or modify it
8 under the terms of the GNU General Public License as published by the
9 Free Software Foundation; either version 3 of the License, or (at your
10 option) any later version.
11 
12 The PPL is distributed in the hope that it will be useful, but WITHOUT
13 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
14 FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
15 for more details.
16 
17 You should have received a copy of the GNU General Public License
18 along with this program; if not, write to the Free Software Foundation,
19 Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02111-1307, USA.
20 
21 For the most up-to-date information see the Parma Polyhedra Library
22 site: http://bugseng.com/products/ppl/ . */
23 
24 #include "ppl-config.h"
25 #include "checked_defs.hh"
26 #include "C_Integer.hh"
27 
28 namespace Parma_Polyhedra_Library {
29 
30 Minus_Infinity MINUS_INFINITY;
31 Plus_Infinity PLUS_INFINITY;
32 Not_A_Number NOT_A_NUMBER;
33 
34 namespace Checked {
35 
36 //! Holds the precision parameter used for irrational calculations.
37 unsigned irrational_precision;
38 
39 struct number_struct {
40   unsigned int base;
41   bool neg_mantissa;
42   bool neg_exponent;
43   std::string mantissa;
44   unsigned int base_for_exponent;
45   unsigned long exponent;
46 };
47 
48 /*! \brief
49   Returns the integer value associated with the ASCII code \p c, in
50   the base \p base positional number system, if there is such an
51   association; returns \f$-1\f$ otherwise.
52 */
53 inline int
get_digit(char c,unsigned int base=10)54 get_digit(char c, unsigned int base = 10) {
55   unsigned int n;
56   switch (c) {
57   case '0': n = 0; break;
58   case '1': n = 1; break;
59   case '2': n = 2; break;
60   case '3': n = 3; break;
61   case '4': n = 4; break;
62   case '5': n = 5; break;
63   case '6': n = 6; break;
64   case '7': n = 7; break;
65   case '8': n = 8; break;
66   case '9': n = 9; break;
67   case 'a': case 'A': n = 10; break;
68   case 'b': case 'B': n = 11; break;
69   case 'c': case 'C': n = 12; break;
70   case 'd': case 'D': n = 13; break;
71   case 'e': case 'E': n = 14; break;
72   case 'f': case 'F': n = 15; break;
73   case 'g': case 'G': n = 16; break;
74   case 'h': case 'H': n = 17; break;
75   case 'i': case 'I': n = 18; break;
76   case 'j': case 'J': n = 19; break;
77   case 'k': case 'K': n = 20; break;
78   case 'l': case 'L': n = 21; break;
79   case 'm': case 'M': n = 22; break;
80   case 'n': case 'N': n = 23; break;
81   case 'o': case 'O': n = 24; break;
82   case 'p': case 'P': n = 25; break;
83   case 'q': case 'Q': n = 26; break;
84   case 'r': case 'R': n = 27; break;
85   case 's': case 'S': n = 28; break;
86   case 't': case 'T': n = 29; break;
87   case 'u': case 'U': n = 30; break;
88   case 'v': case 'V': n = 31; break;
89   case 'w': case 'W': n = 32; break;
90   case 'x': case 'X': n = 33; break;
91   case 'y': case 'Y': n = 34; break;
92   case 'z': case 'Z': n = 35; break;
93   default:
94     return -1;
95   }
96   if (n >= base) {
97     return -1;
98   }
99   return static_cast<int>(n);
100 }
101 
102 /*! \brief
103   Adds the number represented (in the modulus-and-sign representation)
104   by \p b_neg and \p b_mod to the number represented by \p a_neg and
105   \p a_mod, assigning the result to the latter.  Returns
106   <CODE>false</CODE> is the result cannot be represented; returns
107   <CODE>true</CODE> otherwise.
108 */
109 inline bool
sum_sign(bool & a_neg,unsigned long & a_mod,bool b_neg,unsigned long b_mod)110 sum_sign(bool& a_neg, unsigned long& a_mod,
111          bool b_neg, unsigned long b_mod) {
112   if (a_neg == b_neg) {
113     if (a_mod > C_Integer<unsigned long>::max - b_mod) {
114       return false;
115     }
116     a_mod += b_mod;
117   }
118   else if (a_mod >= b_mod) {
119     a_mod -= b_mod;
120   }
121   else {
122     a_neg = !a_neg;
123     a_mod = b_mod - a_mod;
124   }
125   return true;
126 }
127 
128 
129 /*! \brief
130   Helper function for parse_number(): reads the numerator or
131   denominator part of a number from \p is into \p numer, returning the
132   appropriate Result value.
133 */
134 Result
parse_number_part(std::istream & is,number_struct & numer)135 parse_number_part(std::istream& is, number_struct& numer) {
136   enum anonymous_enum { BASE, INTEGER, FRACTIONAL, EXPONENT } state = BASE;
137   PPL_UNINITIALIZED(unsigned long, max_exp_div);
138   PPL_UNINITIALIZED(int, max_exp_rem);
139   bool empty_exponent = true;
140   bool empty_mantissa = true;
141   long exponent_offset = 0;
142   unsigned exponent_offset_scale = 1;
143   numer.base = 10;
144   numer.base_for_exponent = 10;
145   numer.neg_mantissa = false;
146   numer.neg_exponent = false;
147   numer.mantissa.erase();
148   numer.exponent = 0;
149   char c;
150   do {
151     if (!is.get(c)) {
152       return V_CVT_STR_UNK;
153     }
154   } while (is_space(c));
155   switch (c) {
156   case '-':
157     numer.neg_mantissa = true;
158     // Fall through.
159   case '+':
160     if (!is.get(c)) {
161       return V_CVT_STR_UNK;
162     }
163     if (c == 'i' || c == 'I') {
164       goto infinity;
165     }
166     if (c != '.') {
167       break;
168     }
169     // Fall through.
170   case '.':
171     state = FRACTIONAL;
172     if (!is.get(c)) {
173       return V_CVT_STR_UNK;
174     }
175     break;
176   case 'n':
177   case 'N':
178     if (!is.get(c)) {
179       return V_CVT_STR_UNK;
180     }
181     if (c != 'a' && c != 'A') {
182       goto unexpected_char;
183     }
184     if (!is.get(c)) {
185       return V_CVT_STR_UNK;
186     }
187     if (c != 'n' && c != 'N') {
188       goto unexpected_char;
189     }
190     return V_NAN;
191   infinity:
192   case 'i':
193   case 'I':
194     if (!is.get(c)) {
195       return V_CVT_STR_UNK;
196     }
197     if (c != 'n' && c != 'n') {
198       goto unexpected_char;
199     }
200     if (!is.get(c)) {
201       return V_CVT_STR_UNK;
202     }
203     if (c != 'f' && c != 'F') {
204       goto unexpected_char;
205     }
206     return numer.neg_mantissa ? V_EQ_MINUS_INFINITY : V_EQ_PLUS_INFINITY;
207   }
208   if (state != FRACTIONAL) {
209     if (get_digit(c, 10) < 0) {
210       goto unexpected_char;
211     }
212     char d;
213     if (c == '0' && !is.get(d).fail()) {
214       if (d == 'x' || d == 'X') {
215         numer.base = 16;
216         numer.base_for_exponent = 16;
217         state = INTEGER;
218         if (!is.get(c)) {
219           return V_CVT_STR_UNK;
220         }
221       }
222       else {
223         is.unget();
224       }
225     }
226   }
227   do {
228     switch (state) {
229     case BASE:
230       if (get_digit(c, 10) >= 0) {
231         if (c != '0' || !numer.mantissa.empty()) {
232           numer.mantissa += c;
233         }
234         empty_mantissa = false;
235         break;
236       }
237       if (c == '^') {
238         if (!is.get(c)) {
239           return V_CVT_STR_UNK;
240         }
241         if (c != '^') {
242           goto unexpected_char;
243         }
244         numer.base = 0;
245         for (std::string::const_iterator
246                i = numer.mantissa.begin(); i != numer.mantissa.end(); ++i) {
247           numer.base = numer.base * 10 + static_cast<unsigned>(get_digit(*i, 10));
248           if (numer.base > 36) {
249             goto unexpected_char;
250           }
251         }
252         if (numer.base < 2) {
253           goto unexpected_char;
254         }
255         numer.base_for_exponent = numer.base;
256         numer.mantissa.erase();
257         empty_mantissa = true;
258         state = INTEGER;
259         break;
260       }
261       goto integer;
262     case INTEGER:
263       if (get_digit(c, numer.base) >= 0) {
264         if (c != '0' || !numer.mantissa.empty()) {
265           numer.mantissa += c;
266         }
267         empty_mantissa = false;
268         break;
269       }
270     integer:
271       if (c == '.') {
272         state = FRACTIONAL;
273         break;
274       }
275       goto fractional;
276     case FRACTIONAL:
277       if (get_digit(c, numer.base) >= 0) {
278         --exponent_offset;
279         if (c != '0' || !numer.mantissa.empty()) {
280           numer.mantissa += c;
281         }
282         empty_mantissa = false;
283         break;
284       }
285     fractional:
286       if (empty_mantissa) {
287         goto unexpected_char;
288       }
289       if (c == 'e' || c == 'E') {
290         goto exponent;
291       }
292       if (c == 'p' || c == 'P') {
293         if (numer.base == 16) {
294           numer.base_for_exponent = 2;
295           exponent_offset_scale = 4;
296           goto exponent;
297         }
298         else {
299           goto unexpected_char;
300         }
301       }
302       if (c == '*') {
303         if (!is.get(c)) {
304           return V_CVT_STR_UNK;
305         }
306         if (c != '^') {
307           goto unexpected_char;
308         }
309       exponent:
310         state = EXPONENT;
311         PPL_ASSERT(numer.base >= 2);
312         const long l_max = C_Integer<long>::max;
313         max_exp_div = static_cast<unsigned long>(l_max) / numer.base;
314         max_exp_rem = static_cast<int>(l_max % static_cast<long>(numer.base));
315         if (!is.get(c)) {
316           return V_CVT_STR_UNK;
317         }
318         if (c == '-') {
319           numer.neg_exponent = true;
320           break;
321         }
322         if (c == '+') {
323           break;
324         }
325         continue;
326       }
327       is.unget();
328       goto ok;
329     case EXPONENT:
330       const int d = get_digit(c, 10);
331       if (d >= 0) {
332         empty_exponent = false;
333         if (numer.exponent > max_exp_div
334             || (numer.exponent == max_exp_div && d > max_exp_rem)) {
335           return V_CVT_STR_UNK;
336         }
337         numer.exponent = 10 * numer.exponent + static_cast<unsigned long>(d);
338         break;
339       }
340       if (empty_exponent) {
341         goto unexpected_char;
342       }
343       is.unget();
344       goto ok;
345     }
346     is.get(c);
347   } while (!is.fail());
348 
349   if (empty_mantissa || is.bad()) {
350     return V_CVT_STR_UNK;
351   }
352 
353  ok:
354   {
355     std::string::size_type n = numer.mantissa.size();
356     while (n > 0 && numer.mantissa[n - 1] == '0') {
357       --n;
358       ++exponent_offset;
359     }
360     numer.mantissa.erase(n);
361     bool neg;
362     if (exponent_offset < 0) {
363       neg = true;
364       exponent_offset = -exponent_offset;
365     }
366     else {
367       neg = false;
368     }
369     sum_sign(numer.neg_exponent, numer.exponent,
370              neg, static_cast<unsigned long>(exponent_offset) * exponent_offset_scale);
371     return V_EQ;
372   }
373 
374  unexpected_char:
375   is.unget();
376   return V_CVT_STR_UNK;
377 }
378 
379 /*! \brief
380   Reads a number from \p is writing it into \p numer, the numerator,
381   and \p denom, the denominator; the appropriate Result value is
382   returned.
383 */
384 Result
parse_number(std::istream & is,number_struct & numer,number_struct & denom)385 parse_number(std::istream& is, number_struct& numer, number_struct& denom) {
386   // Read the numerator.
387   Result r = parse_number_part(is, numer);
388   if (r != V_EQ) {
389     return r;
390   }
391   char c;
392   is.get(c);
393   if (is.bad()) {
394     return V_CVT_STR_UNK;
395   }
396   if (!is) {
397     denom.base = 0;
398     return r;
399   }
400   if (c != '/') {
401     is.unget();
402     denom.base = 0;
403     return r;
404   }
405   // Read the denominator.
406   r = parse_number_part(is, denom);
407   if (r != V_EQ) {
408     return V_CVT_STR_UNK;
409   }
410   if (numer.base == denom.base
411       && numer.base_for_exponent == denom.base_for_exponent) {
412     if (sum_sign(numer.neg_exponent, numer.exponent,
413                  !denom.neg_exponent, denom.exponent)) {
414       if (numer.neg_exponent) {
415         denom.neg_exponent = false;
416         denom.exponent = numer.exponent;
417         numer.exponent = 0;
418       }
419       else {
420         denom.exponent = 0;
421       }
422     }
423   }
424   return V_EQ;
425 }
426 
427 
428 Result
input_mpq(mpq_class & to,std::istream & is)429 input_mpq(mpq_class& to, std::istream& is) {
430   number_struct numer_struct;
431   number_struct denom_struct;
432   const Result r = parse_number(is, numer_struct, denom_struct);
433   if (r == V_CVT_STR_UNK) {
434     is.setstate(is.failbit);
435     return r;
436   }
437   is.clear(is.rdstate() & ~is.failbit);
438   if (r != V_EQ) {
439     return r;
440   }
441   if (denom_struct.base != 0 && denom_struct.mantissa.empty()) {
442       return V_NAN;
443   }
444   if (numer_struct.mantissa.empty()) {
445     to = 0;
446     return V_EQ;
447   }
448   const mpz_ptr numer = to.get_num().get_mpz_t();
449   const mpz_ptr denom = to.get_den().get_mpz_t();
450   mpz_set_str(numer, numer_struct.mantissa.c_str(),
451               static_cast<int>(numer_struct.base));
452   if (denom_struct.base != 0) {
453     if (numer_struct.neg_mantissa != denom_struct.neg_mantissa) {
454       mpz_neg(numer, numer);
455     }
456     mpz_set_str(denom, denom_struct.mantissa.c_str(),
457                 static_cast<int>(denom_struct.base));
458     if (numer_struct.exponent != 0 || denom_struct.exponent != 0) {
459       // Multiply the exponents into the numerator and denominator.
460       mpz_t z;
461       mpz_init(z);
462       if (numer_struct.exponent != 0) {
463         mpz_ui_pow_ui(z,
464                       numer_struct.base_for_exponent, numer_struct.exponent);
465         if (numer_struct.neg_exponent) {
466           mpz_mul(denom, denom, z);
467         }
468         else {
469           mpz_mul(numer, numer, z);
470         }
471       }
472       if (denom_struct.exponent != 0) {
473         mpz_ui_pow_ui(z,
474                       denom_struct.base_for_exponent, denom_struct.exponent);
475         if (denom_struct.neg_exponent) {
476           mpz_mul(numer, numer, z);
477         }
478         else {
479           mpz_mul(denom, denom, z);
480         }
481       }
482       mpz_clear(z);
483     }
484   }
485   else {
486     if (numer_struct.neg_mantissa) {
487       mpz_neg(numer, numer);
488     }
489     if (numer_struct.exponent != 0) {
490       if (numer_struct.neg_exponent) {
491         // Add the negative exponent as a denominator.
492         mpz_ui_pow_ui(denom,
493                       numer_struct.base_for_exponent, numer_struct.exponent);
494         goto end;
495       }
496       // Multiply the exponent into the numerator.
497       mpz_t z;
498       mpz_init(z);
499       mpz_ui_pow_ui(z,
500                     numer_struct.base_for_exponent, numer_struct.exponent);
501       mpz_mul(numer, numer, z);
502       mpz_clear(z);
503     }
504     mpz_set_ui(denom, 1);
505     return V_EQ;
506   }
507  end:
508   // GMP operators require rationals in canonical form.
509   to.canonicalize();
510   return V_EQ;
511 }
512 
513 /* NOTE: q is overwritten! */
float_mpq_to_string(mpq_class & q)514 std::string float_mpq_to_string(mpq_class& q) {
515   const mpz_ptr n = q.get_num().get_mpz_t();
516   const mpz_ptr d = q.get_den().get_mpz_t();
517   const unsigned long decimals = mpz_sizeinbase(d, 2) - 1;
518   if (decimals != 0) {
519     mpz_ui_pow_ui(d, 5, decimals);
520     mpz_mul(n, n, d);
521   }
522   size_t bufsize = mpz_sizeinbase(n, 10);
523   if (bufsize < decimals) {
524     bufsize = decimals + 4;
525   }
526   else {
527     bufsize += 3;
528   }
529   char buf[bufsize];
530   mpz_get_str(buf, 10, n);
531   if (decimals != 0) {
532     const size_t len = strlen(buf);
533     if (decimals < len) {
534       memmove(&buf[len - decimals + 1], &buf[len - decimals], decimals + 1);
535       buf[len - decimals] = '.';
536     }
537     else {
538       const size_t zeroes = decimals - len;
539       memmove(&buf[2 + zeroes], &buf[0], len + 1);
540       buf[0] = '0';
541       buf[1] = '.';
542       memset(&buf[2], '0', zeroes);
543     }
544   }
545   return buf;
546 }
547 
548 } // namespace Checked
549 
550 } // namespace Parma_Polyhedra_Library
551