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"
30
31namespace operations_research {
32namespace glop {
33
34namespace {
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.
38void 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.
49bool 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
58template <class I, class T>
59double 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
72template <class I, class T>
73double 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.
89template <class I, class T>
90T 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
175ColIndex 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
204ColIndex 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
231void 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
244void 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
260void 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
279const std::vector<ColIndex>& LinearProgram::IntegerVariablesList() const {
280 UpdateAllIntegerVariableLists();
281 return integer_variables_list_;
282}
283
284const std::vector<ColIndex>& LinearProgram::BinaryVariablesList() const {
285 UpdateAllIntegerVariableLists();
286 return binary_variables_list_;
287}
288
289const std::vector<ColIndex>& LinearProgram::NonBinaryVariablesList() const {
290 UpdateAllIntegerVariableLists();
291 return non_binary_variables_list_;
292}
293
295 return variable_types_[col] == VariableType::INTEGER ||
296 variable_types_[col] == VariableType::IMPLIED_INTEGER;
297}
298
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
316void LinearProgram::SetCoefficient(RowIndex row, ColIndex col,
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
359std::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
365std::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
435namespace {
436
437template <typename FractionalValues>
438void 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) {
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
558std::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
645std::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
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
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
753ColIndex 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
948void 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();
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
1029void 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
1063void 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
1135namespace {
1136
1137// Note that we ignore zeros and infinities because they do not matter from a
1138// scaling perspective where this function is used.
1139template <typename FractionalRange>
1140void 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
1151Fractional 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
1162Fractional 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
1174Fractional 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
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) {
1311 return false;
1312 }
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
1332std::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
1425std::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
1458void 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// --------------------------------------------------------
1536std::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(),
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
const std::string & name() const
Definition: lp_data.h:74
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
ColIndex FindOrCreateVariable(const std::string &variable_id)
Definition: lp_data.cc:204
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 DenseColumn & constraint_lower_bounds() const
Definition: lp_data.h:214
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
const DenseRow & variable_upper_bounds() const
Definition: lp_data.h:231
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
const DenseColumn & constraint_upper_bounds() const
Definition: lp_data.h:217
const DenseRow & objective_coefficients() const
Definition: lp_data.h:222
void AddSlackVariablesWhereNecessary(bool detect_integer_constraints)
Definition: lp_data.cc:695
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
const StrictITIVector< ColIndex, VariableType > variable_types() const
Definition: lp_data.h:236
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
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 DenseRow & variable_lower_bounds() const
Definition: lp_data.h:228
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 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
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
SparseColumn * mutable_column(ColIndex col)
Definition: sparse.h:181
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
const SparseColumn & column(ColIndex col) const
Definition: sparse.h:180
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
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
RowIndex ColToRowIndex(ColIndex col)
Definition: lp_types.h:51
std::string StringifyMonomial(const Fractional a, const std::string &x, bool fraction)
bool IsFinite(Fractional value)
Definition: lp_types.h:90
bool IsRightMostSquareMatrixIdentity(const SparseMatrix &matrix)
bool AreBoundsValid(Fractional lower_bound, Fractional upper_bound)
Definition: lp_data.h:683
std::string Stringify(const Fractional x, bool fraction)
Fractional ScalarProduct(const DenseRowOrColumn1 &u, const DenseRowOrColumn2 &v)
const ColIndex kInvalidCol(-1)
StrictITIVector< ColIndex, Fractional > DenseRow
Definition: lp_types.h:299
std::string GetProblemStatusString(ProblemStatus problem_status)
Definition: lp_types.cc:19
std::string GetConstraintStatusString(ConstraintStatus status)
Definition: lp_types.cc:90
Index RowToIntIndex(RowIndex row)
Definition: lp_types.h:57
Index ColToIntIndex(ColIndex col)
Definition: lp_types.h:54
ColIndex RowToColIndex(RowIndex row)
Definition: lp_types.h:48
const RowIndex kInvalidRow(-1)
Fractional PartialScalarProduct(const DenseRowOrColumn &u, const SparseColumn &v, int max_index)
const double kEpsilon
Definition: lp_types.h:86
void ApplyPermutation(const Permutation< IndexType > &perm, const ITIVectorType &b, ITIVectorType *result)
std::string GetVariableStatusString(VariableStatus status)
Definition: lp_types.cc:71
const double kInfinity
Definition: lp_types.h:83
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