cprover
ieee_float.cpp
Go to the documentation of this file.
1 /*******************************************************************\
2 
3 Module:
4 
5 Author: Daniel Kroening, kroening@kroening.com
6 
7 \*******************************************************************/
8 
9 #include "ieee_float.h"
10 
11 #include <cstdint>
12 #include <ostream>
13 #include <cmath>
14 #include <limits>
15 
16 #include "arith_tools.h"
17 #include "invariant.h"
18 #include "std_expr.h"
19 #include "std_types.h"
20 
22 {
23  return power(2, e-1)-1;
24 }
25 
27 {
28  floatbv_typet result;
29  result.set_f(f);
30  result.set_width(width());
31  if(x86_extended)
32  result.set(ID_x86_extended, true);
33  return result;
34 }
35 
37 {
38  return power(2, e)-1;
39 }
40 
42 {
43  return power(2, f)-1;
44 }
45 
47 {
48  std::size_t width=type.get_width();
49  f=type.get_f();
50  DATA_INVARIANT(f != 0, "mantissa must be at least 1 bit");
52  f < width,
53  "mantissa bits must be less than "
54  "originating type width");
55  e=width-f-1;
56  x86_extended=type.get_bool(ID_x86_extended);
57  if(x86_extended)
58  e=e-1; // no hidden bit
59 }
60 
61 void ieee_floatt::print(std::ostream &out) const
62 {
63  out << to_ansi_c_string();
64 }
65 
66 std::string ieee_floatt::format(const format_spect &format_spec) const
67 {
68  std::string result;
69 
70  switch(format_spec.style)
71  {
73  result+=to_string_decimal(format_spec.precision);
74  break;
75 
77  result+=to_string_scientific(format_spec.precision);
78  break;
79 
81  {
82  // On Linux, the man page says:
83  // "Style e is used if the exponent from its conversion
84  // is less than -4 or greater than or equal to the precision."
85  //
86  // On BSD, it's
87  // "The argument is printed in style f (F) or in style e (E)
88  // whichever gives full precision in minimum space."
89 
90  mp_integer _exponent, _fraction;
91  extract_base10(_fraction, _exponent);
92 
93  mp_integer adjusted_exponent=base10_digits(_fraction)+_exponent;
94 
95  if(adjusted_exponent>=format_spec.precision ||
96  adjusted_exponent<-4)
97  {
98  result+=to_string_scientific(format_spec.precision);
99  }
100  else
101  {
102  result+=to_string_decimal(format_spec.precision);
103 
104  // Implementations tested also appear to suppress trailing zeros
105  // and trailing dots.
106 
107  {
108  std::size_t trunc_pos=result.find_last_not_of('0');
109  if(trunc_pos!=std::string::npos)
110  result.resize(trunc_pos+1);
111  }
112 
113  if(!result.empty() && result.back()=='.')
114  result.resize(result.size()-1);
115  }
116  }
117  break;
118  }
119 
120  while(result.size()<format_spec.min_width)
121  result=" "+result;
122 
123  return result;
124 }
125 
127 {
128  PRECONDITION(src >= 0);
129  mp_integer tmp=src;
130  mp_integer result=0;
131  while(tmp!=0) { ++result; tmp/=10; }
132  return result;
133 }
134 
135 std::string ieee_floatt::to_string_decimal(std::size_t precision) const
136 {
137  std::string result;
138 
139  if(sign_flag)
140  result+='-';
141 
142  if((NaN_flag || infinity_flag) && !sign_flag)
143  result+='+';
144 
145  // special cases
146  if(NaN_flag)
147  result+="NaN";
148  else if(infinity_flag)
149  result+="inf";
150  else if(is_zero())
151  {
152  result+='0';
153 
154  // add zeros, if needed
155  if(precision>0)
156  {
157  result+='.';
158  for(std::size_t i=0; i<precision; i++)
159  result+='0';
160  }
161  }
162  else
163  {
164  mp_integer _exponent, _fraction;
165  extract_base2(_fraction, _exponent);
166 
167  // convert to base 10
168  if(_exponent>=0)
169  {
170  result+=integer2string(_fraction*power(2, _exponent));
171 
172  // add dot and zeros, if needed
173  if(precision>0)
174  {
175  result+='.';
176  for(std::size_t i=0; i<precision; i++)
177  result+='0';
178  }
179  }
180  else
181  {
182  #if 1
183  mp_integer position=-_exponent;
184 
185  // 10/2=5 -- this makes it base 10
186  _fraction*=power(5, position);
187 
188  // apply rounding
189  if(position>precision)
190  {
191  mp_integer r=power(10, position-precision);
192  mp_integer remainder=_fraction%r;
193  _fraction/=r;
194  // not sure if this is the right kind of rounding here
195  if(remainder>=r/2)
196  ++_fraction;
197  position=precision;
198  }
199 
200  std::string tmp=integer2string(_fraction);
201 
202  // pad with zeros from the front, if needed
203  while(mp_integer(tmp.size())<=position) tmp="0"+tmp;
204 
205  const std::size_t dot =
206  tmp.size() - numeric_cast_v<std::size_t>(position);
207  result+=std::string(tmp, 0, dot)+'.';
208  result+=std::string(tmp, dot, std::string::npos);
209 
210  // append zeros if needed
211  for(mp_integer i=position; i<precision; ++i)
212  result+='0';
213  #else
214 
215  result+=integer2string(_fraction);
216 
217  if(_exponent!=0)
218  result+="*2^"+integer2string(_exponent);
219 
220  #endif
221  }
222  }
223 
224  return result;
225 }
226 
229 std::string ieee_floatt::to_string_scientific(std::size_t precision) const
230 {
231  std::string result;
232 
233  if(sign_flag)
234  result+='-';
235 
236  if((NaN_flag || infinity_flag) && !sign_flag)
237  result+='+';
238 
239  // special cases
240  if(NaN_flag)
241  result+="NaN";
242  else if(infinity_flag)
243  result+="inf";
244  else if(is_zero())
245  {
246  result+='0';
247 
248  // add zeros, if needed
249  if(precision>0)
250  {
251  result+='.';
252  for(std::size_t i=0; i<precision; i++)
253  result+='0';
254  }
255 
256  result+="e0";
257  }
258  else
259  {
260  mp_integer _exponent, _fraction;
261  extract_base10(_fraction, _exponent);
262 
263  // C99 appears to say that conversion to decimal should
264  // use the currently selected IEEE rounding mode.
265  if(base10_digits(_fraction)>precision+1)
266  {
267  // re-align
268  mp_integer distance=base10_digits(_fraction)-(precision+1);
269  mp_integer p=power(10, distance);
270  mp_integer remainder=_fraction%p;
271  _fraction/=p;
272  _exponent+=distance;
273 
274  if(remainder==p/2)
275  {
276  // need to do rounding mode here
277  ++_fraction;
278  }
279  else if(remainder>p/2)
280  ++_fraction;
281  }
282 
283  std::string decimals=integer2string(_fraction);
284 
285  CHECK_RETURN(!decimals.empty());
286 
287  // First add top digit to result.
288  result+=decimals[0];
289 
290  // Now add dot and further zeros, if needed.
291  if(precision>0)
292  {
293  result+='.';
294 
295  while(decimals.size()<precision+1)
296  decimals+='0';
297 
298  result+=decimals.substr(1, precision);
299  }
300 
301  // add exponent
302  result+='e';
303 
304  std::string exponent_str=
305  integer2string(base10_digits(_fraction)+_exponent-1);
306 
307  if(!exponent_str.empty() && exponent_str[0]!='-')
308  result+='+';
309 
310  result+=exponent_str;
311  }
312 
313  return result;
314 }
315 
317 {
318  PRECONDITION(spec.f != 0);
319  PRECONDITION(spec.e != 0);
320 
321  {
322  mp_integer tmp=i;
323 
324  // split this apart
325  mp_integer pf=power(2, spec.f);
326  fraction=tmp%pf;
327  tmp/=pf;
328 
329  mp_integer pe=power(2, spec.e);
330  exponent=tmp%pe;
331  tmp/=pe;
332 
333  sign_flag=(tmp!=0);
334  }
335 
336  // NaN?
337  if(exponent==spec.max_exponent() && fraction!=0)
338  {
339  make_NaN();
340  }
341  else if(exponent==spec.max_exponent() && fraction==0) // Infinity
342  {
343  NaN_flag=false;
344  infinity_flag=true;
345  }
346  else if(exponent==0 && fraction==0) // zero
347  {
348  NaN_flag=false;
349  infinity_flag=false;
350  }
351  else if(exponent==0) // denormal?
352  {
353  NaN_flag=false;
354  infinity_flag=false;
355  exponent=-spec.bias()+1; // NOT -spec.bias()!
356  }
357  else // normal
358  {
359  NaN_flag=false;
360  infinity_flag=false;
361  fraction+=power(2, spec.f); // hidden bit!
362  exponent-=spec.bias(); // un-bias
363  }
364 }
365 
367 {
368  return fraction>=power(2, spec.f);
369 }
370 
372 {
373  mp_integer result=0;
374 
375  // sign bit
376  if(sign_flag)
377  result+=power(2, spec.e+spec.f);
378 
379  if(NaN_flag)
380  {
381  result+=power(2, spec.f)*spec.max_exponent();
382  result+=1;
383  }
384  else if(infinity_flag)
385  {
386  result+=power(2, spec.f)*spec.max_exponent();
387  }
388  else if(fraction==0 && exponent==0)
389  {
390  // zero
391  }
392  else if(is_normal()) // normal?
393  {
394  // fraction -- need to hide hidden bit
395  result+=fraction-power(2, spec.f); // hidden bit
396 
397  // exponent -- bias!
398  result+=power(2, spec.f)*(exponent+spec.bias());
399  }
400  else // denormal
401  {
402  result+=fraction; // denormal -- no hidden bit
403  // the exponent is zero
404  }
405 
406  return result;
407 }
408 
410  mp_integer &_fraction,
411  mp_integer &_exponent) const
412 {
413  if(is_zero() || is_NaN() || is_infinity())
414  {
415  _fraction=_exponent=0;
416  return;
417  }
418 
419  _exponent=exponent;
420  _fraction=fraction;
421 
422  // adjust exponent
423  _exponent-=spec.f;
424 
425  // try to integer-ize
426  while((_fraction%2)==0)
427  {
428  _fraction/=2;
429  ++_exponent;
430  }
431 }
432 
434  mp_integer &_fraction,
435  mp_integer &_exponent) const
436 {
437  if(is_zero() || is_NaN() || is_infinity())
438  {
439  _fraction=_exponent=0;
440  return;
441  }
442 
443  _exponent=exponent;
444  _fraction=fraction;
445 
446  // adjust exponent
447  _exponent-=spec.f;
448 
449  // now make it base 10
450  if(_exponent>=0)
451  {
452  _fraction*=power(2, _exponent);
453  _exponent=0;
454  }
455  else // _exponent<0
456  {
457  // 10/2=5 -- this makes it base 10
458  _fraction*=power(5, -_exponent);
459  }
460 
461  // try to re-normalize
462  while((_fraction%10)==0)
463  {
464  _fraction/=10;
465  ++_exponent;
466  }
467 }
468 
470  const mp_integer &_fraction,
471  const mp_integer &_exponent)
472 {
473  sign_flag=_fraction<0;
474  fraction=_fraction;
475  if(sign_flag)
477  exponent=_exponent;
478  exponent+=spec.f;
479  align();
480 }
481 
484  const mp_integer &_fraction,
485  const mp_integer &_exponent)
486 {
487  NaN_flag=infinity_flag=false;
488  sign_flag=_fraction<0;
489  fraction=_fraction;
490  if(sign_flag)
492  exponent=spec.f;
493  exponent+=_exponent;
494 
495  if(_exponent<0)
496  {
497  // bring to max. precision
498  mp_integer e_power=power(2, spec.e);
499  fraction*=power(2, e_power);
500  exponent-=e_power;
501  fraction/=power(5, -_exponent);
502  }
503  else if(_exponent>0)
504  {
505  // fix base
506  fraction*=power(5, _exponent);
507  }
508 
509  align();
510 }
511 
513 {
515  exponent=spec.f;
516  fraction=i;
517  align();
518 }
519 
521 {
522  // NaN?
523  if(NaN_flag)
524  {
525  fraction=0;
526  exponent=0;
527  sign_flag=false;
528  return;
529  }
530 
531  // do sign
532  if(fraction<0)
533  {
536  }
537 
538  // zero?
539  if(fraction==0)
540  {
541  exponent=0;
542  return;
543  }
544 
545  // 'usual case'
546  mp_integer f_power=power(2, spec.f);
547  mp_integer f_power_next=power(2, spec.f+1);
548 
549  std::size_t lowPower2=fraction.floorPow2();
550  mp_integer exponent_offset=0;
551 
552  if(lowPower2<spec.f) // too small
553  {
554  exponent_offset-=(spec.f-lowPower2);
555 
556  INVARIANT(
557  fraction * power(2, (spec.f - lowPower2)) >= f_power,
558  "normalisation result must be >= lower bound");
559  INVARIANT(
560  fraction * power(2, (spec.f - lowPower2)) < f_power_next,
561  "normalisation result must be < upper bound");
562  }
563  else if(lowPower2>spec.f) // too large
564  {
565  exponent_offset+=(lowPower2-spec.f);
566 
567  INVARIANT(
568  fraction / power(2, (lowPower2 - spec.f)) >= f_power,
569  "normalisation result must be >= lower bound");
570  INVARIANT(
571  fraction / power(2, (lowPower2 - spec.f)) < f_power_next,
572  "normalisation result must be < upper bound");
573  }
574 
575  mp_integer biased_exponent=exponent+exponent_offset+spec.bias();
576 
577  // exponent too large (infinity)?
578  if(biased_exponent>=spec.max_exponent())
579  {
580  // we need to consider the rounding mode here
581  switch(rounding_mode)
582  {
583  case UNKNOWN:
584  case NONDETERMINISTIC:
585  case ROUND_TO_EVEN:
586  infinity_flag=true;
587  break;
588 
589  case ROUND_TO_MINUS_INF:
590  // the result of the rounding is never larger than the argument
591  if(sign_flag)
592  infinity_flag=true;
593  else
594  make_fltmax();
595  break;
596 
597  case ROUND_TO_PLUS_INF:
598  // the result of the rounding is never smaller than the argument
599  if(sign_flag)
600  {
601  make_fltmax();
602  sign_flag=true; // restore sign
603  }
604  else
605  infinity_flag=true;
606  break;
607 
608  case ROUND_TO_ZERO:
609  if(sign_flag)
610  {
611  make_fltmax();
612  sign_flag=true; // restore sign
613  }
614  else
615  make_fltmax(); // positive
616  break;
617  }
618 
619  return; // done
620  }
621  else if(biased_exponent<=0) // exponent too small?
622  {
623  // produce a denormal (or zero)
624  mp_integer new_exponent=mp_integer(1)-spec.bias();
625  exponent_offset=new_exponent-exponent;
626  }
627 
628  exponent+=exponent_offset;
629 
630  if(exponent_offset>0)
631  {
632  divide_and_round(fraction, power(2, exponent_offset));
633 
634  // rounding might make the fraction too big!
635  if(fraction==f_power_next)
636  {
637  fraction=f_power;
638  ++exponent;
639  }
640  }
641  else if(exponent_offset<0)
642  fraction*=power(2, -exponent_offset);
643 
644  if(fraction==0)
645  exponent=0;
646 }
647 
649  mp_integer &dividend,
650  const mp_integer &divisor)
651 {
652  const mp_integer remainder = dividend % divisor;
653  dividend /= divisor;
654 
655  if(remainder!=0)
656  {
657  switch(rounding_mode)
658  {
659  case ROUND_TO_EVEN:
660  {
661  mp_integer divisor_middle = divisor / 2;
662  if(remainder < divisor_middle)
663  {
664  // crop
665  }
666  else if(remainder > divisor_middle)
667  {
668  ++dividend;
669  }
670  else // exactly in the middle -- go to even
671  {
672  if((dividend % 2) != 0)
673  ++dividend;
674  }
675  }
676  break;
677 
678  case ROUND_TO_ZERO:
679  // this means just crop
680  break;
681 
682  case ROUND_TO_MINUS_INF:
683  if(sign_flag)
684  ++dividend;
685  break;
686 
687  case ROUND_TO_PLUS_INF:
688  if(!sign_flag)
689  ++dividend;
690  break;
691 
692  default:
693  UNREACHABLE;
694  }
695  }
696 }
697 
699 {
701 }
702 
704 {
705  PRECONDITION(other.spec.f == spec.f);
706 
707  // NaN/x = NaN
708  if(NaN_flag)
709  return *this;
710 
711  // x/NaN = NaN
712  if(other.NaN_flag)
713  {
714  make_NaN();
715  return *this;
716  }
717 
718  // 0/0 = NaN
719  if(is_zero() && other.is_zero())
720  {
721  make_NaN();
722  return *this;
723  }
724 
725  // x/0 = +-inf
726  if(other.is_zero())
727  {
728  infinity_flag=true;
729  if(other.sign_flag)
730  negate();
731  return *this;
732  }
733 
734  // x/inf = NaN
735  if(other.infinity_flag)
736  {
737  if(infinity_flag)
738  {
739  make_NaN();
740  return *this;
741  }
742 
743  bool old_sign=sign_flag;
744  make_zero();
745  sign_flag=old_sign;
746 
747  if(other.sign_flag)
748  negate();
749 
750  return *this;
751  } // inf/x = inf
752  else if(infinity_flag)
753  {
754  if(other.sign_flag)
755  negate();
756 
757  return *this;
758  }
759 
760  exponent-=other.exponent;
761  fraction*=power(2, spec.f);
762 
763  // to account for error
764  fraction*=power(2, spec.f);
765  exponent-=spec.f;
766 
767  fraction/=other.fraction;
768 
769  if(other.sign_flag)
770  negate();
771 
772  align();
773 
774  return *this;
775 }
776 
778 {
779  PRECONDITION(other.spec.f == spec.f);
780 
781  if(other.NaN_flag)
782  make_NaN();
783  if(NaN_flag)
784  return *this;
785 
786  if(infinity_flag || other.infinity_flag)
787  {
788  if(is_zero() || other.is_zero())
789  {
790  // special case Inf * 0 is NaN
791  make_NaN();
792  return *this;
793  }
794 
795  if(other.sign_flag)
796  negate();
797  infinity_flag=true;
798  return *this;
799  }
800 
801  exponent+=other.exponent;
802  exponent-=spec.f;
803  fraction*=other.fraction;
804 
805  if(other.sign_flag)
806  negate();
807 
808  align();
809 
810  return *this;
811 }
812 
814 {
815  PRECONDITION(other.spec == spec);
816  ieee_floatt _other=other;
817 
818  if(other.NaN_flag)
819  make_NaN();
820  if(NaN_flag)
821  return *this;
822 
823  if(infinity_flag && other.infinity_flag)
824  {
825  if(sign_flag==other.sign_flag)
826  return *this;
827  make_NaN();
828  return *this;
829  }
830  else if(infinity_flag)
831  return *this;
832  else if(other.infinity_flag)
833  {
834  infinity_flag=true;
835  sign_flag=other.sign_flag;
836  return *this;
837  }
838 
839  // 0 + 0 needs special treatment for the signs
840  if(is_zero() && other.is_zero())
841  {
842  if(get_sign()==other.get_sign())
843  return *this;
844  else
845  {
847  {
848  set_sign(true);
849  return *this;
850  }
851  else
852  {
853  set_sign(false);
854  return *this;
855  }
856  }
857  }
858 
859  // get smaller exponent
860  if(_other.exponent<exponent)
861  {
862  fraction*=power(2, exponent-_other.exponent);
863  exponent=_other.exponent;
864  }
865  else if(exponent<_other.exponent)
866  {
867  _other.fraction*=power(2, _other.exponent-exponent);
868  _other.exponent=exponent;
869  }
870 
871  INVARIANT(
872  exponent == _other.exponent,
873  "prior block equalises the exponents by setting both to the "
874  "minimum of their previous values while adjusting the mantissa");
875 
876  if(sign_flag)
877  fraction.negate();
878  if(_other.sign_flag)
879  _other.fraction.negate();
880 
881  fraction+=_other.fraction;
882 
883  // if the result is zero,
884  // there is some set of rules to get the sign
885  if(fraction==0)
886  {
888  sign_flag=true;
889  else
890  sign_flag=false;
891  }
892  else // fraction!=0
893  {
894  sign_flag=(fraction<0);
895  if(sign_flag)
896  fraction.negate();
897  }
898 
899  align();
900 
901  return *this;
902 }
903 
905 {
906  ieee_floatt _other=other;
907  _other.sign_flag=!_other.sign_flag;
908  return (*this)+=_other;
909 }
910 
911 bool ieee_floatt::operator<(const ieee_floatt &other) const
912 {
913  if(NaN_flag || other.NaN_flag)
914  return false;
915 
916  // check both zero?
917  if(is_zero() && other.is_zero())
918  return false;
919 
920  // one of them zero?
921  if(is_zero())
922  return !other.sign_flag;
923  else if(other.is_zero())
924  return sign_flag;
925 
926  // check sign
927  if(sign_flag!=other.sign_flag)
928  return sign_flag;
929 
930  // handle infinity
931  if(infinity_flag)
932  {
933  if(other.infinity_flag)
934  return false;
935  else
936  return sign_flag;
937  }
938  else if(other.infinity_flag)
939  return !sign_flag;
940 
941  // check exponent
942  if(exponent!=other.exponent)
943  {
944  if(sign_flag) // both negative
945  return exponent>other.exponent;
946  else
947  return exponent<other.exponent;
948  }
949 
950  // check significand
951  if(sign_flag) // both negative
952  return fraction>other.fraction;
953  else
954  return fraction<other.fraction;
955 }
956 
957 bool ieee_floatt::operator<=(const ieee_floatt &other) const
958 {
959  if(NaN_flag || other.NaN_flag)
960  return false;
961 
962  // check zero
963  if(is_zero() && other.is_zero())
964  return true;
965 
966  // handle infinity
967  if(infinity_flag && other.infinity_flag &&
968  sign_flag==other.sign_flag)
969  return true;
970 
971  if(!infinity_flag && !other.infinity_flag &&
972  sign_flag==other.sign_flag &&
973  exponent==other.exponent &&
974  fraction==other.fraction)
975  return true;
976 
977  return *this<other;
978 }
979 
980 bool ieee_floatt::operator>(const ieee_floatt &other) const
981 {
982  return other<*this;
983 }
984 
985 bool ieee_floatt::operator>=(const ieee_floatt &other) const
986 {
987  return other<=*this;
988 }
989 
990 bool ieee_floatt::operator==(const ieee_floatt &other) const
991 {
992  // packed equality!
993  if(NaN_flag && other.NaN_flag)
994  return true;
995  else if(NaN_flag || other.NaN_flag)
996  return false;
997 
998  if(infinity_flag && other.infinity_flag &&
999  sign_flag==other.sign_flag)
1000  return true;
1001  else if(infinity_flag || other.infinity_flag)
1002  return false;
1003 
1004  // if(a.is_zero() && b.is_zero()) return true;
1005 
1006  return
1007  exponent==other.exponent &&
1008  fraction==other.fraction &&
1009  sign_flag==other.sign_flag;
1010 }
1011 
1012 bool ieee_floatt::ieee_equal(const ieee_floatt &other) const
1013 {
1014  if(NaN_flag || other.NaN_flag)
1015  return false;
1016  if(is_zero() && other.is_zero())
1017  return true;
1018  PRECONDITION(spec == other.spec);
1019  return *this==other;
1020 }
1021 
1022 bool ieee_floatt::operator==(int i) const
1023 {
1024  ieee_floatt other(spec);
1025  other.from_integer(i);
1026  return *this==other;
1027 }
1028 
1029 bool ieee_floatt::operator!=(const ieee_floatt &other) const
1030 {
1031  return !(*this==other);
1032 }
1033 
1035 {
1036  if(NaN_flag || other.NaN_flag)
1037  return true; // !!!
1038  if(is_zero() && other.is_zero())
1039  return false;
1040  PRECONDITION(spec == other.spec);
1041  return *this!=other;
1042 }
1043 
1045 {
1046  mp_integer _exponent=exponent-spec.f;
1047  mp_integer _fraction=fraction;
1048 
1049  if(sign_flag)
1050  _fraction.negate();
1051 
1052  spec=dest_spec;
1053 
1054  if(_fraction==0)
1055  {
1056  // We have a zero. It stays a zero.
1057  // Don't call build to preserve sign.
1058  }
1059  else
1060  build(_fraction, _exponent);
1061 }
1062 
1064 {
1066  unpack(bvrep2integer(expr.get_value(), spec.width(), false));
1067 }
1068 
1070 {
1071  if(NaN_flag || infinity_flag || is_zero())
1072  return 0;
1073 
1074  mp_integer result=fraction;
1075 
1076  mp_integer new_exponent=exponent-spec.f;
1077 
1078  // if the exponent is negative, divide
1079  if(new_exponent<0)
1080  result/=power(2, -new_exponent);
1081  else
1082  result*=power(2, new_exponent); // otherwise, multiply
1083 
1084  if(sign_flag)
1085  result.negate();
1086 
1087  return result;
1088 }
1089 
1090 void ieee_floatt::from_double(const double d)
1091 {
1092  spec.f=52;
1093  spec.e=11;
1094  DATA_INVARIANT(spec.width() == 64, "width must be 64 bits");
1095 
1096  static_assert(
1097  std::numeric_limits<double>::is_iec559,
1098  "this requires the double layout is according to the ieee754 "
1099  "standard");
1100  static_assert(
1101  sizeof(double) == sizeof(std::uint64_t), "ensure double has 64 bit width");
1102 
1103  union
1104  {
1105  double d;
1106  std::uint64_t i;
1107  } u;
1108 
1109  u.d=d;
1110 
1111  unpack(u.i);
1112 }
1113 
1114 void ieee_floatt::from_float(const float f)
1115 {
1116  spec.f=23;
1117  spec.e=8;
1118  DATA_INVARIANT(spec.width() == 32, "width must be 32 bits");
1119 
1120  static_assert(
1121  std::numeric_limits<float>::is_iec559,
1122  "this requires the float layout is according to the ieee754 "
1123  "standard");
1124  static_assert(
1125  sizeof(float) == sizeof(std::uint32_t), "ensure float has 32 bit width");
1126 
1127  union
1128  {
1129  float f;
1130  std::uint32_t i;
1131  } u;
1132 
1133  u.f=f;
1134 
1135  unpack(u.i);
1136 }
1137 
1139 {
1140  NaN_flag=true;
1141  // sign=false;
1142  exponent=0;
1143  fraction=0;
1144  infinity_flag=false;
1145 }
1146 
1148 {
1149  mp_integer bit_pattern=
1150  power(2, spec.e+spec.f)-1-power(2, spec.f);
1151  unpack(bit_pattern);
1152 }
1153 
1155 {
1156  unpack(power(2, spec.f));
1157 }
1158 
1160 {
1161  NaN_flag=false;
1162  sign_flag=false;
1163  exponent=0;
1164  fraction=0;
1165  infinity_flag=true;
1166 }
1167 
1169 {
1171  sign_flag=true;
1172 }
1173 
1175 {
1176  return spec.f==52 && spec.e==11;
1177 }
1178 
1180 {
1181  return spec.f==23 && spec.e==8;
1182 }
1183 
1187 {
1189  union { double f; uint64_t i; } a;
1190 
1191  if(infinity_flag)
1192  {
1193  if(sign_flag)
1194  return -std::numeric_limits<double>::infinity();
1195  else
1196  return std::numeric_limits<double>::infinity();
1197  }
1198 
1199  if(NaN_flag)
1200  {
1201  if(sign_flag)
1202  return -std::numeric_limits<double>::quiet_NaN();
1203  else
1204  return std::numeric_limits<double>::quiet_NaN();
1205  }
1206 
1207  mp_integer i=pack();
1208  CHECK_RETURN(i.is_ulong());
1209  CHECK_RETURN(i <= std::numeric_limits<std::uint64_t>::max());
1210 
1211  a.i=i.to_ulong();
1212  return a.f;
1213 }
1214 
1218 {
1219  // if false - ieee_floatt::to_float not supported on this architecture
1220  static_assert(
1221  std::numeric_limits<float>::is_iec559,
1222  "this requires the float layout is according to the IEC 559/IEEE 754 "
1223  "standard");
1224  static_assert(
1225  sizeof(float) == sizeof(uint32_t), "a 32 bit float type is required");
1226 
1227  union { float f; uint32_t i; } a;
1228 
1229  if(infinity_flag)
1230  {
1231  if(sign_flag)
1232  return -std::numeric_limits<float>::infinity();
1233  else
1234  return std::numeric_limits<float>::infinity();
1235  }
1236 
1237  if(NaN_flag)
1238  {
1239  if(sign_flag)
1240  return -std::numeric_limits<float>::quiet_NaN();
1241  else
1242  return std::numeric_limits<float>::quiet_NaN();
1243  }
1244 
1245  mp_integer i=pack();
1246  CHECK_RETURN(i.is_ulong());
1247  CHECK_RETURN(i <= std::numeric_limits<std::uint32_t>::max());
1248 
1249  a.i=i.to_ulong();
1250  return a.f;
1251 }
1252 
1256 {
1257  if(is_NaN())
1258  return;
1259 
1260  bool old_sign=get_sign();
1261 
1262  if(is_zero())
1263  {
1264  unpack(1);
1265  set_sign(!greater);
1266 
1267  return;
1268  }
1269 
1270  if(is_infinity())
1271  {
1272  if(get_sign()==greater)
1273  {
1274  make_fltmax();
1275  set_sign(old_sign);
1276  }
1277  return;
1278  }
1279 
1280  bool dir;
1281  if(greater)
1282  dir=!get_sign();
1283  else
1284  dir=get_sign();
1285 
1286  set_sign(false);
1287 
1288  mp_integer old=pack();
1289  if(dir)
1290  ++old;
1291  else
1292  --old;
1293 
1294  unpack(old);
1295 
1296  // sign change impossible (zero case caught earlier)
1297  set_sign(old_sign);
1298 }
bool get_sign() const
Definition: ieee_float.h:236
std::string to_string_scientific(std::size_t precision) const
format as [-]d.ddde+-d Note that printf always produces at least two digits for the exponent.
Definition: ieee_float.cpp:229
static int8_t r
Definition: irep_hash.h:59
BigInt mp_integer
Definition: mp_arith.h:22
ieee_floatt & operator-=(const ieee_floatt &other)
Definition: ieee_float.cpp:904
bool operator>=(const ieee_floatt &other) const
Definition: ieee_float.cpp:985
bool operator!=(const ieee_floatt &other) const
const std::string integer2string(const mp_integer &n, unsigned base)
Definition: mp_arith.cpp:106
constant_exprt to_expr() const
Definition: ieee_float.cpp:698
Fixed-width bit-vector with IEEE floating-point interpretation.
Definition: std_types.h:1398
ieee_floatt & operator+=(const ieee_floatt &other)
Definition: ieee_float.cpp:813
stylet style
Definition: format_spec.h:28
void from_expr(const constant_exprt &expr)
void build(const mp_integer &exp, const mp_integer &frac)
Definition: ieee_float.cpp:469
bool operator==(const ieee_floatt &other) const
Definition: ieee_float.cpp:990
const irep_idt & get_value() const
Definition: std_expr.h:4403
void from_type(const floatbv_typet &type)
Definition: ieee_float.cpp:46
void next_representable(bool greater)
Sets *this to the next representable number closer to plus infinity (greater = true) or minus infinit...
void make_NaN()
void set_sign(bool _sign)
Definition: ieee_float.h:172
void print(std::ostream &out) const
Definition: ieee_float.cpp:61
irep_idt integer2bvrep(const mp_integer &src, std::size_t width)
convert an integer to bit-vector representation with given width This uses two's complement for negat...
typet & type()
Return the type of the expression.
Definition: expr.h:68
mp_integer max_exponent() const
Definition: ieee_float.cpp:36
bool is_NaN() const
Definition: ieee_float.h:237
A constant literal expression.
Definition: std_expr.h:4384
#define CHECK_RETURN(CONDITION)
Definition: invariant.h:470
bool get_bool(const irep_namet &name) const
Definition: irep.cpp:239
bool is_zero() const
Definition: ieee_float.h:232
mp_integer max_fraction() const
Definition: ieee_float.cpp:41
static ieee_float_spect double_precision()
Definition: ieee_float.h:80
class floatbv_typet to_type() const
Definition: ieee_float.cpp:26
void change_spec(const ieee_float_spect &dest_spec)
void extract_base2(mp_integer &_exponent, mp_integer &_fraction) const
Definition: ieee_float.cpp:409
void dot(const goto_modelt &src, std::ostream &out)
Definition: dot.cpp:352
bool is_infinity() const
Definition: ieee_float.h:238
mp_integer to_integer() const
bool ieee_equal(const ieee_floatt &other) const
void unpack(const mp_integer &i)
Definition: ieee_float.cpp:316
bool is_normal() const
Definition: ieee_float.cpp:366
API to expression classes.
double to_double() const
Note that calling from_double -> to_double can return different bit patterns for NaN.
void make_fltmin()
void from_base10(const mp_integer &exp, const mp_integer &frac)
compute f * (10^e)
Definition: ieee_float.cpp:483
bool operator<=(const ieee_floatt &other) const
Definition: ieee_float.cpp:957
#define PRECONDITION(CONDITION)
Definition: invariant.h:438
void from_float(const float f)
mp_integer fraction
Definition: ieee_float.h:311
mp_integer bvrep2integer(const irep_idt &src, std::size_t width, bool is_signed)
convert a bit-vector representation (possibly signed) to integer
bool sign_flag
Definition: ieee_float.h:309
bool infinity_flag
Definition: ieee_float.h:312
ieee_float_spect spec
Definition: ieee_float.h:134
mp_integer exponent
Definition: ieee_float.h:310
mp_integer bias() const
Definition: ieee_float.cpp:21
std::size_t get_width() const
Definition: std_types.h:1117
std::string format(const format_spect &format_spec) const
Definition: ieee_float.cpp:66
void make_minus_infinity()
void align()
Definition: ieee_float.cpp:520
std::size_t get_f() const
Definition: std_types.cpp:27
unsigned precision
Definition: format_spec.h:19
unsigned min_width
Definition: format_spec.h:18
Pre-defined types.
static mp_integer base10_digits(const mp_integer &src)
Definition: ieee_float.cpp:126
void make_plus_infinity()
const floatbv_typet & to_floatbv_type(const typet &type)
Cast a typet to a floatbv_typet.
Definition: std_types.h:1442
bool NaN_flag
Definition: ieee_float.h:312
void from_double(const double d)
mp_integer pack() const
Definition: ieee_float.cpp:371
void negate()
Definition: ieee_float.h:167
bool is_float() const
bool operator>(const ieee_floatt &other) const
Definition: ieee_float.cpp:980
#define UNREACHABLE
This should be used to mark dead code.
Definition: invariant.h:478
void make_fltmax()
ieee_floatt & operator/=(const ieee_floatt &other)
Definition: ieee_float.cpp:703
std::string to_string_decimal(std::size_t precision) const
Definition: ieee_float.cpp:135
bool is_double() const
void divide_and_round(mp_integer &dividend, const mp_integer &divisor)
Definition: ieee_float.cpp:648
void from_integer(const mp_integer &i)
Definition: ieee_float.cpp:512
void extract_base10(mp_integer &_exponent, mp_integer &_fraction) const
Definition: ieee_float.cpp:433
bool ieee_not_equal(const ieee_floatt &other) const
void set_f(std::size_t b)
Definition: std_types.h:1413
rounding_modet rounding_mode
Definition: ieee_float.h:132
std::size_t width() const
Definition: ieee_float.h:54
ieee_floatt & operator*=(const ieee_floatt &other)
Definition: ieee_float.cpp:777
std::size_t f
Definition: ieee_float.h:30
void make_zero()
Definition: ieee_float.h:175
#define DATA_INVARIANT(CONDITION, REASON)
This condition should be used to document that assumptions that are made on goto_functions,...
Definition: invariant.h:485
std::string to_ansi_c_string() const
Definition: ieee_float.h:269
float to_float() const
Note that calling from_float -> to_float can return different bit patterns for NaN.
std::size_t e
Definition: ieee_float.h:30
mp_integer power(const mp_integer &base, const mp_integer &exponent)
A multi-precision implementation of the power operator.
void set(const irep_namet &name, const irep_idt &value)
Definition: irep.h:286
bool operator<(const ieee_floatt &other) const
Definition: ieee_float.cpp:911
void set_width(std::size_t width)
Definition: std_types.h:1122