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 "Linear_System.defs.hh"
00026 #include "Coefficient.defs.hh"
00027 #include "Row.defs.hh"
00028 #include "Saturation_Matrix.defs.hh"
00029 #include "Scalar_Products.defs.hh"
00030 #include <algorithm>
00031 #include <iostream>
00032 #include <string>
00033 #include <deque>
00034
00035 #include "swapping_sort.icc"
00036
00037 namespace PPL = Parma_Polyhedra_Library;
00038
00039 PPL::dimension_type
00040 PPL::Linear_System::num_lines_or_equalities() const {
00041 assert(num_pending_rows() == 0);
00042 const Linear_System& x = *this;
00043 dimension_type n = 0;
00044 for (dimension_type i = num_rows(); i-- > 0; )
00045 if (x[i].is_line_or_equality())
00046 ++n;
00047 return n;
00048 }
00049
00050 void
00051 PPL::Linear_System::merge_rows_assign(const Linear_System& y) {
00052 assert(row_size >= y.row_size);
00053
00054 assert(check_sorted() && y.check_sorted());
00055 assert(num_pending_rows() == 0 && y.num_pending_rows() == 0);
00056
00057 Linear_System& x = *this;
00058
00059
00060 std::vector<Row> tmp;
00061
00062 tmp.reserve(compute_capacity(x.num_rows() + y.num_rows(), max_num_rows()));
00063
00064 dimension_type xi = 0;
00065 dimension_type x_num_rows = x.num_rows();
00066 dimension_type yi = 0;
00067 dimension_type y_num_rows = y.num_rows();
00068
00069 while (xi < x_num_rows && yi < y_num_rows) {
00070 const int comp = compare(x[xi], y[yi]);
00071 if (comp <= 0) {
00072
00073 std::swap(x[xi++], *tmp.insert(tmp.end(), Linear_Row()));
00074 if (comp == 0)
00075
00076 ++yi;
00077 }
00078 else {
00079
00080 Linear_Row copy(y[yi++], row_size, row_capacity);
00081 std::swap(copy, *tmp.insert(tmp.end(), Linear_Row()));
00082 }
00083 }
00084
00085 if (xi < x_num_rows)
00086 while (xi < x_num_rows)
00087 std::swap(x[xi++], *tmp.insert(tmp.end(), Linear_Row()));
00088 else
00089 while (yi < y_num_rows) {
00090 Linear_Row copy(y[yi++], row_size, row_capacity);
00091 std::swap(copy, *tmp.insert(tmp.end(), Linear_Row()));
00092 }
00093
00094
00095 std::swap(tmp, rows);
00096
00097 unset_pending_rows();
00098 assert(check_sorted());
00099 }
00100
00101 void
00102 PPL::Linear_System::set_rows_topology() {
00103 Linear_System& x = *this;
00104 if (is_necessarily_closed())
00105 for (dimension_type i = num_rows(); i-- > 0; )
00106 x[i].set_necessarily_closed();
00107 else
00108 for (dimension_type i = num_rows(); i-- > 0; )
00109 x[i].set_not_necessarily_closed();
00110 }
00111
00112 void
00113 PPL::Linear_System::ascii_dump(std::ostream& s) const {
00114
00115
00116
00117
00118 const Linear_System& x = *this;
00119 dimension_type x_num_rows = x.num_rows();
00120 dimension_type x_num_columns = x.num_columns();
00121 s << "topology " << (is_necessarily_closed()
00122 ? "NECESSARILY_CLOSED"
00123 : "NOT_NECESSARILY_CLOSED")
00124 << "\n"
00125 << x_num_rows << " x " << x_num_columns
00126 << (x.sorted ? "(sorted)" : "(not_sorted)")
00127 << "\n"
00128 << "index_first_pending " << x.first_pending_row()
00129 << "\n";
00130 for (dimension_type i = 0; i < x_num_rows; ++i)
00131 x[i].ascii_dump(s);
00132 }
00133
00134 PPL_OUTPUT_DEFINITIONS_ASCII_ONLY(Linear_System);
00135
00136 bool
00137 PPL::Linear_System::ascii_load(std::istream& s) {
00138 std::string str;
00139 if (!(s >> str) || str != "topology")
00140 return false;
00141 if (!(s >> str))
00142 return false;
00143 if (str == "NECESSARILY_CLOSED")
00144 set_necessarily_closed();
00145 else {
00146 if (str != "NOT_NECESSARILY_CLOSED")
00147 return false;
00148 set_not_necessarily_closed();
00149 }
00150
00151 dimension_type nrows;
00152 dimension_type ncols;
00153 if (!(s >> nrows))
00154 return false;
00155 if (!(s >> str))
00156 return false;
00157 if (!(s >> ncols))
00158 return false;
00159 resize_no_copy(nrows, ncols);
00160
00161 if (!(s >> str) || (str != "(sorted)" && str != "(not_sorted)"))
00162 return false;
00163 set_sorted(str == "(sorted)");
00164 dimension_type index;
00165 if (!(s >> str) || str != "index_first_pending")
00166 return false;
00167 if (!(s >> index))
00168 return false;
00169 set_index_first_pending_row(index);
00170
00171 Linear_System& x = *this;
00172 for (dimension_type row = 0; row < nrows; ++row)
00173 if (!x[row].ascii_load(s))
00174 return false;
00175
00176
00177 assert(OK(true));
00178 return true;
00179 }
00180
00181 void
00182 PPL::Linear_System::insert(const Linear_Row& r) {
00183
00184
00185 assert(r.check_strong_normalized());
00186 assert(topology() == r.topology());
00187
00188 assert(num_pending_rows() == 0);
00189
00190 const dimension_type old_num_rows = num_rows();
00191 const dimension_type old_num_columns = num_columns();
00192 const dimension_type r_size = r.size();
00193
00194
00195 if (r_size > old_num_columns) {
00196 add_zero_columns(r_size - old_num_columns);
00197 if (!is_necessarily_closed() && old_num_rows != 0)
00198
00199
00200 swap_columns(old_num_columns - 1, r_size - 1);
00201 add_row(r);
00202 }
00203 else if (r_size < old_num_columns)
00204 if (is_necessarily_closed() || old_num_rows == 0)
00205 add_row(Linear_Row(r, old_num_columns, row_capacity));
00206 else {
00207
00208
00209 Linear_Row tmp_row(r, old_num_columns, row_capacity);
00210 std::swap(tmp_row[r_size - 1], tmp_row[old_num_columns - 1]);
00211 add_row(tmp_row);
00212 }
00213 else
00214
00215 add_row(r);
00216
00217
00218 assert(num_pending_rows() == 0);
00219
00220
00221 assert(OK(false));
00222 }
00223
00224 void
00225 PPL::Linear_System::insert_pending(const Linear_Row& r) {
00226
00227
00228 assert(r.check_strong_normalized());
00229 assert(topology() == r.topology());
00230
00231 const dimension_type old_num_rows = num_rows();
00232 const dimension_type old_num_columns = num_columns();
00233 const dimension_type r_size = r.size();
00234
00235
00236 if (r_size > old_num_columns) {
00237 add_zero_columns(r_size - old_num_columns);
00238 if (!is_necessarily_closed() && old_num_rows != 0)
00239
00240
00241 swap_columns(old_num_columns - 1, r_size - 1);
00242 add_pending_row(r);
00243 }
00244 else if (r_size < old_num_columns)
00245 if (is_necessarily_closed() || old_num_rows == 0)
00246 add_pending_row(Linear_Row(r, old_num_columns, row_capacity));
00247 else {
00248
00249
00250 Linear_Row tmp_row(r, old_num_columns, row_capacity);
00251 std::swap(tmp_row[r_size - 1], tmp_row[old_num_columns - 1]);
00252 add_pending_row(tmp_row);
00253 }
00254 else
00255
00256 add_pending_row(r);
00257
00258
00259 assert(num_pending_rows() > 0);
00260
00261
00262 assert(OK(false));
00263 }
00264
00265 void
00266 PPL::Linear_System::add_pending_rows(const Linear_System& y) {
00267 Linear_System& x = *this;
00268 assert(x.row_size == y.row_size);
00269
00270 const dimension_type x_n_rows = x.num_rows();
00271 const dimension_type y_n_rows = y.num_rows();
00272
00273 const bool was_sorted = sorted;
00274 add_zero_rows(y_n_rows, Linear_Row::Flags(row_topology));
00275 sorted = was_sorted;
00276
00277
00278 for (dimension_type i = y_n_rows; i-- > 0; ) {
00279 Row copy(y[i], x.row_size, x.row_capacity);
00280 std::swap(copy, x[x_n_rows+i]);
00281 }
00282
00283
00284 assert(OK(false));
00285 }
00286
00287 void
00288 PPL::Linear_System::add_rows(const Linear_System& y) {
00289 assert(num_pending_rows() == 0);
00290
00291
00292 if (y.num_rows() == 0)
00293 return;
00294
00295
00296 if (is_sorted())
00297 if (!y.is_sorted() || y.num_pending_rows() > 0)
00298 set_sorted(false);
00299 else {
00300
00301 const dimension_type n_rows = num_rows();
00302 if (n_rows > 0)
00303 set_sorted(compare((*this)[n_rows-1], y[0]) <= 0);
00304 }
00305
00306
00307 add_pending_rows(y);
00308
00309 unset_pending_rows();
00310
00311
00312
00313 assert(OK(false));
00314 }
00315
00316 void
00317 PPL::Linear_System::sort_rows() {
00318 const dimension_type num_pending = num_pending_rows();
00319
00320 sort_rows(0, first_pending_row());
00321 set_index_first_pending_row(num_rows() - num_pending);
00322 sorted = true;
00323
00324
00325 assert(OK(false));
00326 }
00327
00328 void
00329 PPL::Linear_System::sort_rows(const dimension_type first_row,
00330 const dimension_type last_row) {
00331 assert(first_row <= last_row && last_row <= num_rows());
00332
00333 assert(first_row >= first_pending_row() || last_row <= first_pending_row());
00334
00335
00336 std::vector<Row>::iterator first = rows.begin() + first_row;
00337 std::vector<Row>::iterator last = rows.begin() + last_row;
00338 swapping_sort(first, last, Row_Less_Than());
00339
00340 std::vector<Row>::iterator new_last = swapping_unique(first, last);
00341
00342 rows.erase(new_last, last);
00343
00344
00345 }
00346
00347 void
00348 PPL::Linear_System::add_row(const Linear_Row& r) {
00349
00350
00351 assert(r.check_strong_normalized());
00352 assert(r.size() == row_size);
00353
00354 assert(num_pending_rows() == 0);
00355
00356 const bool was_sorted = is_sorted();
00357
00358 Matrix::add_row(r);
00359
00360
00361
00362 set_index_first_pending_row(num_rows());
00363
00364 if (was_sorted) {
00365 const dimension_type nrows = num_rows();
00366
00367 if (nrows > 1) {
00368
00369
00370
00371 Linear_System& x = *this;
00372 set_sorted(compare(x[nrows-2], x[nrows-1]) <= 0);
00373 }
00374 else
00375
00376 set_sorted(true);
00377 }
00378
00379 assert(num_pending_rows() == 0);
00380
00381
00382 assert(OK(false));
00383 }
00384
00385 void
00386 PPL::Linear_System::add_pending_row(const Linear_Row& r) {
00387
00388
00389 assert(r.check_strong_normalized());
00390 assert(r.size() == row_size);
00391
00392 const dimension_type new_rows_size = rows.size() + 1;
00393 if (rows.capacity() < new_rows_size) {
00394
00395 std::vector<Row> new_rows;
00396 new_rows.reserve(compute_capacity(new_rows_size, max_num_rows()));
00397 new_rows.insert(new_rows.end(), new_rows_size, Row());
00398
00399 Row new_row(r, row_capacity);
00400 dimension_type i = new_rows_size-1;
00401 std::swap(new_rows[i], new_row);
00402
00403 while (i-- > 0)
00404 new_rows[i].swap(rows[i]);
00405
00406 std::swap(rows, new_rows);
00407 }
00408 else {
00409
00410
00411
00412 Row tmp(r, row_capacity);
00413 std::swap(*rows.insert(rows.end(), Row()), tmp);
00414 }
00415
00416
00417 assert(num_pending_rows() > 0);
00418
00419
00420 assert(OK(false));
00421 }
00422
00423 void
00424 PPL::Linear_System::add_pending_row(const Linear_Row::Flags flags) {
00425 const dimension_type new_rows_size = rows.size() + 1;
00426 if (rows.capacity() < new_rows_size) {
00427
00428 std::vector<Row> new_rows;
00429 new_rows.reserve(compute_capacity(new_rows_size, max_num_rows()));
00430 new_rows.insert(new_rows.end(), new_rows_size, Row());
00431
00432 Linear_Row new_row(row_size, row_capacity, flags);
00433 dimension_type i = new_rows_size-1;
00434 std::swap(new_rows[i], new_row);
00435
00436 while (i-- > 0)
00437 new_rows[i].swap(rows[i]);
00438
00439 std::swap(rows, new_rows);
00440 }
00441 else {
00442
00443
00444
00445 Row& new_row = *rows.insert(rows.end(), Row());
00446 static_cast<Linear_Row&>(new_row).construct(row_size, row_capacity, flags);
00447 }
00448
00449
00450 assert(num_pending_rows() > 0);
00451 }
00452
00453 void
00454 PPL::Linear_System::normalize() {
00455 Linear_System& x = *this;
00456
00457 for (dimension_type i = num_rows(); i-- > 0; )
00458 x[i].normalize();
00459 set_sorted(false);
00460 }
00461
00462 void
00463 PPL::Linear_System::strong_normalize() {
00464 Linear_System& x = *this;
00465
00466 for (dimension_type i = num_rows(); i-- > 0; )
00467 x[i].strong_normalize();
00468 set_sorted(false);
00469 }
00470
00471 void
00472 PPL::Linear_System::sign_normalize() {
00473 Linear_System& x = *this;
00474
00475 for (dimension_type i = num_rows(); i-- > 0; )
00476 x[i].sign_normalize();
00477 set_sorted(false);
00478 }
00479
00481 bool
00482 PPL::operator==(const Linear_System& x, const Linear_System& y) {
00483 if (x.num_columns() != y.num_columns())
00484 return false;
00485 const dimension_type x_num_rows = x.num_rows();
00486 const dimension_type y_num_rows = y.num_rows();
00487 if (x_num_rows != y_num_rows)
00488 return false;
00489 if (x.first_pending_row() != y.first_pending_row())
00490 return false;
00491
00492
00493
00494 for (dimension_type i = x_num_rows; i-- > 0; )
00495 if (x[i] != y[i])
00496 return false;
00497 return true;
00498 }
00499
00500 void
00501 PPL::Linear_System::sort_and_remove_with_sat(Saturation_Matrix& sat) {
00502 Linear_System& sys = *this;
00503
00504 assert(sys.first_pending_row() == sat.num_rows());
00505 if (sys.first_pending_row() <= 1) {
00506 sys.set_sorted(true);
00507 return;
00508 }
00509
00510
00511 With_Saturation_Matrix_iterator first(sys.rows.begin(), sat.rows.begin());
00512 With_Saturation_Matrix_iterator last = first + sat.num_rows();
00513 swapping_sort(first, last, Row_Less_Than());
00514
00515 With_Saturation_Matrix_iterator new_last = swapping_unique(first, last);
00516
00517 const dimension_type num_duplicates = last - new_last;
00518 const dimension_type new_first_pending_row
00519 = sys.first_pending_row() - num_duplicates;
00520
00521 if (sys.num_pending_rows() > 0) {
00522
00523 const dimension_type n_rows = sys.num_rows() - 1;
00524 for (dimension_type i = 0; i < num_duplicates; ++i)
00525 std::swap(sys[new_first_pending_row + i], sys[n_rows - i]);
00526 }
00527
00528 sys.erase_to_end(sys.num_rows() - num_duplicates);
00529 sys.set_index_first_pending_row(new_first_pending_row);
00530
00531 sat.rows_erase_to_end(sat.num_rows() - num_duplicates);
00532 assert(sys.check_sorted());
00533
00534 sys.set_sorted(true);
00535 }
00536
00537 PPL::dimension_type
00538 PPL::Linear_System::gauss(const dimension_type n_lines_or_equalities) {
00539 Linear_System& x = *this;
00540
00541
00542
00543
00544 assert(x.OK(true));
00545 assert(x.num_pending_rows() == 0);
00546 assert(n_lines_or_equalities == x.num_lines_or_equalities());
00547 #ifndef NDEBUG
00548 for (dimension_type i = n_lines_or_equalities; i-- > 0; )
00549 assert(x[i].is_line_or_equality());
00550 #endif
00551
00552 dimension_type rank = 0;
00553
00554 bool changed = false;
00555 for (dimension_type j = x.num_columns(); j-- > 0; )
00556 for (dimension_type i = rank; i < n_lines_or_equalities; ++i) {
00557
00558
00559 if (x[i][j] == 0)
00560 continue;
00561
00562
00563 if (i > rank) {
00564 std::swap(x[i], x[rank]);
00565
00566 changed = true;
00567 }
00568
00569
00570
00571 for (dimension_type k = i + 1; k < n_lines_or_equalities; ++k)
00572 if (x[k][j] != 0) {
00573 x[k].linear_combine(x[rank], j);
00574 changed = true;
00575 }
00576
00577 ++rank;
00578
00579 break;
00580 }
00581 if (changed)
00582 x.set_sorted(false);
00583
00584 assert(x.OK(true));
00585 return rank;
00586 }
00587
00588 void
00589 PPL::Linear_System
00590 ::back_substitute(const dimension_type n_lines_or_equalities) {
00591 Linear_System& x = *this;
00592
00593
00594
00595
00596 assert(x.OK(true));
00597 assert(x.num_pending_rows() == 0);
00598 assert(n_lines_or_equalities <= x.num_lines_or_equalities());
00599 #ifndef NDEBUG
00600 for (dimension_type i = n_lines_or_equalities; i-- > 0; )
00601 assert(x[i].is_line_or_equality());
00602 #endif
00603
00604 const dimension_type nrows = x.num_rows();
00605 const dimension_type ncols = x.num_columns();
00606
00607 bool still_sorted = x.is_sorted();
00608
00609
00610 std::deque<bool> check_for_sortedness;
00611 if (still_sorted)
00612 check_for_sortedness.insert(check_for_sortedness.end(), nrows, false);
00613
00614 for (dimension_type k = n_lines_or_equalities; k-- > 0; ) {
00615
00616
00617
00618 Linear_Row& x_k = x[k];
00619 dimension_type j = ncols - 1;
00620 while (j != 0 && x_k[j] == 0)
00621 --j;
00622
00623
00624 for (dimension_type i = k; i-- > 0; ) {
00625 Linear_Row& x_i = x[i];
00626 if (x_i[j] != 0) {
00627
00628
00629 x_i.linear_combine(x_k, j);
00630 if (still_sorted) {
00631
00632
00633 if (i > 0)
00634 check_for_sortedness[i-1] = true;
00635 check_for_sortedness[i] = true;
00636 }
00637 }
00638 }
00639
00640
00641
00642
00643
00644
00645 const bool have_to_negate = (x_k[j] < 0);
00646 if (have_to_negate)
00647 for (dimension_type h = ncols; h-- > 0; )
00648 PPL::neg_assign(x_k[h]);
00649
00650
00651
00652
00653 for (dimension_type i = n_lines_or_equalities; i < nrows; ++i) {
00654 Linear_Row& x_i = x[i];
00655 if (x_i[j] != 0) {
00656
00657
00658 x_i.linear_combine(x_k, j);
00659 if (still_sorted) {
00660
00661
00662 if (i > n_lines_or_equalities)
00663 check_for_sortedness[i-1] = true;
00664 check_for_sortedness[i] = true;
00665 }
00666 }
00667 }
00668 if (have_to_negate)
00669
00670 for (dimension_type h = ncols; h-- > 0; )
00671 PPL::neg_assign(x_k[h]);
00672 }
00673
00674
00675 for (dimension_type i = 0; still_sorted && i < nrows-1; ++i)
00676 if (check_for_sortedness[i])
00677
00678 still_sorted = (compare(x[i], x[i+1]) <= 0);
00679
00680 x.set_sorted(still_sorted);
00681
00682
00683 assert(x.OK(true));
00684 }
00685
00686 void
00687 PPL::Linear_System::simplify() {
00688 Linear_System& x = *this;
00689
00690
00691 assert(x.OK(true));
00692 assert(x.num_pending_rows() == 0);
00693
00694
00695 dimension_type nrows = x.num_rows();
00696 dimension_type n_lines_or_equalities = 0;
00697 for (dimension_type i = 0; i < nrows; ++i)
00698 if (x[i].is_line_or_equality()) {
00699 if (n_lines_or_equalities < i) {
00700 std::swap(x[i], x[n_lines_or_equalities]);
00701
00702 assert(!x.sorted);
00703 }
00704 ++n_lines_or_equalities;
00705 }
00706
00707 const dimension_type rank = x.gauss(n_lines_or_equalities);
00708
00709 if (rank < n_lines_or_equalities) {
00710 const dimension_type
00711 n_rays_or_points_or_inequalities = nrows - n_lines_or_equalities;
00712 const dimension_type
00713 num_swaps = std::min(n_lines_or_equalities - rank,
00714 n_rays_or_points_or_inequalities);
00715 for (dimension_type i = num_swaps; i-- > 0; )
00716 std::swap(x[--nrows], x[rank + i]);
00717 x.erase_to_end(nrows);
00718 x.unset_pending_rows();
00719 if (n_rays_or_points_or_inequalities > num_swaps)
00720 x.set_sorted(false);
00721 n_lines_or_equalities = rank;
00722 }
00723
00724 x.back_substitute(n_lines_or_equalities);
00725
00726 assert(x.OK(true));
00727 }
00728
00729 void
00730 PPL::Linear_System::add_rows_and_columns(const dimension_type n) {
00731 assert(n > 0);
00732 const bool was_sorted = is_sorted();
00733 const dimension_type old_n_rows = num_rows();
00734 const dimension_type old_n_columns = num_columns();
00735 add_zero_rows_and_columns(n, n, Linear_Row::Flags(row_topology));
00736 Linear_System& x = *this;
00737
00738 for (dimension_type i = old_n_rows; i-- > 0; )
00739 std::swap(x[i], x[i + n]);
00740 for (dimension_type i = n, c = old_n_columns; i-- > 0; ) {
00741
00742
00743
00744 Linear_Row& r = x[i];
00745 r[c++] = 1;
00746 r.set_is_line_or_equality();
00747
00748 }
00749
00750
00751 if (old_n_columns == 0) {
00752 x[n-1].set_is_ray_or_point_or_inequality();
00753
00754
00755 set_sorted(true);
00756 }
00757 else if (was_sorted)
00758 set_sorted(compare(x[n-1], x[n]) <= 0);
00759
00760
00761 assert(OK(true));
00762 }
00763
00764 void
00765 PPL::Linear_System::sort_pending_and_remove_duplicates() {
00766 assert(num_pending_rows() > 0);
00767 assert(is_sorted());
00768 Linear_System& x = *this;
00769
00770
00771
00772 const dimension_type first_pending = x.first_pending_row();
00773 x.sort_rows(first_pending, x.num_rows());
00774
00775
00776 dimension_type num_rows = x.num_rows();
00777
00778 dimension_type k1 = 0;
00779 dimension_type k2 = first_pending;
00780 dimension_type num_duplicates = 0;
00781
00782
00783 while (k1 < first_pending && k2 < num_rows) {
00784 const int cmp = compare(x[k1], x[k2]);
00785 if (cmp == 0) {
00786
00787 ++num_duplicates;
00788 --num_rows;
00789
00790 ++k1;
00791
00792 if (k2 < num_rows)
00793 std::swap(x[k2], x[k2 + num_duplicates]);
00794 }
00795 else if (cmp < 0)
00796
00797 ++k1;
00798 else {
00799
00800
00801
00802 ++k2;
00803 if (num_duplicates > 0 && k2 < num_rows)
00804 std::swap(x[k2], x[k2 + num_duplicates]);
00805 }
00806 }
00807
00808
00809 if (num_duplicates > 0) {
00810 if (k2 < num_rows)
00811 for (++k2; k2 < num_rows; ++k2)
00812 std::swap(x[k2], x[k2 + num_duplicates]);
00813 x.erase_to_end(num_rows);
00814 }
00815
00816
00817 assert(OK(false));
00818 }
00819
00820 bool
00821 PPL::Linear_System::check_sorted() const {
00822 const Linear_System& x = *this;
00823 for (dimension_type i = first_pending_row(); i-- > 1; )
00824 if (compare(x[i], x[i-1]) < 0)
00825 return false;
00826 return true;
00827 }
00828
00829 bool
00830 PPL::Linear_System::OK(const bool check_strong_normalized) const {
00831 #ifndef NDEBUG
00832 using std::endl;
00833 using std::cerr;
00834 #endif
00835
00836
00837 if (first_pending_row() > num_rows()) {
00838 #ifndef NDEBUG
00839 cerr << "Linear_System has a negative number of pending rows!"
00840 << endl;
00841 #endif
00842 return false;
00843 }
00844
00845
00846
00847 if (num_rows() == 0)
00848 if (is_necessarily_closed() || num_columns() != 1)
00849 return true;
00850 else {
00851 #ifndef NDEBUG
00852 cerr << "NNC Linear_System has one column" << endl;
00853 #endif
00854 return false;
00855 }
00856
00857
00858
00859
00860 const dimension_type min_cols = is_necessarily_closed() ? 1 : 2;
00861 if (num_columns() < min_cols) {
00862 #ifndef NDEBUG
00863 cerr << "Linear_System has fewer columns than the minimum "
00864 << "allowed by its topology:"
00865 << endl
00866 << "num_columns is " << num_columns()
00867 << ", minimum is " << min_cols
00868 << endl;
00869 #endif
00870 return false;
00871 }
00872
00873 const Linear_System& x = *this;
00874 const dimension_type n_rows = num_rows();
00875 for (dimension_type i = 0; i < n_rows; ++i) {
00876 if (!x[i].OK(row_size, row_capacity))
00877 return false;
00878
00879 if (x.topology() != x[i].topology()) {
00880 #ifndef NDEBUG
00881 cerr << "Topology mismatch between the system "
00882 << "and one of its rows!"
00883 << endl;
00884 #endif
00885 return false;
00886 }
00887 }
00888
00889 if (check_strong_normalized) {
00890
00891
00892
00893
00894
00895 Linear_System tmp(x, With_Pending());
00896 tmp.strong_normalize();
00897 if (x != tmp) {
00898 #ifndef NDEBUG
00899 cerr << "Linear_System rows are not strongly normalized!"
00900 << endl;
00901 #endif
00902 return false;
00903 }
00904 }
00905
00906 if (sorted && !check_sorted()) {
00907 #ifndef NDEBUG
00908 cerr << "The system declares itself to be sorted but it is not!"
00909 << endl;
00910 #endif
00911 return false;
00912 }
00913
00914
00915 return true;
00916 }