00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023 #ifndef PPL_Polyhedra_Powerset_templates_hh
00024 #define PPL_Polyhedra_Powerset_templates_hh 1
00025
00026 #include "Constraint.defs.hh"
00027 #include "Constraint_System.defs.hh"
00028 #include "Constraint_System.inlines.hh"
00029 #include "C_Polyhedron.defs.hh"
00030 #include "NNC_Polyhedron.defs.hh"
00031 #include <algorithm>
00032 #include <deque>
00033 #include <string>
00034 #include <iostream>
00035 #include <sstream>
00036 #include <stdexcept>
00037
00038 namespace Parma_Polyhedra_Library {
00039
00040 template <typename PH>
00041 void
00042 Polyhedra_Powerset<PH>::add_disjunct(const PH& ph) {
00043 Polyhedra_Powerset& x = *this;
00044 if (x.space_dimension() != ph.space_dimension()) {
00045 std::ostringstream s;
00046 s << "PPL::Polyhedra_Powerset<PH>::add_disjunct(ph):\n"
00047 << "this->space_dimension() == " << x.space_dimension() << ", "
00048 << "ph.space_dimension() == " << ph.space_dimension() << ".";
00049 throw std::invalid_argument(s.str());
00050 }
00051 x.sequence.push_back(Determinate<PH>(ph));
00052 x.reduced = false;
00053 assert(x.OK());
00054 }
00055
00056 template <>
00057 template <typename QH>
00058 Polyhedra_Powerset<NNC_Polyhedron>
00059 ::Polyhedra_Powerset(const Polyhedra_Powerset<QH>& y)
00060 : Base(), space_dim(y.space_dimension()) {
00061 Polyhedra_Powerset& x = *this;
00062 for (typename Polyhedra_Powerset<QH>::const_iterator i = y.begin(),
00063 y_end = y.end(); i != y_end; ++i)
00064 x.sequence.push_back(Determinate<NNC_Polyhedron>(
00065 NNC_Polyhedron(i->element().constraints()))
00066 );
00067 x.reduced = y.reduced;
00068 assert(x.OK());
00069 }
00070
00071 template <>
00072 template <typename QH>
00073 Polyhedra_Powerset<C_Polyhedron>
00074 ::Polyhedra_Powerset(const Polyhedra_Powerset<QH>& y)
00075 : Base(), space_dim(y.space_dimension()) {
00076 Polyhedra_Powerset& x = *this;
00077 for (typename Polyhedra_Powerset<QH>::const_iterator i = y.begin(),
00078 y_end = y.end(); i != y_end; ++i)
00079 x.sequence.push_back(Determinate<C_Polyhedron>(
00080 C_Polyhedron(i->element().constraints()))
00081 );
00082
00083
00084
00085
00086 x.reduced = false;
00087 assert(x.OK());
00088 }
00089
00090 template <typename PH>
00091 void
00092 Polyhedra_Powerset<PH>::concatenate_assign(const Polyhedra_Powerset& y) {
00093 Polyhedra_Powerset& x = *this;
00094
00095 x.omega_reduce();
00096 y.omega_reduce();
00097 Polyhedra_Powerset<PH> new_x(x.space_dim + y.space_dim, EMPTY);
00098 for (const_iterator xi = x.begin(), x_end = x.end(),
00099 y_begin = y.begin(), y_end = y.end(); xi != x_end; ) {
00100 for (const_iterator yi = y_begin; yi != y_end; ++yi) {
00101 CS zi = *xi;
00102 zi.concatenate_assign(*yi);
00103 assert(!zi.is_bottom());
00104 new_x.sequence.push_back(zi);
00105 }
00106 ++xi;
00107 if (abandon_expensive_computations && xi != x_end && y_begin != y_end) {
00108
00109 PH xph = xi->element();
00110 for (++xi; xi != x_end; ++xi)
00111 xph.upper_bound_assign(xi->element());
00112 const_iterator yi = y_begin;
00113 PH yph = yi->element();
00114 for (++yi; yi != y_end; ++yi)
00115 yph.upper_bound_assign(yi->element());
00116 xph.concatenate_assign(yph);
00117 x.swap(new_x);
00118 x.add_disjunct(xph);
00119 assert(x.OK());
00120 return;
00121 }
00122 }
00123 x.swap(new_x);
00124 assert(x.OK());
00125 }
00126
00127 template <typename PH>
00128 void
00129 Polyhedra_Powerset<PH>::add_constraint(const Constraint& c) {
00130 Polyhedra_Powerset& x = *this;
00131 for (Sequence_iterator si = x.sequence.begin(),
00132 s_end = x.sequence.end(); si != s_end; ++si)
00133 si->element().add_constraint(c);
00134 x.reduced = false;
00135 assert(x.OK());
00136 }
00137
00138 template <typename PH>
00139 bool
00140 Polyhedra_Powerset<PH>::add_constraint_and_minimize(const Constraint& c) {
00141 Polyhedra_Powerset& x = *this;
00142 for (Sequence_iterator si = x.sequence.begin(),
00143 s_end = x.sequence.end(); si != s_end; )
00144 if (!si->element().add_constraint_and_minimize(c))
00145 si = x.sequence.erase(si);
00146 else {
00147 x.reduced = false;
00148 ++si;
00149 }
00150 assert(x.OK());
00151 return !x.empty();
00152 }
00153
00154 template <typename PH>
00155 void
00156 Polyhedra_Powerset<PH>::add_constraints(const Constraint_System& cs) {
00157 Polyhedra_Powerset& x = *this;
00158 for (Sequence_iterator si = x.sequence.begin(),
00159 s_end = x.sequence.end(); si != s_end; ++si)
00160 si->element().add_constraints(cs);
00161 x.reduced = false;
00162 assert(x.OK());
00163 }
00164
00165 template <typename PH>
00166 bool
00167 Polyhedra_Powerset<PH>::
00168 add_constraints_and_minimize(const Constraint_System& cs) {
00169 Polyhedra_Powerset& x = *this;
00170 for (Sequence_iterator si = x.sequence.begin(),
00171 s_end = x.sequence.end(); si != s_end; )
00172 if (!si->element().add_constraints_and_minimize(cs))
00173 si = x.sequence.erase(si);
00174 else {
00175 x.reduced = false;
00176 ++si;
00177 }
00178 assert(x.OK());
00179 return !x.empty();
00180 }
00181
00182 template <typename PH>
00183 void
00184 Polyhedra_Powerset<PH>::add_space_dimensions_and_embed(dimension_type m) {
00185 Polyhedra_Powerset& x = *this;
00186 for (Sequence_iterator si = x.sequence.begin(),
00187 s_end = x.sequence.end(); si != s_end; ++si)
00188 si->element().add_space_dimensions_and_embed(m);
00189 x.space_dim += m;
00190 assert(x.OK());
00191 }
00192
00193 template <typename PH>
00194 void
00195 Polyhedra_Powerset<PH>::add_space_dimensions_and_project(dimension_type m) {
00196 Polyhedra_Powerset& x = *this;
00197 for (Sequence_iterator si = x.sequence.begin(),
00198 s_end = x.sequence.end(); si != s_end; ++si)
00199 si->element().add_space_dimensions_and_project(m);
00200 x.space_dim += m;
00201 assert(x.OK());
00202 }
00203
00204 template <typename PH>
00205 void
00206 Polyhedra_Powerset<PH>::
00207 remove_space_dimensions(const Variables_Set& to_be_removed) {
00208 Polyhedra_Powerset& x = *this;
00209 Variables_Set::size_type num_removed = to_be_removed.size();
00210 if (num_removed > 0) {
00211 for (Sequence_iterator si = x.sequence.begin(),
00212 s_end = x.sequence.end(); si != s_end; ++si) {
00213 si->element().remove_space_dimensions(to_be_removed);
00214 x.reduced = false;
00215 }
00216 x.space_dim -= num_removed;
00217 assert(x.OK());
00218 }
00219 }
00220
00221 template <typename PH>
00222 void
00223 Polyhedra_Powerset<PH>::remove_higher_space_dimensions(dimension_type
00224 new_dimension) {
00225 Polyhedra_Powerset& x = *this;
00226 if (new_dimension < x.space_dim) {
00227 for (Sequence_iterator si = x.sequence.begin(),
00228 s_end = x.sequence.end(); si != s_end; ++si) {
00229 si->element().remove_higher_space_dimensions(new_dimension);
00230 x.reduced = false;
00231 }
00232 x.space_dim = new_dimension;
00233 assert(x.OK());
00234 }
00235 }
00236
00237 template <typename PH>
00238 template <typename Partial_Function>
00239 void
00240 Polyhedra_Powerset<PH>::map_space_dimensions(const Partial_Function& pfunc) {
00241 Polyhedra_Powerset& x = *this;
00242 if (x.is_bottom()) {
00243 dimension_type n = 0;
00244 for (dimension_type i = x.space_dim; i-- > 0; ) {
00245 dimension_type new_i;
00246 if (pfunc.maps(i, new_i))
00247 ++n;
00248 }
00249 x.space_dim = n;
00250 }
00251 else {
00252 Sequence_iterator s_begin = x.sequence.begin();
00253 for (Sequence_iterator si = s_begin,
00254 s_end = x.sequence.end(); si != s_end; ++si)
00255 si->element().map_space_dimensions(pfunc);
00256 x.space_dim = s_begin->element().space_dimension();
00257 x.reduced = false;
00258 }
00259 assert(x.OK());
00260 }
00261
00262 template <typename PH>
00263 void
00264 Polyhedra_Powerset<PH>::pairwise_reduce() {
00265 Polyhedra_Powerset& x = *this;
00266
00267 x.omega_reduce();
00268
00269 size_type n = x.size();
00270 size_type deleted;
00271 do {
00272 Polyhedra_Powerset new_x(x.space_dim, EMPTY);
00273 std::deque<bool> marked(n, false);
00274 deleted = 0;
00275 Sequence_iterator s_begin = x.sequence.begin();
00276 Sequence_iterator s_end = x.sequence.end();
00277 unsigned si_index = 0;
00278 for (Sequence_iterator si = s_begin; si != s_end; ++si, ++si_index) {
00279 if (marked[si_index])
00280 continue;
00281 PH& pi = si->element();
00282 Sequence_const_iterator sj = si;
00283 unsigned sj_index = si_index;
00284 for (++sj, ++sj_index; sj != s_end; ++sj, ++sj_index) {
00285 if (marked[sj_index])
00286 continue;
00287 const PH& pj = sj->element();
00288 if (pi.upper_bound_assign_if_exact(pj)) {
00289 marked[si_index] = marked[sj_index] = true;
00290 new_x.add_non_bottom_disjunct(pi);
00291 ++deleted;
00292 goto next;
00293 }
00294 }
00295 next:
00296 ;
00297 }
00298 iterator nx_begin = new_x.begin();
00299 iterator nx_end = new_x.end();
00300 unsigned xi_index = 0;
00301 for (const_iterator xi = x.begin(),
00302 x_end = x.end(); xi != x_end; ++xi, ++xi_index)
00303 if (!marked[xi_index])
00304 nx_begin = new_x.add_non_bottom_disjunct(*xi, nx_begin, nx_end);
00305 std::swap(x.sequence, new_x.sequence);
00306 n -= deleted;
00307 } while (deleted > 0);
00308 assert(x.OK());
00309 }
00310
00311 template <typename PH>
00312 template <typename Widening>
00313 void
00314 Polyhedra_Powerset<PH>::
00315 BGP99_heuristics_assign(const Polyhedra_Powerset& y, Widening wf) {
00316
00317 Polyhedra_Powerset& x = *this;
00318
00319 #ifndef NDEBUG
00320 {
00321
00322 const Polyhedra_Powerset<PH> x_copy = x;
00323 const Polyhedra_Powerset<PH> y_copy = y;
00324 assert(y_copy.definitely_entails(x_copy));
00325 }
00326 #endif
00327
00328 size_type n = x.size();
00329 Polyhedra_Powerset new_x(x.space_dim, EMPTY);
00330 std::deque<bool> marked(n, false);
00331 const_iterator x_begin = x.begin();
00332 const_iterator x_end = x.end();
00333 unsigned i_index = 0;
00334 for (const_iterator i = x_begin,
00335 y_begin = y.begin(), y_end = y.end(); i != x_end; ++i, ++i_index)
00336 for (const_iterator j = y_begin; j != y_end; ++j) {
00337 const PH& pi = i->element();
00338 const PH& pj = j->element();
00339 if (pi.contains(pj)) {
00340 PH pi_copy = pi;
00341 wf(pi_copy, pj);
00342 new_x.add_non_bottom_disjunct(pi_copy);
00343 marked[i_index] = true;
00344 }
00345 }
00346 iterator nx_begin = new_x.begin();
00347 iterator nx_end = new_x.end();
00348 i_index = 0;
00349 for (const_iterator i = x_begin; i != x_end; ++i, ++i_index)
00350 if (!marked[i_index])
00351 nx_begin = new_x.add_non_bottom_disjunct(*i, nx_begin, nx_end);
00352 std::swap(x.sequence, new_x.sequence);
00353 assert(x.OK());
00354 assert(x.is_omega_reduced());
00355 }
00356
00357 template <typename PH>
00358 template <typename Widening>
00359 void
00360 Polyhedra_Powerset<PH>::
00361 BGP99_extrapolation_assign(const Polyhedra_Powerset& y,
00362 Widening wf,
00363 unsigned max_disjuncts) {
00364
00365 Polyhedra_Powerset& x = *this;
00366
00367 #ifndef NDEBUG
00368 {
00369
00370 const Polyhedra_Powerset<PH> x_copy = x;
00371 const Polyhedra_Powerset<PH> y_copy = y;
00372 assert(y_copy.definitely_entails(x_copy));
00373 }
00374 #endif
00375
00376 x.pairwise_reduce();
00377 if (max_disjuncts != 0)
00378 x.collapse(max_disjuncts);
00379 x.BGP99_heuristics_assign(y, wf);
00380 }
00381
00382 template <typename PH>
00383 template <typename Cert>
00384 void
00385 Polyhedra_Powerset<PH>::
00386 collect_certificates(std::map<Cert, size_type,
00387 typename Cert::Compare>& cert_ms) const {
00388 const Polyhedra_Powerset& x = *this;
00389 assert(x.is_omega_reduced());
00390 assert(cert_ms.size() == 0);
00391 for (const_iterator i = x.begin(), end = x.end(); i != end; i++) {
00392 Cert ph_cert(i->element());
00393 ++cert_ms[ph_cert];
00394 }
00395 }
00396
00397 template <typename PH>
00398 template <typename Cert>
00399 bool
00400 Polyhedra_Powerset<PH>::
00401 is_cert_multiset_stabilizing(const std::map<Cert, size_type,
00402 typename Cert::Compare>& y_cert_ms
00403 ) const {
00404 typedef std::map<Cert, size_type, typename Cert::Compare> Cert_Multiset;
00405 Cert_Multiset x_cert_ms;
00406 collect_certificates(x_cert_ms);
00407 typename Cert_Multiset::const_iterator
00408 xi = x_cert_ms.begin(),
00409 xend = x_cert_ms.end(),
00410 yi = y_cert_ms.begin(),
00411 yend = y_cert_ms.end();
00412 while (xi != xend && yi != yend) {
00413 const Cert& xi_cert = xi->first;
00414 const Cert& yi_cert = yi->first;
00415 switch (xi_cert.compare(yi_cert)) {
00416 case 0:
00417
00418 {
00419 const size_type& xi_count = xi->second;
00420 const size_type& yi_count = yi->second;
00421 if (xi_count == yi_count) {
00422
00423 ++xi;
00424 ++yi;
00425 }
00426 else
00427
00428 return xi_count < yi_count;
00429 break;
00430 }
00431 case 1:
00432
00433 return false;
00434
00435 case -1:
00436
00437 return true;
00438 }
00439 }
00440
00441
00442 return yi != yend;
00443 }
00444
00445 template <typename PH>
00446 template <typename Cert, typename Widening>
00447 void
00448 Polyhedra_Powerset<PH>::BHZ03_widening_assign(const Polyhedra_Powerset& y,
00449 Widening wf) {
00450
00451 Polyhedra_Powerset& x = *this;
00452
00453 #ifndef NDEBUG
00454 {
00455
00456 const Polyhedra_Powerset<PH> x_copy = x;
00457 const Polyhedra_Powerset<PH> y_copy = y;
00458 assert(y_copy.definitely_entails(x_copy));
00459 }
00460 #endif
00461
00462
00463
00464
00465 assert(x.size() > 0);
00466 if (y.size() == 0)
00467 return;
00468
00469
00470 PH x_hull(x.space_dim, EMPTY);
00471 for (const_iterator i = x.begin(), x_end = x.end(); i != x_end; ++i)
00472 x_hull.upper_bound_assign(i->element());
00473
00474
00475 PH y_hull(y.space_dim, EMPTY);
00476 for (const_iterator i = y.begin(), y_end = y.end(); i != y_end; ++i)
00477 y_hull.upper_bound_assign(i->element());
00478
00479 const Cert y_hull_cert(y_hull);
00480
00481
00482 int hull_stabilization = y_hull_cert.compare(x_hull);
00483 if (hull_stabilization == 1)
00484 return;
00485
00486
00487 const bool y_is_not_a_singleton = y.size() > 1;
00488
00489
00490
00491 typedef std::map<Cert, size_type, typename Cert::Compare> Cert_Multiset;
00492 Cert_Multiset y_cert_ms;
00493 bool y_cert_ms_computed = false;
00494
00495 if (hull_stabilization == 0 && y_is_not_a_singleton) {
00496
00497 y.collect_certificates(y_cert_ms);
00498 y_cert_ms_computed = true;
00499
00500 if (x.is_cert_multiset_stabilizing(y_cert_ms))
00501 return;
00502 }
00503
00504
00505 Polyhedra_Powerset<PH> bgp99_heuristics = x;
00506 bgp99_heuristics.BGP99_heuristics_assign(y, wf);
00507
00508
00509 PH bgp99_heuristics_hull(x.space_dim, EMPTY);
00510 for (const_iterator i = bgp99_heuristics.begin(),
00511 bh_end = bgp99_heuristics.end(); i != bh_end; ++i)
00512 bgp99_heuristics_hull.upper_bound_assign(i->element());
00513
00514
00515
00516 hull_stabilization = y_hull_cert.compare(bgp99_heuristics_hull);
00517 if (hull_stabilization == 1) {
00518
00519 std::swap(x, bgp99_heuristics);
00520 return;
00521 }
00522 else if (hull_stabilization == 0 && y_is_not_a_singleton) {
00523
00524 if (!y_cert_ms_computed) {
00525 y.collect_certificates(y_cert_ms);
00526 y_cert_ms_computed = true;
00527 }
00528 if (bgp99_heuristics.is_cert_multiset_stabilizing(y_cert_ms)) {
00529 std::swap(x, bgp99_heuristics);
00530 return;
00531 }
00532
00533
00534
00535
00536 Polyhedra_Powerset<PH> reduced_bgp99_heuristics(bgp99_heuristics);
00537 reduced_bgp99_heuristics.pairwise_reduce();
00538 if (reduced_bgp99_heuristics.is_cert_multiset_stabilizing(y_cert_ms)) {
00539 std::swap(x, reduced_bgp99_heuristics);
00540 return;
00541 }
00542 }
00543
00544
00545
00546 if (bgp99_heuristics_hull.strictly_contains(y_hull)) {
00547
00548 PH ph = bgp99_heuristics_hull;
00549 wf(ph, y_hull);
00550
00551 ph.difference_assign(bgp99_heuristics_hull);
00552 x.add_disjunct(ph);
00553 return;
00554 }
00555
00556
00557 Polyhedra_Powerset<PH> x_hull_singleton(x.space_dim, EMPTY);
00558 x_hull_singleton.add_disjunct(x_hull);
00559 std::swap(x, x_hull_singleton);
00560 }
00561
00562 template <typename PH>
00563 void
00564 Polyhedra_Powerset<PH>::ascii_dump(std::ostream& s) const {
00565 const Polyhedra_Powerset& x = *this;
00566 s << "size " << x.size()
00567 << "\nspace_dim " << x.space_dim
00568 << "\n";
00569 for (const_iterator xi = x.begin(), x_end = x.end(); xi != x_end; ++xi)
00570 xi->element().ascii_dump(s);
00571 }
00572
00573 PPL_OUTPUT_TEMPLATE_DEFINITIONS(PH, Polyhedra_Powerset<PH>);
00574
00575 template <typename PH>
00576 bool
00577 Polyhedra_Powerset<PH>::ascii_load(std::istream& s) {
00578 Polyhedra_Powerset& x = *this;
00579 std::string str;
00580
00581 if (!(s >> str) || str != "size")
00582 return false;
00583
00584 size_type sz;
00585
00586 if (!(s >> sz))
00587 return false;
00588
00589 if (!(s >> str) || str != "space_dim")
00590 return false;
00591
00592 if (!(s >> x.space_dim))
00593 return false;
00594
00595 Polyhedra_Powerset new_x(x.space_dim, EMPTY);
00596 while (sz-- > 0) {
00597 PH ph;
00598 if (!ph.ascii_load(s))
00599 return false;
00600 new_x.add_disjunct(ph);
00601 }
00602 x.swap(new_x);
00603
00604
00605 assert(x.OK());
00606 return true;
00607 }
00608
00609 template <typename PH>
00610 bool
00611 Polyhedra_Powerset<PH>::OK() const {
00612 const Polyhedra_Powerset& x = *this;
00613 for (const_iterator xi = x.begin(), x_end = x.end(); xi != x_end; ++xi) {
00614 const PH& pi = xi->element();
00615 if (pi.space_dimension() != x.space_dim) {
00616 #ifndef NDEBUG
00617 std::cerr << "Space dimension mismatch: is " << pi.space_dimension()
00618 << " in an element of the sequence,\nshould be "
00619 << x.space_dim << "."
00620 << std::endl;
00621 #endif
00622 return false;
00623 }
00624 }
00625 return x.Base::OK();
00626 }
00627
00628
00629 namespace Implementation {
00630
00631 namespace Polyhedra_Powersets {
00632
00633 #ifdef PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS
00635
00640 #endif // PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS
00641 template <typename PH>
00642 void
00643 linear_partition_aux(const Constraint& c,
00644 PH& qq,
00645 Polyhedra_Powerset<NNC_Polyhedron>& r) {
00646 Linear_Expression le(c);
00647 Constraint neg_c = c.is_strict_inequality() ? (le <= 0) : (le < 0);
00648 NNC_Polyhedron qqq(qq);
00649 if (qqq.add_constraint_and_minimize(neg_c))
00650 r.add_disjunct(qqq);
00651 qq.add_constraint(c);
00652 }
00653
00654 }
00655
00656 }
00657
00658
00660 template <typename PH>
00661 std::pair<PH, Polyhedra_Powerset<NNC_Polyhedron> >
00662 linear_partition(const PH& p, const PH& q) {
00663 using Implementation::Polyhedra_Powersets::linear_partition_aux;
00664
00665 Polyhedra_Powerset<NNC_Polyhedron> r(p.space_dimension(), EMPTY);
00666 PH qq = q;
00667 const Constraint_System& pcs = p.constraints();
00668 for (Constraint_System::const_iterator i = pcs.begin(),
00669 pcs_end = pcs.end(); i != pcs_end; ++i) {
00670 const Constraint c = *i;
00671 if (c.is_equality()) {
00672 Linear_Expression le(c);
00673 linear_partition_aux(le <= 0, qq, r);
00674 linear_partition_aux(le >= 0, qq, r);
00675 }
00676 else
00677 linear_partition_aux(c, qq, r);
00678 }
00679 return std::pair<PH, Polyhedra_Powerset<NNC_Polyhedron> >(qq, r);
00680 }
00681
00682 }
00683
00684 #endif // !defined(PPL_Polyhedra_Powerset_templates_hh)