1 /* IEC 559 floating point format related functions:
2 non-inline template functions.
3 Copyright (C) 2001-2010 Roberto Bagnara <bagnara@cs.unipr.it>
4 Copyright (C) 2010-2016 BUGSENG srl (http://bugseng.com)
5
6 This file is part of the Parma Polyhedra Library (PPL).
7
8 The PPL is free software; you can redistribute it and/or modify it
9 under the terms of the GNU General Public License as published by the
10 Free Software Foundation; either version 3 of the License, or (at your
11 option) any later version.
12
13 The PPL is distributed in the hope that it will be useful, but WITHOUT
14 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
15 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
16 for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with this program; if not, write to the Free Software Foundation,
20 Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02111-1307, USA.
21
22 For the most up-to-date information see the Parma Polyhedra Library
23 site: http://bugseng.com/products/ppl/ . */
24
25 #ifndef PPL_Float_templates_hh
26 #define PPL_Float_templates_hh 1
27
28 #include "Variable_defs.hh"
29 #include "Linear_Form_defs.hh"
30 #include <cmath>
31
32 namespace Parma_Polyhedra_Library {
33
34 template <typename FP_Interval_Type>
compute_absolute_error(const Floating_Point_Format analyzed_format)35 const FP_Interval_Type& compute_absolute_error(
36 const Floating_Point_Format analyzed_format) {
37 typedef typename FP_Interval_Type::boundary_type analyzer_format;
38
39 // FIXME: check if initializing caches with EMPTY is better.
40 static const FP_Interval_Type ZERO_INTERVAL = FP_Interval_Type(0);
41 // Cached results for each different analyzed format.
42 static FP_Interval_Type ieee754_half_result = ZERO_INTERVAL;
43 static FP_Interval_Type ieee754_single_result = ZERO_INTERVAL;
44 static FP_Interval_Type ieee754_double_result = ZERO_INTERVAL;
45 static FP_Interval_Type ibm_single_result = ZERO_INTERVAL;
46 static FP_Interval_Type ieee754_quad_result = ZERO_INTERVAL;
47 static FP_Interval_Type intel_double_extended_result = ZERO_INTERVAL;
48
49 FP_Interval_Type* to_compute = NULL;
50 // Get the necessary information on the analyzed's format.
51 unsigned int f_base;
52 int f_exponent_bias;
53 unsigned int f_mantissa_bits;
54 switch (analyzed_format) {
55 case IEEE754_HALF:
56 if (ieee754_half_result != ZERO_INTERVAL) {
57 return ieee754_half_result;
58 }
59 to_compute = &ieee754_half_result;
60 f_base = float_ieee754_half::BASE;
61 f_exponent_bias = float_ieee754_half::EXPONENT_BIAS;
62 f_mantissa_bits = float_ieee754_half::MANTISSA_BITS;
63 break;
64 case IEEE754_SINGLE:
65 if (ieee754_single_result != ZERO_INTERVAL) {
66 return ieee754_single_result;
67 }
68
69 to_compute = &ieee754_single_result;
70 f_base = float_ieee754_single::BASE;
71 f_exponent_bias = float_ieee754_single::EXPONENT_BIAS;
72 f_mantissa_bits = float_ieee754_single::MANTISSA_BITS;
73 break;
74 case IEEE754_DOUBLE:
75 if (ieee754_double_result != ZERO_INTERVAL) {
76 return ieee754_double_result;
77 }
78
79 to_compute = &ieee754_double_result;
80 f_base = float_ieee754_double::BASE;
81 f_exponent_bias = float_ieee754_double::EXPONENT_BIAS;
82 f_mantissa_bits = float_ieee754_double::MANTISSA_BITS;
83 break;
84 case IBM_SINGLE:
85 if (ibm_single_result != ZERO_INTERVAL) {
86 return ibm_single_result;
87 }
88
89 to_compute = &ibm_single_result;
90 f_base = float_ibm_single::BASE;
91 f_exponent_bias = float_ibm_single::EXPONENT_BIAS;
92 f_mantissa_bits = float_ibm_single::MANTISSA_BITS;
93 break;
94 case IEEE754_QUAD:
95 if (ieee754_quad_result != ZERO_INTERVAL) {
96 return ieee754_quad_result;
97 }
98
99 to_compute = &ieee754_quad_result;
100 f_base = float_ieee754_quad::BASE;
101 f_exponent_bias = float_ieee754_quad::EXPONENT_BIAS;
102 f_mantissa_bits = float_ieee754_quad::MANTISSA_BITS;
103 break;
104 case INTEL_DOUBLE_EXTENDED:
105 if (intel_double_extended_result != ZERO_INTERVAL) {
106 return intel_double_extended_result;
107 }
108
109 to_compute = &intel_double_extended_result;
110 f_base = float_intel_double_extended::BASE;
111 f_exponent_bias = float_intel_double_extended::EXPONENT_BIAS;
112 f_mantissa_bits = float_intel_double_extended::MANTISSA_BITS;
113 break;
114 default:
115 PPL_UNREACHABLE;
116 break;
117 }
118
119 PPL_ASSERT(to_compute != NULL);
120
121 // We assume that f_base is a power of 2.
122 analyzer_format omega;
123 int power = static_cast<int>(msb_position(f_base))
124 * ((1 - f_exponent_bias) - static_cast<int>(f_mantissa_bits));
125 omega = std::max(static_cast<analyzer_format>(ldexp(1.0, power)),
126 std::numeric_limits<analyzer_format>::denorm_min());
127
128 to_compute->build(i_constraint(GREATER_OR_EQUAL, -omega),
129 i_constraint(LESS_OR_EQUAL, omega));
130 return *to_compute;
131 }
132
133 template <typename FP_Interval_Type>
134 void
discard_occurrences(std::map<dimension_type,Linear_Form<FP_Interval_Type>> & lf_store,Variable var)135 discard_occurrences(std::map<dimension_type,
136 Linear_Form<FP_Interval_Type> >& lf_store,
137 Variable var) {
138 typedef Linear_Form<FP_Interval_Type> FP_Linear_Form;
139 typedef typename std::map<dimension_type, FP_Linear_Form>::iterator Iter;
140 for (Iter i = lf_store.begin(); i != lf_store.end(); ) {
141 if ((i->second).coefficient(var) != 0) {
142 i = lf_store.erase(i);
143 }
144 else {
145 ++i;
146 }
147 }
148 }
149
150 /* FIXME: improve efficiency by adding the list of potentially conflicting
151 variables as an argument. */
152 template <typename FP_Interval_Type>
upper_bound_assign(std::map<dimension_type,Linear_Form<FP_Interval_Type>> & ls1,const std::map<dimension_type,Linear_Form<FP_Interval_Type>> & ls2)153 void upper_bound_assign(std::map<dimension_type,
154 Linear_Form<FP_Interval_Type> >& ls1,
155 const std::map<dimension_type,
156 Linear_Form<FP_Interval_Type> >& ls2) {
157 typedef Linear_Form<FP_Interval_Type> FP_Linear_Form;
158 typedef typename std::map<dimension_type, FP_Linear_Form>::iterator Iter;
159 typedef typename std::map<dimension_type,
160 FP_Linear_Form>::const_iterator Const_Iter;
161
162 Const_Iter i2_end = ls2.end();
163 for (Iter i1 = ls1.begin(), i1_end = ls1.end(); i1 != i1_end; ) {
164 Const_Iter i2 = ls2.find(i1->first);
165 if ((i2 == i2_end) || (i1->second != i2->second)) {
166 i1 = ls1.erase(i1);
167 }
168 else {
169 ++i1;
170 }
171 }
172 }
173
174 } // namespace Parma_Polyhedra_Library
175
176 #endif // !defined(PPL_Float_templates_hh)
177