Google OR-Tools v9.15
a fast and portable software suite for combinatorial optimization
Loading...
Searching...
No Matches
preprocessor.h
Go to the documentation of this file.
1// Copyright 2010-2025 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
14// This file contains the presolving code for a LinearProgram.
15//
16// A classical reference is:
17// E. D. Andersen, K. D. Andersen, "Presolving in linear programming.",
18// Mathematical Programming 71 (1995) 221-245.
19
20#ifndef ORTOOLS_GLOP_PREPROCESSOR_H_
21#define ORTOOLS_GLOP_PREPROCESSOR_H_
22
23#include <deque>
24#include <memory>
25#include <string>
26#include <vector>
27
28#include "absl/strings/string_view.h"
35
36namespace operations_research {
37namespace glop {
38
39// --------------------------------------------------------
40// Preprocessor
41// --------------------------------------------------------
42// This is the base class for preprocessors.
43//
44// TODO(user): On most preprocessors, calling Run() more than once will not work
45// as expected. Fix? or document and crash in debug if this happens.
46class Preprocessor {
47 public:
48 explicit Preprocessor(const GlopParameters* parameters);
49 Preprocessor(const Preprocessor&) = delete;
50 Preprocessor& operator=(const Preprocessor&) = delete;
51 virtual ~Preprocessor();
53 // Runs the preprocessor by modifying the given linear program. Returns true
54 // if a postsolve step will be needed (i.e. RecoverSolution() is not the
55 // identity function). Also updates status_ to something different from
56 // ProblemStatus::INIT if the problem was solved (including bad statuses
57 // like ProblemStatus::ABNORMAL, ProblemStatus::INFEASIBLE, etc.).
58 virtual bool Run(LinearProgram* lp) = 0;
59
60 // Stores the optimal solution of the linear program that was passed to
61 // Run(). The given solution needs to be set to the optimal solution of the
62 // linear program "modified" by Run().
63 virtual void RecoverSolution(ProblemSolution* solution) const = 0;
64
65 // Returns the status of the preprocessor.
66 // A status different from ProblemStatus::INIT means that the problem is
67 // solved and there is not need to call subsequent preprocessors.
68 ProblemStatus status() const { return status_; }
69
70 // Some preprocessors only need minimal changes when used with integer
71 // variables in a MIP context. Setting this to true allows to consider integer
72 // variables as integer in these preprocessors.
73 //
74 // Not all preprocessors handle integer variables correctly, calling this
75 // function on them will cause a LOG(FATAL).
76 virtual void UseInMipContext() { in_mip_context_ = true; }
77
78 void SetTimeLimit(TimeLimit* time_limit) { time_limit_ = time_limit; }
79
80 protected:
81 // Returns true if a is less than b (or slighlty greater than b with a given
82 // tolerance).
84 return ::operations_research::IsSmallerWithinTolerance(
85 a, b, Fractional(parameters_.solution_feasibility_tolerance()));
86 }
88 Fractional b) const {
89 // TODO(user): use an absolute tolerance here to be even more defensive?
90 return ::operations_research::IsSmallerWithinTolerance(
91 a, b, Fractional(parameters_.preprocessor_zero_tolerance()));
92 }
93
97 std::unique_ptr<TimeLimit> infinite_time_limit_;
99};
101// --------------------------------------------------------
102// MainLpPreprocessor
103// --------------------------------------------------------
104// This is the main LP preprocessor responsible for calling all the other
105// preprocessors in this file, possibly more than once.
106class MainLpPreprocessor : public Preprocessor {
107 public:
108 explicit MainLpPreprocessor(const GlopParameters* parameters)
109 : Preprocessor(parameters) {}
112 ~MainLpPreprocessor() override = default;
113
114 bool Run(LinearProgram* lp) final;
117 // Like RecoverSolution but destroys data structures as it goes to reduce peak
118 // RAM use. After calling this the MainLpPreprocessor object may no longer be
119 // used.
121
122 void SetLogger(SolverLogger* logger) { logger_ = logger; }
123
124 private:
125 // Runs the given preprocessor and push it on preprocessors_ for the postsolve
126 // step when needed.
127 void RunAndPushIfRelevant(std::unique_ptr<Preprocessor> preprocessor,
128 absl::string_view name, TimeLimit* time_limit,
129 LinearProgram* lp);
130
131 // Stack of preprocessors currently applied to the lp that needs postsolve.
132 std::vector<std::unique_ptr<Preprocessor>> preprocessors_;
133
134 // Helpers for logging during presolve.
135 SolverLogger default_logger_;
136 SolverLogger* logger_ = &default_logger_;
137
138 // Initial dimension of the lp given to Run(), for displaying purpose.
139 EntryIndex initial_num_entries_;
140 RowIndex initial_num_rows_;
141 ColIndex initial_num_cols_;
142};
143
144// --------------------------------------------------------
145// ColumnDeletionHelper
146// --------------------------------------------------------
147
148// Some preprocessors need to save columns/rows of the matrix for the postsolve.
149// This class helps them do that.
150//
151// Note that we used to simply use a SparseMatrix, which is like a vector of
152// SparseColumn. However on large problem with 10+ millions columns, each empty
153// SparseColumn take 48 bytes, so if we run like 10 presolve step that save as
154// little as 1 columns, we already are at 4GB memory for nothing!
155class ColumnsSaver {
156 public:
157 // Saves a column. The first version CHECKs that it is not already done.
158 void SaveColumn(ColIndex col, const SparseColumn& column);
159 void SaveColumnIfNotAlreadyDone(ColIndex col, const SparseColumn& column);
160
161 // Returns the saved column. The first version CHECKs that it was saved.
162 const SparseColumn& SavedColumn(ColIndex col) const;
163 const SparseColumn& SavedOrEmptyColumn(ColIndex col) const;
164
165 private:
166 SparseColumn empty_column_;
167 absl::flat_hash_map<ColIndex, int> saved_columns_index_;
168
169 // TODO(user): We could optimize further since all these are read only, we
170 // could use a CompactSparseMatrix instead.
171 std::deque<SparseColumn> saved_columns_;
172};
173
174// Help preprocessors deal with column deletion.
176 public:
177 ColumnDeletionHelper() = default;
180
181 // Remember the given column as "deleted" so that it can later be restored
182 // by RestoreDeletedColumns(). Optionally, the caller may indicate the
183 // value and status of the corresponding variable so that it is automatically
184 // restored; if they don't then the restored value and status will be junk
185 // and must be set by the caller.
186 //
187 // The actual deletion is done by LinearProgram::DeleteColumns().
188 void MarkColumnForDeletion(ColIndex col);
189 void MarkColumnForDeletionWithState(ColIndex col, Fractional value,
190 VariableStatus status);
191
192 // From a solution omitting the deleted column, expands it and inserts the
193 // deleted columns. If values and statuses for the corresponding variables
194 // were saved, they'll be restored.
196
197 // Returns whether or not the given column is marked for deletion.
198 bool IsColumnMarked(ColIndex col) const {
199 return col < is_column_deleted_.size() && is_column_deleted_[col];
200 }
201
202 // Returns a Boolean vector of the column to be deleted.
203 const DenseBooleanRow& GetMarkedColumns() const { return is_column_deleted_; }
205 // Returns true if no columns have been marked for deletion.
206 bool IsEmpty() const { return is_column_deleted_.empty(); }
207
208 // Restores the class to its initial state.
209 void Clear();
210
211 // Returns the value that will be restored by
212 // RestoreDeletedColumnInSolution(). Note that only the marked position value
213 // make sense.
214 const DenseRow& GetStoredValue() const { return stored_value_; }
215
216 private:
217 DenseBooleanRow is_column_deleted_;
218
219 // Note that this vector has the same size as is_column_deleted_ and that
220 // the value of the variable corresponding to a deleted column col is stored
221 // at position col. Values of columns not deleted are not used. We use this
222 // data structure so columns can be deleted in any order if needed.
223 DenseRow stored_value_;
224 VariableStatusRow stored_status_;
225};
226
227// --------------------------------------------------------
228// RowDeletionHelper
229// --------------------------------------------------------
230// Help preprocessors deal with row deletion.
231class RowDeletionHelper {
232 public:
233 RowDeletionHelper() = default;
234 RowDeletionHelper(const RowDeletionHelper&) = delete;
236
237 // Returns true if no rows have been marked for deletion.
238 bool IsEmpty() const { return is_row_deleted_.empty(); }
240 // Restores the class to its initial state.
241 void Clear();
243 // Adds a deleted row to the helper.
244 void MarkRowForDeletion(RowIndex row);
245
246 // If the given row was marked for deletion, unmark it.
247 void UnmarkRow(RowIndex row);
248
249 // Returns a Boolean vector of the row to be deleted.
250 const DenseBooleanColumn& GetMarkedRows() const;
251
252 // Returns whether or not the given row is marked for deletion.
253 bool IsRowMarked(RowIndex row) const {
254 return row < is_row_deleted_.size() && is_row_deleted_[row];
255 }
256
257 // From a solution without the deleted rows, expand it by restoring
258 // the deleted rows to a VariableStatus::BASIC status with 0.0 value.
259 // This latter value is important, many preprocessors rely on it.
262 private:
263 DenseBooleanColumn is_row_deleted_;
264};
265
266// --------------------------------------------------------
267// EmptyColumnPreprocessor
268// --------------------------------------------------------
269// Removes the empty columns from the problem.
270class EmptyColumnPreprocessor final : public Preprocessor {
271 public:
272 explicit EmptyColumnPreprocessor(const GlopParameters* parameters)
273 : Preprocessor(parameters) {}
276 ~EmptyColumnPreprocessor() final = default;
277 bool Run(LinearProgram* lp) final;
278 void RecoverSolution(ProblemSolution* solution) const final;
279
280 private:
281 ColumnDeletionHelper column_deletion_helper_;
283
284// --------------------------------------------------------
285// ProportionalColumnPreprocessor
286// --------------------------------------------------------
287// Removes the proportional columns from the problem when possible. Two columns
288// are proportional if one is a non-zero scalar multiple of the other.
289//
290// Note that in the linear programming literature, two proportional columns are
291// usually called duplicates. The notion is the same once the problem has been
292// scaled. However, during presolve the columns can't be assumed to be scaled,
293// so it makes sense to use the more general notion of proportional columns.
294class ProportionalColumnPreprocessor final : public Preprocessor {
295 public:
296 explicit ProportionalColumnPreprocessor(const GlopParameters* parameters)
297 : Preprocessor(parameters) {}
299 delete;
301 const ProportionalColumnPreprocessor&) = delete;
302 ~ProportionalColumnPreprocessor() final = default;
303 bool Run(LinearProgram* lp) final;
304 void RecoverSolution(ProblemSolution* solution) const final;
305 void UseInMipContext() final { LOG(FATAL) << "Not implemented."; }
307 private:
308 // Postsolve information about proportional columns with the same scaled cost
309 // that were merged during presolve.
311 // The proportionality factor of each column. If two columns are proportional
312 // with factor p1 and p2 then p1 times the first column is the same as p2
313 // times the second column.
314 DenseRow column_factors_;
315
316 // If merged_columns_[col] != kInvalidCol, then column col has been merged
317 // into the column merged_columns_[col].
318 ColMapping merged_columns_;
319
320 // The old and new variable bounds.
321 DenseRow lower_bounds_;
322 DenseRow upper_bounds_;
323 DenseRow new_lower_bounds_;
324 DenseRow new_upper_bounds_;
325
326 ColumnDeletionHelper column_deletion_helper_;
327};
328
329// --------------------------------------------------------
330// ProportionalRowPreprocessor
331// --------------------------------------------------------
332// Removes the proportional rows from the problem.
333// The linear programming literature also calls such rows duplicates, see the
334// same remark above for columns in ProportionalColumnPreprocessor.
335class ProportionalRowPreprocessor final : public Preprocessor {
336 public:
337 explicit ProportionalRowPreprocessor(const GlopParameters* parameters)
338 : Preprocessor(parameters) {}
339 ProportionalRowPreprocessor(const ProportionalRowPreprocessor&) = delete;
340 ProportionalRowPreprocessor& operator=(const ProportionalRowPreprocessor&) =
341 delete;
342 ~ProportionalRowPreprocessor() final = default;
343 bool Run(LinearProgram* lp) final;
344 void RecoverSolution(ProblemSolution* solution) const final;
345
346 private:
347 // Informations about proportional rows, only filled for such rows.
348 DenseColumn row_factors_;
349 RowMapping upper_bound_sources_;
350 RowMapping lower_bound_sources_;
352 bool lp_is_maximization_problem_;
353 RowDeletionHelper row_deletion_helper_;
355
356// --------------------------------------------------------
357// SingletonPreprocessor
358// --------------------------------------------------------
359// Removes as many singleton rows and singleton columns as possible from the
360// problem. Note that not all types of singleton columns can be removed. See the
361// comments below on the SingletonPreprocessor functions for more details.
362//
363// TODO(user): Generalize the design used in this preprocessor to a general
364// "propagation" framework in order to apply as many reductions as possible in
365// an efficient manner.
366
367// Holds a triplet (row, col, coefficient).
368struct MatrixEntry {
369 MatrixEntry(RowIndex _row, ColIndex _col, Fractional _coeff)
370 : row(_row), col(_col), coeff(_coeff) {}
371 RowIndex row;
372 ColIndex col;
373 Fractional coeff;
374};
375
376// Stores the information needed to undo a singleton row/column deletion.
377class SingletonUndo {
378 public:
379 // The type of a given operation.
380 typedef enum : uint8_t {
381 ZERO_COST_SINGLETON_COLUMN,
382 SINGLETON_ROW,
383 SINGLETON_COLUMN_IN_EQUALITY,
384 MAKE_CONSTRAINT_AN_EQUALITY,
385 } OperationType;
386
387 // Stores the information, which together with the field deleted_columns_ and
388 // deleted_rows_ of SingletonPreprocessor, are needed to undo an operation
389 // with the given type. Note that all the arguments must refer to the linear
390 // program BEFORE the operation is applied.
391 SingletonUndo(OperationType type, const LinearProgram& lp, MatrixEntry e,
392 ConstraintStatus status);
394 // Undo the operation saved in this class, taking into account the saved
395 // column and row (at the row/col given by Entry()) passed by the calling
396 // instance of SingletonPreprocessor. Note that the operations must be undone
397 // in the reverse order of the one in which they were applied.
398 void Undo(const GlopParameters& parameters, const SparseColumn& saved_column,
399 const SparseColumn& saved_row, ProblemSolution* solution) const;
401 const MatrixEntry& Entry() const { return e_; }
402
403 private:
404 // Actual undo functions for each OperationType.
405 // Undo() just calls the correct one.
406 void SingletonRowUndo(const SparseColumn& saved_column,
408 void ZeroCostSingletonColumnUndo(const GlopParameters& parameters,
409 const SparseColumn& saved_row,
411 void SingletonColumnInEqualityUndo(const GlopParameters& parameters,
412 const SparseColumn& saved_row,
414 void MakeConstraintAnEqualityUndo(ProblemSolution* solution) const;
415
416 // All the information needed during undo.
418 bool is_maximization_;
419 MatrixEntry e_;
420 Fractional cost_;
421
422 // TODO(user): regroup the pair (lower bound, upper bound) in a bound class?
423 Fractional variable_lower_bound_;
424 Fractional variable_upper_bound_;
425 Fractional constraint_lower_bound_;
426 Fractional constraint_upper_bound_;
427
428 // This in only used with MAKE_CONSTRAINT_AN_EQUALITY undo.
429 // TODO(user): Clean that up using many Undo classes and virtual functions.
430 ConstraintStatus constraint_status_;
431};
432
433// Deletes as many singleton rows or singleton columns as possible. Note that
434// each time we delete a row or a column, new singletons may be created.
435class SingletonPreprocessor final : public Preprocessor {
436 public:
437 explicit SingletonPreprocessor(const GlopParameters* parameters)
438 : Preprocessor(parameters) {}
439 SingletonPreprocessor(const SingletonPreprocessor&) = delete;
440 SingletonPreprocessor& operator=(const SingletonPreprocessor&) = delete;
441 ~SingletonPreprocessor() final = default;
442 bool Run(LinearProgram* lp) final;
443 void RecoverSolution(ProblemSolution* solution) const final;
444
445 private:
446 // Returns the MatrixEntry of the given singleton row or column, taking into
447 // account the rows and columns that were already deleted.
448 MatrixEntry GetSingletonColumnMatrixEntry(ColIndex col,
449 const SparseMatrix& matrix);
450 MatrixEntry GetSingletonRowMatrixEntry(RowIndex row,
451 const SparseMatrix& matrix_transpose);
452
453 // A singleton row can always be removed by changing the corresponding
454 // variable bounds to take into account the bounds on this singleton row.
455 void DeleteSingletonRow(MatrixEntry e, LinearProgram* lp);
457 // Internal operation when removing a zero-cost singleton column corresponding
458 // to the given entry. This modifies the constraint bounds to take into
459 // account the bounds of the corresponding variable.
460 void UpdateConstraintBoundsWithVariableBounds(MatrixEntry e,
461 LinearProgram* lp);
462
463 // Checks if all other variables in the constraint are integer and the
464 // coefficients are divisible by the coefficient of the singleton variable.
465 bool IntegerSingletonColumnIsRemovable(const MatrixEntry& matrix_entry,
466 const LinearProgram& lp) const;
467
468 // A singleton column with a cost of zero can always be removed by changing
469 // the corresponding constraint bounds to take into account the bound of this
470 // singleton column.
471 void DeleteZeroCostSingletonColumn(const SparseMatrix& matrix_transpose,
473
474 // Returns true if the constraint associated to the given singleton column was
475 // an equality or could be made one:
476 // If a singleton variable is free in a direction that improves the cost, then
477 // we can always move it as much as possible in this direction. Only the
478 // constraint will stop us, making it an equality. If the constraint doesn't
479 // stop us, then the program is unbounded (provided that there is a feasible
480 // solution).
481 //
482 // Note that this operation does not need any "undo" during the post-solve. At
483 // optimality, the dual value on the constraint row will be of the correct
484 // sign, and relaxing the constraint bound will not impact the dual
485 // feasibility of the solution.
486 //
487 // TODO(user): this operation can be generalized to columns with just one
488 // blocking constraint. Investigate how to use this. The 'reverse' can
489 // probably also be done, relaxing a constraint that is blocking a
490 // unconstrained variable.
491 bool MakeConstraintAnEqualityIfPossible(const SparseMatrix& matrix_transpose,
493
494 // If a singleton column appears in an equality, we can remove its cost by
495 // changing the other variables cost using the constraint. We can then delete
496 // the column like in DeleteZeroCostSingletonColumn().
497 void DeleteSingletonColumnInEquality(const SparseMatrix& matrix_transpose,
499
500 ColumnDeletionHelper column_deletion_helper_;
501 RowDeletionHelper row_deletion_helper_;
502 std::vector<SingletonUndo> undo_stack_;
503
504 // This is used as a "cache" by MakeConstraintAnEqualityIfPossible() to avoid
505 // scanning more than once each row. See the code to see how this is used.
506 util_intops::StrongVector<RowIndex, bool> row_sum_is_cached_;
508 row_lb_sum_;
510 row_ub_sum_;
511
512 // TODO(user): It is annoying that we need to store a part of the matrix that
513 // is not deleted here. This extra memory usage might show the limit of our
514 // presolve architecture that does not require a new matrix factorization on
515 // the original problem to reconstruct the solution.
516 ColumnsSaver columns_saver_;
517 ColumnsSaver rows_saver_;
518};
519
520// --------------------------------------------------------
521// FixedVariablePreprocessor
522// --------------------------------------------------------
523// Removes the fixed variables from the problem.
524class FixedVariablePreprocessor final : public Preprocessor {
525 public:
526 explicit FixedVariablePreprocessor(const GlopParameters* parameters)
527 : Preprocessor(parameters) {}
530 delete;
531 ~FixedVariablePreprocessor() final = default;
532 bool Run(LinearProgram* lp) final;
533 void RecoverSolution(ProblemSolution* solution) const final;
534
535 private:
536 ColumnDeletionHelper column_deletion_helper_;
537};
538
539// --------------------------------------------------------
540// ForcingAndImpliedFreeConstraintPreprocessor
541// --------------------------------------------------------
542// This preprocessor computes for each constraint row the bounds that are
543// implied by the variable bounds and applies one of the following reductions:
545// * If the intersection of the implied bounds and the current constraint bounds
546// is empty (modulo some tolerance), the problem is INFEASIBLE.
548// * If the intersection of the implied bounds and the current constraint bounds
549// is a singleton (modulo some tolerance), then the constraint is said to be
550// forcing and all the variables that appear in it can be fixed to one of their
551// bounds. All these columns and the constraint row is removed.
552//
553// * If the implied bounds are included inside the current constraint bounds
554// (modulo some tolerance) then the constraint is said to be redundant or
555// implied free. Its bounds are relaxed and the constraint will be removed
556// later by the FreeConstraintPreprocessor.
557//
558// * Otherwise, wo do nothing.
560 public:
562 const GlopParameters* parameters)
563 : Preprocessor(parameters) {}
569 bool Run(LinearProgram* lp) final;
570 void RecoverSolution(ProblemSolution* solution) const final;
571
572 private:
573 bool lp_is_maximization_problem_;
574 DenseRow costs_;
575 DenseBooleanColumn is_forcing_up_;
576 ColumnDeletionHelper column_deletion_helper_;
577 RowDeletionHelper row_deletion_helper_;
578 ColumnsSaver columns_saver_;
580
581// --------------------------------------------------------
582// ImpliedFreePreprocessor
583// --------------------------------------------------------
584// It is possible to compute "implied" bounds on a variable from the bounds of
585// all the other variables and the constraints in which this variable take
586// place. If such "implied" bounds are inside the variable bounds, then the
587// variable bounds can be relaxed and the variable is said to be "implied free".
589// This preprocessor detects the implied free variables and make as many as
590// possible free with a priority towards low-degree columns. This transformation
591// will make the simplex algorithm more efficient later, but will also make it
592// possible to reduce the problem by applying subsequent transformations:
593//
594// * The SingletonPreprocessor already deals with implied free singleton
595// variables and removes the columns and the rows in which they appear.
596//
597// * Any multiple of the column of a free variable can be added to any other
598// column without changing the linear program solution. This is the dual
599// counterpart of the fact that any multiple of an equality row can be added to
600// any row.
601//
602// TODO(user): Only process doubleton columns so we have more chance in the
603// later passes to create more doubleton columns? Such columns lead to a smaller
604// problem thanks to the DoubletonFreeColumnPreprocessor.
605class ImpliedFreePreprocessor final : public Preprocessor {
606 public:
607 explicit ImpliedFreePreprocessor(const GlopParameters* parameters)
608 : Preprocessor(parameters) {}
611 ~ImpliedFreePreprocessor() final = default;
612 bool Run(LinearProgram* lp) final;
613 void RecoverSolution(ProblemSolution* solution) const final;
614
615 private:
616 // This preprocessor adds fixed offsets to some variables. We remember those
617 // here to un-offset them in RecoverSolution().
618 DenseRow variable_offsets_;
619
620 // This preprocessor causes some variables who would normally be
621 // AT_{LOWER,UPPER}_BOUND to be VariableStatus::FREE. We store the restore
622 // value of these variables; which will only be used (eg. restored) if the
623 // variable actually turns out to be VariableStatus::FREE.
624 VariableStatusRow postsolve_status_of_free_variables_;
625};
626
627// --------------------------------------------------------
628// DoubletonFreeColumnPreprocessor
629// --------------------------------------------------------
630// This preprocessor removes one of the two rows in which a doubleton column of
631// a free variable appears. Since we can add any multiple of such a column to
632// any other column, the way this works is that we can always remove all the
633// entries on one row.
634//
635// Actually, we can remove all the entries except the one of the free column.
636// But we will be left with a singleton row that we can delete in the same way
637// as what is done in SingletonPreprocessor. That is by reporting the constraint
638// bounds into the one of the originally free variable. After this operation,
639// the doubleton free column will become a singleton and may or may not be
640// removed later by the SingletonPreprocessor.
641//
642// Note that this preprocessor can be seen as the dual of the
643// DoubletonEqualityRowPreprocessor since when taking the dual, an equality row
644// becomes a free variable and vice versa.
645//
646// Note(user): As far as I know, this doubleton free column procedure is more
647// general than what can be found in the research papers or in any of the linear
648// solver open source codes as of July 2013. All of them only process such
649// columns if one of the two rows is also an equality which is not actually
650// required. Most probably, commercial solvers do use it though.
652 public:
653 explicit DoubletonFreeColumnPreprocessor(const GlopParameters* parameters)
654 : Preprocessor(parameters) {}
656 delete;
658 const DoubletonFreeColumnPreprocessor&) = delete;
659 ~DoubletonFreeColumnPreprocessor() final = default;
660 bool Run(LinearProgram* lp) final;
661 void RecoverSolution(ProblemSolution* solution) const final;
662
663 private:
664 enum RowChoice : int {
665 DELETED = 0,
666 MODIFIED = 1,
667 // This is just a constant for the number of rows in a doubleton column.
668 // That is 2, one will be DELETED, the other MODIFIED.
669 NUM_ROWS = 2,
670 };
671 struct RestoreInfo {
672 // The index of the original free doubleton column and its objective.
673 ColIndex col;
674 Fractional objective_coefficient;
676 // The row indices of the two involved rows and their coefficients on
677 // column col.
678 RowIndex row[NUM_ROWS];
679 Fractional coeff[NUM_ROWS];
680
681 // The deleted row as a column.
682 SparseColumn deleted_row_as_column;
683 };
684
685 std::vector<RestoreInfo> restore_stack_;
686 RowDeletionHelper row_deletion_helper_;
687};
688
689// --------------------------------------------------------
690// UnconstrainedVariablePreprocessor
691// --------------------------------------------------------
692// If for a given variable, none of the constraints block it in one direction
693// and this direction improves the objective, then this variable can be fixed to
694// its bound in this direction. If this bound is infinite and the variable cost
695// is non-zero, then the problem is unbounded.
696//
697// More generally, by using the constraints and the variables that are unbounded
698// on one side, one can derive bounds on the dual values. These can be
699// translated into bounds on the reduced costs or the columns, which may force
700// variables to their bounds. This is called forcing and dominated columns in
701// the Andersen & Andersen paper.
702class UnconstrainedVariablePreprocessor final : public Preprocessor {
703 public:
704 explicit UnconstrainedVariablePreprocessor(const GlopParameters* parameters)
705 : Preprocessor(parameters) {}
706 UnconstrainedVariablePreprocessor(const UnconstrainedVariablePreprocessor&) =
707 delete;
708 UnconstrainedVariablePreprocessor& operator=(
709 const UnconstrainedVariablePreprocessor&) = delete;
710 ~UnconstrainedVariablePreprocessor() final = default;
711 bool Run(LinearProgram* lp) final;
712 void RecoverSolution(ProblemSolution* solution) const final;
713
714 // Removes the given variable and all the rows in which it appears: If a
715 // variable is unconstrained with a zero cost, then all the constraints in
716 // which it appears can be made free! More precisely, during postsolve, if
717 // such a variable is unconstrained towards +kInfinity, for any activity value
718 // of the involved constraints, an M exists such that for each value of the
719 // variable >= M the problem will be feasible.
720 //
721 // The algorithm during postsolve is to find a feasible value for all such
722 // variables while trying to keep their magnitudes small (for better numerical
723 // behavior). target_bound should take only two possible values: +/-kInfinity.
724 void RemoveZeroCostUnconstrainedVariable(ColIndex col,
725 Fractional target_bound,
726 LinearProgram* lp);
727
728 private:
729 // Lower/upper bounds on the feasible dual value. We use constraints and
730 // variables unbounded in one direction to derive these bounds. We use these
731 // bounds to compute bounds on the reduced costs of the problem variables.
732 // Note that any finite bounds on a reduced cost means that the variable
733 // (ignoring its domain) can move freely in one direction.
734 DenseColumn dual_lb_;
735 DenseColumn dual_ub_;
737 // Indicates if a given column may have participated in the current lb/ub
738 // on the reduced cost of the same column.
739 DenseBooleanRow may_have_participated_ub_;
740 DenseBooleanRow may_have_participated_lb_;
741
742 ColumnDeletionHelper column_deletion_helper_;
743 RowDeletionHelper row_deletion_helper_;
744 ColumnsSaver rows_saver_;
745 DenseColumn rhs_;
746 DenseColumn activity_sign_correction_;
747 DenseBooleanRow is_unbounded_;
748};
749
750// --------------------------------------------------------
751// FreeConstraintPreprocessor
752// --------------------------------------------------------
753// Removes the constraints with no bounds from the problem.
754class FreeConstraintPreprocessor final : public Preprocessor {
755 public:
756 explicit FreeConstraintPreprocessor(const GlopParameters* parameters)
757 : Preprocessor(parameters) {}
760 delete;
761 ~FreeConstraintPreprocessor() final = default;
762 bool Run(LinearProgram* lp) final;
763 void RecoverSolution(ProblemSolution* solution) const final;
764
765 private:
766 RowDeletionHelper row_deletion_helper_;
767};
768
769// --------------------------------------------------------
770// EmptyConstraintPreprocessor
771// --------------------------------------------------------
772// Removes the constraints with no coefficients from the problem.
773class EmptyConstraintPreprocessor final : public Preprocessor {
774 public:
775 explicit EmptyConstraintPreprocessor(const GlopParameters* parameters)
776 : Preprocessor(parameters) {}
777 EmptyConstraintPreprocessor(const EmptyConstraintPreprocessor&) = delete;
778 EmptyConstraintPreprocessor& operator=(const EmptyConstraintPreprocessor&) =
779 delete;
780 ~EmptyConstraintPreprocessor() final = default;
781 bool Run(LinearProgram* lp) final;
783
784 private:
785 RowDeletionHelper row_deletion_helper_;
788// --------------------------------------------------------
789// SingletonColumnSignPreprocessor
790// --------------------------------------------------------
791// Make sure that the only coefficient of all singleton columns (i.e. column
792// with only one entry) is positive. This is because this way the column will
793// be transformed in an identity column by the scaling. This will lead to more
794// efficient solve when this column is involved.
796 public:
797 explicit SingletonColumnSignPreprocessor(const GlopParameters* parameters)
798 : Preprocessor(parameters) {}
800 delete;
802 const SingletonColumnSignPreprocessor&) = delete;
804 bool Run(LinearProgram* lp) final;
806
807 private:
808 std::vector<ColIndex> changed_columns_;
809};
811// --------------------------------------------------------
812// DoubletonEqualityRowPreprocessor
813// --------------------------------------------------------
814// Reduce equality constraints involving two variables (i.e. aX + bY = c),
815// by substitution (and thus removal) of one of the variables by the other
816// in all the constraints that it is involved in.
818 public:
819 explicit DoubletonEqualityRowPreprocessor(const GlopParameters* parameters)
820 : Preprocessor(parameters) {}
822 delete;
824 const DoubletonEqualityRowPreprocessor&) = delete;
825 ~DoubletonEqualityRowPreprocessor() final = default;
826 bool Run(LinearProgram* lp) final;
828
829 private:
830 enum ColChoice : int {
831 DELETED = 0,
832 MODIFIED = 1,
833 // For `for()` loops iterating over the ColChoice values and/or arrays.
834 NUM_DOUBLETON_COLS = 2,
835 };
836 static ColChoice OtherColChoice(ColChoice x) {
837 return x == DELETED ? MODIFIED : DELETED;
838 }
839
840 ColumnDeletionHelper column_deletion_helper_;
841 RowDeletionHelper row_deletion_helper_;
842
843 struct RestoreInfo {
844 // The row index of the doubleton equality constraint, and its constant.
845 RowIndex row;
846 Fractional rhs; // The constant c in the equality aX + bY = c.
847
848 // The indices and the data of the two columns that we touched, exactly
849 // as they were beforehand.
850 ColIndex col[NUM_DOUBLETON_COLS];
851 Fractional coeff[NUM_DOUBLETON_COLS];
852 Fractional lb[NUM_DOUBLETON_COLS];
853 Fractional ub[NUM_DOUBLETON_COLS];
854 Fractional objective_coefficient[NUM_DOUBLETON_COLS];
856 // If the modified variable has status AT_[LOWER,UPPER]_BOUND, then we'll
857 // set one of the two original variables to one of its bounds, and set the
858 // other to VariableStatus::BASIC. We store this information (which variable
859 // will be set to one of its bounds, and which bound) for each possible
860 // outcome.
861 struct ColChoiceAndStatus {
862 ColChoice col_choice;
864 Fractional value;
865 ColChoiceAndStatus() : col_choice(), status(), value(0.0) {}
866 ColChoiceAndStatus(ColChoice c, VariableStatus s, Fractional v)
867 : col_choice(c), status(s), value(v) {}
868 };
869 ColChoiceAndStatus bound_backtracking_at_lower_bound;
870 ColChoiceAndStatus bound_backtracking_at_upper_bound;
871 };
872 void SwapDeletedAndModifiedVariableRestoreInfo(RestoreInfo* r);
873
874 std::vector<RestoreInfo> restore_stack_;
875 DenseColumn saved_row_lower_bounds_;
876 DenseColumn saved_row_upper_bounds_;
877
878 ColumnsSaver columns_saver_;
879 DenseRow saved_objective_;
880};
881
882// Because of numerical imprecision, a preprocessor like
883// DoubletonEqualityRowPreprocessor can transform a constraint/variable domain
884// like [1, 1+1e-7] to a fixed domain (for ex by multiplying the above domain by
885// 1e9). This causes an issue because at postsolve, a FIXED_VALUE status now
886// needs to be transformed to a AT_LOWER_BOUND/AT_UPPER_BOUND status. This is
887// what this function is doing for the constraint statuses only.
888//
889// TODO(user): A better solution would simply be to get rid of the FIXED status
890// altogether, it is better to simply use AT_LOWER_BOUND/AT_UPPER_BOUND
891// depending on the constraining bound in the optimal solution. Note that we can
892// always at the end transform any variable/constraint with a fixed domain to
893// FIXED_VALUE if needed to keep the same external API.
894void FixConstraintWithFixedStatuses(const DenseColumn& row_lower_bounds,
895 const DenseColumn& row_upper_bounds,
898// --------------------------------------------------------
899// DualizerPreprocessor
900// --------------------------------------------------------
901// DualizerPreprocessor may change the given program to its dual depending
902// on the value of the parameter solve_dual_problem.
903//
904// IMPORTANT: FreeConstraintPreprocessor() must be called first since this
905// preprocessor does not deal correctly with free constraints.
906class DualizerPreprocessor final : public Preprocessor {
907 public:
908 explicit DualizerPreprocessor(const GlopParameters* parameters)
909 : Preprocessor(parameters) {}
912 ~DualizerPreprocessor() final = default;
913 bool Run(LinearProgram* lp) final;
914 void RecoverSolution(ProblemSolution* solution) const final;
915 void UseInMipContext() final {
916 LOG(FATAL) << "In the presence of integer variables, "
917 << "there is no notion of a dual problem.";
918 }
919
920 // Convert the given problem status to the one of its dual.
921 ProblemStatus ChangeStatusToDualStatus(ProblemStatus status) const;
922
923 private:
924 DenseRow variable_lower_bounds_;
925 DenseRow variable_upper_bounds_;
926
927 RowIndex primal_num_rows_;
928 ColIndex primal_num_cols_;
929 bool primal_is_maximization_problem_;
930 RowToColMapping duplicated_rows_;
931
932 // For postsolving the variable/constraint statuses.
933 VariableStatusRow dual_status_correspondence_;
934 ColMapping slack_or_surplus_mapping_;
935};
936
937// --------------------------------------------------------
938// ShiftVariableBoundsPreprocessor
939// --------------------------------------------------------
940// For each variable, inspects its bounds and "shift" them if necessary, so that
941// its domain contains zero. A variable that was shifted will always have at
942// least one of its bounds to zero. Doing it all at once allows to have a better
943// precision when modifying the constraint bounds by using an accurate summation
944// algorithm.
945//
946// Example:
947// - A variable with bound [1e10, infinity] will be shifted to [0, infinity].
948// - A variable with domain [-1e10, 1e10] will not be shifted. Note that
949// compared to the first case, doing so here may introduce unnecessary
950// numerical errors if the variable value in the final solution is close to
951// zero.
952//
953// The expected impact of this is:
954// - Better behavior of the scaling.
955// - Better precision and numerical accuracy of the simplex method.
956// - Slightly improved speed (because adding a column with a variable value of
957// zero takes no work later).
958//
959// TODO(user): Having for each variable one of their bounds at zero is a
960// requirement for the DualizerPreprocessor and for the implied free column in
961// the ImpliedFreePreprocessor. However, shifting a variable with a domain like
962// [-1e10, 1e10] may introduce numerical issues. Relax the definition of
963// a free variable so that only having a domain containing 0.0 is enough?
965 public:
966 explicit ShiftVariableBoundsPreprocessor(const GlopParameters* parameters)
967 : Preprocessor(parameters) {}
969 delete;
971 const ShiftVariableBoundsPreprocessor&) = delete;
972 ~ShiftVariableBoundsPreprocessor() final = default;
973 bool Run(LinearProgram* lp) final;
974 void RecoverSolution(ProblemSolution* solution) const final;
975
976 const DenseRow& offsets() const { return offsets_; }
977
978 private:
979 // Contains for each variable by how much its bounds where shifted during
980 // presolve. Note that the shift was negative (new bound = initial bound -
981 // offset).
982 DenseRow offsets_;
983 // Contains the initial problem bounds. They are needed to get the perfect
984 // numerical accuracy for variables at their bound after postsolve.
985 DenseRow variable_initial_lbs_;
986 DenseRow variable_initial_ubs_;
987};
988
989// --------------------------------------------------------
990// ScalingPreprocessor
991// --------------------------------------------------------
992// Scales the SparseMatrix of the linear program using a SparseMatrixScaler.
993// This is only applied if the parameter use_scaling is true.
994class ScalingPreprocessor final : public Preprocessor {
995 public:
996 explicit ScalingPreprocessor(const GlopParameters* parameters)
997 : Preprocessor(parameters) {}
998 ScalingPreprocessor(const ScalingPreprocessor&) = delete;
999 ScalingPreprocessor& operator=(const ScalingPreprocessor&) = delete;
1000 ~ScalingPreprocessor() final = default;
1001 bool Run(LinearProgram* lp) final;
1003 void UseInMipContext() final { LOG(FATAL) << "Not implemented."; }
1005 private:
1006 DenseRow variable_lower_bounds_;
1007 DenseRow variable_upper_bounds_;
1008 Fractional cost_scaling_factor_;
1009 Fractional bound_scaling_factor_;
1011};
1012
1013// --------------------------------------------------------
1014// ToMinimizationPreprocessor
1015// --------------------------------------------------------
1016// Changes the problem from maximization to minimization (if applicable).
1017class ToMinimizationPreprocessor final : public Preprocessor {
1018 public:
1019 explicit ToMinimizationPreprocessor(const GlopParameters* parameters)
1020 : Preprocessor(parameters) {}
1023 delete;
1024 ~ToMinimizationPreprocessor() final = default;
1025 bool Run(LinearProgram* lp) final;
1026 void RecoverSolution(ProblemSolution* solution) const final;
1027};
1028
1029// --------------------------------------------------------
1030// AddSlackVariablesPreprocessor
1031// --------------------------------------------------------
1032// Transforms the linear program to the equation form
1033// min c.x, s.t. A.x = 0. This is done by:
1034// 1. Introducing slack variables for all constraints; all these variables are
1035// introduced with coefficient 1.0, and their bounds are set to be negative
1036// bounds of the corresponding constraint.
1037// 2. Changing the bounds of all constraints to (0, 0) to make them an equality.
1039// As a consequence, the matrix of the linear program always has full row rank
1040// after this preprocessor. Note that the slack variables are always added last,
1041// so that the rightmost square sub-matrix is always the identity matrix.
1042//
1043// TODO(user): Do not require this step to talk to the revised simplex. On large
1044// LPs like supportcase11.mps, this step alone can add 1.5 GB to the solver peak
1045// memory for no good reason. The internal matrix representation used in glop is
1046// a lot more efficient, and there is no point keeping the slacks in
1047// LinearProgram. It is also bad for incrementaly modifying the LP.
1048class AddSlackVariablesPreprocessor final : public Preprocessor {
1049 public:
1050 explicit AddSlackVariablesPreprocessor(const GlopParameters* parameters)
1051 : Preprocessor(parameters) {}
1054 const AddSlackVariablesPreprocessor&) = delete;
1055 ~AddSlackVariablesPreprocessor() final = default;
1056 bool Run(LinearProgram* lp) final;
1057 void RecoverSolution(ProblemSolution* solution) const final;
1058
1059 private:
1060 ColIndex first_slack_col_;
1062
1063} // namespace glop
1064} // namespace operations_research
1065
1066#endif // ORTOOLS_GLOP_PREPROCESSOR_H_
ColumnDeletionHelper & operator=(const ColumnDeletionHelper &)=delete
const DenseBooleanRow & GetMarkedColumns() const
void RestoreDeletedColumns(ProblemSolution *solution) const
void MarkColumnForDeletionWithState(ColIndex col, Fractional value, VariableStatus status)
const SparseColumn & SavedOrEmptyColumn(ColIndex col) const
void SaveColumn(ColIndex col, const SparseColumn &column)
const SparseColumn & SavedColumn(ColIndex col) const
void SaveColumnIfNotAlreadyDone(ColIndex col, const SparseColumn &column)
void RecoverSolution(ProblemSolution *solution) const final
DoubletonEqualityRowPreprocessor & operator=(const DoubletonEqualityRowPreprocessor &)=delete
DualizerPreprocessor & operator=(const DualizerPreprocessor &)=delete
void RecoverSolution(ProblemSolution *solution) const final
EmptyColumnPreprocessor & operator=(const EmptyColumnPreprocessor &)=delete
EmptyColumnPreprocessor(const GlopParameters *parameters)
void RecoverSolution(ProblemSolution *solution) const final
void RecoverSolution(ProblemSolution *solution) const final
EmptyConstraintPreprocessor & operator=(const EmptyConstraintPreprocessor &)=delete
FixedVariablePreprocessor & operator=(const FixedVariablePreprocessor &)=delete
void RecoverSolution(ProblemSolution *solution) const final
ForcingAndImpliedFreeConstraintPreprocessor & operator=(const ForcingAndImpliedFreeConstraintPreprocessor &)=delete
FreeConstraintPreprocessor & operator=(const FreeConstraintPreprocessor &)=delete
void RecoverSolution(ProblemSolution *solution) const final
void RecoverSolution(ProblemSolution *solution) const final
ImpliedFreePreprocessor & operator=(const ImpliedFreePreprocessor &)=delete
MainLpPreprocessor & operator=(const MainLpPreprocessor &)=delete
void DestructiveRecoverSolution(ProblemSolution *solution)
void RecoverSolution(ProblemSolution *solution) const override
MainLpPreprocessor(const GlopParameters *parameters)
bool IsSmallerWithinFeasibilityTolerance(Fractional a, Fractional b) const
virtual void RecoverSolution(ProblemSolution *solution) const =0
Preprocessor(const GlopParameters *parameters)
std::unique_ptr< TimeLimit > infinite_time_limit_
void SetTimeLimit(TimeLimit *time_limit)
bool IsSmallerWithinPreprocessorZeroTolerance(Fractional a, Fractional b) const
Preprocessor & operator=(const Preprocessor &)=delete
RowDeletionHelper & operator=(const RowDeletionHelper &)=delete
const DenseBooleanColumn & GetMarkedRows() const
void RestoreDeletedRows(ProblemSolution *solution) const
ScalingPreprocessor & operator=(const ScalingPreprocessor &)=delete
void RecoverSolution(ProblemSolution *solution) const final
ShiftVariableBoundsPreprocessor & operator=(const ShiftVariableBoundsPreprocessor &)=delete
void RecoverSolution(ProblemSolution *solution) const final
void RecoverSolution(ProblemSolution *solution) const final
SingletonPreprocessor & operator=(const SingletonPreprocessor &)=delete
void RecoverSolution(ProblemSolution *solution) const final
void Undo(const GlopParameters &parameters, const SparseColumn &saved_column, const SparseColumn &saved_row, ProblemSolution *solution) const
const MatrixEntry & Entry() const
UnconstrainedVariablePreprocessor & operator=(const UnconstrainedVariablePreprocessor &)=delete
void RecoverSolution(ProblemSolution *solution) const final
double Fractional
Definition lp_types.h:81
StrictITIVector< RowIndex, ColIndex > RowToColMapping
Definition lp_types.h:394
StrictITIVector< RowIndex, bool > DenseBooleanColumn
Definition lp_types.h:383
StrictITIVector< ColIndex, ColIndex > ColMapping
Definition lp_types.h:357
SumWithOneMissing< true > SumWithPositiveInfiniteAndOneMissing
Definition lp_utils.h:395
StrictITIVector< RowIndex, RowIndex > RowMapping
Definition lp_types.h:389
StrictITIVector< ColIndex, VariableStatus > VariableStatusRow
Definition lp_types.h:372
StrictITIVector< RowIndex, Fractional > DenseColumn
Definition lp_types.h:380
StrictITIVector< ColIndex, Fractional > DenseRow
Definition lp_types.h:351
StrictITIVector< ColIndex, bool > DenseBooleanRow
Definition lp_types.h:354
SumWithOneMissing< false > SumWithNegativeInfiniteAndOneMissing
Definition lp_utils.h:396
OR-Tools root namespace.
Select next search node to expand Select next item_i to add this new search node to the search Generate a new search node where item_i is not in the knapsack Check validity of this new partial solution(using propagators) - If valid
STL namespace.
MatrixEntry(RowIndex _row, ColIndex _col, Fractional _coeff)