1 /* Grid class implementation (non-inline public 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 #include "ppl-config.h"
25 #include "Grid_defs.hh"
26 #include "Topology_types.hh"
27 #include "Scalar_Products_defs.hh"
28 #include "Scalar_Products_inlines.hh"
29 #include "Polyhedron_defs.hh"
30 #include "assertions.hh"
31 #include <iostream>
32
33 namespace PPL = Parma_Polyhedra_Library;
34
35 // TODO: In the Grid constructors adapt and use the given system if it
36 // is modifiable, instead of using a copy.
37
Grid(const Grid & y,Complexity_Class)38 PPL::Grid::Grid(const Grid& y, Complexity_Class)
39 : con_sys(),
40 gen_sys(),
41 status(y.status),
42 space_dim(y.space_dim),
43 dim_kinds(y.dim_kinds) {
44 if (space_dim == 0) {
45 con_sys = y.con_sys;
46 gen_sys = y.gen_sys;
47 }
48 else {
49 if (y.congruences_are_up_to_date()) {
50 con_sys = y.con_sys;
51 }
52 else {
53 con_sys.set_space_dimension(space_dim);
54 }
55 if (y.generators_are_up_to_date()) {
56 gen_sys = y.gen_sys;
57 }
58 else {
59 gen_sys = Grid_Generator_System(y.space_dim);
60 }
61 }
62 }
63
Grid(const Constraint_System & cs)64 PPL::Grid::Grid(const Constraint_System& cs)
65 : con_sys(check_space_dimension_overflow(cs.space_dimension(),
66 max_space_dimension(),
67 "PPL::Grid::",
68 "Grid(cs)",
69 "the space dimension of cs "
70 "exceeds the maximum allowed "
71 "space dimension")),
72 gen_sys(cs.space_dimension()) {
73 space_dim = cs.space_dimension();
74
75 if (space_dim == 0) {
76 // See if an inconsistent constraint has been passed.
77 for (Constraint_System::const_iterator i = cs.begin(),
78 cs_end = cs.end(); i != cs_end; ++i) {
79 if (i->is_inconsistent()) {
80 // Inconsistent constraint found: the grid is empty.
81 status.set_empty();
82 // Insert the zero dim false congruence system into `con_sys'.
83 // `gen_sys' is already in empty form.
84 con_sys.insert(Congruence::zero_dim_false());
85 PPL_ASSERT(OK());
86 return;
87 }
88 }
89 set_zero_dim_univ();
90 PPL_ASSERT(OK());
91 return;
92 }
93
94 Congruence_System cgs(cs.space_dimension());
95 for (Constraint_System::const_iterator i = cs.begin(),
96 cs_end = cs.end(); i != cs_end; ++i) {
97 if (i->is_equality()) {
98 cgs.insert(*i);
99 }
100 else {
101 throw_invalid_constraints("Grid(cs)", "cs");
102 }
103 }
104 construct(cgs);
105 }
106
Grid(Constraint_System & cs,Recycle_Input)107 PPL::Grid::Grid(Constraint_System& cs, Recycle_Input)
108 : con_sys(check_space_dimension_overflow(cs.space_dimension(),
109 max_space_dimension(),
110 "PPL::Grid::",
111 "Grid(cs, recycle)",
112 "the space dimension of cs "
113 "exceeds the maximum allowed "
114 "space dimension")),
115 gen_sys(cs.space_dimension()) {
116 space_dim = cs.space_dimension();
117
118 if (space_dim == 0) {
119 // See if an inconsistent constraint has been passed.
120 for (Constraint_System::const_iterator i = cs.begin(),
121 cs_end = cs.end(); i != cs_end; ++i) {
122 if (i->is_inconsistent()) {
123 // Inconsistent constraint found: the grid is empty.
124 status.set_empty();
125 // Insert the zero dim false congruence system into `con_sys'.
126 // `gen_sys' is already in empty form.
127 con_sys.insert(Congruence::zero_dim_false());
128 PPL_ASSERT(OK());
129 return;
130 }
131 }
132 set_zero_dim_univ();
133 PPL_ASSERT(OK());
134 return;
135 }
136
137 Congruence_System cgs(space_dim);
138 for (Constraint_System::const_iterator i = cs.begin(),
139 cs_end = cs.end(); i != cs_end; ++i) {
140 if (i->is_equality()) {
141 cgs.insert(*i);
142 }
143 else {
144 throw_invalid_constraint("Grid(cs)", "cs");
145 }
146 }
147 construct(cgs);
148 }
149
Grid(const Polyhedron & ph,Complexity_Class complexity)150 PPL::Grid::Grid(const Polyhedron& ph,
151 Complexity_Class complexity)
152 : con_sys(check_space_dimension_overflow(ph.space_dimension(),
153 max_space_dimension(),
154 "PPL::Grid::",
155 "Grid(ph)",
156 "the space dimension of ph "
157 "exceeds the maximum allowed "
158 "space dimension")),
159 gen_sys(ph.space_dimension()) {
160 space_dim = ph.space_dimension();
161
162 // A zero-dim polyhedron causes no complexity problems.
163 if (space_dim == 0) {
164 if (ph.is_empty()) {
165 set_empty();
166 }
167 else {
168 set_zero_dim_univ();
169 }
170 return;
171 }
172
173 // A polyhedron known to be empty causes no complexity problems.
174 if (ph.marked_empty()) {
175 set_empty();
176 return;
177 }
178
179 const bool use_constraints = ph.constraints_are_minimized()
180 || !ph.generators_are_up_to_date();
181
182 // Minimize the constraint description if it is needed and
183 // the complexity allows it.
184 if (use_constraints && complexity == ANY_COMPLEXITY) {
185 if (!ph.minimize()) {
186 set_empty();
187 return;
188 }
189 }
190 if (use_constraints) {
191 // Only the equality constraints need be used.
192 PPL_ASSERT(ph.constraints_are_up_to_date());
193 const Constraint_System& cs = ph.constraints();
194 Congruence_System cgs(space_dim);
195 for (Constraint_System::const_iterator i = cs.begin(),
196 cs_end = cs.end(); i != cs_end; ++i) {
197 if (i->is_equality()) {
198 cgs.insert(*i);
199 }
200 }
201 construct(cgs);
202 }
203 else {
204 // First find a point or closure point and convert it to a
205 // grid point and add to the (initially empty) set of grid generators.
206 PPL_ASSERT(ph.generators_are_up_to_date());
207 const Generator_System& gs = ph.generators();
208 Grid_Generator_System ggs(space_dim);
209 Linear_Expression point_expr;
210 point_expr.set_space_dimension(space_dim);
211 PPL_DIRTY_TEMP_COEFFICIENT(point_divisor);
212 for (Generator_System::const_iterator g = gs.begin(),
213 gs_end = gs.end(); g != gs_end; ++g) {
214 if (g->is_point() || g->is_closure_point()) {
215 point_expr.linear_combine(g->expr, Coefficient_one(), Coefficient_one(),
216 1, space_dim + 1);
217 point_divisor = g->divisor();
218 ggs.insert(grid_point(point_expr, point_divisor));
219 break;
220 }
221 }
222 // Add grid lines for all the other generators.
223 // If the polyhedron's generator is a (closure) point, the grid line must
224 // have the direction given by a line that joins the grid point already
225 // inserted and the new point.
226 for (Generator_System::const_iterator g = gs.begin(),
227 gs_end = gs.end(); g != gs_end; ++g) {
228 Linear_Expression e;
229 e.set_space_dimension(space_dim);
230 if (g->is_point() || g->is_closure_point()) {
231 e.linear_combine(point_expr, Coefficient_one(), g->divisor(),
232 1, space_dim + 1);
233 e.linear_combine(g->expr, Coefficient_one(), -point_divisor,
234 1, space_dim + 1);
235 if (e.all_homogeneous_terms_are_zero()) {
236 continue;
237 }
238 }
239 else {
240 e.linear_combine(g->expr, Coefficient_one(), Coefficient_one(),
241 1, space_dim + 1);
242 }
243 ggs.insert(grid_line(e));
244 }
245 construct(ggs);
246 }
247 PPL_ASSERT(OK());
248 }
249
250 PPL::Grid&
operator =(const Grid & y)251 PPL::Grid::operator=(const Grid& y) {
252 space_dim = y.space_dim;
253 dim_kinds = y.dim_kinds;
254 if (y.marked_empty()) {
255 set_empty();
256 }
257 else if (space_dim == 0) {
258 set_zero_dim_univ();
259 }
260 else {
261 status = y.status;
262 if (y.congruences_are_up_to_date()) {
263 con_sys = y.con_sys;
264 }
265 if (y.generators_are_up_to_date()) {
266 gen_sys = y.gen_sys;
267 }
268 }
269 return *this;
270 }
271
272 PPL::dimension_type
affine_dimension() const273 PPL::Grid::affine_dimension() const {
274 if (space_dim == 0 || is_empty()) {
275 return 0;
276 }
277
278 if (generators_are_up_to_date()) {
279 if (generators_are_minimized()) {
280 return gen_sys.num_rows() - 1;
281 }
282 if (!(congruences_are_up_to_date() && congruences_are_minimized())) {
283 return minimized_grid_generators().num_rows() - 1;
284 }
285 }
286 else {
287 minimized_congruences();
288 }
289 PPL_ASSERT(congruences_are_minimized());
290 dimension_type d = space_dim;
291 for (dimension_type i = con_sys.num_rows(); i-- > 0; ) {
292 if (con_sys[i].is_equality()) {
293 --d;
294 }
295 }
296 return d;
297 }
298
299 const PPL::Congruence_System&
congruences() const300 PPL::Grid::congruences() const {
301 if (marked_empty()) {
302 return con_sys;
303 }
304
305 if (space_dim == 0) {
306 // Zero-dimensional universe.
307 PPL_ASSERT(con_sys.num_rows() == 0 && con_sys.space_dimension() == 0);
308 return con_sys;
309 }
310
311 if (!congruences_are_up_to_date()) {
312 update_congruences();
313 }
314
315 return con_sys;
316 }
317
318 const PPL::Congruence_System&
minimized_congruences() const319 PPL::Grid::minimized_congruences() const {
320 if (congruences_are_up_to_date() && !congruences_are_minimized()) {
321 // Minimize the congruences.
322 Grid& gr = const_cast<Grid&>(*this);
323 if (gr.simplify(gr.con_sys, gr.dim_kinds)) {
324 gr.set_empty();
325 }
326 else {
327 gr.set_congruences_minimized();
328 }
329 }
330 return congruences();
331 }
332
333 const PPL::Grid_Generator_System&
grid_generators() const334 PPL::Grid::grid_generators() const {
335 if (space_dim == 0) {
336 PPL_ASSERT(gen_sys.space_dimension() == 0
337 && gen_sys.num_rows() == (marked_empty() ? 0U : 1U));
338 return gen_sys;
339 }
340
341 if (marked_empty()) {
342 PPL_ASSERT(gen_sys.has_no_rows());
343 return gen_sys;
344 }
345
346 if (!generators_are_up_to_date() && !update_generators()) {
347 // Updating found the grid empty.
348 const_cast<Grid&>(*this).set_empty();
349 return gen_sys;
350 }
351
352 return gen_sys;
353 }
354
355 const PPL::Grid_Generator_System&
minimized_grid_generators() const356 PPL::Grid::minimized_grid_generators() const {
357 if (space_dim == 0) {
358 PPL_ASSERT(gen_sys.space_dimension() == 0
359 && gen_sys.num_rows() == (marked_empty() ? 0U : 1U));
360 return gen_sys;
361 }
362
363 if (marked_empty()) {
364 PPL_ASSERT(gen_sys.has_no_rows());
365 return gen_sys;
366 }
367
368 if (generators_are_up_to_date()) {
369 if (!generators_are_minimized()) {
370 // Minimize the generators.
371 Grid& gr = const_cast<Grid&>(*this);
372 gr.simplify(gr.gen_sys, gr.dim_kinds);
373 gr.set_generators_minimized();
374 }
375 }
376 else if (!update_generators()) {
377 // Updating found the grid empty.
378 const_cast<Grid&>(*this).set_empty();
379 return gen_sys;
380 }
381
382 return gen_sys;
383 }
384
385 PPL::Poly_Con_Relation
relation_with(const Congruence & cg) const386 PPL::Grid::relation_with(const Congruence& cg) const {
387 // Dimension-compatibility check.
388 if (space_dim < cg.space_dimension()) {
389 throw_dimension_incompatible("relation_with(cg)", "cg", cg);
390 }
391
392 if (marked_empty()) {
393 return Poly_Con_Relation::saturates()
394 && Poly_Con_Relation::is_included()
395 && Poly_Con_Relation::is_disjoint();
396 }
397
398 if (space_dim == 0) {
399 if (cg.is_inconsistent()) {
400 return Poly_Con_Relation::is_disjoint();
401 }
402 else if (cg.is_equality()) {
403 return Poly_Con_Relation::saturates()
404 && Poly_Con_Relation::is_included();
405 }
406 else if (cg.inhomogeneous_term() % cg.modulus() == 0) {
407 return Poly_Con_Relation::saturates()
408 && Poly_Con_Relation::is_included();
409 }
410 }
411
412 if (!generators_are_up_to_date() && !update_generators()) {
413 // Updating found the grid empty.
414 return Poly_Con_Relation::saturates()
415 && Poly_Con_Relation::is_included()
416 && Poly_Con_Relation::is_disjoint();
417 }
418
419 // Return one of the relations
420 // 'strictly_intersects' a strict subset of the grid points satisfy cg
421 // 'is_included' every grid point satisfies cg
422 // 'is_disjoint' cg and the grid occupy separate spaces.
423 // There is always a point.
424 // Scalar product of the congruence and the first point that
425 // satisfies the congruence.
426 PPL_DIRTY_TEMP_COEFFICIENT(point_sp);
427 point_sp = 0;
428
429 PPL_DIRTY_TEMP_COEFFICIENT(div);
430 div = cg.modulus();
431
432 PPL_DIRTY_TEMP_COEFFICIENT(sp);
433
434 bool known_to_intersect = false;
435
436 for (Grid_Generator_System::const_iterator i = gen_sys.begin(),
437 i_end = gen_sys.end(); i != i_end; ++i) {
438 const Grid_Generator& g = *i;
439 Scalar_Products::assign(sp, cg, g);
440
441 switch (g.type()) {
442
443 case Grid_Generator::POINT:
444 if (cg.is_proper_congruence()) {
445 sp %= div;
446 }
447 if (sp == 0) {
448 // The point satisfies the congruence.
449 if (point_sp == 0) {
450 // Any previous points satisfied the congruence.
451 known_to_intersect = true;
452 }
453 else {
454 return Poly_Con_Relation::strictly_intersects();
455 }
456 }
457 else {
458 if (point_sp == 0) {
459 if (known_to_intersect) {
460 return Poly_Con_Relation::strictly_intersects();
461 }
462 // Assign `sp' to `point_sp' as `sp' is the scalar product
463 // of cg and a point g and is non-zero.
464 point_sp = sp;
465 }
466 else {
467 // A previously considered point p failed to satisfy cg such that
468 // `point_sp' = `scalar_prod(p, cg)'
469 // so, if we consider the parameter g-p instead of g, we have
470 // scalar_prod(g-p, cg) = scalar_prod(g, cg) - scalar_prod(p, cg)
471 // = sp - point_sp
472 sp -= point_sp;
473
474 if (sp != 0) {
475 // Find the GCD between sp and the previous GCD.
476 gcd_assign(div, div, sp);
477 if (point_sp % div == 0) {
478 // There is a point in the grid satisfying cg.
479 return Poly_Con_Relation::strictly_intersects();
480 }
481 }
482 }
483 }
484 break;
485
486 case Grid_Generator::PARAMETER:
487 if (cg.is_proper_congruence()) {
488 sp %= (div * g.divisor());
489 }
490 if (sp == 0) {
491 // Parameter g satisfies the cg so the relation depends
492 // entirely on the other generators.
493 break;
494 }
495 if (known_to_intersect) {
496 // At least one point satisfies cg. However, the sum of such
497 // a point and the parameter g fails to satisfy cg (due to g).
498 return Poly_Con_Relation::strictly_intersects();
499 }
500 // Find the GCD between sp and the previous GCD.
501 gcd_assign(div, div, sp);
502 if (point_sp != 0) {
503 // At least one of any previously encountered points fails to
504 // satisfy cg.
505 if (point_sp % div == 0) {
506 // There is also a grid point that satisfies cg.
507 return Poly_Con_Relation::strictly_intersects();
508 }
509 }
510 break;
511
512 case Grid_Generator::LINE:
513 if (sp == 0) {
514 // Line g satisfies the cg so the relation depends entirely on
515 // the other generators.
516 break;
517 }
518
519 // Line g intersects the congruence.
520 //
521 // There is a point p in the grid. Suppose <p*cg> = p_sp. Then
522 // (-p_sp/sp)*g + p is a point that satisfies cg: <((-p_sp/sp)*g
523 // + p).cg> = -(p_sp/sp)*sp + p_sp) = 0. If p does not satisfy
524 // `cg' and hence is not in the grid defined by `cg', the grid
525 // `*this' strictly intersects the `cg' grid. On the other
526 // hand, if `p' is in the grid defined by `cg' so that p_sp = 0,
527 // then <p+g.cg> = p_sp + sp != 0; thus `p+g' is a point in
528 // *this that does not satisfy `cg' and hence `p+g' is a point
529 // in *this not in the grid defined by `cg'; therefore `*this'
530 // strictly intersects the `cg' grid.
531 return Poly_Con_Relation::strictly_intersects();
532 }
533 }
534
535 if (point_sp == 0) {
536 if (cg.is_equality()) {
537 // Every generator satisfied the cg.
538 return Poly_Con_Relation::is_included()
539 && Poly_Con_Relation::saturates();
540 }
541 else {
542 // Every generator satisfied the cg.
543 return Poly_Con_Relation::is_included();
544 }
545 }
546
547 PPL_ASSERT(!known_to_intersect);
548 return Poly_Con_Relation::is_disjoint();
549 }
550
551 PPL::Poly_Gen_Relation
relation_with(const Grid_Generator & g) const552 PPL::Grid::relation_with(const Grid_Generator& g) const {
553 // Dimension-compatibility check.
554 if (space_dim < g.space_dimension()) {
555 throw_dimension_incompatible("relation_with(g)", "g", g);
556 }
557
558 // The empty grid cannot subsume a generator.
559 if (marked_empty()) {
560 return Poly_Gen_Relation::nothing();
561 }
562
563 // A universe grid in a zero-dimensional space subsumes all the
564 // generators of a zero-dimensional space.
565 if (space_dim == 0) {
566 return Poly_Gen_Relation::subsumes();
567 }
568
569 if (!congruences_are_up_to_date()) {
570 update_congruences();
571 }
572
573 return
574 con_sys.satisfies_all_congruences(g)
575 ? Poly_Gen_Relation::subsumes()
576 : Poly_Gen_Relation::nothing();
577 }
578
579 PPL::Poly_Gen_Relation
relation_with(const Generator & g) const580 PPL::Grid::relation_with(const Generator& g) const {
581 const dimension_type g_space_dim = g.space_dimension();
582
583 // Dimension-compatibility check.
584 if (space_dim < g_space_dim) {
585 throw_dimension_incompatible("relation_with(g)", "g", g);
586 }
587
588 // The empty grid cannot subsume a generator.
589 if (marked_empty()) {
590 return Poly_Gen_Relation::nothing();
591 }
592
593 // A universe grid in a zero-dimensional space subsumes all the
594 // generators of a zero-dimensional space.
595 if (space_dim == 0) {
596 return Poly_Gen_Relation::subsumes();
597 }
598
599 if (!congruences_are_up_to_date()) {
600 update_congruences();
601 }
602
603 const Linear_Expression expr(g.expression());
604 Grid_Generator gg(grid_point());
605 if (g.is_point() || g.is_closure_point()) {
606 // Points and closure points are converted to grid points.
607 gg = grid_point(expr, g.divisor());
608 }
609 else {
610 // The generator is a ray or line.
611 // In both cases, we convert it to a grid line
612 gg = grid_line(expr);
613 }
614
615 return
616 con_sys.satisfies_all_congruences(gg)
617 ? Poly_Gen_Relation::subsumes()
618 : Poly_Gen_Relation::nothing();
619 }
620
621 PPL::Poly_Con_Relation
relation_with(const Constraint & c) const622 PPL::Grid::relation_with(const Constraint& c) const {
623 // Dimension-compatibility check.
624 if (space_dim < c.space_dimension()) {
625 throw_dimension_incompatible("relation_with(c)", "c", c);
626 }
627
628 if (c.is_equality()) {
629 const Congruence cg(c);
630 return relation_with(cg);
631 }
632
633 if (marked_empty()) {
634 return Poly_Con_Relation::saturates()
635 && Poly_Con_Relation::is_included()
636 && Poly_Con_Relation::is_disjoint();
637 }
638
639 if (space_dim == 0) {
640 if (c.is_inconsistent()) {
641 if (c.is_strict_inequality() && c.inhomogeneous_term() == 0) {
642 // The constraint 0 > 0 implicitly defines the hyperplane 0 = 0;
643 // thus, the zero-dimensional point also saturates it.
644 return Poly_Con_Relation::saturates()
645 && Poly_Con_Relation::is_disjoint();
646 }
647 else {
648 return Poly_Con_Relation::is_disjoint();
649 }
650 }
651 else if (c.inhomogeneous_term() == 0) {
652 return Poly_Con_Relation::saturates()
653 && Poly_Con_Relation::is_included();
654 }
655 else {
656 // The zero-dimensional point saturates
657 // neither the positivity constraint 1 >= 0,
658 // nor the strict positivity constraint 1 > 0.
659 return Poly_Con_Relation::is_included();
660 }
661 }
662
663 if (!generators_are_up_to_date() && !update_generators()) {
664 // Updating found the grid empty.
665 return Poly_Con_Relation::saturates()
666 && Poly_Con_Relation::is_included()
667 && Poly_Con_Relation::is_disjoint();
668 }
669
670 // Return one of the relations
671 // 'strictly_intersects' a strict subset of the grid points satisfy c
672 // 'is_included' every grid point satisfies c
673 // 'is_disjoint' c and the grid occupy separate spaces.
674 // There is always a point.
675
676 bool point_is_included = false;
677 bool point_saturates = false;
678 const Grid_Generator* first_point = 0;
679
680 for (Grid_Generator_System::const_iterator i = gen_sys.begin(),
681 i_end = gen_sys.end(); i != i_end; ++i) {
682 const Grid_Generator& g = *i;
683 switch (g.type()) {
684 case Grid_Generator::POINT:
685 {
686 if (first_point == 0) {
687 first_point = &g;
688 const int sign = Scalar_Products::sign(c, g);
689 if (sign == 0) {
690 point_saturates = !c.is_strict_inequality();
691 }
692 else if (sign > 0) {
693 point_is_included = !c.is_equality();
694 }
695 break;
696 }
697 // Not the first point: convert `g' to be a parameter
698 // and fall through into the parameter case.
699 Grid_Generator& gen = const_cast<Grid_Generator&>(g);
700 const Grid_Generator& point = *first_point;
701 const Coefficient& p_div = point.divisor();
702 const Coefficient& g_div = gen.divisor();
703 gen.expr.linear_combine(point.expr, p_div, -g_div,
704 1, gen.expr.space_dimension());
705 gen.expr.set_inhomogeneous_term(g_div * p_div);
706 gen.strong_normalize();
707 gen.set_is_parameter();
708 PPL_ASSERT(gen.OK());
709 }
710 // Intentionally fall through.
711
712 case Grid_Generator::PARAMETER:
713 case Grid_Generator::LINE:
714 {
715 const int sign = c.is_strict_inequality()
716 ? Scalar_Products::reduced_sign(c.expr, g.expr)
717 : Scalar_Products::sign(c.expr, g.expr);
718 if (sign != 0) {
719 return Poly_Con_Relation::strictly_intersects();
720 }
721 }
722 break;
723 } // switch
724 }
725
726 // If this program point is reached, then all lines and parameters
727 // saturate the constraint. Hence, the result is determined by
728 // the previosly computed relation with the point.
729 if (point_saturates) {
730 return Poly_Con_Relation::saturates()
731 && Poly_Con_Relation::is_included();
732 }
733
734 if (point_is_included) {
735 return Poly_Con_Relation::is_included();
736 }
737
738 return Poly_Con_Relation::is_disjoint();
739 }
740
741 bool
is_empty() const742 PPL::Grid::is_empty() const {
743 if (marked_empty()) {
744 return true;
745 }
746 // Try a fast-fail test: if generators are up-to-date then the
747 // generator system (since it is well formed) contains a point.
748 if (generators_are_up_to_date()) {
749 return false;
750 }
751 if (space_dim == 0) {
752 return false;
753 }
754 if (congruences_are_minimized()) {
755 // If the grid was empty it would be marked empty.
756 return false;
757 }
758 // Minimize the congruences to check if the grid is empty.
759 Grid& gr = const_cast<Grid&>(*this);
760 if (gr.simplify(gr.con_sys, gr.dim_kinds)) {
761 gr.set_empty();
762 return true;
763 }
764 gr.set_congruences_minimized();
765 return false;
766 }
767
768 bool
is_universe() const769 PPL::Grid::is_universe() const {
770 if (marked_empty()) {
771 return false;
772 }
773
774 if (space_dim == 0) {
775 return true;
776 }
777
778 if (congruences_are_up_to_date()) {
779 if (congruences_are_minimized()) {
780 // The minimized universe congruence system has only one row,
781 // the integrality congruence.
782 return con_sys.num_rows() == 1 && con_sys[0].is_tautological();
783 }
784 }
785 else {
786 update_congruences();
787 return con_sys.num_rows() == 1 && con_sys[0].is_tautological();
788 }
789
790 // Test con_sys's inclusion in a universe generator system.
791
792 // The zero dimension cases are handled above.
793 for (dimension_type i = space_dim; i-- > 0; ) {
794 Linear_Expression expr;
795 expr.set_space_dimension(space_dim);
796 expr += Variable(i);
797 if (!con_sys.satisfies_all_congruences(grid_line(expr))) {
798 return false;
799 }
800 }
801 #ifndef NDEBUG
802 Linear_Expression expr;
803 expr.set_space_dimension(space_dim);
804 PPL_ASSERT(con_sys.satisfies_all_congruences(grid_point(expr)));
805 #endif
806 return true;
807 }
808
809 bool
is_bounded() const810 PPL::Grid::is_bounded() const {
811 // A zero-dimensional or empty grid is bounded.
812 if (space_dim == 0
813 || marked_empty()
814 || (!generators_are_up_to_date() && !update_generators())) {
815 return true;
816 }
817 // TODO: Consider using con_sys when gen_sys is out of date.
818
819 if (gen_sys.num_rows() > 1) {
820 // Check if all generators are the same point.
821 const Grid_Generator& first_point = gen_sys[0];
822 if (first_point.is_line_or_parameter()) {
823 return false;
824 }
825 for (dimension_type row = gen_sys.num_rows(); row-- > 0; ) {
826 const Grid_Generator& gen = gen_sys[row];
827 if (gen.is_line_or_parameter() || gen != first_point) {
828 return false;
829 }
830 }
831 }
832 return true;
833 }
834
835 bool
is_discrete() const836 PPL::Grid::is_discrete() const {
837 // A zero-dimensional or empty grid is discrete.
838 if (space_dim == 0
839 || marked_empty()
840 || (!generators_are_up_to_date() && !update_generators())) {
841 return true;
842 }
843 // Search for lines in the generator system.
844 for (dimension_type row = gen_sys.num_rows(); row-- > 1; ) {
845 if (gen_sys[row].is_line()) {
846 return false;
847 }
848 }
849
850 // The system of generators is composed only by
851 // points and parameters: the grid is discrete.
852 return true;
853 }
854
855 bool
is_topologically_closed() const856 PPL::Grid::is_topologically_closed() const {
857 return true;
858 }
859
860 bool
contains_integer_point() const861 PPL::Grid::contains_integer_point() const {
862 // Empty grids have no points.
863 if (marked_empty()) {
864 return false;
865 }
866
867 // A zero-dimensional, universe grid has, by convention, an
868 // integer point.
869 if (space_dim == 0) {
870 return true;
871 }
872
873 // A grid has an integer point if its intersection with the integer
874 // grid is non-empty.
875 Congruence_System cgs;
876 for (dimension_type var_index = space_dim; var_index-- > 0; ) {
877 cgs.insert(Variable(var_index) %= 0);
878 }
879
880 Grid gr = *this;
881 gr.add_recycled_congruences(cgs);
882 return !gr.is_empty();
883 }
884
885 bool
constrains(const Variable var) const886 PPL::Grid::constrains(const Variable var) const {
887 // `var' should be one of the dimensions of the grid.
888 const dimension_type var_space_dim = var.space_dimension();
889 if (space_dim < var_space_dim) {
890 throw_dimension_incompatible("constrains(v)", "v", var);
891 }
892
893 // An empty grid constrains all variables.
894 if (marked_empty()) {
895 return true;
896 }
897
898 if (generators_are_up_to_date()) {
899 // Since generators are up-to-date, the generator system (since it is
900 // well formed) contains a point. Hence the grid is not empty.
901 if (congruences_are_up_to_date()) {
902 // Here a variable is constrained if and only if it is
903 // syntactically constrained.
904 goto syntactic_check;
905 }
906
907 if (generators_are_minimized()) {
908 // Try a quick, incomplete check for the universe grid:
909 // a universe grid constrains no variable.
910 // Count the number of lines (they are linearly independent).
911 dimension_type num_lines = 0;
912 for (dimension_type i = gen_sys.num_rows(); i-- > 0; ) {
913 if (gen_sys[i].is_line()) {
914 ++num_lines;
915 }
916 }
917
918 if (num_lines == space_dim) {
919 return false;
920 }
921 }
922
923 // Scan generators: perhaps we will find line(var).
924 for (dimension_type i = gen_sys.num_rows(); i-- > 0; ) {
925 const Grid_Generator& g_i = gen_sys[i];
926 if (!g_i.is_line()) {
927 continue;
928 }
929 if (sgn(g_i.coefficient(var)) != 0) {
930 if (g_i.expression().all_zeroes(1, var.space_dimension())
931 && g_i.expression().all_zeroes(var.space_dimension() + 1, space_dim + 1)) {
932 // The only nonzero coefficient in g_i is the one of var.
933 return true;
934 }
935 }
936 }
937
938 // We are still here: at least we know that the grid is not empty.
939 update_congruences();
940 goto syntactic_check;
941 }
942
943 // We must minimize to detect emptiness and obtain constraints.
944 if (!minimize()) {
945 return true;
946 }
947
948 syntactic_check:
949 for (dimension_type i = con_sys.num_rows(); i-- > 0; ) {
950 if (con_sys[i].coefficient(var) != 0) {
951 return true;
952 }
953 }
954 return false;
955 }
956
957 bool
OK(bool check_not_empty) const958 PPL::Grid::OK(bool check_not_empty) const {
959 #ifndef NDEBUG
960 using std::endl;
961 using std::cerr;
962 #endif
963
964 // Check whether the status information is legal.
965 if (!status.OK()) {
966 goto fail;
967 }
968
969 if (marked_empty()) {
970 if (check_not_empty) {
971 // The caller does not want the grid to be empty.
972 #ifndef NDEBUG
973 cerr << "Empty grid!" << endl;
974 #endif
975 goto fail;
976 }
977
978 if (con_sys.space_dimension() != space_dim) {
979 #ifndef NDEBUG
980 cerr << "The grid is in a space of dimension " << space_dim
981 << " while the system of congruences is in a space of dimension "
982 << con_sys.space_dimension()
983 << endl;
984 #endif
985 goto fail;
986 }
987 return true;
988 }
989
990 // A zero-dimensional universe grid is legal only if the system of
991 // congruences `con_sys' is empty, and the generator system contains
992 // one point.
993 if (space_dim == 0) {
994 if (con_sys.has_no_rows()) {
995 if (gen_sys.num_rows() == 1 && gen_sys[0].is_point()) {
996 return true;
997 }
998 }
999 #ifndef NDEBUG
1000 cerr << "Zero-dimensional grid should have an empty congruence" << endl
1001 << "system and a generator system of a single point." << endl;
1002 #endif
1003 goto fail;
1004 }
1005
1006 // A grid is defined by a system of congruences or a system of
1007 // generators. At least one of them must be up to date.
1008 if (!congruences_are_up_to_date() && !generators_are_up_to_date()) {
1009 #ifndef NDEBUG
1010 cerr << "Grid not empty, not zero-dimensional" << endl
1011 << "and with neither congruences nor generators up-to-date!"
1012 << endl;
1013 #endif
1014 goto fail;
1015 }
1016
1017 {
1018 // The expected number of columns in the congruence and generator
1019 // systems, if they are not empty.
1020 const dimension_type num_columns = space_dim + 1;
1021
1022 // Here we check if the size of the matrices is consistent.
1023 // Let us suppose that all the matrices are up-to-date; this means:
1024 // `con_sys' : number of congruences x poly_num_columns
1025 // `gen_sys' : number of generators x poly_num_columns
1026 if (congruences_are_up_to_date()) {
1027 if (con_sys.space_dimension() != space_dim) {
1028 #ifndef NDEBUG
1029 cerr << "Incompatible size! (con_sys and space_dim)"
1030 << endl;
1031 #endif
1032 goto fail;
1033 }
1034 }
1035 if (generators_are_up_to_date()) {
1036 if (gen_sys.space_dimension() != space_dim) {
1037 #ifndef NDEBUG
1038 cerr << "Incompatible size! (gen_sys and space_dim)"
1039 << endl;
1040 #endif
1041 goto fail;
1042 }
1043
1044 // A non-empty system of generators describing a grid is valid
1045 // if and only if it contains a point.
1046 if (!gen_sys.has_no_rows() && !gen_sys.has_points()) {
1047 #ifndef NDEBUG
1048 cerr << "Non-empty generator system declared up-to-date "
1049 << "has no points!"
1050 << endl;
1051 #endif
1052 goto fail;
1053 }
1054
1055 if (generators_are_minimized()) {
1056 Grid_Generator_System gs = gen_sys;
1057
1058 if (dim_kinds.size() != num_columns) {
1059 #ifndef NDEBUG
1060 cerr << "Size of dim_kinds should equal the number of columns."
1061 << endl;
1062 #endif
1063 goto fail;
1064 }
1065
1066 if (!upper_triangular(gs, dim_kinds)) {
1067 #ifndef NDEBUG
1068 cerr << "Reduced generators should be upper triangular."
1069 << endl;
1070 #endif
1071 goto fail;
1072 }
1073
1074 // Check that dim_kinds corresponds to the row kinds in gen_sys.
1075 for (dimension_type dim = space_dim,
1076 row = gen_sys.num_rows(); dim > 0; --dim) {
1077 if (dim_kinds[dim] == GEN_VIRTUAL) {
1078 goto ok;
1079 }
1080 if (gen_sys[--row].is_parameter_or_point()
1081 && dim_kinds[dim] == PARAMETER) {
1082 goto ok;
1083 }
1084 PPL_ASSERT(gen_sys[row].is_line());
1085 if (dim_kinds[dim] == LINE) {
1086 goto ok;
1087 }
1088 #ifndef NDEBUG
1089 cerr << "Kinds in dim_kinds should match those in gen_sys."
1090 << endl;
1091 #endif
1092 goto fail;
1093 ok:
1094 PPL_ASSERT(row <= dim);
1095 }
1096
1097 // A reduced generator system must be the same as a temporary
1098 // reduced copy.
1099 Dimension_Kinds dim_kinds_copy = dim_kinds;
1100 // `gs' is minimized and marked_empty returned false, so `gs'
1101 // should contain rows.
1102 PPL_ASSERT(!gs.has_no_rows());
1103 simplify(gs, dim_kinds_copy);
1104 // gs contained rows before being reduced, so it should
1105 // contain at least a single point afterward.
1106 PPL_ASSERT(!gs.has_no_rows());
1107 for (dimension_type row = gen_sys.num_rows(); row-- > 0; ) {
1108 const Grid_Generator& g = gs[row];
1109 const Grid_Generator& g_copy = gen_sys[row];
1110 if (g.is_equal_to(g_copy)) {
1111 continue;
1112 }
1113 #ifndef NDEBUG
1114 cerr << "Generators are declared minimized,"
1115 " but they change under reduction.\n"
1116 << "Here is the generator system:\n";
1117 gen_sys.ascii_dump(cerr);
1118 cerr << "and here is the minimized form of the temporary copy:\n";
1119 gs.ascii_dump(cerr);
1120 #endif
1121 goto fail;
1122 }
1123 }
1124
1125 } // if (congruences_are_up_to_date())
1126 }
1127
1128 if (congruences_are_up_to_date()) {
1129 // Check if the system of congruences is well-formed.
1130 if (!con_sys.OK()) {
1131 goto fail;
1132 }
1133
1134 Grid tmp_gr = *this;
1135 // Make a copy here, before changing tmp_gr, to check later.
1136 const Congruence_System cs_copy = tmp_gr.con_sys;
1137
1138 // Clear the generators in tmp_gr.
1139 Grid_Generator_System gs(space_dim);
1140 tmp_gr.gen_sys.m_swap(gs);
1141 tmp_gr.clear_generators_up_to_date();
1142
1143 if (!tmp_gr.update_generators()) {
1144 if (check_not_empty) {
1145 // Want to know the satisfiability of the congruences.
1146 #ifndef NDEBUG
1147 cerr << "Unsatisfiable system of congruences!"
1148 << endl;
1149 #endif
1150 goto fail;
1151 }
1152 // The grid is empty, all checks are done.
1153 return true;
1154 }
1155
1156 if (congruences_are_minimized()) {
1157 // A reduced congruence system must be lower triangular.
1158 if (!lower_triangular(con_sys, dim_kinds)) {
1159 #ifndef NDEBUG
1160 cerr << "Reduced congruences should be lower triangular." << endl;
1161 #endif
1162 goto fail;
1163 }
1164
1165 // If the congruences are minimized, all the elements in the
1166 // congruence system must be the same as those in the temporary,
1167 // minimized system `cs_copy'.
1168 if (!con_sys.is_equal_to(cs_copy)) {
1169 #ifndef NDEBUG
1170 cerr << "Congruences are declared minimized, but they change under reduction!"
1171 << endl
1172 << "Here is the minimized form of the congruence system:"
1173 << endl;
1174 cs_copy.ascii_dump(cerr);
1175 cerr << endl;
1176 #endif
1177 goto fail;
1178 }
1179
1180 if (dim_kinds.size() != con_sys.space_dimension() + 1 /* inhomogeneous term */) {
1181 #ifndef NDEBUG
1182 cerr << "Size of dim_kinds should equal the number of columns."
1183 << endl;
1184 #endif
1185 goto fail;
1186 }
1187
1188 // Check that dim_kinds corresponds to the row kinds in con_sys.
1189 for (dimension_type dim = space_dim, row = 0; dim > 0; --dim) {
1190 if (dim_kinds[dim] == CON_VIRTUAL) {
1191 continue;
1192 }
1193 if (con_sys[row++].is_proper_congruence()
1194 && dim_kinds[dim] == PROPER_CONGRUENCE) {
1195 continue;
1196 }
1197 PPL_ASSERT(con_sys[row-1].is_equality());
1198 if (dim_kinds[dim] == EQUALITY) {
1199 continue;
1200 }
1201 #ifndef NDEBUG
1202 cerr << "Kinds in dim_kinds should match those in con_sys." << endl;
1203 #endif
1204 goto fail;
1205 }
1206 }
1207 }
1208
1209 return true;
1210
1211 fail:
1212 #ifndef NDEBUG
1213 cerr << "Here is the grid under check:" << endl;
1214 ascii_dump(cerr);
1215 #endif
1216 return false;
1217 }
1218
1219 void
add_constraints(const Constraint_System & cs)1220 PPL::Grid::add_constraints(const Constraint_System& cs) {
1221 // The dimension of `cs' must be at most `space_dim'.
1222 if (space_dim < cs.space_dimension()) {
1223 throw_dimension_incompatible("add_constraints(cs)", "cs", cs);
1224 }
1225 if (marked_empty()) {
1226 return;
1227 }
1228
1229 for (Constraint_System::const_iterator i = cs.begin(),
1230 cs_end = cs.end(); i != cs_end; ++i) {
1231 add_constraint_no_check(*i);
1232 if (marked_empty()) {
1233 return;
1234 }
1235 }
1236 }
1237
1238 void
add_grid_generator(const Grid_Generator & g)1239 PPL::Grid::add_grid_generator(const Grid_Generator& g) {
1240 // The dimension of `g' must be at most space_dim.
1241 const dimension_type g_space_dim = g.space_dimension();
1242 if (space_dim < g_space_dim) {
1243 throw_dimension_incompatible("add_grid_generator(g)", "g", g);
1244 }
1245
1246 // Deal with zero-dimension case first.
1247 if (space_dim == 0) {
1248 // Points and parameters are the only zero-dimension generators
1249 // that can be created.
1250 if (marked_empty()) {
1251 if (g.is_parameter()) {
1252 throw_invalid_generator("add_grid_generator(g)", "g");
1253 }
1254 set_zero_dim_univ();
1255 }
1256 PPL_ASSERT(OK());
1257 return;
1258 }
1259
1260 if (marked_empty()
1261 || (!generators_are_up_to_date() && !update_generators())) {
1262 // Here the grid is empty: the specification says we can only
1263 // insert a point.
1264 if (g.is_line_or_parameter()) {
1265 throw_invalid_generator("add_grid_generator(g)", "g");
1266 }
1267 gen_sys.insert(g);
1268 clear_empty();
1269 }
1270 else {
1271 PPL_ASSERT(generators_are_up_to_date());
1272 gen_sys.insert(g);
1273 if (g.is_parameter_or_point()) {
1274 normalize_divisors(gen_sys);
1275 }
1276 }
1277
1278 // With the added generator, congruences are out of date.
1279 clear_congruences_up_to_date();
1280
1281 clear_generators_minimized();
1282 set_generators_up_to_date();
1283 PPL_ASSERT(OK());
1284 }
1285
1286 void
add_recycled_congruences(Congruence_System & cgs)1287 PPL::Grid::add_recycled_congruences(Congruence_System& cgs) {
1288 // Dimension-compatibility check.
1289 const dimension_type cgs_space_dim = cgs.space_dimension();
1290 if (space_dim < cgs_space_dim) {
1291 throw_dimension_incompatible("add_recycled_congruences(cgs)", "cgs", cgs);
1292 }
1293
1294 if (cgs.has_no_rows()) {
1295 return;
1296 }
1297
1298 if (marked_empty()) {
1299 return;
1300 }
1301
1302 if (space_dim == 0) {
1303 // In a 0-dimensional space the congruences are trivial (e.g., 0
1304 // == 0 or 1 %= 0) or false (e.g., 1 == 0). In a system of
1305 // congruences `begin()' and `end()' are equal if and only if the
1306 // system contains only trivial congruences.
1307 if (cgs.begin() != cgs.end()) {
1308 // There is a congruence, it must be false, the grid becomes empty.
1309 set_empty();
1310 }
1311 return;
1312 }
1313
1314 // The congruences are required.
1315 if (!congruences_are_up_to_date()) {
1316 update_congruences();
1317 }
1318
1319 // Swap (instead of copying) the coefficients of `cgs' (which is
1320 // writable).
1321 con_sys.insert(cgs, Recycle_Input());
1322
1323 // Congruences may not be minimized and generators are out of date.
1324 clear_congruences_minimized();
1325 clear_generators_up_to_date();
1326 // Note: the congruence system may have become unsatisfiable, thus
1327 // we do not check for satisfiability.
1328 PPL_ASSERT(OK());
1329 }
1330
1331 void
add_recycled_grid_generators(Grid_Generator_System & gs)1332 PPL::Grid::add_recycled_grid_generators(Grid_Generator_System& gs) {
1333 // Dimension-compatibility check.
1334 const dimension_type gs_space_dim = gs.space_dimension();
1335 if (space_dim < gs_space_dim) {
1336 throw_dimension_incompatible("add_recycled_grid_generators(gs)", "gs", gs);
1337 }
1338
1339 // Adding no generators leaves the grid the same.
1340 if (gs.has_no_rows()) {
1341 return;
1342 }
1343
1344 // Adding valid generators to a zero-dimensional grid transforms it
1345 // to the zero-dimensional universe grid.
1346 if (space_dim == 0) {
1347 if (marked_empty()) {
1348 set_zero_dim_univ();
1349 }
1350 else {
1351 PPL_ASSERT(gs.has_points());
1352 }
1353 PPL_ASSERT(OK(true));
1354 return;
1355 }
1356
1357 if (!marked_empty()) {
1358 // The grid contains at least one point.
1359
1360 if (!generators_are_up_to_date()) {
1361 update_generators();
1362 }
1363 normalize_divisors(gs, gen_sys);
1364
1365 gen_sys.insert(gs, Recycle_Input());
1366
1367 // Congruences are out of date and generators are not minimized.
1368 clear_congruences_up_to_date();
1369 clear_generators_minimized();
1370
1371 PPL_ASSERT(OK(true));
1372 return;
1373 }
1374
1375 // The grid is empty.
1376
1377 // `gs' must contain at least one point.
1378 if (!gs.has_points()) {
1379 throw_invalid_generators("add_recycled_grid_generators(gs)", "gs");
1380 }
1381 // Adjust `gs' to the right dimension.
1382 gs.set_space_dimension(space_dim);
1383
1384 gen_sys.m_swap(gs);
1385
1386 normalize_divisors(gen_sys);
1387
1388 // The grid is no longer empty and generators are up-to-date.
1389 set_generators_up_to_date();
1390 clear_empty();
1391
1392 PPL_ASSERT(OK());
1393 }
1394
1395 void
add_grid_generators(const Grid_Generator_System & gs)1396 PPL::Grid::add_grid_generators(const Grid_Generator_System& gs) {
1397 // TODO: this is just an executable specification.
1398 Grid_Generator_System gs_copy = gs;
1399 add_recycled_grid_generators(gs_copy);
1400 }
1401
1402 void
refine_with_constraint(const Constraint & c)1403 PPL::Grid::refine_with_constraint(const Constraint& c) {
1404 // The dimension of `c' must be at most `space_dim'.
1405 if (space_dim < c.space_dimension()) {
1406 throw_dimension_incompatible("refine_with_constraint(c)", "c", c);
1407 }
1408 if (marked_empty()) {
1409 return;
1410 }
1411 refine_no_check(c);
1412 }
1413
1414 void
refine_with_constraints(const Constraint_System & cs)1415 PPL::Grid::refine_with_constraints(const Constraint_System& cs) {
1416 // The dimension of `cs' must be at most `space_dim'.
1417 if (space_dim < cs.space_dimension()) {
1418 throw_dimension_incompatible("refine_with_constraints(cs)", "cs", cs);
1419 }
1420
1421 for (Constraint_System::const_iterator i = cs.begin(),
1422 cs_end = cs.end(); !marked_empty() && i != cs_end; ++i) {
1423 refine_no_check(*i);
1424 }
1425 }
1426
1427 void
unconstrain(const Variable var)1428 PPL::Grid::unconstrain(const Variable var) {
1429 // Dimension-compatibility check.
1430 if (space_dim < var.space_dimension()) {
1431 throw_dimension_incompatible("unconstrain(var)", var.space_dimension());
1432 }
1433
1434 // Do something only if the grid is non-empty.
1435 if (marked_empty()
1436 || (!generators_are_up_to_date() && !update_generators())) {
1437 // Empty: do nothing.
1438 return;
1439 }
1440 PPL_ASSERT(generators_are_up_to_date());
1441 Grid_Generator l = grid_line(var);
1442 gen_sys.insert(l, Recycle_Input());
1443 // With the added generator, congruences are out of date.
1444 clear_congruences_up_to_date();
1445 clear_generators_minimized();
1446 PPL_ASSERT(OK());
1447 }
1448
1449 void
unconstrain(const Variables_Set & vars)1450 PPL::Grid::unconstrain(const Variables_Set& vars) {
1451 // The cylindrification with respect to no dimensions is a no-op.
1452 // This case also captures the only legal cylindrification
1453 // of a grid in a 0-dim space.
1454 if (vars.empty()) {
1455 return;
1456 }
1457 // Dimension-compatibility check.
1458 const dimension_type min_space_dim = vars.space_dimension();
1459 if (space_dim < min_space_dim) {
1460 throw_dimension_incompatible("unconstrain(vs)", min_space_dim);
1461 }
1462
1463 // Do something only if the grid is non-empty.
1464 if (marked_empty()
1465 || (!generators_are_up_to_date() && !update_generators())) {
1466 // Empty: do nothing.
1467 return;
1468 }
1469
1470 PPL_ASSERT(generators_are_up_to_date());
1471 // Since `gen_sys' is not empty, the space dimension of the inserted
1472 // generators are automatically adjusted.
1473 for (Variables_Set::const_iterator vsi = vars.begin(),
1474 vsi_end = vars.end(); vsi != vsi_end; ++vsi) {
1475 Grid_Generator l = grid_line(Variable(*vsi));
1476 gen_sys.insert(l, Recycle_Input());
1477 }
1478 // Constraints are no longer up-to-date.
1479 clear_generators_minimized();
1480 clear_congruences_up_to_date();
1481 PPL_ASSERT(OK());
1482 }
1483
1484 void
intersection_assign(const Grid & y)1485 PPL::Grid::intersection_assign(const Grid& y) {
1486 Grid& x = *this;
1487 // Dimension-compatibility check.
1488 if (x.space_dim != y.space_dim) {
1489 throw_dimension_incompatible("intersection_assign(y)", "y", y);
1490 }
1491 // If one of the two grids is empty, the intersection is empty.
1492 if (x.marked_empty()) {
1493 return;
1494 }
1495 if (y.marked_empty()) {
1496 x.set_empty();
1497 return;
1498 }
1499
1500 // If both grids are zero-dimensional, then at this point they are
1501 // necessarily universe, so the intersection is also universe.
1502 if (x.space_dim == 0) {
1503 return;
1504 }
1505
1506 // The congruences must be up-to-date.
1507 if (!x.congruences_are_up_to_date()) {
1508 x.update_congruences();
1509 }
1510 if (!y.congruences_are_up_to_date()) {
1511 y.update_congruences();
1512 }
1513
1514 if (!y.con_sys.has_no_rows()) {
1515 x.con_sys.insert(y.con_sys);
1516 // Grid_Generators may be out of date and congruences may have changed
1517 // from minimal form.
1518 x.clear_generators_up_to_date();
1519 x.clear_congruences_minimized();
1520 }
1521
1522 PPL_ASSERT(x.OK() && y.OK());
1523 }
1524
1525 void
upper_bound_assign(const Grid & y)1526 PPL::Grid::upper_bound_assign(const Grid& y) {
1527 Grid& x = *this;
1528 // Dimension-compatibility check.
1529 if (x.space_dim != y.space_dim) {
1530 throw_dimension_incompatible("upper_bound_assign(y)", "y", y);
1531 }
1532 // The join of a grid `gr' with an empty grid is `gr'.
1533 if (y.marked_empty()) {
1534 return;
1535 }
1536 if (x.marked_empty()) {
1537 x = y;
1538 return;
1539 }
1540
1541 // If both grids are zero-dimensional, then they are necessarily
1542 // universe grids, and so is their join.
1543 if (x.space_dim == 0) {
1544 return;
1545 }
1546
1547 // The generators must be up-to-date.
1548 if (!x.generators_are_up_to_date() && !x.update_generators()) {
1549 // Discovered `x' empty when updating generators.
1550 x = y;
1551 return;
1552 }
1553 if (!y.generators_are_up_to_date() && !y.update_generators()) {
1554 // Discovered `y' empty when updating generators.
1555 return;
1556 }
1557 // Match the divisors of the x and y generator systems.
1558 Grid_Generator_System gs(y.gen_sys);
1559 normalize_divisors(x.gen_sys, gs);
1560 x.gen_sys.insert(gs, Recycle_Input());
1561 // Congruences may be out of date and generators may have lost
1562 // minimal form.
1563 x.clear_congruences_up_to_date();
1564 x.clear_generators_minimized();
1565
1566 // At this point both `x' and `y' are not empty.
1567 PPL_ASSERT(x.OK(true) && y.OK(true));
1568 }
1569
1570 bool
upper_bound_assign_if_exact(const Grid & y)1571 PPL::Grid::upper_bound_assign_if_exact(const Grid& y) {
1572 const Grid& x = *this;
1573
1574 // Dimension-compatibility check.
1575 if (x.space_dim != y.space_dim) {
1576 throw_dimension_incompatible("upper_bound_assign_if_exact(y)", "y", y);
1577 }
1578
1579 if (x.marked_empty()
1580 || y.marked_empty()
1581 || x.space_dim == 0
1582 || x.is_included_in(y)
1583 || y.is_included_in(x)) {
1584 upper_bound_assign(y);
1585 return true;
1586 }
1587
1588 // The above test 'x.is_included_in(y)' will ensure the generators of x
1589 // are up to date.
1590 PPL_ASSERT(generators_are_up_to_date());
1591
1592 Grid x_copy = x;
1593 x_copy.upper_bound_assign(y);
1594 x_copy.difference_assign(y);
1595 if (x_copy.is_included_in(x)) {
1596 upper_bound_assign(y);
1597 return true;
1598 }
1599
1600 return false;
1601 }
1602
1603 void
difference_assign(const Grid & y)1604 PPL::Grid::difference_assign(const Grid& y) {
1605 Grid& x = *this;
1606 // Dimension-compatibility check.
1607 if (x.space_dim != y.space_dim) {
1608 throw_dimension_incompatible("difference_assign(y)", "y", y);
1609 }
1610
1611 if (y.marked_empty() || x.marked_empty()) {
1612 return;
1613 }
1614
1615 // If both grids are zero-dimensional, then they are necessarily
1616 // universe grids, so the result is empty.
1617 if (x.space_dim == 0) {
1618 x.set_empty();
1619 return;
1620 }
1621
1622 if (y.contains(x)) {
1623 x.set_empty();
1624 return;
1625 }
1626
1627 Grid new_grid(x.space_dim, EMPTY);
1628
1629 const Congruence_System& y_cgs = y.congruences();
1630 for (Congruence_System::const_iterator i = y_cgs.begin(),
1631 y_cgs_end = y_cgs.end(); i != y_cgs_end; ++i) {
1632 const Congruence& cg = *i;
1633
1634 // The 2-complement cg2 of cg = ((e %= 0) / m) is the congruence
1635 // defining the sets of points exactly half-way between successive
1636 // hyperplanes e = km and e = (k+1)m, for any integer k; that is,
1637 // the hyperplanes defined by 2e = (2k + 1)m, for any integer k.
1638 // Thus `cg2' is the congruence ((2e %= m) / 2m).
1639
1640 // As the grid difference must be a grid, only add the
1641 // 2-complement congruence to x if the resulting grid includes all
1642 // the points in x that did not satisfy `cg'.
1643
1644 // The 2-complement of cg can be included in the result only if x
1645 // holds points other than those in cg.
1646 if (x.relation_with(cg).implies(Poly_Con_Relation::is_included())) {
1647 continue;
1648 }
1649
1650 if (cg.is_proper_congruence()) {
1651 const Linear_Expression e(cg.expression());
1652 // Congruence cg is ((e %= 0) / m).
1653 const Coefficient& m = cg.modulus();
1654 // If x is included in the grid defined by the congruences cg
1655 // and its 2-complement (i.e. the grid defined by the congruence
1656 // (2e %= 0) / m) then add the 2-complement to the potential
1657 // result.
1658 if (x.relation_with((2*e %= 0) / m)
1659 .implies(Poly_Con_Relation::is_included())) {
1660 Grid z = x;
1661 z.add_congruence_no_check((2*e %= m) / (2*m));
1662 new_grid.upper_bound_assign(z);
1663 continue;
1664 }
1665 }
1666 return;
1667 }
1668
1669 *this = new_grid;
1670
1671 PPL_ASSERT(OK());
1672 }
1673
1674 namespace {
1675
1676 struct Ruled_Out_Pair {
1677 PPL::dimension_type congruence_index;
1678 PPL::dimension_type num_ruled_out;
1679 };
1680
1681 struct Ruled_Out_Less_Than {
operator ()__anon85d1942c0111::Ruled_Out_Less_Than1682 bool operator()(const Ruled_Out_Pair& x,
1683 const Ruled_Out_Pair& y) const {
1684 return x.num_ruled_out > y.num_ruled_out;
1685 }
1686 };
1687
1688 } // namespace
1689
1690 bool
simplify_using_context_assign(const Grid & y)1691 PPL::Grid::simplify_using_context_assign(const Grid& y) {
1692 Grid& x = *this;
1693 // Dimension-compatibility check.
1694 if (x.space_dim != y.space_dim) {
1695 throw_dimension_incompatible("simplify_using_context_assign(y)", "y", y);
1696 }
1697
1698 // Filter away the zero-dimensional case.
1699 if (x.space_dim == 0) {
1700 if (y.is_empty()) {
1701 set_zero_dim_univ();
1702 PPL_ASSERT(OK());
1703 return false;
1704 }
1705 else {
1706 return !x.is_empty();
1707 }
1708 }
1709
1710 // If `y' is empty, the biggest enlargement for `x' is the universe.
1711 if (!y.minimize()) {
1712 Grid gr(x.space_dim, UNIVERSE);
1713 m_swap(gr);
1714 return false;
1715 }
1716
1717 // If `x' is empty, the intersection is empty.
1718 if (!x.minimize()) {
1719 // Search for a congruence of `y' that is not a tautology.
1720 PPL_ASSERT(y.congruences_are_up_to_date());
1721 Grid gr(x.space_dim, UNIVERSE);
1722 for (dimension_type i = y.con_sys.num_rows(); i-- > 0; ) {
1723 const Congruence& y_con_sys_i = y.con_sys[i];
1724 if (!y_con_sys_i.is_tautological()) {
1725 // Found: we obtain a congruence `c' contradicting the one we
1726 // found, and assign to `x' the grid `gr' with `c' as
1727 // the only congruence.
1728 const Linear_Expression le(y_con_sys_i.expression());
1729 if (y_con_sys_i.is_equality()) {
1730 gr.refine_no_check(le == 1);
1731 break;
1732 }
1733 else {
1734 const Coefficient& y_modulus_i = y_con_sys_i.modulus();
1735 if (y_modulus_i > 1) {
1736 gr.refine_no_check(le == 1);
1737 }
1738 else {
1739 Linear_Expression le2 = le;
1740 le2 *= 2;
1741 gr.refine_no_check(le2 == y_modulus_i);
1742 }
1743 break;
1744 }
1745 }
1746 }
1747 m_swap(gr);
1748 PPL_ASSERT(OK());
1749 return false;
1750 }
1751
1752 PPL_ASSERT(x.congruences_are_minimized()
1753 && y.generators_are_minimized());
1754
1755 const Congruence_System& x_cs = x.con_sys;
1756 const dimension_type x_cs_num_rows = x_cs.num_rows();
1757
1758 // Record into `redundant_by_y' the info about which congruences of
1759 // `x' are redundant in the context `y'. Count the number of
1760 // redundancies found.
1761 std::vector<bool> redundant_by_y(x_cs_num_rows, false);
1762 dimension_type num_redundant_by_y = 0;
1763 for (dimension_type i = 0; i < x_cs_num_rows; ++i) {
1764 if (y.relation_with(x_cs[i]).implies(Poly_Con_Relation::is_included())) {
1765 redundant_by_y[i] = true;
1766 ++num_redundant_by_y;
1767 }
1768 }
1769
1770 if (num_redundant_by_y < x_cs_num_rows) {
1771
1772 // Some congruences were not identified as redundant.
1773
1774 Congruence_System result_cs;
1775 const Congruence_System& y_cs = y.con_sys;
1776 const dimension_type y_cs_num_rows = y_cs.num_rows();
1777 // Compute into `z' the intersection of `x' and `y'.
1778 const bool x_first = (x_cs_num_rows > y_cs_num_rows);
1779 Grid z(x_first ? x : y);
1780 if (x_first) {
1781 z.add_congruences(y_cs);
1782 }
1783 else {
1784 // Only copy (and then recycle) the non-redundant congruences.
1785 Congruence_System tmp_cs;
1786 for (dimension_type i = 0; i < x_cs_num_rows; ++i) {
1787 if (!redundant_by_y[i]) {
1788 tmp_cs.insert(x_cs[i]);
1789 }
1790 }
1791 z.add_recycled_congruences(tmp_cs);
1792 }
1793
1794 // Congruences are added to `w' until it equals `z'.
1795 // We do not care about minimization or maximization, since
1796 // we are only interested in satisfiability.
1797 Grid w;
1798 w.add_space_dimensions_and_embed(x.space_dim);
1799 // First add the congruences for `y'.
1800 w.add_congruences(y_cs);
1801
1802 // We apply the following heuristics here: congruences of `x' that
1803 // are not made redundant by `y' are added to `w' depending on
1804 // the number of generators of `y' they rule out (the more generators)
1805 // (they rule out, the sooner they are added). Of course, as soon
1806 // as `w' becomes empty, we stop adding.
1807 std::vector<Ruled_Out_Pair>
1808 ruled_out_vec(x_cs_num_rows - num_redundant_by_y);
1809
1810 PPL_DIRTY_TEMP_COEFFICIENT(sp);
1811 PPL_DIRTY_TEMP_COEFFICIENT(div);
1812
1813 for (dimension_type i = 0, j = 0; i < x_cs_num_rows; ++i) {
1814 if (!redundant_by_y[i]) {
1815 const Congruence& c = x_cs[i];
1816 const Coefficient& modulus = c.modulus();
1817 div = modulus;
1818
1819 const Grid_Generator_System& y_gs = y.gen_sys;
1820 dimension_type num_ruled_out_generators = 0;
1821 for (Grid_Generator_System::const_iterator k = y_gs.begin(),
1822 y_gs_end = y_gs.end(); k != y_gs_end; ++k) {
1823 const Grid_Generator& g = *k;
1824 // If the generator is not to be ruled out,
1825 // it must saturate the congruence.
1826 Scalar_Products::assign(sp, c, g);
1827 // If `c' is a proper congruence the scalar product must be
1828 // reduced modulo a (possibly scaled) modulus.
1829 if (c.is_proper_congruence()) {
1830 // If `g' is a parameter the congruence modulus must be scaled
1831 // up by the divisor of the generator.
1832 if (g.is_parameter()) {
1833 sp %= (div * g.divisor());
1834 }
1835 else {
1836 if (g.is_point()) {
1837 sp %= div;
1838 }
1839 }
1840 }
1841 if (sp == 0) {
1842 continue;
1843 }
1844 ++num_ruled_out_generators;
1845 }
1846 ruled_out_vec[j].congruence_index = i;
1847 ruled_out_vec[j].num_ruled_out = num_ruled_out_generators;
1848 ++j;
1849 }
1850 }
1851 std::sort(ruled_out_vec.begin(), ruled_out_vec.end(),
1852 Ruled_Out_Less_Than());
1853
1854 const bool empty_intersection = (!z.minimize());
1855
1856 // Add the congruences in the "ruled out" order to `w'
1857 // until the result is the intersection.
1858 for (std::vector<Ruled_Out_Pair>::const_iterator
1859 j = ruled_out_vec.begin(), ruled_out_vec_end = ruled_out_vec.end();
1860 j != ruled_out_vec_end;
1861 ++j) {
1862 const Congruence& c = x_cs[j->congruence_index];
1863 result_cs.insert(c);
1864 w.add_congruence(c);
1865 if ((empty_intersection && w.is_empty())
1866 || (!empty_intersection && w.is_included_in(z))) {
1867 Grid result_gr(x.space_dim, UNIVERSE);
1868 result_gr.add_congruences(result_cs);
1869 x.m_swap(result_gr);
1870 PPL_ASSERT(x.OK());
1871 return !empty_intersection;
1872 }
1873 }
1874 // Cannot exit from here.
1875 PPL_UNREACHABLE;
1876 }
1877
1878 // All the congruences are redundant so that the simplified grid
1879 // is the universe.
1880 Grid result_gr(x.space_dim, UNIVERSE);
1881 x.m_swap(result_gr);
1882 PPL_ASSERT(x.OK());
1883 return true;
1884 }
1885
1886 void
affine_image(const Variable var,const Linear_Expression & expr,Coefficient_traits::const_reference denominator)1887 PPL::Grid::affine_image(const Variable var,
1888 const Linear_Expression& expr,
1889 Coefficient_traits::const_reference denominator) {
1890 // The denominator cannot be zero.
1891 if (denominator == 0) {
1892 throw_invalid_argument("affine_image(v, e, d)", "d == 0");
1893 }
1894
1895 // Dimension-compatibility checks.
1896 // The dimension of `expr' must be at most the dimension of `*this'.
1897 const dimension_type expr_space_dim = expr.space_dimension();
1898 if (space_dim < expr_space_dim) {
1899 throw_dimension_incompatible("affine_image(v, e, d)", "e", expr);
1900 }
1901 // `var' must be one of the dimensions of the grid.
1902 const dimension_type var_space_dim = var.space_dimension();
1903 if (space_dim < var_space_dim) {
1904 throw_dimension_incompatible("affine_image(v, e, d)", "v", var);
1905 }
1906
1907 if (marked_empty()) {
1908 return;
1909 }
1910
1911 Coefficient_traits::const_reference expr_var = expr.coefficient(var);
1912
1913 if (var_space_dim <= expr_space_dim
1914 && expr_var != 0) {
1915 // The transformation is invertible.
1916 if (generators_are_up_to_date()) {
1917 // Grid_Generator_System::affine_image() requires the third argument
1918 // to be a positive Coefficient.
1919 if (denominator > 0) {
1920 gen_sys.affine_image(var, expr, denominator);
1921 }
1922 else {
1923 gen_sys.affine_image(var, -expr, -denominator);
1924 }
1925 clear_generators_minimized();
1926 // Strong normalization in gs::affine_image may have modified
1927 // divisors.
1928 normalize_divisors(gen_sys);
1929 }
1930 if (congruences_are_up_to_date()) {
1931 // To build the inverse transformation,
1932 // after copying and negating `expr',
1933 // we exchange the roles of `expr[var_space_dim]' and `denominator'.
1934 Linear_Expression inverse;
1935 if (expr_var > 0) {
1936 inverse = -expr;
1937 inverse.set_coefficient(var, denominator);
1938 con_sys.affine_preimage(var, inverse, expr_var);
1939 }
1940 else {
1941 // The new denominator is negative: we negate everything once
1942 // more, as Congruence_System::affine_preimage() requires the
1943 // third argument to be positive.
1944 inverse = expr;
1945 inverse.set_coefficient(var, -denominator);
1946 con_sys.affine_preimage(var, inverse, -expr_var);
1947 }
1948 clear_congruences_minimized();
1949 }
1950 }
1951 else {
1952 // The transformation is not invertible.
1953 // We need an up-to-date system of generators.
1954 if (!generators_are_up_to_date()) {
1955 minimize();
1956 }
1957 if (!marked_empty()) {
1958 // Grid_Generator_System::affine_image() requires the third argument
1959 // to be a positive Coefficient.
1960 if (denominator > 0) {
1961 gen_sys.affine_image(var, expr, denominator);
1962 }
1963 else {
1964 gen_sys.affine_image(var, -expr, -denominator);
1965 }
1966
1967 clear_congruences_up_to_date();
1968 clear_generators_minimized();
1969 // Strong normalization in gs::affine_image may have modified
1970 // divisors.
1971 normalize_divisors(gen_sys);
1972 }
1973 }
1974 PPL_ASSERT(OK());
1975 }
1976
1977 void
1978 PPL::Grid::
affine_preimage(const Variable var,const Linear_Expression & expr,Coefficient_traits::const_reference denominator)1979 affine_preimage(const Variable var,
1980 const Linear_Expression& expr,
1981 Coefficient_traits::const_reference denominator) {
1982 // The denominator cannot be zero.
1983 if (denominator == 0) {
1984 throw_invalid_argument("affine_preimage(v, e, d)", "d == 0");
1985 }
1986
1987 // Dimension-compatibility checks.
1988 // The dimension of `expr' should not be greater than the dimension
1989 // of `*this'.
1990 const dimension_type expr_space_dim = expr.space_dimension();
1991 if (space_dim < expr_space_dim) {
1992 throw_dimension_incompatible("affine_preimage(v, e, d)", "e", expr);
1993 }
1994 // `var' should be one of the dimensions of the grid.
1995 const dimension_type var_space_dim = var.space_dimension();
1996 if (space_dim < var_space_dim) {
1997 throw_dimension_incompatible("affine_preimage(v, e, d)", "v", var);
1998 }
1999
2000 if (marked_empty()) {
2001 return;
2002 }
2003
2004 Coefficient_traits::const_reference expr_var = expr.coefficient(var);
2005
2006 if (var_space_dim <= expr_space_dim && expr_var != 0) {
2007 // The transformation is invertible.
2008 if (congruences_are_up_to_date()) {
2009 // Congruence_System::affine_preimage() requires the third argument
2010 // to be a positive Coefficient.
2011 if (denominator > 0) {
2012 con_sys.affine_preimage(var, expr, denominator);
2013 }
2014 else {
2015 con_sys.affine_preimage(var, -expr, -denominator);
2016 }
2017 clear_congruences_minimized();
2018 }
2019 if (generators_are_up_to_date()) {
2020 // To build the inverse transformation,
2021 // after copying and negating `expr',
2022 // we exchange the roles of `expr[var_space_dim]' and `denominator'.
2023 Linear_Expression inverse;
2024 if (expr_var > 0) {
2025 inverse = -expr;
2026 inverse.set_coefficient(var, denominator);
2027 gen_sys.affine_image(var, inverse, expr_var);
2028 }
2029 else {
2030 // The new denominator is negative: we negate everything once
2031 // more, as Grid_Generator_System::affine_image() requires the
2032 // third argument to be positive.
2033 inverse = expr;
2034 inverse.set_coefficient(var, -denominator);
2035 gen_sys.affine_image(var, inverse, -expr_var);
2036 }
2037 clear_generators_minimized();
2038 }
2039 }
2040 else {
2041 // The transformation is not invertible.
2042 // We need an up-to-date system of congruences.
2043 if (!congruences_are_up_to_date()) {
2044 minimize();
2045 }
2046 // Congruence_System::affine_preimage() requires the third argument
2047 // to be a positive Coefficient.
2048 if (denominator > 0) {
2049 con_sys.affine_preimage(var, expr, denominator);
2050 }
2051 else {
2052 con_sys.affine_preimage(var, -expr, -denominator);
2053 }
2054
2055 clear_generators_up_to_date();
2056 clear_congruences_minimized();
2057 }
2058 PPL_ASSERT(OK());
2059 }
2060
2061 void
2062 PPL::Grid::
generalized_affine_image(const Variable var,const Relation_Symbol relsym,const Linear_Expression & expr,Coefficient_traits::const_reference denominator,Coefficient_traits::const_reference modulus)2063 generalized_affine_image(const Variable var,
2064 const Relation_Symbol relsym,
2065 const Linear_Expression& expr,
2066 Coefficient_traits::const_reference denominator,
2067 Coefficient_traits::const_reference modulus) {
2068 // The denominator cannot be zero.
2069 if (denominator == 0) {
2070 throw_invalid_argument("generalized_affine_image(v, r, e, d, m)",
2071 "d == 0");
2072 }
2073
2074 // Dimension-compatibility checks.
2075 // The dimension of `expr' should not be greater than the dimension
2076 // of `*this'.
2077 const dimension_type expr_space_dim = expr.space_dimension();
2078 if (space_dim < expr_space_dim) {
2079 throw_dimension_incompatible("generalized_affine_image(v, r, e, d, m)",
2080 "e", expr);
2081 }
2082 // `var' should be one of the dimensions of the grid.
2083 const dimension_type var_space_dim = var.space_dimension();
2084 if (space_dim < var_space_dim) {
2085 throw_dimension_incompatible("generalized_affine_image(v, r, e, d, m)",
2086 "v", var);
2087 }
2088 // The relation symbol cannot be a disequality.
2089 if (relsym == NOT_EQUAL) {
2090 throw_invalid_argument("generalized_affine_image(v, r, e, d, m)",
2091 "r is the disequality relation symbol");
2092 }
2093
2094 // Any image of an empty grid is empty.
2095 if (marked_empty()) {
2096 return;
2097 }
2098
2099 // If relsym is not EQUAL, then we return a safe approximation
2100 // by adding a line in the direction of var.
2101 if (relsym != EQUAL) {
2102
2103 if (modulus != 0) {
2104 throw_invalid_argument("generalized_affine_image(v, r, e, d, m)",
2105 "r != EQUAL && m != 0");
2106 }
2107
2108 if (!generators_are_up_to_date()) {
2109 minimize();
2110 }
2111
2112 // Any image of an empty grid is empty.
2113 if (marked_empty()) {
2114 return;
2115 }
2116
2117 add_grid_generator(grid_line(var));
2118
2119 PPL_ASSERT(OK());
2120 return;
2121 }
2122
2123 PPL_ASSERT(relsym == EQUAL);
2124
2125 affine_image(var, expr, denominator);
2126
2127 if (modulus == 0) {
2128 return;
2129 }
2130
2131 // Modulate dimension `var' according to `modulus'.
2132
2133 if (!generators_are_up_to_date()) {
2134 minimize();
2135 }
2136
2137 // Test if minimization, possibly in affine_image, found an empty
2138 // grid.
2139 if (marked_empty()) {
2140 return;
2141 }
2142
2143 if (modulus < 0) {
2144 gen_sys.insert(parameter(-modulus * var));
2145 }
2146 else {
2147 gen_sys.insert(parameter(modulus * var));
2148 }
2149
2150 normalize_divisors(gen_sys);
2151
2152 clear_generators_minimized();
2153 clear_congruences_up_to_date();
2154
2155 PPL_ASSERT(OK());
2156 }
2157
2158 void
2159 PPL::Grid::
generalized_affine_preimage(const Variable var,const Relation_Symbol relsym,const Linear_Expression & expr,Coefficient_traits::const_reference denominator,Coefficient_traits::const_reference modulus)2160 generalized_affine_preimage(const Variable var,
2161 const Relation_Symbol relsym,
2162 const Linear_Expression& expr,
2163 Coefficient_traits::const_reference denominator,
2164 Coefficient_traits::const_reference modulus) {
2165 // The denominator cannot be zero.
2166 if (denominator == 0) {
2167 throw_invalid_argument("generalized_affine_preimage(v, r, e, d, m)",
2168 "d == 0");
2169 }
2170
2171 // The dimension of `expr' should be at most the dimension of
2172 // `*this'.
2173 const dimension_type expr_space_dim = expr.space_dimension();
2174 if (space_dim < expr_space_dim) {
2175 throw_dimension_incompatible("generalized_affine_preimage(v, r, e, d, m)",
2176 "e", expr);
2177 }
2178 // `var' should be one of the dimensions of the grid.
2179 const dimension_type var_space_dim = var.space_dimension();
2180 if (space_dim < var_space_dim) {
2181 throw_dimension_incompatible("generalized_affine_preimage(v, r, e, d, m)",
2182 "v", var);
2183 }
2184 // The relation symbol cannot be a disequality.
2185 if (relsym == NOT_EQUAL) {
2186 throw_invalid_argument("generalized_affine_preimage(v, r, e, d, m)",
2187 "r is the disequality relation symbol");
2188 }
2189
2190 // If relsym is not EQUAL, then we return a safe approximation
2191 // by adding a line in the direction of var.
2192 if (relsym != EQUAL) {
2193
2194 if (modulus != 0) {
2195 throw_invalid_argument("generalized_affine_preimage(v, r, e, d, m)",
2196 "r != EQUAL && m != 0");
2197 }
2198
2199 if (!generators_are_up_to_date()) {
2200 minimize();
2201 }
2202
2203 // Any image of an empty grid is empty.
2204 if (marked_empty()) {
2205 return;
2206 }
2207
2208 add_grid_generator(grid_line(var));
2209
2210 PPL_ASSERT(OK());
2211 return;
2212 }
2213
2214 PPL_ASSERT(relsym == EQUAL);
2215 // Any image of an empty grid is empty.
2216 if (marked_empty()) {
2217 return;
2218 }
2219 // Check whether the affine relation is an affine function.
2220 if (modulus == 0) {
2221 affine_preimage(var, expr, denominator);
2222 return;
2223 }
2224
2225 // Check whether the preimage of this affine relation can be easily
2226 // computed as the image of its inverse relation.
2227 const Coefficient& var_coefficient = expr.coefficient(var);
2228 if (var_space_dim <= expr_space_dim && var_coefficient != 0) {
2229 const Linear_Expression inverse_expr
2230 = expr - (denominator + var_coefficient) * var;
2231 PPL_DIRTY_TEMP_COEFFICIENT(inverse_denominator);
2232 neg_assign(inverse_denominator, var_coefficient);
2233 if (modulus < 0) {
2234 generalized_affine_image(var, EQUAL, inverse_expr, inverse_denominator,
2235 - modulus);
2236 }
2237 else {
2238 generalized_affine_image(var, EQUAL, inverse_expr, inverse_denominator,
2239 modulus);
2240 }
2241 return;
2242 }
2243
2244 // Here `var_coefficient == 0', so that the preimage cannot be
2245 // easily computed by inverting the affine relation. Add the
2246 // congruence induced by the affine relation.
2247 {
2248 Congruence cg((denominator*var %= expr) / denominator);
2249 if (modulus < 0) {
2250 cg /= -modulus;
2251 }
2252 else {
2253 cg /= modulus;
2254 }
2255 add_congruence_no_check(cg);
2256 }
2257
2258 // If the resulting grid is empty, its preimage is empty too.
2259 // Note: DO check for emptiness here, as we will later add a line.
2260 if (is_empty()) {
2261 return;
2262 }
2263 add_grid_generator(grid_line(var));
2264 PPL_ASSERT(OK());
2265 }
2266
2267 void
2268 PPL::Grid::
generalized_affine_image(const Linear_Expression & lhs,const Relation_Symbol relsym,const Linear_Expression & rhs,Coefficient_traits::const_reference modulus)2269 generalized_affine_image(const Linear_Expression& lhs,
2270 const Relation_Symbol relsym,
2271 const Linear_Expression& rhs,
2272 Coefficient_traits::const_reference modulus) {
2273 // Dimension-compatibility checks.
2274 // The dimension of `lhs' should be at most the dimension of
2275 // `*this'.
2276 dimension_type lhs_space_dim = lhs.space_dimension();
2277 if (space_dim < lhs_space_dim) {
2278 throw_dimension_incompatible("generalized_affine_image(e1, r, e2, m)",
2279 "e1", lhs);
2280 }
2281 // The dimension of `rhs' should be at most the dimension of
2282 // `*this'.
2283 const dimension_type rhs_space_dim = rhs.space_dimension();
2284 if (space_dim < rhs_space_dim) {
2285 throw_dimension_incompatible("generalized_affine_image(e1, r, e2, m)",
2286 "e2", rhs);
2287 }
2288 // The relation symbol cannot be a disequality.
2289 if (relsym == NOT_EQUAL) {
2290 throw_invalid_argument("generalized_affine_image(e1, r, e2, m)",
2291 "r is the disequality relation symbol");
2292 }
2293
2294 // Any image of an empty grid is empty.
2295 if (marked_empty()) {
2296 return;
2297 }
2298
2299 // If relsym is not EQUAL, then we return a safe approximation
2300 // by adding a line in the direction of var.
2301 if (relsym != EQUAL) {
2302
2303 if (modulus != 0) {
2304 throw_invalid_argument("generalized_affine_image(e1, r, e2, m)",
2305 "r != EQUAL && m != 0");
2306 }
2307
2308 if (!generators_are_up_to_date()) {
2309 minimize();
2310 }
2311
2312 // Any image of an empty grid is empty.
2313 if (marked_empty()) {
2314 return;
2315 }
2316
2317 for (Linear_Expression::const_iterator i = lhs.begin(), i_end = lhs.end();
2318 i != i_end; ++i) {
2319 add_grid_generator(grid_line(i.variable()));
2320 }
2321
2322 PPL_ASSERT(OK());
2323 return;
2324 }
2325
2326 PPL_ASSERT(relsym == EQUAL);
2327
2328 PPL_DIRTY_TEMP_COEFFICIENT(tmp_modulus);
2329 tmp_modulus = modulus;
2330 if (tmp_modulus < 0) {
2331 neg_assign(tmp_modulus);
2332 }
2333
2334 // Compute the actual space dimension of `lhs',
2335 // i.e., the highest dimension having a non-zero coefficient in `lhs'.
2336
2337 lhs_space_dim = lhs.last_nonzero();
2338 if (lhs_space_dim == 0) {
2339 // All variables have zero coefficients, so `lhs' is a constant.
2340 add_congruence_no_check((lhs %= rhs) / tmp_modulus);
2341 return;
2342 }
2343
2344 // Gather in `new_lines' the collections of all the lines having the
2345 // direction of variables occurring in `lhs'.
2346 Grid_Generator_System new_lines;
2347 for (Linear_Expression::const_iterator i = lhs.begin(),
2348 i_end = lhs.lower_bound(Variable(lhs_space_dim)); i != i_end; ++i) {
2349 new_lines.insert(grid_line(i.variable()));
2350 }
2351
2352 const dimension_type num_common_dims = std::min(lhs_space_dim, rhs_space_dim);
2353 if (lhs.have_a_common_variable(rhs, Variable(0), Variable(num_common_dims))) {
2354 // Some variables in `lhs' also occur in `rhs'.
2355 // To ease the computation, add an additional dimension.
2356 const Variable new_var(space_dim);
2357 add_space_dimensions_and_embed(1);
2358
2359 // Constrain the new dimension to be equal to the right hand side.
2360 // TODO: Use add_congruence() when it has been updated.
2361 Congruence_System new_cgs1(new_var == rhs);
2362 add_recycled_congruences(new_cgs1);
2363 if (!is_empty()) {
2364 // The grid still contains points.
2365
2366 // Existentially quantify all the variables occurring in the left
2367 // hand side expression.
2368
2369 // Adjust `new_lines' to the right dimension.
2370 new_lines.set_space_dimension(space_dim);
2371 // Add the lines to `gen_sys' (first make sure they are up-to-date).
2372 update_generators();
2373 gen_sys.insert(new_lines, Recycle_Input());
2374 normalize_divisors(gen_sys);
2375 // Update the flags.
2376 clear_congruences_up_to_date();
2377 clear_generators_minimized();
2378
2379 // Constrain the new dimension so that it is congruent to the left
2380 // hand side expression modulo `modulus'.
2381 // TODO: Use add_congruence() when it has been updated.
2382 Congruence_System new_cgs2((lhs %= new_var) / tmp_modulus);
2383 add_recycled_congruences(new_cgs2);
2384 }
2385
2386 // Remove the temporarily added dimension.
2387 remove_higher_space_dimensions(space_dim-1);
2388 }
2389 else {
2390 // `lhs' and `rhs' variables are disjoint:
2391 // there is no need to add a further dimension.
2392
2393 // Only add the lines and congruence if there are points.
2394 if (is_empty()) {
2395 return;
2396 }
2397
2398 // Existentially quantify all the variables occurring in the left hand
2399 // side expression.
2400 add_recycled_grid_generators(new_lines);
2401
2402 // Constrain the left hand side expression so that it is congruent to
2403 // the right hand side expression modulo `modulus'.
2404 add_congruence_no_check((lhs %= rhs) / tmp_modulus);
2405 }
2406
2407 PPL_ASSERT(OK());
2408 }
2409
2410 void
2411 PPL::Grid::
generalized_affine_preimage(const Linear_Expression & lhs,const Relation_Symbol relsym,const Linear_Expression & rhs,Coefficient_traits::const_reference modulus)2412 generalized_affine_preimage(const Linear_Expression& lhs,
2413 const Relation_Symbol relsym,
2414 const Linear_Expression& rhs,
2415 Coefficient_traits::const_reference modulus) {
2416 // The dimension of `lhs' must be at most the dimension of `*this'.
2417 dimension_type lhs_space_dim = lhs.space_dimension();
2418 if (space_dim < lhs_space_dim) {
2419 throw_dimension_incompatible("generalized_affine_preimage(e1, r, e2, m)",
2420 "lhs", lhs);
2421 }
2422 // The dimension of `rhs' must be at most the dimension of `*this'.
2423 const dimension_type rhs_space_dim = rhs.space_dimension();
2424 if (space_dim < rhs_space_dim) {
2425 throw_dimension_incompatible("generalized_affine_preimage(e1, r, e2, m)",
2426 "e2", rhs);
2427 }
2428 // The relation symbol cannot be a disequality.
2429 if (relsym == NOT_EQUAL) {
2430 throw_invalid_argument("generalized_affine_preimage(e1, r, e2, m)",
2431 "r is the disequality relation symbol");
2432 }
2433
2434 // Any preimage of an empty grid is empty.
2435 if (marked_empty()) {
2436 return;
2437 }
2438
2439 // If relsym is not EQUAL, then we return a safe approximation
2440 // by adding a line in the direction of var.
2441 if (relsym != EQUAL) {
2442
2443 if (modulus != 0) {
2444 throw_invalid_argument("generalized_affine_preimage(e1, r, e2, m)",
2445 "r != EQUAL && m != 0");
2446 }
2447
2448 if (!generators_are_up_to_date()) {
2449 minimize();
2450 }
2451 // Any image of an empty grid is empty.
2452 if (marked_empty()) {
2453 return;
2454 }
2455
2456 for (Linear_Expression::const_iterator i = lhs.begin(), i_end = lhs.end();
2457 i != i_end; ++i) {
2458 add_grid_generator(grid_line(i.variable()));
2459 }
2460
2461 PPL_ASSERT(OK());
2462 return;
2463 }
2464
2465 PPL_ASSERT(relsym == EQUAL);
2466
2467 PPL_DIRTY_TEMP_COEFFICIENT(tmp_modulus);
2468 tmp_modulus = modulus;
2469 if (tmp_modulus < 0) {
2470 neg_assign(tmp_modulus);
2471 }
2472
2473 // Compute the actual space dimension of `lhs',
2474 // i.e., the highest dimension having a non-zero coefficient in `lhs'.
2475 lhs_space_dim = lhs.last_nonzero();
2476 if (lhs_space_dim == 0) {
2477 // All variables have zero coefficients, so `lhs' is a constant.
2478 // In this case, preimage and image happen to be the same.
2479 add_congruence_no_check((lhs %= rhs) / tmp_modulus);
2480 return;
2481 }
2482
2483 // Gather in `new_lines' the collections of all the lines having
2484 // the direction of variables occurring in `lhs'.
2485 Grid_Generator_System new_lines;
2486 for (Linear_Expression::const_iterator i = lhs.begin(),
2487 i_end = lhs.lower_bound(Variable(lhs_space_dim)); i != i_end; ++i) {
2488 new_lines.insert(grid_line(i.variable()));
2489 }
2490
2491 const dimension_type num_common_dims
2492 = std::min(lhs_space_dim, rhs_space_dim);
2493 if (lhs.have_a_common_variable(rhs, Variable(0), Variable(num_common_dims))) {
2494 // Some variables in `lhs' also occur in `rhs'.
2495 // To ease the computation, add an additional dimension.
2496 const Variable new_var(space_dim);
2497 add_space_dimensions_and_embed(1);
2498
2499 // Constrain the new dimension to be equal to `lhs'
2500 // TODO: Use add_congruence() when it has been updated.
2501 Congruence_System new_cgs1(new_var == lhs);
2502 add_recycled_congruences(new_cgs1);
2503 if (!is_empty()) {
2504 // The grid still contains points.
2505
2506 // Existentially quantify all the variables occurring in the left
2507 // hand side
2508
2509 // Adjust `new_lines' to the right dimension.
2510 new_lines.set_space_dimension(space_dim);
2511 // Add the lines to `gen_sys' (first make sure they are up-to-date).
2512 update_generators();
2513 gen_sys.insert(new_lines, Recycle_Input());
2514 normalize_divisors(gen_sys);
2515 // Update the flags.
2516 clear_congruences_up_to_date();
2517 clear_generators_minimized();
2518
2519 // Constrain the new dimension so that it is related to
2520 // the right hand side modulo `modulus'.
2521 // TODO: Use add_congruence() when it has been updated.
2522 Congruence_System new_cgs2((rhs %= new_var) / tmp_modulus);
2523 add_recycled_congruences(new_cgs2);
2524 }
2525
2526 // Remove the temporarily added dimension.
2527 remove_higher_space_dimensions(space_dim-1);
2528 }
2529 else {
2530 // `lhs' and `rhs' variables are disjoint:
2531 // there is no need to add a further dimension.
2532
2533 // Constrain the left hand side expression so that it is congruent to
2534 // the right hand side expression modulo `mod'.
2535 add_congruence_no_check((lhs %= rhs) / tmp_modulus);
2536
2537 // Any image of an empty grid is empty.
2538 if (is_empty()) {
2539 return;
2540 }
2541
2542 // Existentially quantify all the variables occurring in `lhs'.
2543 add_recycled_grid_generators(new_lines);
2544 }
2545 PPL_ASSERT(OK());
2546 }
2547
2548 void
2549 PPL::Grid::
bounded_affine_image(const Variable var,const Linear_Expression & lb_expr,const Linear_Expression & ub_expr,Coefficient_traits::const_reference denominator)2550 bounded_affine_image(const Variable var,
2551 const Linear_Expression& lb_expr,
2552 const Linear_Expression& ub_expr,
2553 Coefficient_traits::const_reference denominator) {
2554
2555 // The denominator cannot be zero.
2556 if (denominator == 0) {
2557 throw_invalid_argument("bounded_affine_image(v, lb, ub, d)", "d == 0");
2558 }
2559 // Dimension-compatibility checks.
2560 // `var' should be one of the dimensions of the grid.
2561 const dimension_type var_space_dim = var.space_dimension();
2562 if (space_dim < var_space_dim) {
2563 throw_dimension_incompatible("bounded_affine_image(v, lb, ub, d)",
2564 "v", var);
2565 }
2566 // The dimension of `lb_expr' and `ub_expr' should not be
2567 // greater than the dimension of `*this'.
2568 const dimension_type lb_space_dim = lb_expr.space_dimension();
2569 if (space_dim < lb_space_dim) {
2570 throw_dimension_incompatible("bounded_affine_image(v, lb, ub, d)",
2571 "lb", lb_expr);
2572 }
2573 const dimension_type ub_space_dim = ub_expr.space_dimension();
2574 if (space_dim < ub_space_dim) {
2575 throw_dimension_incompatible("bounded_affine_image(v, lb, ub, d)",
2576 "ub", ub_expr);
2577 }
2578
2579 // Any image of an empty grid is empty.
2580 if (marked_empty()) {
2581 return;
2582 }
2583 // In all other cases, generalized_affine_preimage() must
2584 // just add a line in the direction of var.
2585 generalized_affine_image(var,
2586 LESS_OR_EQUAL,
2587 ub_expr,
2588 denominator);
2589
2590 PPL_ASSERT(OK());
2591 }
2592
2593
2594 void
2595 PPL::Grid::
bounded_affine_preimage(const Variable var,const Linear_Expression & lb_expr,const Linear_Expression & ub_expr,Coefficient_traits::const_reference denominator)2596 bounded_affine_preimage(const Variable var,
2597 const Linear_Expression& lb_expr,
2598 const Linear_Expression& ub_expr,
2599 Coefficient_traits::const_reference denominator) {
2600
2601 // The denominator cannot be zero.
2602 if (denominator == 0) {
2603 throw_invalid_argument("bounded_affine_preimage(v, lb, ub, d)", "d == 0");
2604 }
2605 // Dimension-compatibility checks.
2606 // `var' should be one of the dimensions of the grid.
2607 const dimension_type var_space_dim = var.space_dimension();
2608 if (space_dim < var_space_dim) {
2609 throw_dimension_incompatible("bounded_affine_preimage(v, lb, ub, d)",
2610 "v", var);
2611 }
2612 // The dimension of `lb_expr' and `ub_expr' should not be
2613 // greater than the dimension of `*this'.
2614 const dimension_type lb_space_dim = lb_expr.space_dimension();
2615 if (space_dim < lb_space_dim) {
2616 throw_dimension_incompatible("bounded_affine_preimage(v, lb, ub, d)",
2617 "lb", lb_expr);
2618 }
2619 const dimension_type ub_space_dim = ub_expr.space_dimension();
2620 if (space_dim < ub_space_dim) {
2621 throw_dimension_incompatible("bounded_affine_preimage(v, lb, ub, d)",
2622 "ub", ub_expr);
2623 }
2624
2625 // Any preimage of an empty grid is empty.
2626 if (marked_empty()) {
2627 return;
2628 }
2629
2630 // In all other cases, generalized_affine_preimage() must
2631 // just add a line in the direction of var.
2632 generalized_affine_preimage(var,
2633 LESS_OR_EQUAL,
2634 ub_expr,
2635 denominator);
2636
2637 PPL_ASSERT(OK());
2638 }
2639
2640 void
time_elapse_assign(const Grid & y)2641 PPL::Grid::time_elapse_assign(const Grid& y) {
2642 Grid& x = *this;
2643 // Check dimension-compatibility.
2644 if (x.space_dim != y.space_dim) {
2645 throw_dimension_incompatible("time_elapse_assign(y)", "y", y);
2646 }
2647
2648 // Deal with the zero-dimensional case.
2649 if (x.space_dim == 0) {
2650 if (y.marked_empty()) {
2651 x.set_empty();
2652 }
2653 return;
2654 }
2655
2656 // If either one of `x' or `y' is empty, the result is empty too.
2657 if (x.marked_empty()) {
2658 return;
2659 }
2660 if (y.marked_empty()
2661 || (!x.generators_are_up_to_date() && !x.update_generators())
2662 || (!y.generators_are_up_to_date() && !y.update_generators())) {
2663 x.set_empty();
2664 return;
2665 }
2666
2667 // At this point both generator systems are up-to-date.
2668 Grid_Generator_System gs = y.gen_sys;
2669 const dimension_type gs_num_rows = gs.num_rows();
2670
2671 normalize_divisors(gs, gen_sys);
2672
2673 for (dimension_type i = gs_num_rows; i-- > 0; ) {
2674 Grid_Generator& g = gs.sys.rows[i];
2675 if (g.is_point()) {
2676 // Transform the point into a parameter.
2677 g.set_is_parameter();
2678 }
2679 }
2680
2681 PPL_ASSERT(gs.sys.OK());
2682
2683 if (gs_num_rows == 0) {
2684 // `y' was the grid containing a single point at the origin, so
2685 // the result is `x'.
2686 return;
2687 }
2688
2689 // Append `gs' to the generators of `x'.
2690
2691 gen_sys.insert(gs, Recycle_Input());
2692
2693 x.clear_congruences_up_to_date();
2694 x.clear_generators_minimized();
2695
2696 PPL_ASSERT(x.OK(true) && y.OK(true));
2697 }
2698
2699 bool
frequency(const Linear_Expression & expr,Coefficient & freq_n,Coefficient & freq_d,Coefficient & val_n,Coefficient & val_d) const2700 PPL::Grid::frequency(const Linear_Expression& expr,
2701 Coefficient& freq_n, Coefficient& freq_d,
2702 Coefficient& val_n, Coefficient& val_d) const {
2703 // The dimension of `expr' must be at most the dimension of *this.
2704 if (space_dim < expr.space_dimension()) {
2705 throw_dimension_incompatible("frequency(e, ...)", "e", expr);
2706 }
2707
2708 // Space dimension is 0: if empty, then return false;
2709 // otherwise the frequency is 1 and the value is 0.
2710 if (space_dim == 0) {
2711 if (is_empty()) {
2712 return false;
2713 }
2714 freq_n = 0;
2715 freq_d = 1;
2716 val_n = 0;
2717 val_d = 1;
2718 return true;
2719 }
2720 if (!generators_are_minimized() && !minimize()) {
2721 // Minimizing found `this' empty.
2722 return false;
2723 }
2724
2725 return frequency_no_check(expr, freq_n, freq_d, val_n, val_d);
2726 }
2727
2728 /*! \relates Parma_Polyhedra_Library::Grid */
2729 bool
operator ==(const Grid & x,const Grid & y)2730 PPL::operator==(const Grid& x, const Grid& y) {
2731 if (x.space_dim != y.space_dim) {
2732 return false;
2733 }
2734 if (x.marked_empty()) {
2735 return y.is_empty();
2736 }
2737 if (y.marked_empty()) {
2738 return x.is_empty();
2739 }
2740 if (x.space_dim == 0) {
2741 return true;
2742 }
2743
2744 switch (x.quick_equivalence_test(y)) {
2745 case Grid::TVB_TRUE:
2746 return true;
2747
2748 case Grid::TVB_FALSE:
2749 return false;
2750
2751 default:
2752 if (x.is_included_in(y)) {
2753 if (x.marked_empty()) {
2754 return y.is_empty();
2755 }
2756 return y.is_included_in(x);
2757 }
2758 return false;
2759 }
2760 }
2761
2762 bool
contains(const Grid & y) const2763 PPL::Grid::contains(const Grid& y) const {
2764 const Grid& x = *this;
2765
2766 // Dimension-compatibility check.
2767 if (x.space_dim != y.space_dim) {
2768 throw_dimension_incompatible("contains(y)", "y", y);
2769 }
2770
2771 if (y.marked_empty()) {
2772 return true;
2773 }
2774 if (x.marked_empty()) {
2775 return y.is_empty();
2776 }
2777 if (y.space_dim == 0) {
2778 return true;
2779 }
2780 if (x.quick_equivalence_test(y) == Grid::TVB_TRUE) {
2781 return true;
2782 }
2783 return y.is_included_in(x);
2784 }
2785
2786 bool
is_disjoint_from(const Grid & y) const2787 PPL::Grid::is_disjoint_from(const Grid& y) const {
2788 // Dimension-compatibility check.
2789 if (space_dim != y.space_dim) {
2790 throw_dimension_incompatible("is_disjoint_from(y)", "y", y);
2791 }
2792 Grid z = *this;
2793 z.intersection_assign(y);
2794 return z.is_empty();
2795 }
2796
2797 void
ascii_dump(std::ostream & s) const2798 PPL::Grid::ascii_dump(std::ostream& s) const {
2799 using std::endl;
2800
2801 s << "space_dim "
2802 << space_dim
2803 << endl;
2804 status.ascii_dump(s);
2805 s << "con_sys ("
2806 << (congruences_are_up_to_date() ? "" : "not_")
2807 << "up-to-date)"
2808 << endl;
2809 con_sys.ascii_dump(s);
2810 s << "gen_sys ("
2811 << (generators_are_up_to_date() ? "" : "not_")
2812 << "up-to-date)"
2813 << endl;
2814 gen_sys.ascii_dump(s);
2815 s << "dimension_kinds";
2816 if ((generators_are_up_to_date() && generators_are_minimized())
2817 || (congruences_are_up_to_date() && congruences_are_minimized())) {
2818 for (Dimension_Kinds::const_iterator i = dim_kinds.begin();
2819 i != dim_kinds.end();
2820 ++i) {
2821 s << " " << *i;
2822 }
2823 }
2824 s << endl;
2825 }
2826
PPL_OUTPUT_DEFINITIONS(Grid)2827 PPL_OUTPUT_DEFINITIONS(Grid)
2828
2829 bool
2830 PPL::Grid::ascii_load(std::istream& s) {
2831 std::string str;
2832
2833 if (!(s >> str) || str != "space_dim") {
2834 return false;
2835 }
2836
2837 if (!(s >> space_dim)) {
2838 return false;
2839 }
2840
2841 if (!status.ascii_load(s)) {
2842 return false;
2843 }
2844
2845 if (!(s >> str) || str != "con_sys") {
2846 return false;
2847 }
2848
2849 if (!(s >> str)) {
2850 return false;
2851 }
2852 if (str == "(up-to-date)") {
2853 set_congruences_up_to_date();
2854 }
2855 else if (str != "(not_up-to-date)") {
2856 return false;
2857 }
2858
2859 if (!con_sys.ascii_load(s)) {
2860 return false;
2861 }
2862
2863 if (!(s >> str) || str != "gen_sys") {
2864 return false;
2865 }
2866
2867 if (!(s >> str)) {
2868 return false;
2869 }
2870 if (str == "(up-to-date)") {
2871 set_generators_up_to_date();
2872 }
2873 else if (str != "(not_up-to-date)") {
2874 return false;
2875 }
2876
2877 if (!gen_sys.ascii_load(s)) {
2878 return false;
2879 }
2880
2881 if (!(s >> str) || str != "dimension_kinds") {
2882 return false;
2883 }
2884
2885 if (!marked_empty()
2886 && ((generators_are_up_to_date() && generators_are_minimized())
2887 || (congruences_are_up_to_date() && congruences_are_minimized()))) {
2888 dim_kinds.resize(space_dim + 1);
2889 for (Dimension_Kinds::size_type dim = 0; dim <= space_dim; ++dim) {
2890 short unsigned int dim_kind;
2891 if (!(s >> dim_kind)) {
2892 return false;
2893 }
2894 switch (dim_kind) {
2895 case 0:
2896 dim_kinds[dim] = PARAMETER;
2897 break;
2898 case 1:
2899 dim_kinds[dim] = LINE;
2900 break;
2901 case 2:
2902 dim_kinds[dim] = GEN_VIRTUAL;
2903 break;
2904 default:
2905 return false;
2906 }
2907 }
2908 }
2909
2910 // Check invariants.
2911 PPL_ASSERT(OK());
2912 return true;
2913 }
2914
2915 PPL::memory_size_type
external_memory_in_bytes() const2916 PPL::Grid::external_memory_in_bytes() const {
2917 return
2918 con_sys.external_memory_in_bytes()
2919 + gen_sys.external_memory_in_bytes();
2920 }
2921
2922 void
wrap_assign(const Variables_Set & vars,Bounded_Integer_Type_Width w,Bounded_Integer_Type_Representation r,Bounded_Integer_Type_Overflow o,const Constraint_System * cs_p,unsigned,bool)2923 PPL::Grid::wrap_assign(const Variables_Set& vars,
2924 Bounded_Integer_Type_Width w,
2925 Bounded_Integer_Type_Representation r,
2926 Bounded_Integer_Type_Overflow o,
2927 const Constraint_System* cs_p,
2928 unsigned /* complexity_threshold */,
2929 bool /* wrap_individually */) {
2930
2931 // Dimension-compatibility check of `*cs_p', if any.
2932 if (cs_p != 0) {
2933 const dimension_type cs_p_space_dim = cs_p->space_dimension();
2934 if (cs_p->space_dimension() > space_dim) {
2935 throw_dimension_incompatible("wrap_assign(vs, ...)", cs_p_space_dim);
2936 }
2937 }
2938
2939 // Wrapping no variable is a no-op.
2940 if (vars.empty()) {
2941 return;
2942 }
2943
2944 // Dimension-compatibility check of `vars'.
2945 const dimension_type min_space_dim = vars.space_dimension();
2946 if (space_dim < min_space_dim) {
2947 throw_dimension_incompatible("wrap_assign(vs, ...)", min_space_dim);
2948 }
2949
2950 // Wrapping an empty grid is a no-op.
2951 if (marked_empty()) {
2952 return;
2953 }
2954 if (!generators_are_minimized() && !minimize()) {
2955 // Minimizing found `this' empty.
2956 return;
2957 }
2958
2959 // Set the wrap frequency for variables of width `w'.
2960 // This is independent of the signedness `s'.
2961 PPL_DIRTY_TEMP_COEFFICIENT(wrap_frequency);
2962 mul_2exp_assign(wrap_frequency, Coefficient_one(), w);
2963 // Set `min_value' and `max_value' to the minimum and maximum values
2964 // a variable of width `w' and signedness `s' can take.
2965 PPL_DIRTY_TEMP_COEFFICIENT(min_value);
2966 PPL_DIRTY_TEMP_COEFFICIENT(max_value);
2967 if (r == UNSIGNED) {
2968 min_value = 0;
2969 mul_2exp_assign(max_value, Coefficient_one(), w);
2970 --max_value;
2971 }
2972 else {
2973 PPL_ASSERT(r == SIGNED_2_COMPLEMENT);
2974 mul_2exp_assign(max_value, Coefficient_one(), w-1);
2975 neg_assign(min_value, max_value);
2976 --max_value;
2977 }
2978
2979 // Generators are up-to-date and minimized.
2980 const Grid gr = *this;
2981
2982 // Overflow is impossible or wraps.
2983 if (o == OVERFLOW_IMPOSSIBLE || o == OVERFLOW_WRAPS) {
2984 PPL_DIRTY_TEMP_COEFFICIENT(f_n);
2985 PPL_DIRTY_TEMP_COEFFICIENT(f_d);
2986 PPL_DIRTY_TEMP_COEFFICIENT(v_n);
2987 PPL_DIRTY_TEMP_COEFFICIENT(v_d);
2988 for (Variables_Set::const_iterator i = vars.begin(),
2989 vars_end = vars.end(); i != vars_end; ++i) {
2990 const Variable x(*i);
2991 // Find the frequency and a value for `x' in `gr'.
2992 if (!gr.frequency_no_check(x, f_n, f_d, v_n, v_d)) {
2993 continue;
2994 }
2995 if (f_n == 0) {
2996
2997 // `x' is a constant in `gr'.
2998 if (v_d != 1) {
2999 // The value for `x' is not integral (`frequency_no_check()'
3000 // ensures that `v_n' and `v_d' have no common divisors).
3001 set_empty();
3002 return;
3003 }
3004
3005 // `x' is a constant and has an integral value.
3006 if ((v_n > max_value) || (v_n < min_value)) {
3007 // The value is outside the range of the bounded integer type.
3008 if (o == OVERFLOW_IMPOSSIBLE) {
3009 // Then `x' has no possible value and hence set empty.
3010 set_empty();
3011 return;
3012 }
3013 PPL_ASSERT(o == OVERFLOW_WRAPS);
3014 // The value v_n for `x' is wrapped modulo the 'wrap_frequency'.
3015 v_n %= wrap_frequency;
3016 // `v_n' is the value closest to 0 and may be negative.
3017 if (r == UNSIGNED && v_n < 0) {
3018 v_n += wrap_frequency;
3019 }
3020 unconstrain(x);
3021 add_constraint(x == v_n);
3022 }
3023 continue;
3024 }
3025
3026 // `x' is not a constant in `gr'.
3027 PPL_ASSERT(f_n != 0);
3028
3029 if (f_d % v_d != 0) {
3030 // Then `x' has no integral value and hence `gr' is set empty.
3031 set_empty();
3032 return;
3033 }
3034 if (f_d != 1) {
3035 // `x' has non-integral values, so add the integrality
3036 // congruence for `x'.
3037 add_congruence((x %= 0) / 1);
3038 }
3039 if (o == OVERFLOW_WRAPS && f_n != wrap_frequency) {
3040 // We know that `x' is not a constant, so, if overflow wraps,
3041 // `x' may wrap to a value modulo the `wrap_frequency'.
3042 add_grid_generator(parameter(wrap_frequency * x));
3043 }
3044 else if ((o == OVERFLOW_IMPOSSIBLE && 2*f_n >= wrap_frequency)
3045 || (f_n == wrap_frequency)) {
3046 // In these cases, `x' can only take a unique (ie constant)
3047 // value.
3048 if (r == UNSIGNED && v_n < 0) {
3049 // `v_n' is the value closest to 0 and may be negative.
3050 v_n += f_n;
3051 }
3052 unconstrain(x);
3053 add_constraint(x == v_n);
3054 }
3055 else {
3056 // If overflow is impossible but the grid frequency is less than
3057 // half the wrap frequency, then there is more than one possible
3058 // value for `x' in the range of the bounded integer type,
3059 // so the grid is unchanged.
3060 PPL_ASSERT(o == OVERFLOW_IMPOSSIBLE && 2*f_n < wrap_frequency);
3061 }
3062 }
3063 return;
3064 }
3065
3066 PPL_ASSERT(o == OVERFLOW_UNDEFINED);
3067 // If overflow is undefined, then all we know is that the variable
3068 // may take any integer within the range of the bounded integer type.
3069 const Grid_Generator& point = gr.gen_sys[0];
3070 const Coefficient& div = point.divisor();
3071 max_value *= div;
3072 min_value *= div;
3073 for (Variables_Set::const_iterator i = vars.begin(),
3074 vars_end = vars.end(); i != vars_end; ++i) {
3075 const Variable x(*i);
3076 if (!gr.bounds_no_check(x)) {
3077 // `x' is not a constant in `gr'.
3078 // We know that `x' is not a constant, so `x' may wrap to any
3079 // value `x + z' where z is an integer.
3080 if (point.coefficient(x) % div != 0) {
3081 // We know that `x' can take non-integral values, so enforce
3082 // integrality.
3083 unconstrain(x);
3084 add_congruence(x %= 0);
3085 }
3086 else {
3087 // `x' has at least one integral value.
3088 // `x' may also take other non-integral values,
3089 // but checking could be costly so we ignore this.
3090 add_grid_generator(parameter(x));
3091 }
3092 }
3093 else {
3094 // `x' is a constant `v' in `gr'.
3095 const Coefficient& coeff_x = point.coefficient(x);
3096 // `x' should be integral.
3097 if (coeff_x % div != 0) {
3098 set_empty();
3099 return;
3100 }
3101 // If the value `v' for `x' is not within the range for the
3102 // bounded integer type, then `x' may wrap to any value `v + z'
3103 // where `z' is an integer; otherwise `x' is unchanged.
3104 if (coeff_x > max_value || coeff_x < min_value) {
3105 add_grid_generator(parameter(x));
3106 }
3107 }
3108 }
3109 }
3110
3111 void
drop_some_non_integer_points(Complexity_Class)3112 PPL::Grid::drop_some_non_integer_points(Complexity_Class) {
3113 if (marked_empty() || space_dim == 0) {
3114 return;
3115 }
3116
3117 // By adding integral congruences for each dimension to the
3118 // congruence system, defining \p *this, the grid will keep only
3119 // those points that have integral coordinates. All points in \p
3120 // *this with non-integral coordinates are removed.
3121 for (dimension_type i = space_dim; i-- > 0; ) {
3122 add_congruence(Variable(i) %= 0);
3123 }
3124
3125 PPL_ASSERT(OK());
3126 }
3127
3128 void
drop_some_non_integer_points(const Variables_Set & vars,Complexity_Class)3129 PPL::Grid::drop_some_non_integer_points(const Variables_Set& vars,
3130 Complexity_Class) {
3131 // Dimension-compatibility check.
3132 const dimension_type min_space_dim = vars.space_dimension();
3133 if (space_dimension() < min_space_dim) {
3134 throw_dimension_incompatible("drop_some_non_integer_points(vs, cmpl)",
3135 min_space_dim);
3136 }
3137
3138 if (marked_empty() || min_space_dim == 0) {
3139 return;
3140 }
3141
3142 // By adding the integral congruences for each dimension in vars to
3143 // the congruence system defining \p *this, the grid will keep only
3144 // those points that have integer coordinates for all the dimensions
3145 // in vars. All points in \p *this with non-integral coordinates for
3146 // the dimensions in vars are removed.
3147 for (Variables_Set::const_iterator i = vars.begin(),
3148 vars_end = vars.end(); i != vars_end; ++i) {
3149 add_congruence(Variable(*i) %= 0);
3150 }
3151
3152 PPL_ASSERT(OK());
3153 }
3154
3155 /*! \relates Parma_Polyhedra_Library::Grid */
3156 std::ostream&
operator <<(std::ostream & s,const Grid & gr)3157 PPL::IO_Operators::operator<<(std::ostream& s, const Grid& gr) {
3158 if (gr.is_empty()) {
3159 s << "false";
3160 }
3161 else if (gr.is_universe()) {
3162 s << "true";
3163 }
3164 else {
3165 s << gr.minimized_congruences();
3166 }
3167 return s;
3168 }
3169