00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023 #include <config.h>
00024
00025 #include "Constraint_System.defs.hh"
00026 #include "Constraint_System.inlines.hh"
00027 #include "Generator.defs.hh"
00028 #include "Scalar_Products.defs.hh"
00029 #include <cassert>
00030 #include <string>
00031 #include <vector>
00032 #include <iostream>
00033 #include <stdexcept>
00034
00035 namespace PPL = Parma_Polyhedra_Library;
00036
00037 bool
00038 PPL::Constraint_System::
00039 adjust_topology_and_space_dimension(const Topology new_topology,
00040 const dimension_type new_space_dim) {
00041 assert(space_dimension() <= new_space_dim);
00042
00043 const dimension_type old_space_dim = space_dimension();
00044 const Topology old_topology = topology();
00045 dimension_type cols_to_be_added = new_space_dim - old_space_dim;
00046
00047
00048 if (num_rows() == 0) {
00049 if (num_columns() == 0)
00050 if (new_topology == NECESSARILY_CLOSED) {
00051 add_zero_columns(++cols_to_be_added);
00052 set_necessarily_closed();
00053 }
00054 else {
00055 cols_to_be_added += 2;
00056 add_zero_columns(cols_to_be_added);
00057 set_not_necessarily_closed();
00058 }
00059 else
00060
00061 if (old_topology != new_topology)
00062 if (new_topology == NECESSARILY_CLOSED) {
00063 switch (cols_to_be_added) {
00064 case 0:
00065 remove_trailing_columns(1);
00066 break;
00067 case 1:
00068
00069 break;
00070 default:
00071 add_zero_columns(--cols_to_be_added);
00072 }
00073 set_necessarily_closed();
00074 }
00075 else {
00076
00077
00078 add_zero_columns(++cols_to_be_added);
00079 set_not_necessarily_closed();
00080 }
00081 else {
00082
00083 if (cols_to_be_added > 0)
00084 add_zero_columns(cols_to_be_added);
00085 }
00086 assert(OK());
00087 return true;
00088 }
00089
00090
00091 if (cols_to_be_added > 0)
00092 if (old_topology != new_topology)
00093 if (new_topology == NECESSARILY_CLOSED) {
00094
00095
00096
00097 if (has_strict_inequalities())
00098 return false;
00099
00100
00101
00102
00103
00104
00105 Constraint_System& cs = *this;
00106 const dimension_type eps_index = old_space_dim + 1;
00107 dimension_type cs_num_rows = cs.num_rows();
00108 bool was_sorted = cs.is_sorted();
00109 if (was_sorted)
00110 cs.set_sorted(false);
00111
00112
00113
00114 if (cs.num_pending_rows() == 0) {
00115 for (dimension_type i = cs_num_rows; i-- > 0; )
00116 if (cs[i][eps_index] != 0) {
00117 --cs_num_rows;
00118 std::swap(cs[i], cs[cs_num_rows]);
00119 }
00120 cs.erase_to_end(cs_num_rows);
00121 cs.unset_pending_rows();
00122 }
00123 else {
00124
00125
00126
00127
00128 const dimension_type old_first_pending = cs.first_pending_row();
00129 dimension_type new_first_pending = old_first_pending;
00130 for (dimension_type i = new_first_pending; i-- > 0; )
00131 if (cs[i][eps_index] != 0) {
00132 --new_first_pending;
00133 std::swap(cs[i], cs[new_first_pending]);
00134 }
00135 const dimension_type num_swaps
00136 = old_first_pending - new_first_pending;
00137 cs.set_index_first_pending_row(new_first_pending);
00138
00139 for (dimension_type i = num_swaps; i-- > 0; )
00140 std::swap(cs[old_first_pending - i], cs[cs_num_rows - i]);
00141 cs_num_rows -= num_swaps;
00142
00143 for (dimension_type i = cs_num_rows; i-- > new_first_pending; )
00144 if (cs[i][eps_index] != 0) {
00145 --cs_num_rows;
00146 std::swap(cs[i], cs[cs_num_rows]);
00147 }
00148 cs.erase_to_end(cs_num_rows);
00149 }
00150
00151
00152 if (was_sorted)
00153 cs.sort_rows();
00154 if (--cols_to_be_added > 0)
00155 add_zero_columns(cols_to_be_added);
00156 set_necessarily_closed();
00157 }
00158 else {
00159
00160
00161
00162 add_zero_columns(++cols_to_be_added);
00163 set_not_necessarily_closed();
00164 }
00165 else {
00166
00167 add_zero_columns(cols_to_be_added);
00168
00169
00170 if (old_topology == NOT_NECESSARILY_CLOSED)
00171 swap_columns(old_space_dim + 1, new_space_dim + 1);
00172 }
00173 else
00174
00175 if (old_topology != new_topology)
00176 if (new_topology == NECESSARILY_CLOSED) {
00177
00178
00179
00180 if (has_strict_inequalities())
00181 return false;
00182
00183 remove_trailing_columns(1);
00184 set_necessarily_closed();
00185 }
00186 else {
00187
00188 add_zero_columns(1);
00189 set_not_necessarily_closed();
00190 }
00191
00192 assert(OK());
00193 return true;
00194 }
00195
00196 bool
00197 PPL::Constraint_System::has_strict_inequalities() const {
00198 if (is_necessarily_closed())
00199 return false;
00200 const Constraint_System& cs = *this;
00201 const dimension_type eps_index = cs.num_columns() - 1;
00202
00203
00204 for (dimension_type i = cs.num_rows(); i-- > 0; ) {
00205 const Constraint& c = cs[i];
00206
00207
00208
00209
00210 if (c[eps_index] < 0 && !c.is_tautological())
00211 return true;
00212 }
00213 return false;
00214 }
00215
00216 void
00217 PPL::Constraint_System::insert(const Constraint& c) {
00218
00219
00220 assert(num_pending_rows() == 0);
00221 if (topology() == c.topology())
00222 Linear_System::insert(c);
00223 else
00224
00225 if (is_necessarily_closed()) {
00226
00227
00228 add_zero_columns(1);
00229 set_not_necessarily_closed();
00230 Linear_System::insert(c);
00231 }
00232 else {
00233
00234
00235
00236
00237
00238 const dimension_type new_size = 2 + std::max(c.space_dimension(),
00239 space_dimension());
00240 Constraint tmp_c(c, new_size);
00241 tmp_c.set_not_necessarily_closed();
00242 Linear_System::insert(tmp_c);
00243 }
00244 assert(OK());
00245 }
00246
00247 void
00248 PPL::Constraint_System::insert_pending(const Constraint& c) {
00249 if (topology() == c.topology())
00250 Linear_System::insert_pending(c);
00251 else
00252
00253 if (is_necessarily_closed()) {
00254
00255
00256 add_zero_columns(1);
00257 set_not_necessarily_closed();
00258 Linear_System::insert_pending(c);
00259 }
00260 else {
00261
00262
00263
00264 const dimension_type new_size = 2 + std::max(c.space_dimension(),
00265 space_dimension());
00266 Constraint tmp_c(c, new_size);
00267 tmp_c.set_not_necessarily_closed();
00268 Linear_System::insert_pending(tmp_c);
00269 }
00270 assert(OK());
00271 }
00272
00273 PPL::dimension_type
00274 PPL::Constraint_System::num_inequalities() const {
00275
00276
00277 assert(num_pending_rows() == 0);
00278 const Constraint_System& cs = *this;
00279 dimension_type n = 0;
00280
00281
00282 if (is_sorted())
00283 for (dimension_type i = num_rows(); i > 0 && cs[--i].is_inequality(); )
00284 ++n;
00285 else
00286 for (dimension_type i = num_rows(); i-- > 0 ; )
00287 if (cs[i].is_inequality())
00288 ++n;
00289 return n;
00290 }
00291
00292 PPL::dimension_type
00293 PPL::Constraint_System::num_equalities() const {
00294
00295
00296 assert(num_pending_rows() == 0);
00297 return num_rows() - num_inequalities();
00298 }
00299
00300 void
00301 PPL::Constraint_System::const_iterator::skip_forward() {
00302 const Linear_System::const_iterator csp_end = csp->end();
00303 while (i != csp_end && (*this)->is_tautological())
00304 ++i;
00305 }
00306
00307 bool
00308 PPL::Constraint_System::satisfies_all_constraints(const Generator& g) const {
00309 assert(g.space_dimension() <= space_dimension());
00310
00311
00312
00313
00314 Topology_Adjusted_Scalar_Product_Sign sps(g);
00315
00316 const Constraint_System& cs = *this;
00317 if (cs.is_necessarily_closed()) {
00318 if (g.is_line()) {
00319
00320 for (dimension_type i = cs.num_rows(); i-- > 0; )
00321 if (sps(g, cs[i]) != 0)
00322 return false;
00323 }
00324 else
00325
00326 for (dimension_type i = cs.num_rows(); i-- > 0; ) {
00327 const Constraint& c = cs[i];
00328 const int sp_sign = sps(g, c);
00329 if (c.is_inequality()) {
00330
00331
00332 if (sp_sign < 0)
00333 return false;
00334 }
00335 else
00336
00337 if (sp_sign != 0)
00338 return false;
00339 }
00340 }
00341 else
00342
00343 switch (g.type()) {
00344
00345 case Generator::LINE:
00346
00347 for (dimension_type i = cs.num_rows(); i-- > 0; )
00348 if (sps(g, cs[i]) != 0)
00349 return false;
00350 break;
00351
00352 case Generator::POINT:
00353
00354
00355 for (dimension_type i = cs.num_rows(); i-- > 0; ) {
00356 const Constraint& c = cs[i];
00357 const int sp_sign = sps(g, c);
00358 switch (c.type()) {
00359 case Constraint::EQUALITY:
00360 if (sp_sign != 0)
00361 return false;
00362 break;
00363 case Constraint::NONSTRICT_INEQUALITY:
00364 if (sp_sign < 0)
00365 return false;
00366 break;
00367 case Constraint::STRICT_INEQUALITY:
00368 if (sp_sign <= 0)
00369 return false;
00370 break;
00371 }
00372 }
00373 break;
00374
00375 case Generator::RAY:
00376
00377 case Generator::CLOSURE_POINT:
00378 for (dimension_type i = cs.num_rows(); i-- > 0; ) {
00379 const Constraint& c = cs[i];
00380 const int sp_sign = sps(g, c);
00381 if (c.is_inequality()) {
00382
00383 if (sp_sign < 0)
00384 return false;
00385 }
00386 else
00387
00388 if (sp_sign != 0)
00389 return false;
00390 }
00391 break;
00392 }
00393
00394
00395 return true;
00396 }
00397
00398
00399 void
00400 PPL::Constraint_System
00401 ::affine_preimage(const dimension_type v,
00402 const Linear_Expression& expr,
00403 Coefficient_traits::const_reference denominator) {
00404 Constraint_System& x = *this;
00405
00406
00407
00408 assert(v > 0 && v <= x.space_dimension());
00409 assert(expr.space_dimension() <= x.space_dimension());
00410 assert(denominator > 0);
00411
00412 const dimension_type n_columns = x.num_columns();
00413 const dimension_type n_rows = x.num_rows();
00414 const dimension_type expr_size = expr.size();
00415 const bool not_invertible = (v >= expr_size || expr[v] == 0);
00416
00417 if (denominator != 1)
00418 for (dimension_type i = n_rows; i-- > 0; ) {
00419 Constraint& row = x[i];
00420 Coefficient& row_v = row[v];
00421 if (row_v != 0) {
00422 for (dimension_type j = n_columns; j-- > 0; )
00423 if (j != v) {
00424 Coefficient& row_j = row[j];
00425 row_j *= denominator;
00426 if (j < expr_size)
00427 add_mul_assign(row_j, row_v, expr[j]);
00428 }
00429 if (not_invertible)
00430 row_v = 0;
00431 else
00432 row_v *= expr[v];
00433 }
00434 }
00435 else
00436
00437
00438 for (dimension_type i = n_rows; i-- > 0; ) {
00439 Constraint& row = x[i];
00440 Coefficient& row_v = row[v];
00441 if (row_v != 0) {
00442 for (dimension_type j = expr_size; j-- > 0; )
00443 if (j != v)
00444 add_mul_assign(row[j], row_v, expr[j]);
00445 if (not_invertible)
00446 row_v = 0;
00447 else
00448 row_v *= expr[v];
00449 }
00450 }
00451
00452 x.strong_normalize();
00453 }
00454
00455 void
00456 PPL::Constraint_System::ascii_dump(std::ostream& s) const {
00457 const Constraint_System& x = *this;
00458 const dimension_type x_num_rows = x.num_rows();
00459 const dimension_type x_num_columns = x.num_columns();
00460 s << "topology " << (is_necessarily_closed()
00461 ? "NECESSARILY_CLOSED"
00462 : "NOT_NECESSARILY_CLOSED")
00463 << "\n"
00464 << x_num_rows << " x " << x_num_columns << ' '
00465 << (x.is_sorted() ? "(sorted)" : "(not_sorted)")
00466 << "\n"
00467 << "index_first_pending " << x.first_pending_row()
00468 << "\n";
00469 for (dimension_type i = 0; i < x_num_rows; ++i) {
00470 const Constraint& c = x[i];
00471 for (dimension_type j = 0; j < x_num_columns; ++j)
00472 s << c[j] << ' ';
00473 switch (c.type()) {
00474 case Constraint::EQUALITY:
00475 s << "=";
00476 break;
00477 case Constraint::NONSTRICT_INEQUALITY:
00478 s << ">=";
00479 break;
00480 case Constraint::STRICT_INEQUALITY:
00481 s << ">";
00482 break;
00483 }
00484 s << "\n";
00485 }
00486 }
00487
00488 PPL_OUTPUT_DEFINITIONS(Constraint_System);
00489
00490 bool
00491 PPL::Constraint_System::ascii_load(std::istream& s) {
00492 std::string str;
00493 if (!(s >> str) || str != "topology")
00494 return false;
00495 if (!(s >> str))
00496 return false;
00497 if (str == "NECESSARILY_CLOSED")
00498 set_necessarily_closed();
00499 else {
00500 if (str != "NOT_NECESSARILY_CLOSED")
00501 return false;
00502 set_not_necessarily_closed();
00503 }
00504
00505 dimension_type nrows;
00506 dimension_type ncols;
00507 if (!(s >> nrows))
00508 return false;
00509 if (!(s >> str))
00510 return false;
00511 if (!(s >> ncols))
00512 return false;
00513 resize_no_copy(nrows, ncols);
00514
00515 if (!(s >> str) || (str != "(sorted)" && str != "(not_sorted)"))
00516 return false;
00517 set_sorted(str == "(sorted)");
00518 dimension_type index;
00519 if (!(s >> str) || str != "index_first_pending")
00520 return false;
00521 if (!(s >> index))
00522 return false;
00523 set_index_first_pending_row(index);
00524
00525 Constraint_System& x = *this;
00526 for (dimension_type i = 0; i < x.num_rows(); ++i) {
00527 for (dimension_type j = 0; j < x.num_columns(); ++j)
00528 if (!(s >> x[i][j]))
00529 return false;
00530
00531 if (!(s >> str))
00532 return false;
00533 if (str == "=")
00534 x[i].set_is_equality();
00535 else
00536 x[i].set_is_inequality();
00537
00538
00539 switch (x[i].type()) {
00540 case Constraint::EQUALITY:
00541 if (str == "=")
00542 continue;
00543 break;
00544 case Constraint::NONSTRICT_INEQUALITY:
00545 if (str == ">=")
00546 continue;
00547 break;
00548 case Constraint::STRICT_INEQUALITY:
00549 if (str == ">")
00550 continue;
00551 break;
00552 }
00553
00554 return false;
00555 }
00556
00557 assert(OK());
00558 return true;
00559 }
00560
00561 bool
00562 PPL::Constraint_System::OK() const {
00563
00564
00565
00566 if (!Linear_System::OK(false))
00567 return false;
00568
00569
00570 const Constraint_System& x = *this;
00571 for (dimension_type i = num_rows(); i-- > 0; )
00572 if (!x[i].OK())
00573 return false;
00574
00575
00576 return true;
00577 }
00578
00580 std::ostream&
00581 PPL::IO_Operators::operator<<(std::ostream& s, const Constraint_System& cs) {
00582 Constraint_System::const_iterator i = cs.begin();
00583 const Constraint_System::const_iterator cs_end = cs.end();
00584 if (i == cs_end)
00585 s << "true";
00586 else {
00587 while (i != cs_end) {
00588 s << *i++;
00589 if (i != cs_end)
00590 s << ", ";
00591 }
00592 }
00593 return s;
00594 }