OR-Tools  8.2
lp_data.cc
Go to the documentation of this file.
1 // Copyright 2010-2018 Google LLC
2 // Licensed under the Apache License, Version 2.0 (the "License");
3 // you may not use this file except in compliance with the License.
4 // You may obtain a copy of the License at
5 //
6 // http://www.apache.org/licenses/LICENSE-2.0
7 //
8 // Unless required by applicable law or agreed to in writing, software
9 // distributed under the License is distributed on an "AS IS" BASIS,
10 // WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
11 // See the License for the specific language governing permissions and
12 // limitations under the License.
13 
15 
16 #include <algorithm>
17 #include <cmath>
18 #include <string>
19 #include <utility>
20 
21 #include "absl/container/flat_hash_map.h"
22 #include "absl/strings/str_cat.h"
23 #include "absl/strings/str_format.h"
24 #include "ortools/base/logging.h"
30 
31 namespace operations_research {
32 namespace glop {
33 
34 namespace {
35 
36 // This should be the same as DCHECK(AreBoundsValid()), but the DCHECK() are
37 // split to give more meaningful information to the user in case of failure.
38 void DebugCheckBoundsValid(Fractional lower_bound, Fractional upper_bound) {
39  DCHECK(!std::isnan(lower_bound));
40  DCHECK(!std::isnan(upper_bound));
41  DCHECK(!(lower_bound == kInfinity && upper_bound == kInfinity));
42  DCHECK(!(lower_bound == -kInfinity && upper_bound == -kInfinity));
43  DCHECK_LE(lower_bound, upper_bound);
44  DCHECK(AreBoundsValid(lower_bound, upper_bound));
45 }
46 
47 // Returns true if the bounds are the ones of a free or boxed row. Note that
48 // a fixed row is not counted as boxed.
49 bool AreBoundsFreeOrBoxed(Fractional lower_bound, Fractional upper_bound) {
50  if (lower_bound == -kInfinity && upper_bound == kInfinity) return true;
51  if (lower_bound != -kInfinity && upper_bound != kInfinity &&
52  lower_bound != upper_bound) {
53  return true;
54  }
55  return false;
56 }
57 
58 template <class I, class T>
59 double Average(const absl::StrongVector<I, T>& v) {
60  const size_t size = v.size();
61  DCHECK_LT(0, size);
62  double sum = 0.0;
63  double n = 0.0; // n is used in a calculation involving doubles.
64  for (I i(0); i < size; ++i) {
65  if (v[i] == 0.0) continue;
66  ++n;
67  sum += static_cast<double>(v[i].value());
68  }
69  return n == 0.0 ? 0.0 : sum / n;
70 }
71 
72 template <class I, class T>
73 double StandardDeviation(const absl::StrongVector<I, T>& v) {
74  const size_t size = v.size();
75  double n = 0.0; // n is used in a calculation involving doubles.
76  double sigma_square = 0.0;
77  double sigma = 0.0;
78  for (I i(0); i < size; ++i) {
79  double sample = static_cast<double>(v[i].value());
80  if (sample == 0.0) continue;
81  sigma_square += sample * sample;
82  sigma += sample;
83  ++n;
84  }
85  return n == 0.0 ? 0.0 : sqrt((sigma_square - sigma * sigma / n) / n);
86 }
87 
88 // Returns 0 when the vector is empty.
89 template <class I, class T>
90 T GetMaxElement(const absl::StrongVector<I, T>& v) {
91  const size_t size = v.size();
92  if (size == 0) {
93  return T(0);
94  }
95 
96  T max_index = v[I(0)];
97  for (I i(1); i < size; ++i) {
98  if (max_index < v[i]) {
99  max_index = v[i];
100  }
101  }
102  return max_index;
103 }
104 
105 } // anonymous namespace
106 
107 // --------------------------------------------------------
108 // LinearProgram
109 // --------------------------------------------------------
111  : matrix_(),
112  transpose_matrix_(),
113  constraint_lower_bounds_(),
114  constraint_upper_bounds_(),
115  constraint_names_(),
116  objective_coefficients_(),
117  variable_lower_bounds_(),
118  variable_upper_bounds_(),
119  variable_names_(),
120  variable_types_(),
121  integer_variables_list_(),
122  variable_table_(),
123  constraint_table_(),
124  objective_offset_(0.0),
125  objective_scaling_factor_(1.0),
126  maximize_(false),
127  columns_are_known_to_be_clean_(true),
128  transpose_matrix_is_consistent_(true),
129  integer_variables_list_is_consistent_(true),
130  name_(),
131  first_slack_variable_(kInvalidCol) {}
132 
134  matrix_.Clear();
135  transpose_matrix_.Clear();
136 
137  constraint_lower_bounds_.clear();
138  constraint_upper_bounds_.clear();
139  constraint_names_.clear();
140 
141  objective_coefficients_.clear();
142  variable_lower_bounds_.clear();
143  variable_upper_bounds_.clear();
144  variable_types_.clear();
145  integer_variables_list_.clear();
146  variable_names_.clear();
147 
148  constraint_table_.clear();
149  variable_table_.clear();
150 
151  maximize_ = false;
152  objective_offset_ = 0.0;
153  objective_scaling_factor_ = 1.0;
154  columns_are_known_to_be_clean_ = true;
155  transpose_matrix_is_consistent_ = true;
156  integer_variables_list_is_consistent_ = true;
157  name_.clear();
158  first_slack_variable_ = kInvalidCol;
159 }
160 
162  DCHECK_EQ(kInvalidCol, first_slack_variable_)
163  << "New variables can't be added to programs that already have slack "
164  "variables. Consider calling LinearProgram::DeleteSlackVariables() "
165  "before adding new variables to the problem.";
166  objective_coefficients_.push_back(0.0);
167  variable_lower_bounds_.push_back(0);
168  variable_upper_bounds_.push_back(kInfinity);
169  variable_types_.push_back(VariableType::CONTINUOUS);
170  variable_names_.push_back("");
171  transpose_matrix_is_consistent_ = false;
172  return matrix_.AppendEmptyColumn();
173 }
174 
175 ColIndex LinearProgram::CreateNewSlackVariable(bool is_integer_slack_variable,
176  Fractional lower_bound,
177  Fractional upper_bound,
178  const std::string& name) {
179  objective_coefficients_.push_back(0.0);
180  variable_lower_bounds_.push_back(lower_bound);
181  variable_upper_bounds_.push_back(upper_bound);
182  variable_types_.push_back(is_integer_slack_variable
185  variable_names_.push_back(name);
186  transpose_matrix_is_consistent_ = false;
187  return matrix_.AppendEmptyColumn();
188 }
189 
191  DCHECK_EQ(kInvalidCol, first_slack_variable_)
192  << "New constraints can't be added to programs that already have slack "
193  "variables. Consider calling LinearProgram::DeleteSlackVariables() "
194  "before adding new variables to the problem.";
195  const RowIndex row(constraint_names_.size());
196  matrix_.SetNumRows(row + 1);
197  constraint_lower_bounds_.push_back(Fractional(0.0));
198  constraint_upper_bounds_.push_back(Fractional(0.0));
199  constraint_names_.push_back("");
200  transpose_matrix_is_consistent_ = false;
201  return row;
202 }
203 
204 ColIndex LinearProgram::FindOrCreateVariable(const std::string& variable_id) {
205  const absl::flat_hash_map<std::string, ColIndex>::iterator it =
206  variable_table_.find(variable_id);
207  if (it != variable_table_.end()) {
208  return it->second;
209  } else {
210  const ColIndex col = CreateNewVariable();
211  variable_names_[col] = variable_id;
212  variable_table_[variable_id] = col;
213  return col;
214  }
215 }
216 
218  const std::string& constraint_id) {
219  const absl::flat_hash_map<std::string, RowIndex>::iterator it =
220  constraint_table_.find(constraint_id);
221  if (it != constraint_table_.end()) {
222  return it->second;
223  } else {
224  const RowIndex row = CreateNewConstraint();
225  constraint_names_[row] = constraint_id;
226  constraint_table_[constraint_id] = row;
227  return row;
228  }
229 }
230 
231 void LinearProgram::SetVariableName(ColIndex col, absl::string_view name) {
232  variable_names_[col] = std::string(name);
233 }
234 
236  const bool var_was_integer = IsVariableInteger(col);
237  variable_types_[col] = type;
238  const bool var_is_integer = IsVariableInteger(col);
239  if (var_is_integer != var_was_integer) {
240  integer_variables_list_is_consistent_ = false;
241  }
242 }
243 
244 void LinearProgram::SetConstraintName(RowIndex row, absl::string_view name) {
245  constraint_names_[row] = std::string(name);
246 }
247 
249  Fractional upper_bound) {
250  if (dcheck_bounds_) DebugCheckBoundsValid(lower_bound, upper_bound);
251  const bool var_was_binary = IsVariableBinary(col);
252  variable_lower_bounds_[col] = lower_bound;
253  variable_upper_bounds_[col] = upper_bound;
254  const bool var_is_binary = IsVariableBinary(col);
255  if (var_is_binary != var_was_binary) {
256  integer_variables_list_is_consistent_ = false;
257  }
258 }
259 
260 void LinearProgram::UpdateAllIntegerVariableLists() const {
261  if (integer_variables_list_is_consistent_) return;
262  integer_variables_list_.clear();
263  binary_variables_list_.clear();
264  non_binary_variables_list_.clear();
265  const ColIndex num_cols = num_variables();
266  for (ColIndex col(0); col < num_cols; ++col) {
267  if (IsVariableInteger(col)) {
268  integer_variables_list_.push_back(col);
269  if (IsVariableBinary(col)) {
270  binary_variables_list_.push_back(col);
271  } else {
272  non_binary_variables_list_.push_back(col);
273  }
274  }
275  }
276  integer_variables_list_is_consistent_ = true;
277 }
278 
279 const std::vector<ColIndex>& LinearProgram::IntegerVariablesList() const {
280  UpdateAllIntegerVariableLists();
281  return integer_variables_list_;
282 }
283 
284 const std::vector<ColIndex>& LinearProgram::BinaryVariablesList() const {
285  UpdateAllIntegerVariableLists();
286  return binary_variables_list_;
287 }
288 
289 const std::vector<ColIndex>& LinearProgram::NonBinaryVariablesList() const {
290  UpdateAllIntegerVariableLists();
291  return non_binary_variables_list_;
292 }
293 
294 bool LinearProgram::IsVariableInteger(ColIndex col) const {
295  return variable_types_[col] == VariableType::INTEGER ||
296  variable_types_[col] == VariableType::IMPLIED_INTEGER;
297 }
298 
299 bool LinearProgram::IsVariableBinary(ColIndex col) const {
300  // TODO(user,user): bounds of binary variables (and of integer ones) should
301  // be integer. Add a preprocessor for that.
302  return IsVariableInteger(col) && (variable_lower_bounds_[col] < kEpsilon) &&
303  (variable_lower_bounds_[col] > Fractional(-1)) &&
304  (variable_upper_bounds_[col] > Fractional(1) - kEpsilon) &&
305  (variable_upper_bounds_[col] < 2);
306 }
307 
309  Fractional upper_bound) {
310  if (dcheck_bounds_) DebugCheckBoundsValid(lower_bound, upper_bound);
311  ResizeRowsIfNeeded(row);
312  constraint_lower_bounds_[row] = lower_bound;
313  constraint_upper_bounds_[row] = upper_bound;
314 }
315 
316 void LinearProgram::SetCoefficient(RowIndex row, ColIndex col,
317  Fractional value) {
319  ResizeRowsIfNeeded(row);
320  columns_are_known_to_be_clean_ = false;
321  transpose_matrix_is_consistent_ = false;
323 }
324 
327  objective_coefficients_[col] = value;
328 }
329 
332  objective_offset_ = objective_offset;
333 }
334 
336  Fractional objective_scaling_factor) {
339  objective_scaling_factor_ = objective_scaling_factor;
340 }
341 
343  maximize_ = maximize;
344 }
345 
347  if (columns_are_known_to_be_clean_) return;
348  matrix_.CleanUp();
349  columns_are_known_to_be_clean_ = true;
350  transpose_matrix_is_consistent_ = false;
351 }
352 
354  if (columns_are_known_to_be_clean_) return true;
355  columns_are_known_to_be_clean_ = matrix_.IsCleanedUp();
356  return columns_are_known_to_be_clean_;
357 }
358 
359 std::string LinearProgram::GetVariableName(ColIndex col) const {
360  return col >= variable_names_.size() || variable_names_[col].empty()
361  ? absl::StrFormat("c%d", col.value())
362  : variable_names_[col];
363 }
364 
365 std::string LinearProgram::GetConstraintName(RowIndex row) const {
366  return row >= constraint_names_.size() || constraint_names_[row].empty()
367  ? absl::StrFormat("r%d", row.value())
368  : constraint_names_[row];
369 }
370 
372  return variable_types_[col];
373 }
374 
376  if (!transpose_matrix_is_consistent_) {
377  transpose_matrix_.PopulateFromTranspose(matrix_);
378  transpose_matrix_is_consistent_ = true;
379  }
380  DCHECK_EQ(transpose_matrix_.num_rows().value(), matrix_.num_cols().value());
381  DCHECK_EQ(transpose_matrix_.num_cols().value(), matrix_.num_rows().value());
382  return transpose_matrix_;
383 }
384 
386  if (!transpose_matrix_is_consistent_) {
387  transpose_matrix_.PopulateFromTranspose(matrix_);
388  }
389  // This enables a client to start modifying the matrix and then abort and not
390  // call UseTransposeMatrixAsReference(). Then, the other client of
391  // GetTransposeSparseMatrix() will still see the correct matrix.
392  transpose_matrix_is_consistent_ = false;
393  return &transpose_matrix_;
394 }
395 
397  DCHECK_EQ(transpose_matrix_.num_rows().value(), matrix_.num_cols().value());
398  DCHECK_EQ(transpose_matrix_.num_cols().value(), matrix_.num_rows().value());
399  matrix_.PopulateFromTranspose(transpose_matrix_);
400  transpose_matrix_is_consistent_ = true;
401 }
402 
404  transpose_matrix_.Clear();
405  transpose_matrix_is_consistent_ = false;
406 }
407 
409  return matrix_.column(col);
410 }
411 
413  columns_are_known_to_be_clean_ = false;
414  transpose_matrix_is_consistent_ = false;
415  return matrix_.mutable_column(col);
416 }
417 
419  ColIndex col) const {
420  return maximize_ ? -objective_coefficients()[col]
422 }
423 
425  Fractional min_magnitude = 0.0;
426  Fractional max_magnitude = 0.0;
427  matrix_.ComputeMinAndMaxMagnitudes(&min_magnitude, &max_magnitude);
428  return absl::StrFormat(
429  "%d rows, %d columns, %d entries with magnitude in [%e, %e]",
431  // static_cast<int64> is needed because the Android port uses int32.
432  static_cast<int64>(num_entries().value()), min_magnitude, max_magnitude);
433 }
434 
435 namespace {
436 
437 template <typename FractionalValues>
438 void UpdateStats(const FractionalValues& values, int64* num_non_zeros,
439  Fractional* min_value, Fractional* max_value) {
440  for (const Fractional v : values) {
441  if (v == 0 || v == kInfinity || v == -kInfinity) continue;
442  *min_value = std::min(*min_value, v);
443  *max_value = std::max(*max_value, v);
444  ++(*num_non_zeros);
445  }
446 }
447 
448 } // namespace
449 
451  int64 num_non_zeros = 0;
452  Fractional min_value = +kInfinity;
453  Fractional max_value = -kInfinity;
454  UpdateStats(objective_coefficients_, &num_non_zeros, &min_value, &max_value);
455  if (num_non_zeros == 0) {
456  return "No objective term. This is a pure feasibility problem.";
457  } else {
458  return absl::StrFormat("%d non-zeros, range [%e, %e]", num_non_zeros,
459  min_value, max_value);
460  }
461 }
462 
464  int64 num_non_zeros = 0;
465  Fractional min_value = +kInfinity;
466  Fractional max_value = -kInfinity;
467  UpdateStats(variable_lower_bounds_, &num_non_zeros, &min_value, &max_value);
468  UpdateStats(variable_upper_bounds_, &num_non_zeros, &min_value, &max_value);
469  UpdateStats(constraint_lower_bounds_, &num_non_zeros, &min_value, &max_value);
470  UpdateStats(constraint_upper_bounds_, &num_non_zeros, &min_value, &max_value);
471  if (num_non_zeros == 0) {
472  return "All variables/constraints bounds are zero or +/- infinity.";
473  } else {
474  return absl::StrFormat("%d non-zeros, range [%e, %e]", num_non_zeros,
475  min_value, max_value);
476  }
477 }
478 
480  const DenseRow& solution, Fractional absolute_tolerance) const {
481  DCHECK_EQ(solution.size(), num_variables());
482  if (solution.size() != num_variables()) return false;
483  const ColIndex num_cols = num_variables();
484  for (ColIndex col = ColIndex(0); col < num_cols; ++col) {
485  if (!IsFinite(solution[col])) return false;
486  const Fractional lb_error = variable_lower_bounds()[col] - solution[col];
487  const Fractional ub_error = solution[col] - variable_upper_bounds()[col];
488  if (lb_error > absolute_tolerance || ub_error > absolute_tolerance) {
489  return false;
490  }
491  }
492  return true;
493 }
494 
496  Fractional absolute_tolerance) const {
497  if (!SolutionIsWithinVariableBounds(solution, absolute_tolerance)) {
498  return false;
499  }
500  const SparseMatrix& transpose = GetTransposeSparseMatrix();
501  const RowIndex num_rows = num_constraints();
502  for (RowIndex row = RowIndex(0); row < num_rows; ++row) {
503  const Fractional sum =
504  ScalarProduct(solution, transpose.column(RowToColIndex(row)));
505  if (!IsFinite(sum)) return false;
506  const Fractional lb_error = constraint_lower_bounds()[row] - sum;
507  const Fractional ub_error = sum - constraint_upper_bounds()[row];
508  if (lb_error > absolute_tolerance || ub_error > absolute_tolerance) {
509  return false;
510  }
511  }
512  return true;
513 }
514 
516  Fractional absolute_tolerance) const {
517  DCHECK_EQ(solution.size(), num_variables());
518  if (solution.size() != num_variables()) return false;
519  for (ColIndex col : IntegerVariablesList()) {
520  if (!IsFinite(solution[col])) return false;
521  const Fractional fractionality = fabs(solution[col] - round(solution[col]));
522  if (fractionality > absolute_tolerance) return false;
523  }
524  return true;
525 }
526 
528  Fractional absolute_tolerance) const {
529  return SolutionIsLPFeasible(solution, absolute_tolerance) &&
530  SolutionIsInteger(solution, absolute_tolerance);
531 }
532 
534  CHECK(solution != nullptr);
535  const ColIndex num_cols = GetFirstSlackVariable();
536  const SparseMatrix& transpose = GetTransposeSparseMatrix();
537  const RowIndex num_rows = num_constraints();
538  CHECK_EQ(solution->size(), num_variables());
539  for (RowIndex row = RowIndex(0); row < num_rows; ++row) {
540  const Fractional sum = PartialScalarProduct(
541  *solution, transpose.column(RowToColIndex(row)), num_cols.value());
542  const ColIndex slack_variable = GetSlackVariable(row);
543  CHECK_NE(slack_variable, kInvalidCol);
544  (*solution)[slack_variable] = -sum;
545  }
546 }
547 
549  Fractional value) const {
551 }
552 
554  Fractional value) const {
556 }
557 
558 std::string LinearProgram::Dump() const {
559  // Objective line.
560  std::string output = maximize_ ? "max:" : "min:";
561  if (objective_offset_ != 0.0) {
562  output += Stringify(objective_offset_);
563  }
564  const ColIndex num_cols = num_variables();
565  for (ColIndex col(0); col < num_cols; ++col) {
566  const Fractional coeff = objective_coefficients()[col];
567  if (coeff != 0.0) {
568  output += StringifyMonomial(coeff, GetVariableName(col), false);
569  }
570  }
571  output.append(";\n");
572 
573  // Constraints.
574  const RowIndex num_rows = num_constraints();
575  for (RowIndex row(0); row < num_rows; ++row) {
576  const Fractional lower_bound = constraint_lower_bounds()[row];
577  const Fractional upper_bound = constraint_upper_bounds()[row];
578  output += GetConstraintName(row);
579  output += ":";
580  if (AreBoundsFreeOrBoxed(lower_bound, upper_bound)) {
581  output += " ";
582  output += Stringify(lower_bound);
583  output += " <=";
584  }
585  for (ColIndex col(0); col < num_cols; ++col) {
586  const Fractional coeff = matrix_.LookUpValue(row, col);
587  output += StringifyMonomial(coeff, GetVariableName(col), false);
588  }
589  if (AreBoundsFreeOrBoxed(lower_bound, upper_bound)) {
590  output += " <= ";
591  output += Stringify(upper_bound);
592  } else if (lower_bound == upper_bound) {
593  output += " = ";
594  output += Stringify(upper_bound);
595  } else if (lower_bound != -kInfinity) {
596  output += " >= ";
597  output += Stringify(lower_bound);
598  } else if (lower_bound != kInfinity) {
599  output += " <= ";
600  output += Stringify(upper_bound);
601  }
602  output += ";\n";
603  }
604 
605  // Variables.
606  for (ColIndex col(0); col < num_cols; ++col) {
607  const Fractional lower_bound = variable_lower_bounds()[col];
608  const Fractional upper_bound = variable_upper_bounds()[col];
609  if (AreBoundsFreeOrBoxed(lower_bound, upper_bound)) {
610  output += Stringify(lower_bound);
611  output += " <= ";
612  }
613  output += GetVariableName(col);
614  if (AreBoundsFreeOrBoxed(lower_bound, upper_bound)) {
615  output += " <= ";
616  output += Stringify(upper_bound);
617  } else if (lower_bound == upper_bound) {
618  output += " = ";
619  output += Stringify(upper_bound);
620  } else if (lower_bound != -kInfinity) {
621  output += " >= ";
622  output += Stringify(lower_bound);
623  } else if (lower_bound != kInfinity) {
624  output += " <= ";
625  output += Stringify(upper_bound);
626  }
627  output += ";\n";
628  }
629 
630  // Integer variables.
631  // TODO(user): if needed provide similar output for binary variables.
632  const std::vector<ColIndex>& integer_variables = IntegerVariablesList();
633  if (!integer_variables.empty()) {
634  output += "int";
635  for (ColIndex col : integer_variables) {
636  output += " ";
637  output += GetVariableName(col);
638  }
639  output += ";\n";
640  }
641 
642  return output;
643 }
644 
645 std::string LinearProgram::DumpSolution(const DenseRow& variable_values) const {
646  DCHECK_EQ(variable_values.size(), num_variables());
647  std::string output;
648  for (ColIndex col(0); col < variable_values.size(); ++col) {
649  if (!output.empty()) absl::StrAppend(&output, ", ");
650  absl::StrAppend(&output, GetVariableName(col), " = ",
651  (variable_values[col]));
652  }
653  return output;
654 }
655 
656 std::string LinearProgram::GetProblemStats() const {
657  return ProblemStatFormatter(
658  "%d,%d,%d,%d,%d,%d,%d,%d,%d,%d,%d,%d,%d,%d,"
659  "%d,%d,%d,%d");
660 }
661 
663  return ProblemStatFormatter(
664  "Number of rows : %d\n"
665  "Number of variables in file : %d\n"
666  "Number of entries (non-zeros) : %d\n"
667  "Number of entries in the objective : %d\n"
668  "Number of entries in the right-hand side : %d\n"
669  "Number of <= constraints : %d\n"
670  "Number of >= constraints : %d\n"
671  "Number of = constraints : %d\n"
672  "Number of range constraints : %d\n"
673  "Number of non-negative variables : %d\n"
674  "Number of boxed variables : %d\n"
675  "Number of free variables : %d\n"
676  "Number of fixed variables : %d\n"
677  "Number of other variables : %d\n"
678  "Number of integer variables : %d\n"
679  "Number of binary variables : %d\n"
680  "Number of non-binary integer variables : %d\n"
681  "Number of continuous variables : %d\n");
682 }
683 
684 std::string LinearProgram::GetNonZeroStats() const {
685  return NonZeroStatFormatter("%.2f%%,%d,%.2f,%.2f,%d,%.2f,%.2f");
686 }
687 
689  return NonZeroStatFormatter(
690  "Fill rate : %.2f%%\n"
691  "Entries in row (Max / average / std. dev.) : %d / %.2f / %.2f\n"
692  "Entries in column (Max / average / std. dev.): %d / %.2f / %.2f\n");
693 }
694 
696  bool detect_integer_constraints) {
697  // Clean up the matrix. We're going to add entries, but we'll only be adding
698  // them to new columns, and only one entry per column, which does not
699  // invalidate the "cleanness" of the matrix.
700  CleanUp();
701 
702  // Detect which constraints produce an integer slack variable. A constraint
703  // has an integer slack variable, if it contains only integer variables with
704  // integer coefficients. We do not check the bounds of the constraints,
705  // because in such case, they will be tightened to integer values by the
706  // preprocessors.
707  //
708  // We don't use the transpose, because it might not be valid and it would be
709  // inefficient to update it and invalidate it again at the end of this
710  // preprocessor.
711  DenseBooleanColumn has_integer_slack_variable(num_constraints(),
712  detect_integer_constraints);
713  if (detect_integer_constraints) {
714  for (ColIndex col(0); col < num_variables(); ++col) {
715  const SparseColumn& column = matrix_.column(col);
716  const bool is_integer_variable = IsVariableInteger(col);
717  for (const SparseColumn::Entry& entry : column) {
718  const RowIndex row = entry.row();
719  has_integer_slack_variable[row] =
720  has_integer_slack_variable[row] && is_integer_variable &&
721  round(entry.coefficient()) == entry.coefficient();
722  }
723  }
724  }
725 
726  // Extend the matrix of the problem with an identity matrix.
727  const ColIndex original_num_variables = num_variables();
728  for (RowIndex row(0); row < num_constraints(); ++row) {
729  ColIndex slack_variable_index = GetSlackVariable(row);
730  if (slack_variable_index != kInvalidCol &&
731  slack_variable_index < original_num_variables) {
732  // Slack variable is already present in this constraint.
733  continue;
734  }
735  const ColIndex slack_col = CreateNewSlackVariable(
736  has_integer_slack_variable[row], -constraint_upper_bounds_[row],
737  -constraint_lower_bounds_[row], absl::StrCat("s", row.value()));
738  SetCoefficient(row, slack_col, 1.0);
739  SetConstraintBounds(row, 0.0, 0.0);
740  }
741 
742  columns_are_known_to_be_clean_ = true;
743  transpose_matrix_is_consistent_ = false;
744  if (first_slack_variable_ == kInvalidCol) {
745  first_slack_variable_ = original_num_variables;
746  }
747 }
748 
750  return first_slack_variable_;
751 }
752 
753 ColIndex LinearProgram::GetSlackVariable(RowIndex row) const {
754  DCHECK_GE(row, RowIndex(0));
756  if (first_slack_variable_ == kInvalidCol) {
757  return kInvalidCol;
758  }
759  return first_slack_variable_ + RowToColIndex(row);
760 }
761 
763  RowToColMapping* duplicated_rows) {
764  const ColIndex dual_num_variables = dual.num_variables();
765  const RowIndex dual_num_constraints = dual.num_constraints();
766  Clear();
767 
768  // We always take the dual in its minimization form thanks to the
769  // GetObjectiveCoefficientForMinimizationVersion() below, so this will always
770  // be a maximization problem.
772 
773  // Taking the dual does not change the offset nor the objective scaling
774  // factor.
777 
778  // Create the dual variables y, with bounds depending on the type
779  // of constraints in the primal.
780  for (RowIndex dual_row(0); dual_row < dual_num_constraints; ++dual_row) {
782  const ColIndex col = RowToColIndex(dual_row);
783  const Fractional lower_bound = dual.constraint_lower_bounds()[dual_row];
784  const Fractional upper_bound = dual.constraint_upper_bounds()[dual_row];
785  if (lower_bound == upper_bound) {
787  SetObjectiveCoefficient(col, lower_bound);
788  } else if (upper_bound != kInfinity) {
789  // Note that for a ranged constraint, the loop will be continued here.
790  // This is wanted because we want to deal with the lower_bound afterwards.
792  SetObjectiveCoefficient(col, upper_bound);
793  } else if (lower_bound != -kInfinity) {
795  SetObjectiveCoefficient(col, lower_bound);
796  } else {
797  // This code does not support free rows in linear_program.
798  LOG(DFATAL) << "PopulateFromDual() was called with a program "
799  << "containing free constraints.";
800  }
801  }
802  // Create the dual slack variables v.
803  for (ColIndex dual_col(0); dual_col < dual_num_variables; ++dual_col) {
804  const Fractional lower_bound = dual.variable_lower_bounds()[dual_col];
805  if (lower_bound != -kInfinity) {
806  const ColIndex col = CreateNewVariable();
807  SetObjectiveCoefficient(col, lower_bound);
809  const RowIndex row = ColToRowIndex(dual_col);
811  }
812  }
813  // Create the dual surplus variables w.
814  for (ColIndex dual_col(0); dual_col < dual_num_variables; ++dual_col) {
815  const Fractional upper_bound = dual.variable_upper_bounds()[dual_col];
816  if (upper_bound != kInfinity) {
817  const ColIndex col = CreateNewVariable();
818  SetObjectiveCoefficient(col, upper_bound);
820  const RowIndex row = ColToRowIndex(dual_col);
822  }
823  }
824  // Store the transpose of the matrix.
825  for (ColIndex dual_col(0); dual_col < dual_num_variables; ++dual_col) {
826  const RowIndex row = ColToRowIndex(dual_col);
827  const Fractional row_bound =
829  SetConstraintBounds(row, row_bound, row_bound);
830  for (const SparseColumn::Entry e : dual.GetSparseColumn(dual_col)) {
831  const RowIndex dual_row = e.row();
832  const ColIndex col = RowToColIndex(dual_row);
833  SetCoefficient(row, col, e.coefficient());
834  }
835  }
836 
837  // Take care of ranged constraints.
838  duplicated_rows->assign(dual_num_constraints, kInvalidCol);
839  for (RowIndex dual_row(0); dual_row < dual_num_constraints; ++dual_row) {
840  const Fractional lower_bound = dual.constraint_lower_bounds()[dual_row];
841  const Fractional upper_bound = dual.constraint_upper_bounds()[dual_row];
842  if (AreBoundsFreeOrBoxed(lower_bound, upper_bound)) {
843  DCHECK(upper_bound != kInfinity || lower_bound != -kInfinity);
844 
845  // upper_bound was done in a loop above, now do the lower_bound.
846  const ColIndex col = CreateNewVariable();
848  SetObjectiveCoefficient(col, lower_bound);
850  matrix_.column(RowToColIndex(dual_row)));
851  (*duplicated_rows)[dual_row] = col;
852  }
853  }
854 
855  // We know that the columns are ordered by rows.
856  columns_are_known_to_be_clean_ = true;
857  transpose_matrix_is_consistent_ = false;
858 }
859 
861  const LinearProgram& linear_program) {
862  matrix_.PopulateFromSparseMatrix(linear_program.matrix_);
863  if (linear_program.transpose_matrix_is_consistent_) {
864  transpose_matrix_is_consistent_ = true;
865  transpose_matrix_.PopulateFromSparseMatrix(
866  linear_program.transpose_matrix_);
867  } else {
868  transpose_matrix_is_consistent_ = false;
869  transpose_matrix_.Clear();
870  }
871 
872  constraint_lower_bounds_ = linear_program.constraint_lower_bounds_;
873  constraint_upper_bounds_ = linear_program.constraint_upper_bounds_;
874  constraint_names_ = linear_program.constraint_names_;
875  constraint_table_.clear();
876 
877  PopulateNameObjectiveAndVariablesFromLinearProgram(linear_program);
878  first_slack_variable_ = linear_program.first_slack_variable_;
879 }
880 
882  const LinearProgram& lp, const RowPermutation& row_permutation,
883  const ColumnPermutation& col_permutation) {
884  DCHECK(lp.IsCleanedUp());
885  DCHECK_EQ(row_permutation.size(), lp.num_constraints());
886  DCHECK_EQ(col_permutation.size(), lp.num_variables());
888  Clear();
889 
890  // Populate matrix coefficients.
891  ColumnPermutation inverse_col_permutation;
892  inverse_col_permutation.PopulateFromInverse(col_permutation);
893  matrix_.PopulateFromPermutedMatrix(lp.matrix_, row_permutation,
894  inverse_col_permutation);
896 
897  // Populate constraints.
898  ApplyPermutation(row_permutation, lp.constraint_lower_bounds(),
899  &constraint_lower_bounds_);
900  ApplyPermutation(row_permutation, lp.constraint_upper_bounds(),
901  &constraint_upper_bounds_);
902 
903  // Populate variables.
904  ApplyPermutation(col_permutation, lp.objective_coefficients(),
905  &objective_coefficients_);
906  ApplyPermutation(col_permutation, lp.variable_lower_bounds(),
907  &variable_lower_bounds_);
908  ApplyPermutation(col_permutation, lp.variable_upper_bounds(),
909  &variable_upper_bounds_);
910  ApplyPermutation(col_permutation, lp.variable_types(), &variable_types_);
911  integer_variables_list_is_consistent_ = false;
912 
913  // There is no vector based accessor to names, because they may be created
914  // on the fly.
915  constraint_names_.resize(lp.num_constraints());
916  for (RowIndex old_row(0); old_row < lp.num_constraints(); ++old_row) {
917  const RowIndex new_row = row_permutation[old_row];
918  constraint_names_[new_row] = lp.constraint_names_[old_row];
919  }
920  variable_names_.resize(lp.num_variables());
921  for (ColIndex old_col(0); old_col < lp.num_variables(); ++old_col) {
922  const ColIndex new_col = col_permutation[old_col];
923  variable_names_[new_col] = lp.variable_names_[old_col];
924  }
925 
926  // Populate singular fields.
927  maximize_ = lp.maximize_;
928  objective_offset_ = lp.objective_offset_;
929  objective_scaling_factor_ = lp.objective_scaling_factor_;
930  name_ = lp.name_;
931 }
932 
934  const LinearProgram& linear_program) {
935  matrix_.PopulateFromZero(RowIndex(0), linear_program.num_variables());
936  first_slack_variable_ = kInvalidCol;
937  transpose_matrix_is_consistent_ = false;
938  transpose_matrix_.Clear();
939 
940  constraint_lower_bounds_.clear();
941  constraint_upper_bounds_.clear();
942  constraint_names_.clear();
943  constraint_table_.clear();
944 
945  PopulateNameObjectiveAndVariablesFromLinearProgram(linear_program);
946 }
947 
948 void LinearProgram::PopulateNameObjectiveAndVariablesFromLinearProgram(
949  const LinearProgram& linear_program) {
950  objective_coefficients_ = linear_program.objective_coefficients_;
951  variable_lower_bounds_ = linear_program.variable_lower_bounds_;
952  variable_upper_bounds_ = linear_program.variable_upper_bounds_;
953  variable_names_ = linear_program.variable_names_;
954  variable_types_ = linear_program.variable_types_;
955  integer_variables_list_is_consistent_ =
956  linear_program.integer_variables_list_is_consistent_;
957  integer_variables_list_ = linear_program.integer_variables_list_;
958  binary_variables_list_ = linear_program.binary_variables_list_;
959  non_binary_variables_list_ = linear_program.non_binary_variables_list_;
960  variable_table_.clear();
961 
962  maximize_ = linear_program.maximize_;
963  objective_offset_ = linear_program.objective_offset_;
964  objective_scaling_factor_ = linear_program.objective_scaling_factor_;
965  columns_are_known_to_be_clean_ =
966  linear_program.columns_are_known_to_be_clean_;
967  name_ = linear_program.name_;
968 }
969 
971  const SparseMatrix& coefficients, const DenseColumn& left_hand_sides,
972  const DenseColumn& right_hand_sides,
974  const RowIndex num_new_constraints = coefficients.num_rows();
975  DCHECK_EQ(num_variables(), coefficients.num_cols());
976  DCHECK_EQ(num_new_constraints, left_hand_sides.size());
977  DCHECK_EQ(num_new_constraints, right_hand_sides.size());
978  DCHECK_EQ(num_new_constraints, names.size());
979 
981  transpose_matrix_is_consistent_ = false;
982  transpose_matrix_.Clear();
983  columns_are_known_to_be_clean_ = false;
984 
985  // Copy constraint bounds and names from linear_program.
986  constraint_lower_bounds_.insert(constraint_lower_bounds_.end(),
987  left_hand_sides.begin(),
988  left_hand_sides.end());
989  constraint_upper_bounds_.insert(constraint_upper_bounds_.end(),
990  right_hand_sides.begin(),
991  right_hand_sides.end());
992  constraint_names_.insert(constraint_names_.end(), names.begin(), names.end());
993 }
994 
996  const SparseMatrix& coefficients, const DenseColumn& left_hand_sides,
997  const DenseColumn& right_hand_sides,
999  bool detect_integer_constraints_for_slack) {
1000  AddConstraints(coefficients, left_hand_sides, right_hand_sides, names);
1001  AddSlackVariablesWhereNecessary(detect_integer_constraints_for_slack);
1002 }
1003 
1005  const DenseRow& variable_lower_bounds,
1006  const DenseRow& variable_upper_bounds) {
1007  const ColIndex num_vars = num_variables();
1008  DCHECK_EQ(variable_lower_bounds.size(), num_vars);
1009  DCHECK_EQ(variable_upper_bounds.size(), num_vars);
1010 
1011  DenseRow new_lower_bounds(num_vars, 0);
1012  DenseRow new_upper_bounds(num_vars, 0);
1013  for (ColIndex i(0); i < num_vars; ++i) {
1014  const Fractional new_lower_bound =
1015  std::max(variable_lower_bounds[i], variable_lower_bounds_[i]);
1016  const Fractional new_upper_bound =
1017  std::min(variable_upper_bounds[i], variable_upper_bounds_[i]);
1018  if (new_lower_bound > new_upper_bound) {
1019  return false;
1020  }
1021  new_lower_bounds[i] = new_lower_bound;
1022  new_upper_bounds[i] = new_upper_bound;
1023  }
1024  variable_lower_bounds_.swap(new_lower_bounds);
1025  variable_upper_bounds_.swap(new_upper_bounds);
1026  return true;
1027 }
1028 
1029 void LinearProgram::Swap(LinearProgram* linear_program) {
1030  matrix_.Swap(&linear_program->matrix_);
1031  transpose_matrix_.Swap(&linear_program->transpose_matrix_);
1032 
1033  constraint_lower_bounds_.swap(linear_program->constraint_lower_bounds_);
1034  constraint_upper_bounds_.swap(linear_program->constraint_upper_bounds_);
1035  constraint_names_.swap(linear_program->constraint_names_);
1036 
1037  objective_coefficients_.swap(linear_program->objective_coefficients_);
1038  variable_lower_bounds_.swap(linear_program->variable_lower_bounds_);
1039  variable_upper_bounds_.swap(linear_program->variable_upper_bounds_);
1040  variable_names_.swap(linear_program->variable_names_);
1041  variable_types_.swap(linear_program->variable_types_);
1042  integer_variables_list_.swap(linear_program->integer_variables_list_);
1043  binary_variables_list_.swap(linear_program->binary_variables_list_);
1044  non_binary_variables_list_.swap(linear_program->non_binary_variables_list_);
1045 
1046  variable_table_.swap(linear_program->variable_table_);
1047  constraint_table_.swap(linear_program->constraint_table_);
1048 
1049  std::swap(maximize_, linear_program->maximize_);
1050  std::swap(objective_offset_, linear_program->objective_offset_);
1051  std::swap(objective_scaling_factor_,
1052  linear_program->objective_scaling_factor_);
1053  std::swap(columns_are_known_to_be_clean_,
1054  linear_program->columns_are_known_to_be_clean_);
1055  std::swap(transpose_matrix_is_consistent_,
1056  linear_program->transpose_matrix_is_consistent_);
1057  std::swap(integer_variables_list_is_consistent_,
1058  linear_program->integer_variables_list_is_consistent_);
1059  name_.swap(linear_program->name_);
1060  std::swap(first_slack_variable_, linear_program->first_slack_variable_);
1061 }
1062 
1063 void LinearProgram::DeleteColumns(const DenseBooleanRow& columns_to_delete) {
1064  if (columns_to_delete.empty()) return;
1065  integer_variables_list_is_consistent_ = false;
1066  const ColIndex num_cols = num_variables();
1067  ColumnPermutation permutation(num_cols);
1068  ColIndex new_index(0);
1069  for (ColIndex col(0); col < num_cols; ++col) {
1070  permutation[col] = new_index;
1071  if (col >= columns_to_delete.size() || !columns_to_delete[col]) {
1072  objective_coefficients_[new_index] = objective_coefficients_[col];
1073  variable_lower_bounds_[new_index] = variable_lower_bounds_[col];
1074  variable_upper_bounds_[new_index] = variable_upper_bounds_[col];
1075  variable_names_[new_index] = variable_names_[col];
1076  variable_types_[new_index] = variable_types_[col];
1077  ++new_index;
1078  } else {
1079  permutation[col] = kInvalidCol;
1080  }
1081  }
1082 
1083  matrix_.DeleteColumns(columns_to_delete);
1084  objective_coefficients_.resize(new_index, 0.0);
1085  variable_lower_bounds_.resize(new_index, 0.0);
1086  variable_upper_bounds_.resize(new_index, 0.0);
1087  variable_types_.resize(new_index, VariableType::CONTINUOUS);
1088  variable_names_.resize(new_index, "");
1089 
1090  // Remove the id of the deleted columns and adjust the index of the other.
1091  absl::flat_hash_map<std::string, ColIndex>::iterator it =
1092  variable_table_.begin();
1093  while (it != variable_table_.end()) {
1094  const ColIndex col = it->second;
1095  if (col >= columns_to_delete.size() || !columns_to_delete[col]) {
1096  it->second = permutation[col];
1097  ++it;
1098  } else {
1099  // This safely deletes the entry and moves the iterator one step ahead.
1100  variable_table_.erase(it++);
1101  }
1102  }
1103 
1104  // Eventually update transpose_matrix_.
1105  if (transpose_matrix_is_consistent_) {
1106  transpose_matrix_.DeleteRows(
1107  ColToRowIndex(new_index),
1108  reinterpret_cast<const RowPermutation&>(permutation));
1109  }
1110 }
1111 
1113  DCHECK_NE(first_slack_variable_, kInvalidCol);
1114  DenseBooleanRow slack_variables(matrix_.num_cols(), false);
1115  // Restore the bounds on the constraints corresponding to the slack variables.
1116  for (ColIndex slack_variable = first_slack_variable_;
1117  slack_variable < matrix_.num_cols(); ++slack_variable) {
1118  const SparseColumn& column = matrix_.column(slack_variable);
1119  // Slack variables appear only in the constraints for which they were
1120  // created. We can find this constraint by looking at the (only) entry in
1121  // the columnm of the slack variable.
1122  DCHECK_EQ(column.num_entries(), 1);
1123  const RowIndex row = column.EntryRow(EntryIndex(0));
1124  DCHECK_EQ(constraint_lower_bounds_[row], 0.0);
1125  DCHECK_EQ(constraint_upper_bounds_[row], 0.0);
1126  SetConstraintBounds(row, -variable_upper_bounds_[slack_variable],
1127  -variable_lower_bounds_[slack_variable]);
1128  slack_variables[slack_variable] = true;
1129  }
1130 
1131  DeleteColumns(slack_variables);
1132  first_slack_variable_ = kInvalidCol;
1133 }
1134 
1135 namespace {
1136 
1137 // Note that we ignore zeros and infinities because they do not matter from a
1138 // scaling perspective where this function is used.
1139 template <typename FractionalRange>
1140 void UpdateMinAndMaxMagnitude(const FractionalRange& range,
1141  Fractional* min_magnitude,
1142  Fractional* max_magnitude) {
1143  for (const Fractional value : range) {
1144  const Fractional magnitude = std::abs(value);
1145  if (magnitude == 0 || magnitude == kInfinity) continue;
1146  *min_magnitude = std::min(*min_magnitude, magnitude);
1147  *max_magnitude = std::max(*max_magnitude, magnitude);
1148  }
1149 }
1150 
1151 Fractional GetMedianScalingFactor(const DenseRow& range) {
1152  std::vector<Fractional> median;
1153  for (const Fractional value : range) {
1154  if (value == 0.0) continue;
1155  median.push_back(std::abs(value));
1156  }
1157  if (median.empty()) return 1.0;
1158  std::sort(median.begin(), median.end());
1159  return median[median.size() / 2];
1160 }
1161 
1162 Fractional GetMeanScalingFactor(const DenseRow& range) {
1163  Fractional mean = 0.0;
1164  int num_non_zeros = 0;
1165  for (const Fractional value : range) {
1166  if (value == 0.0) continue;
1167  ++num_non_zeros;
1168  mean += std::abs(value);
1169  }
1170  if (num_non_zeros == 0.0) return 1.0;
1171  return mean / static_cast<Fractional>(num_non_zeros);
1172 }
1173 
1174 Fractional ComputeDivisorSoThatRangeContainsOne(Fractional min_magnitude,
1175  Fractional max_magnitude) {
1176  if (min_magnitude > 1.0 && min_magnitude < kInfinity) {
1177  return min_magnitude;
1178  } else if (max_magnitude > 0.0 && max_magnitude < 1.0) {
1179  return max_magnitude;
1180  }
1181  return 1.0;
1182 }
1183 
1184 } // namespace
1185 
1187  GlopParameters::CostScalingAlgorithm method) {
1188  Fractional min_magnitude = kInfinity;
1189  Fractional max_magnitude = 0.0;
1190  UpdateMinAndMaxMagnitude(objective_coefficients(), &min_magnitude,
1191  &max_magnitude);
1192  Fractional cost_scaling_factor = 1.0;
1193  switch (method) {
1194  case GlopParameters::NO_COST_SCALING:
1195  break;
1196  case GlopParameters::CONTAIN_ONE_COST_SCALING:
1197  cost_scaling_factor =
1198  ComputeDivisorSoThatRangeContainsOne(min_magnitude, max_magnitude);
1199  break;
1200  case GlopParameters::MEAN_COST_SCALING:
1201  cost_scaling_factor = GetMeanScalingFactor(objective_coefficients());
1202  break;
1203  case GlopParameters::MEDIAN_COST_SCALING:
1204  cost_scaling_factor = GetMedianScalingFactor(objective_coefficients());
1205  break;
1206  }
1207  if (cost_scaling_factor != 1.0) {
1208  for (ColIndex col(0); col < num_variables(); ++col) {
1209  if (objective_coefficients()[col] == 0.0) continue;
1211  col, objective_coefficients()[col] / cost_scaling_factor);
1212  }
1213  SetObjectiveScalingFactor(objective_scaling_factor() * cost_scaling_factor);
1214  SetObjectiveOffset(objective_offset() / cost_scaling_factor);
1215  }
1216  VLOG(1) << "Objective magnitude range is [" << min_magnitude << ", "
1217  << max_magnitude << "] (dividing by " << cost_scaling_factor << ").";
1218  return cost_scaling_factor;
1219 }
1220 
1222  Fractional min_magnitude = kInfinity;
1223  Fractional max_magnitude = 0.0;
1224  UpdateMinAndMaxMagnitude(variable_lower_bounds(), &min_magnitude,
1225  &max_magnitude);
1226  UpdateMinAndMaxMagnitude(variable_upper_bounds(), &min_magnitude,
1227  &max_magnitude);
1228  UpdateMinAndMaxMagnitude(constraint_lower_bounds(), &min_magnitude,
1229  &max_magnitude);
1230  UpdateMinAndMaxMagnitude(constraint_upper_bounds(), &min_magnitude,
1231  &max_magnitude);
1232  const Fractional bound_scaling_factor =
1233  ComputeDivisorSoThatRangeContainsOne(min_magnitude, max_magnitude);
1234  if (bound_scaling_factor != 1.0) {
1236  bound_scaling_factor);
1237  SetObjectiveOffset(objective_offset() / bound_scaling_factor);
1238  for (ColIndex col(0); col < num_variables(); ++col) {
1240  variable_lower_bounds()[col] / bound_scaling_factor,
1241  variable_upper_bounds()[col] / bound_scaling_factor);
1242  }
1243  for (RowIndex row(0); row < num_constraints(); ++row) {
1245  row, constraint_lower_bounds()[row] / bound_scaling_factor,
1246  constraint_upper_bounds()[row] / bound_scaling_factor);
1247  }
1248  }
1249 
1250  VLOG(1) << "Bounds magnitude range is [" << min_magnitude << ", "
1251  << max_magnitude << "] (dividing bounds by " << bound_scaling_factor
1252  << ").";
1253  return bound_scaling_factor;
1254 }
1255 
1256 void LinearProgram::DeleteRows(const DenseBooleanColumn& rows_to_delete) {
1257  if (rows_to_delete.empty()) return;
1258 
1259  // Deal with row-indexed data and construct the row mapping that will need to
1260  // be applied to every column entry.
1261  const RowIndex num_rows = num_constraints();
1262  RowPermutation permutation(num_rows);
1263  RowIndex new_index(0);
1264  for (RowIndex row(0); row < num_rows; ++row) {
1265  if (row >= rows_to_delete.size() || !rows_to_delete[row]) {
1266  constraint_lower_bounds_[new_index] = constraint_lower_bounds_[row];
1267  constraint_upper_bounds_[new_index] = constraint_upper_bounds_[row];
1268  constraint_names_[new_index].swap(constraint_names_[row]);
1269  permutation[row] = new_index;
1270  ++new_index;
1271  } else {
1272  permutation[row] = kInvalidRow;
1273  }
1274  }
1275  constraint_lower_bounds_.resize(new_index, 0.0);
1276  constraint_upper_bounds_.resize(new_index, 0.0);
1277  constraint_names_.resize(new_index, "");
1278 
1279  // Remove the rows from the matrix.
1280  matrix_.DeleteRows(new_index, permutation);
1281 
1282  // Remove the id of the deleted rows and adjust the index of the other.
1283  absl::flat_hash_map<std::string, RowIndex>::iterator it =
1284  constraint_table_.begin();
1285  while (it != constraint_table_.end()) {
1286  const RowIndex row = it->second;
1287  if (permutation[row] != kInvalidRow) {
1288  it->second = permutation[row];
1289  ++it;
1290  } else {
1291  // This safely deletes the entry and moves the iterator one step ahead.
1292  constraint_table_.erase(it++);
1293  }
1294  }
1295 
1296  // Eventually update transpose_matrix_.
1297  if (transpose_matrix_is_consistent_) {
1298  transpose_matrix_.DeleteColumns(
1299  reinterpret_cast<const DenseBooleanRow&>(rows_to_delete));
1300  }
1301 }
1302 
1304  if (!IsFinite(objective_offset_)) return false;
1305  if (!IsFinite(objective_scaling_factor_)) return false;
1306  if (objective_scaling_factor_ == 0.0) return false;
1307  const ColIndex num_cols = num_variables();
1308  for (ColIndex col(0); col < num_cols; ++col) {
1310  variable_upper_bounds()[col])) {
1311  return false;
1312  }
1313  if (!IsFinite(objective_coefficients()[col])) {
1314  return false;
1315  }
1316  for (const SparseColumn::Entry e : GetSparseColumn(col)) {
1317  if (!IsFinite(e.coefficient())) return false;
1318  }
1319  }
1320  if (constraint_upper_bounds_.size() != constraint_lower_bounds_.size()) {
1321  return false;
1322  }
1323  for (RowIndex row(0); row < constraint_lower_bounds_.size(); ++row) {
1326  return false;
1327  }
1328  }
1329  return true;
1330 }
1331 
1332 std::string LinearProgram::ProblemStatFormatter(
1333  const absl::string_view format) const {
1334  int num_objective_non_zeros = 0;
1335  int num_non_negative_variables = 0;
1336  int num_boxed_variables = 0;
1337  int num_free_variables = 0;
1338  int num_fixed_variables = 0;
1339  int num_other_variables = 0;
1340  const ColIndex num_cols = num_variables();
1341  for (ColIndex col(0); col < num_cols; ++col) {
1342  if (objective_coefficients()[col] != 0.0) {
1343  ++num_objective_non_zeros;
1344  }
1345 
1346  const Fractional lower_bound = variable_lower_bounds()[col];
1347  const Fractional upper_bound = variable_upper_bounds()[col];
1348  const bool lower_bounded = (lower_bound != -kInfinity);
1349  const bool upper_bounded = (upper_bound != kInfinity);
1350 
1351  if (!lower_bounded && !upper_bounded) {
1352  ++num_free_variables;
1353  } else if (lower_bound == 0.0 && !upper_bounded) {
1354  ++num_non_negative_variables;
1355  } else if (!upper_bounded || !lower_bounded) {
1356  ++num_other_variables;
1357  } else if (lower_bound == upper_bound) {
1358  ++num_fixed_variables;
1359  } else {
1360  ++num_boxed_variables;
1361  }
1362  }
1363 
1364  int num_range_constraints = 0;
1365  int num_less_than_constraints = 0;
1366  int num_greater_than_constraints = 0;
1367  int num_equal_constraints = 0;
1368  int num_rhs_non_zeros = 0;
1369  const RowIndex num_rows = num_constraints();
1370  for (RowIndex row(0); row < num_rows; ++row) {
1371  const Fractional lower_bound = constraint_lower_bounds()[row];
1372  const Fractional upper_bound = constraint_upper_bounds()[row];
1373  if (AreBoundsFreeOrBoxed(lower_bound, upper_bound)) {
1374  // TODO(user): we currently count a free row as a range constraint.
1375  // Add a new category?
1376  ++num_range_constraints;
1377  continue;
1378  }
1379  if (lower_bound == upper_bound) {
1380  ++num_equal_constraints;
1381  if (lower_bound != 0) {
1382  ++num_rhs_non_zeros;
1383  }
1384  continue;
1385  }
1386  if (lower_bound == -kInfinity) {
1387  ++num_less_than_constraints;
1388  if (upper_bound != 0) {
1389  ++num_rhs_non_zeros;
1390  }
1391  continue;
1392  }
1393  if (upper_bound == kInfinity) {
1394  ++num_greater_than_constraints;
1395  if (lower_bound != 0) {
1396  ++num_rhs_non_zeros;
1397  }
1398  continue;
1399  }
1400  LOG(DFATAL) << "There is a bug since all possible cases for the row bounds "
1401  "should have been accounted for. row="
1402  << row;
1403  }
1404 
1405  const int num_integer_variables = IntegerVariablesList().size();
1406  const int num_binary_variables = BinaryVariablesList().size();
1407  const int num_non_binary_variables = NonBinaryVariablesList().size();
1408  const int num_continuous_variables =
1409  ColToIntIndex(num_variables()) - num_integer_variables;
1410  auto format_runtime =
1411  absl::ParsedFormat<'d', 'd', 'd', 'd', 'd', 'd', 'd', 'd', 'd', 'd', 'd',
1412  'd', 'd', 'd', 'd', 'd', 'd', 'd'>::New(format);
1413  CHECK(format_runtime);
1414  return absl::StrFormat(
1415  *format_runtime, RowToIntIndex(num_constraints()),
1416  ColToIntIndex(num_variables()), matrix_.num_entries().value(),
1417  num_objective_non_zeros, num_rhs_non_zeros, num_less_than_constraints,
1418  num_greater_than_constraints, num_equal_constraints,
1419  num_range_constraints, num_non_negative_variables, num_boxed_variables,
1420  num_free_variables, num_fixed_variables, num_other_variables,
1421  num_integer_variables, num_binary_variables, num_non_binary_variables,
1422  num_continuous_variables);
1423 }
1424 
1425 std::string LinearProgram::NonZeroStatFormatter(
1426  const absl::string_view format) const {
1427  StrictITIVector<RowIndex, EntryIndex> num_entries_in_row(num_constraints(),
1428  EntryIndex(0));
1429  StrictITIVector<ColIndex, EntryIndex> num_entries_in_column(num_variables(),
1430  EntryIndex(0));
1431  EntryIndex num_entries(0);
1432  const ColIndex num_cols = num_variables();
1433  for (ColIndex col(0); col < num_cols; ++col) {
1434  const SparseColumn& sparse_column = GetSparseColumn(col);
1435  num_entries += sparse_column.num_entries();
1436  num_entries_in_column[col] = sparse_column.num_entries();
1437  for (const SparseColumn::Entry e : sparse_column) {
1438  ++num_entries_in_row[e.row()];
1439  }
1440  }
1441 
1442  // To avoid division by 0 if there are no columns or no rows, we set
1443  // height and width to be at least one.
1444  const int64 height = std::max(RowToIntIndex(num_constraints()), 1);
1445  const int64 width = std::max(ColToIntIndex(num_variables()), 1);
1446  const double fill_rate = 100.0 * static_cast<double>(num_entries.value()) /
1447  static_cast<double>(height * width);
1448 
1449  auto format_runtime =
1450  absl::ParsedFormat<'f', 'd', 'f', 'f', 'd', 'f', 'f'>::New(format);
1451  return absl::StrFormat(
1452  *format_runtime, fill_rate, GetMaxElement(num_entries_in_row).value(),
1453  Average(num_entries_in_row), StandardDeviation(num_entries_in_row),
1454  GetMaxElement(num_entries_in_column).value(),
1455  Average(num_entries_in_column), StandardDeviation(num_entries_in_column));
1456 }
1457 
1458 void LinearProgram::ResizeRowsIfNeeded(RowIndex row) {
1459  DCHECK_GE(row, 0);
1460  if (row >= num_constraints()) {
1461  transpose_matrix_is_consistent_ = false;
1462  matrix_.SetNumRows(row + 1);
1463  constraint_lower_bounds_.resize(row + 1, Fractional(0.0));
1464  constraint_upper_bounds_.resize(row + 1, Fractional(0.0));
1465  constraint_names_.resize(row + 1, "");
1466  }
1467 }
1468 
1470  for (RowIndex constraint(0); constraint < num_constraints(); ++constraint) {
1471  if (constraint_lower_bounds_[constraint] != 0.0 ||
1472  constraint_upper_bounds_[constraint] != 0.0) {
1473  return false;
1474  }
1475  }
1476  const ColIndex num_slack_variables =
1478  return num_constraints().value() == num_slack_variables.value() &&
1480 }
1481 
1483  Fractional tolerance) const {
1484  for (const ColIndex col : IntegerVariablesList()) {
1485  if ((IsFinite(variable_lower_bounds_[col]) &&
1486  !IsIntegerWithinTolerance(variable_lower_bounds_[col], tolerance)) ||
1487  (IsFinite(variable_upper_bounds_[col]) &&
1488  !IsIntegerWithinTolerance(variable_upper_bounds_[col], tolerance))) {
1489  VLOG(1) << "Bounds of variable " << col.value() << " are non-integer ("
1490  << variable_lower_bounds_[col] << ", "
1491  << variable_upper_bounds_[col] << ").";
1492  return false;
1493  }
1494  }
1495  return true;
1496 }
1497 
1499  Fractional tolerance) const {
1500  // Using transpose for this is faster (complexity = O(number of non zeros in
1501  // matrix)) than directly iterating through entries (complexity = O(number of
1502  // constraints * number of variables)).
1503  const SparseMatrix& transpose = GetTransposeSparseMatrix();
1504  for (RowIndex row = RowIndex(0); row < num_constraints(); ++row) {
1505  bool integer_constraint = true;
1506  for (const SparseColumn::Entry var : transpose.column(RowToColIndex(row))) {
1507  if (!IsVariableInteger(RowToColIndex(var.row()))) {
1508  integer_constraint = false;
1509  break;
1510  }
1511  if (!IsIntegerWithinTolerance(var.coefficient(), tolerance)) {
1512  integer_constraint = false;
1513  break;
1514  }
1515  }
1516  if (integer_constraint) {
1517  if ((IsFinite(constraint_lower_bounds_[row]) &&
1518  !IsIntegerWithinTolerance(constraint_lower_bounds_[row],
1519  tolerance)) ||
1520  (IsFinite(constraint_upper_bounds_[row]) &&
1521  !IsIntegerWithinTolerance(constraint_upper_bounds_[row],
1522  tolerance))) {
1523  VLOG(1) << "Bounds of constraint " << row.value()
1524  << " are non-integer (" << constraint_lower_bounds_[row] << ", "
1525  << constraint_upper_bounds_[row] << ").";
1526  return false;
1527  }
1528  }
1529  }
1530  return true;
1531 }
1532 
1533 // --------------------------------------------------------
1534 // ProblemSolution
1535 // --------------------------------------------------------
1536 std::string ProblemSolution::DebugString() const {
1537  std::string s = "Problem status: " + GetProblemStatusString(status);
1538  for (ColIndex col(0); col < primal_values.size(); ++col) {
1539  absl::StrAppendFormat(&s, "\n Var #%d: %s %g", col.value(),
1541  primal_values[col]);
1542  }
1543  s += "\n------------------------------";
1544  for (RowIndex row(0); row < dual_values.size(); ++row) {
1545  absl::StrAppendFormat(&s, "\n Constraint #%d: %s %g", row.value(),
1547  dual_values[row]);
1548  }
1549  return s;
1550 }
1551 
1552 } // namespace glop
1553 } // namespace operations_research
int64 min
Definition: alldiff_cst.cc:138
int64 max
Definition: alldiff_cst.cc:139
#define CHECK(condition)
Definition: base/logging.h:495
#define DCHECK_LE(val1, val2)
Definition: base/logging.h:887
#define DCHECK_NE(val1, val2)
Definition: base/logging.h:886
#define CHECK_EQ(val1, val2)
Definition: base/logging.h:697
#define DCHECK_GE(val1, val2)
Definition: base/logging.h:889
#define CHECK_NE(val1, val2)
Definition: base/logging.h:698
#define DCHECK_LT(val1, val2)
Definition: base/logging.h:888
#define LOG(severity)
Definition: base/logging.h:420
#define DCHECK(condition)
Definition: base/logging.h:884
#define DCHECK_EQ(val1, val2)
Definition: base/logging.h:885
#define VLOG(verboselevel)
Definition: base/logging.h:978
iterator insert(const_iterator pos, const value_type &x)
size_type size() const
bool empty() const
void push_back(const value_type &x)
void swap(StrongVector &x)
SparseMatrix * GetMutableTransposeSparseMatrix()
Definition: lp_data.cc:385
std::string GetObjectiveStatsString() const
Definition: lp_data.cc:450
void SetObjectiveScalingFactor(Fractional objective_scaling_factor)
Definition: lp_data.cc:335
void PopulateFromPermutedLinearProgram(const LinearProgram &lp, const RowPermutation &row_permutation, const ColumnPermutation &col_permutation)
Definition: lp_data.cc:881
void SetVariableBounds(ColIndex col, Fractional lower_bound, Fractional upper_bound)
Definition: lp_data.cc:248
std::string GetVariableName(ColIndex col) const
Definition: lp_data.cc:359
void SetConstraintName(RowIndex row, absl::string_view name)
Definition: lp_data.cc:244
const SparseMatrix & GetTransposeSparseMatrix() const
Definition: lp_data.cc:375
bool SolutionIsWithinVariableBounds(const DenseRow &solution, Fractional absolute_tolerance) const
Definition: lp_data.cc:479
bool BoundsOfIntegerConstraintsAreInteger(Fractional tolerance) const
Definition: lp_data.cc:1498
void SetObjectiveOffset(Fractional objective_offset)
Definition: lp_data.cc:330
void PopulateFromLinearProgram(const LinearProgram &linear_program)
Definition: lp_data.cc:860
std::string GetPrettyProblemStats() const
Definition: lp_data.cc:662
bool SolutionIsMIPFeasible(const DenseRow &solution, Fractional absolute_tolerance) const
Definition: lp_data.cc:527
void SetCoefficient(RowIndex row, ColIndex col, Fractional value)
Definition: lp_data.cc:316
bool BoundsOfIntegerVariablesAreInteger(Fractional tolerance) const
Definition: lp_data.cc:1482
void SetVariableName(ColIndex col, absl::string_view name)
Definition: lp_data.cc:231
std::string DumpSolution(const DenseRow &variable_values) const
Definition: lp_data.cc:645
ColIndex GetSlackVariable(RowIndex row) const
Definition: lp_data.cc:753
const DenseRow & variable_lower_bounds() const
Definition: lp_data.h:228
ColIndex FindOrCreateVariable(const std::string &variable_id)
Definition: lp_data.cc:204
const DenseColumn & constraint_lower_bounds() const
Definition: lp_data.h:214
std::string GetBoundsStatsString() const
Definition: lp_data.cc:463
Fractional ScaleObjective(GlopParameters::CostScalingAlgorithm method)
Definition: lp_data.cc:1186
const std::vector< ColIndex > & BinaryVariablesList() const
Definition: lp_data.cc:284
const DenseRow & objective_coefficients() const
Definition: lp_data.h:222
Fractional RemoveObjectiveScalingAndOffset(Fractional value) const
Definition: lp_data.cc:553
const std::vector< ColIndex > & IntegerVariablesList() const
Definition: lp_data.cc:279
Fractional GetObjectiveCoefficientForMinimizationVersion(ColIndex col) const
Definition: lp_data.cc:418
void SetConstraintBounds(RowIndex row, Fractional lower_bound, Fractional upper_bound)
Definition: lp_data.cc:308
ColIndex CreateNewSlackVariable(bool is_integer_slack_variable, Fractional lower_bound, Fractional upper_bound, const std::string &name)
Definition: lp_data.cc:175
VariableType GetVariableType(ColIndex col) const
Definition: lp_data.cc:371
RowIndex FindOrCreateConstraint(const std::string &constraint_id)
Definition: lp_data.cc:217
void Swap(LinearProgram *linear_program)
Definition: lp_data.cc:1029
void AddConstraints(const SparseMatrix &coefficients, const DenseColumn &left_hand_sides, const DenseColumn &right_hand_sides, const StrictITIVector< RowIndex, std::string > &names)
Definition: lp_data.cc:970
std::string GetPrettyNonZeroStats() const
Definition: lp_data.cc:688
void SetVariableType(ColIndex col, VariableType type)
Definition: lp_data.cc:235
const std::vector< ColIndex > & NonBinaryVariablesList() const
Definition: lp_data.cc:289
bool SolutionIsInteger(const DenseRow &solution, Fractional absolute_tolerance) const
Definition: lp_data.cc:515
SparseColumn * GetMutableSparseColumn(ColIndex col)
Definition: lp_data.cc:412
std::string GetConstraintName(RowIndex row) const
Definition: lp_data.cc:365
void AddSlackVariablesWhereNecessary(bool detect_integer_constraints)
Definition: lp_data.cc:695
const DenseColumn & constraint_upper_bounds() const
Definition: lp_data.h:217
void ComputeSlackVariableValues(DenseRow *solution) const
Definition: lp_data.cc:533
bool SolutionIsLPFeasible(const DenseRow &solution, Fractional absolute_tolerance) const
Definition: lp_data.cc:495
bool IsVariableInteger(ColIndex col) const
Definition: lp_data.cc:294
void SetObjectiveCoefficient(ColIndex col, Fractional value)
Definition: lp_data.cc:325
bool IsVariableBinary(ColIndex col) const
Definition: lp_data.cc:299
Fractional ApplyObjectiveScalingAndOffset(Fractional value) const
Definition: lp_data.cc:548
void DeleteRows(const DenseBooleanColumn &rows_to_delete)
Definition: lp_data.cc:1256
void DeleteColumns(const DenseBooleanRow &columns_to_delete)
Definition: lp_data.cc:1063
const DenseRow & variable_upper_bounds() const
Definition: lp_data.h:231
bool UpdateVariableBoundsToIntersection(const DenseRow &variable_lower_bounds, const DenseRow &variable_upper_bounds)
Definition: lp_data.cc:1004
void PopulateFromDual(const LinearProgram &dual, RowToColMapping *duplicated_rows)
Definition: lp_data.cc:762
const std::string & name() const
Definition: lp_data.h:74
void PopulateFromLinearProgramVariables(const LinearProgram &linear_program)
Definition: lp_data.cc:933
std::string GetDimensionString() const
Definition: lp_data.cc:424
Fractional objective_scaling_factor() const
Definition: lp_data.h:260
void SetMaximizationProblem(bool maximize)
Definition: lp_data.cc:342
void AddConstraintsWithSlackVariables(const SparseMatrix &coefficients, const DenseColumn &left_hand_sides, const DenseColumn &right_hand_sides, const StrictITIVector< RowIndex, std::string > &names, bool detect_integer_constraints_for_slack)
Definition: lp_data.cc:995
const StrictITIVector< ColIndex, VariableType > variable_types() const
Definition: lp_data.h:236
const SparseColumn & GetSparseColumn(ColIndex col) const
Definition: lp_data.cc:408
void PopulateFromInverse(const Permutation &inverse)
RowIndex EntryRow(EntryIndex i) const
Definition: sparse_column.h:51
SparseColumn * mutable_column(ColIndex col)
Definition: sparse.h:181
void PopulateFromPermutedMatrix(const Matrix &a, const RowPermutation &row_perm, const ColumnPermutation &inverse_col_perm)
Definition: sparse.cc:212
void PopulateFromTranspose(const Matrix &input)
Definition: sparse.cc:181
void SetNumRows(RowIndex num_rows)
Definition: sparse.cc:143
Fractional LookUpValue(RowIndex row, ColIndex col) const
Definition: sparse.cc:323
void Swap(SparseMatrix *matrix)
Definition: sparse.cc:158
void ComputeMinAndMaxMagnitudes(Fractional *min_magnitude, Fractional *max_magnitude) const
Definition: sparse.cc:369
void DeleteRows(RowIndex num_rows, const RowPermutation &permutation)
Definition: sparse.cc:289
bool AppendRowsFromSparseMatrix(const SparseMatrix &matrix)
Definition: sparse.cc:302
void DeleteColumns(const DenseBooleanRow &columns_to_delete)
Definition: sparse.cc:276
void PopulateFromSparseMatrix(const SparseMatrix &matrix)
Definition: sparse.cc:206
const SparseColumn & column(ColIndex col) const
Definition: sparse.h:180
void PopulateFromZero(RowIndex num_rows, ColIndex num_cols)
Definition: sparse.cc:164
void SetCoefficient(Index index, Fractional value)
void PopulateFromSparseVector(const SparseVector &sparse_vector)
void assign(IntType size, const T &v)
Definition: lp_types.h:274
const std::string name
int64 value
IntVar * var
Definition: expr_array.cc:1858
int64_t int64
ColIndex col
Definition: markowitz.cc:176
RowIndex row
Definition: markowitz.cc:175
std::string StringifyMonomial(const Fractional a, const std::string &x, bool fraction)
bool IsRightMostSquareMatrixIdentity(const SparseMatrix &matrix)
bool AreBoundsValid(Fractional lower_bound, Fractional upper_bound)
Definition: lp_data.h:683
const RowIndex kInvalidRow(-1)
std::string Stringify(const Fractional x, bool fraction)
Fractional ScalarProduct(const DenseRowOrColumn1 &u, const DenseRowOrColumn2 &v)
StrictITIVector< ColIndex, Fractional > DenseRow
Definition: lp_types.h:299
std::string GetProblemStatusString(ProblemStatus problem_status)
Definition: lp_types.cc:19
Index ColToIntIndex(ColIndex col)
Definition: lp_types.h:54
std::string GetConstraintStatusString(ConstraintStatus status)
Definition: lp_types.cc:90
ColIndex RowToColIndex(RowIndex row)
Definition: lp_types.h:48
bool IsFinite(Fractional value)
Definition: lp_types.h:90
RowIndex ColToRowIndex(ColIndex col)
Definition: lp_types.h:51
const double kEpsilon
Definition: lp_types.h:86
void ApplyPermutation(const Permutation< IndexType > &perm, const ITIVectorType &b, ITIVectorType *result)
Fractional PartialScalarProduct(const DenseRowOrColumn &u, const SparseColumn &v, int max_index)
std::string GetVariableStatusString(VariableStatus status)
Definition: lp_types.cc:71
Index RowToIntIndex(RowIndex row)
Definition: lp_types.h:57
const double kInfinity
Definition: lp_types.h:83
const ColIndex kInvalidCol(-1)
The vehicle routing library lets one model and solve generic vehicle routing problems ranging from th...
bool IsIntegerWithinTolerance(FloatType x, FloatType tolerance)
Definition: fp_utils.h:161
std::vector< double > coefficients
const bool maximize_
Definition: search.cc:2499
ConstraintStatusColumn constraint_statuses
Definition: lp_data.h:676