1 /* Grid class implementation: inline functions.
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 #ifndef PPL_Grid_templates_hh
25 #define PPL_Grid_templates_hh 1
26
27 #include "Grid_Generator_defs.hh"
28 #include "Grid_Generator_System_defs.hh"
29 #include "Grid_Generator_System_inlines.hh"
30 #include <algorithm>
31 #include <deque>
32
33 namespace Parma_Polyhedra_Library {
34
35 template <typename Interval>
Grid(const Box<Interval> & box,Complexity_Class)36 Grid::Grid(const Box<Interval>& box, Complexity_Class)
37 : con_sys(),
38 gen_sys() {
39 space_dim = check_space_dimension_overflow(box.space_dimension(),
40 max_space_dimension(),
41 "PPL::Grid::",
42 "Grid(box, from_bounding_box)",
43 "the space dimension of box "
44 "exceeds the maximum allowed "
45 "space dimension");
46
47 if (box.is_empty()) {
48 // Empty grid.
49 set_empty();
50 PPL_ASSERT(OK());
51 return;
52 }
53
54 if (space_dim == 0) {
55 set_zero_dim_univ();
56 }
57 else {
58 // Initialize the space dimension as indicated by the box.
59 con_sys.set_space_dimension(space_dim);
60 gen_sys.set_space_dimension(space_dim);
61 // Add congruences and generators according to `box'.
62 PPL_DIRTY_TEMP_COEFFICIENT(l_n);
63 PPL_DIRTY_TEMP_COEFFICIENT(l_d);
64 PPL_DIRTY_TEMP_COEFFICIENT(u_n);
65 PPL_DIRTY_TEMP_COEFFICIENT(u_d);
66 gen_sys.insert(grid_point());
67 for (dimension_type k = space_dim; k-- > 0; ) {
68 const Variable v_k = Variable(k);
69 bool closed = false;
70 // TODO: Consider producing the system(s) in minimized form.
71 if (box.has_lower_bound(v_k, l_n, l_d, closed)) {
72 if (box.has_upper_bound(v_k, u_n, u_d, closed)) {
73 if (l_n * u_d == u_n * l_d) {
74 // A point interval sets dimension k of every point to a
75 // single value.
76 con_sys.insert(l_d * v_k == l_n);
77
78 // This is declared here because it may be invalidated
79 // by the call to gen_sys.insert() at the end of the loop.
80 Grid_Generator& point = gen_sys.sys.rows[0];
81
82 // Scale the point to use as divisor the lcm of the
83 // divisors of the existing point and the lower bound.
84 const Coefficient& point_divisor = point.divisor();
85 gcd_assign(u_n, l_d, point_divisor);
86 // `u_n' now holds the gcd.
87 exact_div_assign(u_n, point_divisor, u_n);
88 if (l_d < 0) {
89 neg_assign(u_n);
90 }
91 // l_d * u_n == abs(l_d * (point_divisor / gcd(l_d, point_divisor)))
92 point.scale_to_divisor(l_d * u_n);
93 // Set dimension k of the point to the lower bound.
94 if (l_d < 0) {
95 neg_assign(u_n);
96 }
97 // point[k + 1] = l_n * point_divisor / gcd(l_d, point_divisor)
98 point.expr.set(Variable(k), l_n * u_n);
99 PPL_ASSERT(point.OK());
100
101 PPL_ASSERT(gen_sys.sys.OK());
102
103 continue;
104 }
105 }
106 }
107 // A universe interval allows any value in dimension k.
108 gen_sys.insert(grid_line(v_k));
109 }
110 set_congruences_up_to_date();
111 set_generators_up_to_date();
112 }
113
114 PPL_ASSERT(OK());
115 }
116
117 template <typename Partial_Function>
118 void
map_space_dimensions(const Partial_Function & pfunc)119 Grid::map_space_dimensions(const Partial_Function& pfunc) {
120 if (space_dim == 0) {
121 return;
122 }
123
124 if (pfunc.has_empty_codomain()) {
125 // All dimensions vanish: the grid becomes zero_dimensional.
126 if (marked_empty()
127 || (!generators_are_up_to_date() && !update_generators())) {
128 // Removing all dimensions from the empty grid.
129 space_dim = 0;
130 set_empty();
131 }
132 else {
133 // Removing all dimensions from a non-empty grid.
134 set_zero_dim_univ();
135 }
136
137 PPL_ASSERT(OK());
138 return;
139 }
140
141 dimension_type new_space_dimension = pfunc.max_in_codomain() + 1;
142
143 if (new_space_dimension == space_dim) {
144 // The partial function `pfunc' is indeed total and thus specifies
145 // a permutation, that is, a renaming of the dimensions. For
146 // maximum efficiency, we will simply permute the columns of the
147 // constraint system and/or the generator system.
148
149 std::vector<Variable> cycle;
150 cycle.reserve(space_dim);
151
152 // Used to mark elements as soon as they are inserted in a cycle.
153 std::deque<bool> visited(space_dim);
154
155 for (dimension_type i = space_dim; i-- > 0; ) {
156 if (!visited[i]) {
157 dimension_type j = i;
158 do {
159 visited[j] = true;
160 // The following initialization is only to make the compiler happy.
161 dimension_type k = 0;
162 if (!pfunc.maps(j, k)) {
163 throw_invalid_argument("map_space_dimensions(pfunc)",
164 " pfunc is inconsistent");
165 }
166 if (k == j) {
167 break;
168 }
169
170 cycle.push_back(Variable(j));
171 // Go along the cycle.
172 j = k;
173 } while (!visited[j]);
174
175 // End of cycle.
176
177 // Avoid calling clear_*_minimized() if cycle.size() is less than 2,
178 // to improve efficiency.
179 if (cycle.size() >= 2) {
180 // Permute all that is up-to-date.
181 if (congruences_are_up_to_date()) {
182 con_sys.permute_space_dimensions(cycle);
183 clear_congruences_minimized();
184 }
185
186 if (generators_are_up_to_date()) {
187 gen_sys.permute_space_dimensions(cycle);
188 clear_generators_minimized();
189 }
190 }
191
192 cycle.clear();
193 }
194 }
195
196 PPL_ASSERT(OK());
197 return;
198 }
199
200 // If control gets here, then `pfunc' is not a permutation and some
201 // dimensions must be projected away.
202
203 const Grid_Generator_System& old_gensys = grid_generators();
204
205 if (old_gensys.has_no_rows()) {
206 // The grid is empty.
207 Grid new_grid(new_space_dimension, EMPTY);
208 m_swap(new_grid);
209 PPL_ASSERT(OK());
210 return;
211 }
212
213 // Make a local copy of the partial function.
214 std::vector<dimension_type> pfunc_maps(space_dim, not_a_dimension());
215 for (dimension_type j = space_dim; j-- > 0; ) {
216 dimension_type pfunc_j;
217 if (pfunc.maps(j, pfunc_j)) {
218 pfunc_maps[j] = pfunc_j;
219 }
220 }
221
222 Grid_Generator_System new_gensys;
223 // Set sortedness, for the assertion met via gs::insert.
224 new_gensys.set_sorted(false);
225 // Get the divisor of the first point.
226 Grid_Generator_System::const_iterator i;
227 Grid_Generator_System::const_iterator old_gensys_end = old_gensys.end();
228 for (i = old_gensys.begin(); i != old_gensys_end; ++i) {
229 if (i->is_point()) {
230 break;
231 }
232 }
233 PPL_ASSERT(i != old_gensys_end);
234 const Coefficient& system_divisor = i->divisor();
235 for (i = old_gensys.begin(); i != old_gensys_end; ++i) {
236 const Grid_Generator& old_g = *i;
237 const Grid_Generator::expr_type old_g_e = old_g.expression();
238 Linear_Expression expr;
239 expr.set_space_dimension(new_space_dimension);
240 bool all_zeroes = true;
241 for (Grid_Generator::expr_type::const_iterator j = old_g_e.begin(),
242 j_end = old_g_e.end(); j != j_end; ++j) {
243 const dimension_type mapped_id = pfunc_maps[j.variable().id()];
244 if (mapped_id != not_a_dimension()) {
245 add_mul_assign(expr, *j, Variable(mapped_id));
246 all_zeroes = false;
247 }
248 }
249 switch (old_g.type()) {
250 case Grid_Generator::LINE:
251 if (!all_zeroes) {
252 new_gensys.insert(grid_line(expr));
253 }
254 break;
255 case Grid_Generator::PARAMETER:
256 if (!all_zeroes) {
257 new_gensys.insert(parameter(expr, system_divisor));
258 }
259 break;
260 case Grid_Generator::POINT:
261 new_gensys.insert(grid_point(expr, old_g.divisor()));
262 break;
263 }
264 }
265
266 Grid new_grid(new_gensys);
267 m_swap(new_grid);
268
269 PPL_ASSERT(OK(true));
270 }
271
272 // Needed for converting the congruence or grid_generator system
273 // to "strong minimal form".
274 template <typename M>
275 void
reduce_reduced(Swapping_Vector<typename M::row_type> & rows,const dimension_type dim,const dimension_type pivot_index,const dimension_type start,const dimension_type end,const Dimension_Kinds & sys_dim_kinds,const bool generators)276 Grid::reduce_reduced(Swapping_Vector<typename M::row_type>& rows,
277 const dimension_type dim,
278 const dimension_type pivot_index,
279 const dimension_type start,
280 const dimension_type end,
281 const Dimension_Kinds& sys_dim_kinds,
282 const bool generators) {
283 // TODO: Remove this.
284 typedef typename M::row_type M_row_type;
285
286 const M_row_type& pivot = rows[pivot_index];
287 const Coefficient& pivot_dim = pivot.expr.get(dim);
288
289 if (pivot_dim == 0) {
290 return;
291 }
292
293 PPL_DIRTY_TEMP_COEFFICIENT(pivot_dim_half);
294 pivot_dim_half = (pivot_dim + 1) / 2;
295 const Dimension_Kind row_kind = sys_dim_kinds[dim];
296 const bool row_is_line_or_equality
297 = (row_kind == (generators ? LINE : EQUALITY));
298
299 PPL_DIRTY_TEMP_COEFFICIENT(num_rows_to_subtract);
300 PPL_DIRTY_TEMP_COEFFICIENT(row_dim_remainder);
301 for (dimension_type kinds_index = dim,
302 row_index = pivot_index; row_index-- > 0; ) {
303 if (generators) {
304 --kinds_index;
305 // Move over any virtual rows.
306 while (sys_dim_kinds[kinds_index] == GEN_VIRTUAL) {
307 --kinds_index;
308 }
309 }
310 else {
311 ++kinds_index;
312 // Move over any virtual rows.
313 while (sys_dim_kinds[kinds_index] == CON_VIRTUAL) {
314 ++kinds_index;
315 }
316 }
317
318 // row_kind CONGRUENCE is included as PARAMETER
319 if (row_is_line_or_equality
320 || (row_kind == PARAMETER
321 && sys_dim_kinds[kinds_index] == PARAMETER)) {
322 M_row_type& row = rows[row_index];
323
324 const Coefficient& row_dim = row.expr.get(dim);
325 // num_rows_to_subtract may be positive or negative.
326 num_rows_to_subtract = row_dim / pivot_dim;
327
328 // Ensure that after subtracting num_rows_to_subtract * r_dim
329 // from row_dim, -pivot_dim_half < row_dim <= pivot_dim_half.
330 // E.g., if pivot[dim] = 9, then after this reduction
331 // -5 < row_dim <= 5.
332 row_dim_remainder = row_dim % pivot_dim;
333 if (row_dim_remainder < 0) {
334 if (row_dim_remainder <= -pivot_dim_half) {
335 --num_rows_to_subtract;
336 }
337 }
338 else if (row_dim_remainder > 0 && row_dim_remainder > pivot_dim_half) {
339 ++num_rows_to_subtract;
340 }
341 // Subtract num_rows_to_subtract copies of pivot from row i. Only the
342 // entries from dim need to be subtracted, as the preceding
343 // entries are all zero.
344 // If num_rows_to_subtract is negative, these copies of pivot are
345 // added to row i.
346 if (num_rows_to_subtract != 0) {
347 row.expr.linear_combine(pivot.expr,
348 Coefficient_one(), -num_rows_to_subtract,
349 start, end + 1);
350 }
351 }
352 }
353 }
354
355 } // namespace Parma_Polyhedra_Library
356
357 #endif // !defined(PPL_Grid_templates_hh)
358