1 // Copyright 2010-2021 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 /**
15 * \file
16 * A C++ wrapper that provides a simple and unified interface to
17 * several linear programming and mixed integer programming solvers:
18 * GLOP, GLPK, CLP, CBC, and SCIP. The wrapper can also be used in Java, C#,
19 * and Python via SWIG.
20 *
21 * What is Linear Programming?
22 *
23 * In mathematics, linear programming (LP) is a technique for optimization of
24 * a linear objective function, subject to linear equality and linear
25 * inequality constraints. Informally, linear programming determines the way
26 * to achieve the best outcome (such as maximum profit or lowest cost) in a
27 * given mathematical model and given some list of requirements represented
28 * as linear equations.
29 *
30 * The most widely used technique for solving a linear program is the Simplex
31 * algorithm, devised by George Dantzig in 1947. It performs very well on
32 * most instances, for which its running time is polynomial. A lot of effort
33 * has been put into improving the algorithm and its implementation. As a
34 * byproduct, it has however been shown that one can always construct
35 * problems that take exponential time for the Simplex algorithm to solve.
36 * Research has thus focused on trying to find a polynomial algorithm for
37 * linear programming, or to prove that linear programming is indeed
38 * polynomial.
39 *
40 * Leonid Khachiyan first exhibited in 1979 a weakly polynomial algorithm for
41 * linear programming. "Weakly polynomial" means that the running time of the
42 * algorithm is in O(P(n) * 2^p) where P(n) is a polynomial of the size of the
43 * problem, and p is the precision of computations expressed in number of
44 * bits. With a fixed-precision, floating-point-based implementation, a weakly
45 * polynomial algorithm will thus run in polynomial time. No implementation
46 * of Khachiyan's algorithm has proved efficient, but a larger breakthrough in
47 * the field came in 1984 when Narendra Karmarkar introduced a new interior
48 * point method for solving linear programming problems. Interior point
49 * algorithms have proved efficient on very large linear programs.
50 *
51 * Check Wikipedia for more detail:
52 * http://en.wikipedia.org/wiki/Linear_programming
53 *
54 * -----------------------------------
55 *
56 * Example of a Linear Program
57 *
58 * maximize:
59 * 3x + y
60 * subject to:
61 * 1.5 x + 2 y <= 12
62 * 0 <= x <= 3
63 * 0 <= y <= 5
64 *
65 * A linear program has:
66 * 1) a linear objective function
67 * 2) linear constraints that can be equalities or inequalities
68 * 3) bounds on variables that can be positive, negative, finite or
69 * infinite.
70 *
71 * -----------------------------------
72 *
73 * What is Mixed Integer Programming?
74 *
75 * Here, the constraints and the objective are still linear but
76 * there are additional integrality requirements for variables. If
77 * all variables are required to take integer values, then the
78 * problem is called an integer program (IP). In most cases, only
79 * some variables are required to be integer and the rest of the
80 * variables are continuous: this is called a mixed integer program
81 * (MIP). IPs and MIPs are generally NP-hard.
82 *
83 * Integer variables can be used to model discrete decisions (build a
84 * datacenter in city A or city B), logical relationships (only
85 * place machines in datacenter A if we have decided to build
86 * datacenter A) and approximate non-linear functions with piecewise
87 * linear functions (for example, the cost of machines as a function
88 * of how many machines are bought, or the latency of a server as a
89 * function of its load).
90 *
91 * -----------------------------------
92 *
93 * How to use the wrapper
94 *
95 * The user builds the model and solves it through the MPSolver class,
96 * then queries the solution through the MPSolver, MPVariable and
97 * MPConstraint classes. To be able to query a solution, you need the
98 * following:
99 * - A solution exists: MPSolver::Solve has been called and a solution
100 * has been found.
101 * - The model has not been modified since the last time
102 * MPSolver::Solve was called. Otherwise, the solution obtained
103 * before the model modification may not longer be feasible or
104 * optimal.
105 *
106 * @see ../examples/linear_programming.cc for a simple LP example.
107 *
108 * @see ../examples/integer_programming.cc for a simple MIP example.
109 *
110 * All methods cannot be called successfully in all cases. For
111 * example: you cannot query a solution when no solution exists, you
112 * cannot query a reduced cost value (which makes sense only on
113 * continuous problems) on a discrete problem. When a method is
114 * called in an unsuitable context, it aborts with a
115 * LOG(FATAL).
116 * TODO(user): handle failures gracefully.
117 *
118 * -----------------------------------
119 *
120 * For developers: How the wrapper works
121 *
122 * MPSolver stores a representation of the model (variables,
123 * constraints and objective) in its own data structures and a
124 * pointer to a MPSolverInterface that wraps the underlying solver
125 * (GLOP, CBC, CLP, GLPK, or SCIP) that does the actual work. The
126 * underlying solver also keeps a representation of the model in its
127 * own data structures. The model representations in MPSolver and in
128 * the underlying solver are kept in sync by the 'extraction'
129 * mechanism: synchronously for some changes and asynchronously
130 * (when MPSolver::Solve is called) for others. Synchronicity
131 * depends on the modification applied and on the underlying solver.
132 */
133
134 #ifndef OR_TOOLS_LINEAR_SOLVER_LINEAR_SOLVER_H_
135 #define OR_TOOLS_LINEAR_SOLVER_LINEAR_SOLVER_H_
136
137 #include <atomic>
138 #include <cstdint>
139 #include <functional>
140 #include <limits>
141 #include <map>
142 #include <memory>
143 #include <string>
144 #include <utility>
145 #include <vector>
146
147 #include "absl/base/port.h"
148 #include "absl/container/flat_hash_map.h"
149 #include "absl/flags/parse.h"
150 #include "absl/flags/usage.h"
151 #include "absl/status/status.h"
152 #include "absl/strings/match.h"
153 #include "absl/strings/str_format.h"
154 #include "absl/types/optional.h"
155 #include "ortools/base/integral_types.h"
156 #include "ortools/base/logging.h"
157 #include "ortools/base/macros.h"
158 #include "ortools/base/timer.h"
159 #include "ortools/linear_solver/linear_expr.h"
160 #include "ortools/linear_solver/linear_solver.pb.h"
161 #include "ortools/linear_solver/linear_solver_callback.h"
162 #include "ortools/port/proto_utils.h"
163
164 ABSL_DECLARE_FLAG(bool, linear_solver_enable_verbose_output);
165
166 namespace operations_research {
167
168 constexpr double kDefaultPrimalTolerance = 1e-07;
169
170 class MPConstraint;
171 class MPObjective;
172 class MPSolverInterface;
173 class MPSolverParameters;
174 class MPVariable;
175
176 // There is a homonymous version taking a MPSolver::OptimizationProblemType.
177 bool SolverTypeIsMip(MPModelRequest::SolverType solver_type);
178
179 /**
180 * This mathematical programming (MP) solver class is the main class
181 * though which users build and solve problems.
182 */
183 class MPSolver {
184 public:
185 /**
186 * The type of problems (LP or MIP) that will be solved and the underlying
187 * solver (GLOP, GLPK, CLP, CBC or SCIP) that will solve them. This must
188 * remain consistent with MPModelRequest::OptimizationProblemType
189 * (take particular care of the open-source version).
190 */
191 enum OptimizationProblemType {
192 // Linear programming problems.
193 // ----------------------------
194 CLP_LINEAR_PROGRAMMING = 0,
195 GLPK_LINEAR_PROGRAMMING = 1,
196 GLOP_LINEAR_PROGRAMMING = 2, // Recommended default value. Made in Google.
197
198 // Integer programming problems.
199 // -----------------------------
200 SCIP_MIXED_INTEGER_PROGRAMMING = 3, // Recommended default value.
201 GLPK_MIXED_INTEGER_PROGRAMMING = 4,
202 CBC_MIXED_INTEGER_PROGRAMMING = 5,
203
204 // Commercial software (need license).
205 GUROBI_LINEAR_PROGRAMMING = 6,
206 GUROBI_MIXED_INTEGER_PROGRAMMING = 7,
207 CPLEX_LINEAR_PROGRAMMING = 10,
208 CPLEX_MIXED_INTEGER_PROGRAMMING = 11,
209 XPRESS_LINEAR_PROGRAMMING = 101,
210 XPRESS_MIXED_INTEGER_PROGRAMMING = 102,
211
212 // Boolean optimization problem (requires only integer variables and works
213 // best with only Boolean variables).
214 BOP_INTEGER_PROGRAMMING = 12,
215
216 // SAT based solver (requires only integer and Boolean variables).
217 // If you pass it mixed integer problems, it will scale coefficients to
218 // integer values, and solver continuous variables as integral variables.
219 SAT_INTEGER_PROGRAMMING = 14,
220
221 // Dedicated knapsack solvers.
222 KNAPSACK_MIXED_INTEGER_PROGRAMMING = 13,
223 };
224
225 /// Create a solver with the given name and underlying solver backend.
226 MPSolver(const std::string& name, OptimizationProblemType problem_type);
227 virtual ~MPSolver();
228
229 /**
230 * Recommended factory method to create a MPSolver instance, especially in
231 * non C++ languages.
232 *
233 * It returns a newly created solver instance if successful, or a nullptr
234 * otherwise. This can occur if the relevant interface is not linked in, or if
235 * a needed license is not accessible for commercial solvers.
236 *
237 * Ownership of the solver is passed on to the caller of this method.
238 * It will accept both string names of the OptimizationProblemType enum, as
239 * well as a short version (i.e. "SCIP_MIXED_INTEGER_PROGRAMMING" or "SCIP").
240 *
241 * solver_id is case insensitive, and the following names are supported:
242 * - CLP_LINEAR_PROGRAMMING or CLP
243 * - CBC_MIXED_INTEGER_PROGRAMMING or CBC
244 * - GLOP_LINEAR_PROGRAMMING or GLOP
245 * - BOP_INTEGER_PROGRAMMING or BOP
246 * - SAT_INTEGER_PROGRAMMING or SAT or CP_SAT
247 * - SCIP_MIXED_INTEGER_PROGRAMMING or SCIP
248 * - GUROBI_LINEAR_PROGRAMMING or GUROBI_LP
249 * - GUROBI_MIXED_INTEGER_PROGRAMMING or GUROBI or GUROBI_MIP
250 * - CPLEX_LINEAR_PROGRAMMING or CPLEX_LP
251 * - CPLEX_MIXED_INTEGER_PROGRAMMING or CPLEX or CPLEX_MIP
252 * - XPRESS_LINEAR_PROGRAMMING or XPRESS_LP
253 * - XPRESS_MIXED_INTEGER_PROGRAMMING or XPRESS or XPRESS_MIP
254 * - GLPK_LINEAR_PROGRAMMING or GLPK_LP
255 * - GLPK_MIXED_INTEGER_PROGRAMMING or GLPK or GLPK_MIP
256 */
257 static MPSolver* CreateSolver(const std::string& solver_id);
258
259 /**
260 * Whether the given problem type is supported (this will depend on the
261 * targets that you linked).
262 */
263 static bool SupportsProblemType(OptimizationProblemType problem_type);
264
265 /**
266 * Parses the name of the solver. Returns true if the solver type is
267 * successfully parsed as one of the OptimizationProblemType.
268 * See the documentation of CreateSolver() for the list of supported names.
269 */
270 static bool ParseSolverType(absl::string_view solver_id,
271 OptimizationProblemType* type);
272
273 /**
274 * Parses the name of the solver and returns the correct optimization type or
275 * dies. Invariant: ParseSolverTypeOrDie(ToString(type)) = type.
276 */
277 static OptimizationProblemType ParseSolverTypeOrDie(
278 const std::string& solver_id);
279
280 bool IsMIP() const;
281
282 /// Returns the name of the model set at construction.
Name()283 const std::string& Name() const {
284 return name_; // Set at construction.
285 }
286
287 /// Returns the optimization problem type set at construction.
ProblemType()288 virtual OptimizationProblemType ProblemType() const {
289 return problem_type_; // Set at construction.
290 }
291
292 /**
293 * Clears the objective (including the optimization direction), all variables
294 * and constraints. All the other properties of the MPSolver (like the time
295 * limit) are kept untouched.
296 */
297 void Clear();
298
299 /// Returns the number of variables.
NumVariables()300 int NumVariables() const { return variables_.size(); }
301
302 /**
303 * Returns the array of variables handled by the MPSolver. (They are listed in
304 * the order in which they were created.)
305 */
variables()306 const std::vector<MPVariable*>& variables() const { return variables_; }
307
308 /**
309 * Returns the variable at position index.
310 */
variable(int index)311 MPVariable* variable(int index) const { return variables_[index]; }
312
313 /**
314 * Looks up a variable by name, and returns nullptr if it does not exist. The
315 * first call has a O(n) complexity, as the variable name index is lazily
316 * created upon first use. Will crash if variable names are not unique.
317 */
318 MPVariable* LookupVariableOrNull(const std::string& var_name) const;
319
320 /**
321 * Creates a variable with the given bounds, integrality requirement and
322 * name. Bounds can be finite or +/- MPSolver::infinity(). The MPSolver owns
323 * the variable (i.e. the returned pointer is borrowed). Variable names are
324 * optional. If you give an empty name, name() will auto-generate one for you
325 * upon request.
326 */
327 MPVariable* MakeVar(double lb, double ub, bool integer,
328 const std::string& name);
329
330 /// Creates a continuous variable.
331 MPVariable* MakeNumVar(double lb, double ub, const std::string& name);
332
333 /// Creates an integer variable.
334 MPVariable* MakeIntVar(double lb, double ub, const std::string& name);
335
336 /// Creates a boolean variable.
337 MPVariable* MakeBoolVar(const std::string& name);
338
339 /**
340 * Creates an array of variables. All variables created have the same bounds
341 * and integrality requirement. If nb <= 0, no variables are created, the
342 * function crashes in non-opt mode.
343 *
344 * @param nb the number of variables to create.
345 * @param lb the lower bound of created variables
346 * @param ub the upper bound of created variables
347 * @param integer controls whether the created variables are continuous or
348 * integral.
349 * @param name_prefix the prefix of the variable names. Variables are named
350 * name_prefix0, name_prefix1, ...
351 * @param[out] vars the vector of variables to fill with variables.
352 */
353 void MakeVarArray(int nb, double lb, double ub, bool integer,
354 const std::string& name_prefix,
355 std::vector<MPVariable*>* vars);
356
357 /// Creates an array of continuous variables.
358 void MakeNumVarArray(int nb, double lb, double ub, const std::string& name,
359 std::vector<MPVariable*>* vars);
360
361 /// Creates an array of integer variables.
362 void MakeIntVarArray(int nb, double lb, double ub, const std::string& name,
363 std::vector<MPVariable*>* vars);
364
365 /// Creates an array of boolean variables.
366 void MakeBoolVarArray(int nb, const std::string& name,
367 std::vector<MPVariable*>* vars);
368
369 /// Returns the number of constraints.
NumConstraints()370 int NumConstraints() const { return constraints_.size(); }
371
372 /**
373 * Returns the array of constraints handled by the MPSolver.
374 *
375 * They are listed in the order in which they were created.
376 */
constraints()377 const std::vector<MPConstraint*>& constraints() const { return constraints_; }
378
379 /** Returns the constraint at the given index. */
constraint(int index)380 MPConstraint* constraint(int index) const { return constraints_[index]; }
381
382 /**
383 * Looks up a constraint by name, and returns nullptr if it does not exist.
384 *
385 * The first call has a O(n) complexity, as the constraint name index is
386 * lazily created upon first use. Will crash if constraint names are not
387 * unique.
388 */
389 MPConstraint* LookupConstraintOrNull(
390 const std::string& constraint_name) const;
391
392 /**
393 * Creates a linear constraint with given bounds.
394 *
395 * Bounds can be finite or +/- MPSolver::infinity(). The MPSolver class
396 * assumes ownership of the constraint.
397 *
398 * @return a pointer to the newly created constraint.
399 */
400 MPConstraint* MakeRowConstraint(double lb, double ub);
401
402 /// Creates a constraint with -infinity and +infinity bounds.
403 MPConstraint* MakeRowConstraint();
404
405 /// Creates a named constraint with given bounds.
406 MPConstraint* MakeRowConstraint(double lb, double ub,
407 const std::string& name);
408
409 /// Creates a named constraint with -infinity and +infinity bounds.
410 MPConstraint* MakeRowConstraint(const std::string& name);
411
412 /**
413 * Creates a constraint owned by MPSolver enforcing:
414 * range.lower_bound() <= range.linear_expr() <= range.upper_bound()
415 */
416 MPConstraint* MakeRowConstraint(const LinearRange& range);
417
418 /// As above, but also names the constraint.
419 MPConstraint* MakeRowConstraint(const LinearRange& range,
420 const std::string& name);
421
422 /**
423 * Returns the objective object.
424 *
425 * Note that the objective is owned by the solver, and is initialized to its
426 * default value (see the MPObjective class below) at construction.
427 */
Objective()428 const MPObjective& Objective() const { return *objective_; }
429
430 /// Returns the mutable objective object.
MutableObjective()431 MPObjective* MutableObjective() { return objective_.get(); }
432
433 /**
434 * The status of solving the problem. The straightforward translation to
435 * homonymous enum values of MPSolverResponseStatus (see
436 * ./linear_solver.proto) is guaranteed by ./enum_consistency_test.cc, you may
437 * rely on it.
438 */
439 enum ResultStatus {
440 /// optimal.
441 OPTIMAL,
442 /// feasible, or stopped by limit.
443 FEASIBLE,
444 /// proven infeasible.
445 INFEASIBLE,
446 /// proven unbounded.
447 UNBOUNDED,
448 /// abnormal, i.e., error of some kind.
449 ABNORMAL,
450 /// the model is trivially invalid (NaN coefficients, etc).
451 MODEL_INVALID,
452 /// not been solved yet.
453 NOT_SOLVED = 6
454 };
455
456 /// Solves the problem using the default parameter values.
457 ResultStatus Solve();
458
459 /// Solves the problem using the specified parameter values.
460 ResultStatus Solve(const MPSolverParameters& param);
461
462 /**
463 * Writes the model using the solver internal write function. Currently only
464 * available for Gurobi.
465 */
466 void Write(const std::string& file_name);
467
468 /**
469 * Advanced usage: compute the "activities" of all constraints, which are the
470 * sums of their linear terms. The activities are returned in the same order
471 * as constraints(), which is the order in which constraints were added; but
472 * you can also use MPConstraint::index() to get a constraint's index.
473 */
474 std::vector<double> ComputeConstraintActivities() const;
475
476 /**
477 * Advanced usage: Verifies the *correctness* of the solution.
478 *
479 * It verifies that all variables must be within their domains, all
480 * constraints must be satisfied, and the reported objective value must be
481 * accurate.
482 *
483 * Usage:
484 * - This can only be called after Solve() was called.
485 * - "tolerance" is interpreted as an absolute error threshold.
486 * - For the objective value only, if the absolute error is too large,
487 * the tolerance is interpreted as a relative error threshold instead.
488 * - If "log_errors" is true, every single violation will be logged.
489 * - If "tolerance" is negative, it will be set to infinity().
490 *
491 * Most users should just set the --verify_solution flag and not bother using
492 * this method directly.
493 */
494 bool VerifySolution(double tolerance, bool log_errors) const;
495
496 /**
497 * Advanced usage: resets extracted model to solve from scratch.
498 *
499 * This won't reset the parameters that were set with
500 * SetSolverSpecificParametersAsString() or set_time_limit() or even clear the
501 * linear program. It will just make sure that next Solve() will be as if
502 * everything was reconstructed from scratch.
503 */
504 void Reset();
505
506 /** Interrupts the Solve() execution to terminate processing if possible.
507 *
508 * If the underlying interface supports interruption; it does that and returns
509 * true regardless of whether there's an ongoing Solve() or not. The Solve()
510 * call may still linger for a while depending on the conditions. If
511 * interruption is not supported; returns false and does nothing.
512 * MPSolver::SolverTypeSupportsInterruption can be used to check if
513 * interruption is supported for a given solver type.
514 */
515 bool InterruptSolve();
516
517 /**
518 * Loads model from protocol buffer.
519 *
520 * Returns MPSOLVER_MODEL_IS_VALID if the model is valid, and another status
521 * otherwise (currently only MPSOLVER_MODEL_INVALID and MPSOLVER_INFEASIBLE).
522 * If the model isn't valid, populates "error_message".
523 */
524 MPSolverResponseStatus LoadModelFromProto(const MPModelProto& input_model,
525 std::string* error_message);
526 /**
527 * Loads model from protocol buffer.
528 *
529 * The same as above, except that the loading keeps original variable and
530 * constraint names. Caller should make sure that all variable names and
531 * constraint names are unique, respectively.
532 */
533 MPSolverResponseStatus LoadModelFromProtoWithUniqueNamesOrDie(
534 const MPModelProto& input_model, std::string* error_message);
535
536 /// Encodes the current solution in a solution response protocol buffer.
537 void FillSolutionResponseProto(MPSolutionResponse* response) const;
538
539 /**
540 * Solves the model encoded by a MPModelRequest protocol buffer and fills the
541 * solution encoded as a MPSolutionResponse. The solve is stopped prematurely
542 * if interrupt is non-null at set to true during (or before) solving.
543 * Interruption is only supported if SolverTypeSupportsInterruption() returns
544 * true for the requested solver. Passing a non-null interruption with any
545 * other solver type immediately returns an MPSOLVER_INCOMPATIBLE_OPTIONS
546 * error.
547 *
548 * Note(user): This attempts to first use `DirectlySolveProto()` (if
549 * implemented). Consequently, this most likely does *not* override any of
550 * the default parameters of the underlying solver. This behavior *differs*
551 * from `MPSolver::Solve()` which by default sets the feasibility tolerance
552 * and the gap limit (as of 2020/02/11, to 1e-7 and 0.0001, respectively).
553 */
554 static void SolveWithProto(const MPModelRequest& model_request,
555 MPSolutionResponse* response,
556 // `interrupt` is non-const because the internal
557 // solver may set it to true itself, in some cases.
558 std::atomic<bool>* interrupt = nullptr);
559
SolverTypeSupportsInterruption(const MPModelRequest::SolverType solver)560 static bool SolverTypeSupportsInterruption(
561 const MPModelRequest::SolverType solver) {
562 // Interruption requires that MPSolver::InterruptSolve is supported for the
563 // underlying solver. Interrupting requests using SCIP is also not supported
564 // as of 2021/08/23, since InterruptSolve is not go/thread-safe
565 // for SCIP (see e.g. cl/350545631 for details).
566 return solver == MPModelRequest::GLOP_LINEAR_PROGRAMMING ||
567 solver == MPModelRequest::GUROBI_LINEAR_PROGRAMMING ||
568 solver == MPModelRequest::GUROBI_MIXED_INTEGER_PROGRAMMING ||
569 solver == MPModelRequest::SAT_INTEGER_PROGRAMMING;
570 }
571
572 /// Exports model to protocol buffer.
573 void ExportModelToProto(MPModelProto* output_model) const;
574
575 /**
576 * Load a solution encoded in a protocol buffer onto this solver for easy
577 access via the MPSolver interface.
578 *
579 * IMPORTANT: This may only be used in conjunction with ExportModel(),
580 following this example:
581 *
582 \code
583 MPSolver my_solver;
584 ... add variables and constraints ...
585 MPModelProto model_proto;
586 my_solver.ExportModelToProto(&model_proto);
587 MPSolutionResponse solver_response;
588 MPSolver::SolveWithProto(model_proto, &solver_response);
589 if (solver_response.result_status() == MPSolutionResponse::OPTIMAL) {
590 CHECK_OK(my_solver.LoadSolutionFromProto(solver_response));
591 ... inspect the solution using the usual API: solution_value(), etc...
592 }
593 \endcode
594 *
595 * The response must be in OPTIMAL or FEASIBLE status.
596 *
597 * Returns a non-OK status if a problem arised (typically, if it wasn't used
598 * like it should be):
599 * - loading a solution whose variables don't correspond to the solver's
600 * current variables
601 * - loading a dual solution whose constraints don't correspond to the
602 * solver's current constraints
603 * - loading a solution with a status other than OPTIMAL / FEASIBLE.
604 *
605 * Note: the objective value isn't checked. You can use VerifySolution() for
606 * that.
607 */
608 absl::Status LoadSolutionFromProto(
609 const MPSolutionResponse& response,
610 double tolerance = std::numeric_limits<double>::infinity());
611
612 /**
613 * Resets values of out of bound variables to the corresponding bound and
614 * returns an error if any of the variables have NaN value.
615 */
616 absl::Status ClampSolutionWithinBounds();
617
618 /**
619 * Shortcuts to the homonymous MPModelProtoExporter methods, via exporting to
620 * a MPModelProto with ExportModelToProto() (see above).
621 *
622 * Produces empty std::string on portable platforms (e.g. android, ios).
623 */
624 bool ExportModelAsLpFormat(bool obfuscate, std::string* model_str) const;
625 bool ExportModelAsMpsFormat(bool fixed_format, bool obfuscate,
626 std::string* model_str) const;
627
628 /**
629 * Sets the number of threads to use by the underlying solver.
630 *
631 * Returns OkStatus if the operation was successful. num_threads must be equal
632 * to or greater than 1. Note that the behaviour of this call depends on the
633 * underlying solver. E.g., it may set the exact number of threads or the max
634 * number of threads (check the solver's interface implementation for
635 * details). Also, some solvers may not (yet) support this function, but still
636 * enable multi-threading via SetSolverSpecificParametersAsString().
637 */
638 absl::Status SetNumThreads(int num_threads);
639
640 /// Returns the number of threads to be used during solve.
GetNumThreads()641 int GetNumThreads() const { return num_threads_; }
642
643 /**
644 * Advanced usage: pass solver specific parameters in text format.
645 *
646 * The format is solver-specific and is the same as the corresponding solver
647 * configuration file format. Returns true if the operation was successful.
648 */
649 bool SetSolverSpecificParametersAsString(const std::string& parameters);
GetSolverSpecificParametersAsString()650 std::string GetSolverSpecificParametersAsString() const {
651 return solver_specific_parameter_string_;
652 }
653
654 /**
655 * Sets a hint for solution.
656 *
657 * If a feasible or almost-feasible solution to the problem is already known,
658 * it may be helpful to pass it to the solver so that it can be used. A solver
659 * that supports this feature will try to use this information to create its
660 * initial feasible solution.
661 *
662 * Note: It may not always be faster to give a hint like this to the
663 * solver. There is also no guarantee that the solver will use this hint or
664 * try to return a solution "close" to this assignment in case of multiple
665 * optimal solutions.
666 */
667 void SetHint(std::vector<std::pair<const MPVariable*, double> > hint);
668
669 /**
670 * Advanced usage: possible basis status values for a variable and the slack
671 * variable of a linear constraint.
672 */
673 enum BasisStatus {
674 FREE = 0,
675 AT_LOWER_BOUND,
676 AT_UPPER_BOUND,
677 FIXED_VALUE,
678 BASIC
679 };
680
681 /**
682 * Advanced usage: Incrementality.
683 *
684 * This function takes a starting basis to be used in the next LP Solve()
685 * call. The statuses of a current solution can be retrieved via the
686 * basis_status() function of a MPVariable or a MPConstraint.
687 *
688 * WARNING: With Glop, you should disable presolve when using this because
689 * this information will not be modified in sync with the presolve and will
690 * likely not mean much on the presolved problem.
691 */
692 void SetStartingLpBasis(
693 const std::vector<MPSolver::BasisStatus>& variable_statuses,
694 const std::vector<MPSolver::BasisStatus>& constraint_statuses);
695
696 /**
697 * Infinity.
698 *
699 * You can use -MPSolver::infinity() for negative infinity.
700 */
infinity()701 static double infinity() { return std::numeric_limits<double>::infinity(); }
702
703 /**
704 * Controls (or queries) the amount of output produced by the underlying
705 * solver. The output can surface to LOGs, or to stdout or stderr, depending
706 * on the implementation. The amount of output will greatly vary with each
707 * implementation and each problem.
708 *
709 * Output is suppressed by default.
710 */
711 bool OutputIsEnabled() const;
712
713 /// Enables solver logging.
714 void EnableOutput();
715
716 /// Suppresses solver logging.
717 void SuppressOutput();
718
TimeLimit()719 absl::Duration TimeLimit() const { return time_limit_; }
SetTimeLimit(absl::Duration time_limit)720 void SetTimeLimit(absl::Duration time_limit) {
721 DCHECK_GE(time_limit, absl::ZeroDuration());
722 time_limit_ = time_limit;
723 }
724
DurationSinceConstruction()725 absl::Duration DurationSinceConstruction() const {
726 return absl::Now() - construction_time_;
727 }
728
729 /// Returns the number of simplex iterations.
730 int64_t iterations() const;
731
732 /**
733 * Returns the number of branch-and-bound nodes evaluated during the solve.
734 *
735 * Only available for discrete problems.
736 */
737 int64_t nodes() const;
738
739 /// Returns a string describing the underlying solver and its version.
740 std::string SolverVersion() const;
741
742 /**
743 * Advanced usage: returns the underlying solver.
744 *
745 * Returns the underlying solver so that the user can use solver-specific
746 * features or features that are not exposed in the simple API of MPSolver.
747 * This method is for advanced users, use at your own risk! In particular, if
748 * you modify the model or the solution by accessing the underlying solver
749 * directly, then the underlying solver will be out of sync with the
750 * information kept in the wrapper (MPSolver, MPVariable, MPConstraint,
751 * MPObjective). You need to cast the void* returned back to its original type
752 * that depends on the interface (CBC: OsiClpSolverInterface*, CLP:
753 * ClpSimplex*, GLPK: glp_prob*, SCIP: SCIP*).
754 */
755 void* underlying_solver();
756
757 /** Advanced usage: computes the exact condition number of the current scaled
758 * basis: L1norm(B) * L1norm(inverse(B)), where B is the scaled basis.
759 *
760 * This method requires that a basis exists: it should be called after Solve.
761 * It is only available for continuous problems. It is implemented for GLPK
762 * but not CLP because CLP does not provide the API for doing it.
763 *
764 * The condition number measures how well the constraint matrix is conditioned
765 * and can be used to predict whether numerical issues will arise during the
766 * solve: the model is declared infeasible whereas it is feasible (or
767 * vice-versa), the solution obtained is not optimal or violates some
768 * constraints, the resolution is slow because of repeated singularities.
769 *
770 * The rule of thumb to interpret the condition number kappa is:
771 * - o kappa <= 1e7: virtually no chance of numerical issues
772 * - o 1e7 < kappa <= 1e10: small chance of numerical issues
773 * - o 1e10 < kappa <= 1e13: medium chance of numerical issues
774 * - o kappa > 1e13: high chance of numerical issues
775 *
776 * The computation of the condition number depends on the quality of the LU
777 * decomposition, so it is not very accurate when the matrix is ill
778 * conditioned.
779 */
780 double ComputeExactConditionNumber() const;
781
782 /**
783 * Some solvers (MIP only, not LP) can produce multiple solutions to the
784 * problem. Returns true when another solution is available, and updates the
785 * MPVariable* objects to make the new solution queryable. Call only after
786 * calling solve.
787 *
788 * The optimality properties of the additional solutions found, and whether or
789 * not the solver computes them ahead of time or when NextSolution() is called
790 * is solver specific.
791 *
792 * As of 2020-02-10, only Gurobi and SCIP support NextSolution(), see
793 * linear_solver_interfaces_test for an example of how to configure these
794 * solvers for multiple solutions. Other solvers return false unconditionally.
795 */
796 ABSL_MUST_USE_RESULT bool NextSolution();
797
798 // Does not take ownership of "mp_callback".
799 //
800 // As of 2019-10-22, only SCIP and Gurobi support Callbacks.
801 // SCIP does not support suggesting a heuristic solution in the callback.
802 //
803 // See go/mpsolver-callbacks for additional documentation.
804 void SetCallback(MPCallback* mp_callback);
805 bool SupportsCallbacks() const;
806
807 // Global counters of variables and constraints ever created across all
808 // MPSolver instances. Those are only updated after the destruction
809 // (or Clear()) of each MPSolver instance.
810 static int64_t global_num_variables();
811 static int64_t global_num_constraints();
812
813 // DEPRECATED: Use TimeLimit() and SetTimeLimit(absl::Duration) instead.
814 // NOTE: These deprecated functions used the convention time_limit = 0 to mean
815 // "no limit", which now corresponds to time_limit_ = InfiniteDuration().
time_limit()816 int64_t time_limit() const {
817 return time_limit_ == absl::InfiniteDuration()
818 ? 0
819 : absl::ToInt64Milliseconds(time_limit_);
820 }
set_time_limit(int64_t time_limit_milliseconds)821 void set_time_limit(int64_t time_limit_milliseconds) {
822 SetTimeLimit(time_limit_milliseconds == 0
823 ? absl::InfiniteDuration()
824 : absl::Milliseconds(time_limit_milliseconds));
825 }
time_limit_in_secs()826 double time_limit_in_secs() const {
827 return static_cast<double>(time_limit()) / 1000.0;
828 }
829
830 // DEPRECATED: Use DurationSinceConstruction() instead.
wall_time()831 int64_t wall_time() const {
832 return absl::ToInt64Milliseconds(DurationSinceConstruction());
833 }
834
835 friend class GLPKInterface;
836 friend class CLPInterface;
837 friend class CBCInterface;
838 friend class SCIPInterface;
839 friend class GurobiInterface;
840 friend class CplexInterface;
841 friend class XpressInterface;
842 friend class SLMInterface;
843 friend class MPSolverInterface;
844 friend class GLOPInterface;
845 friend class BopInterface;
846 friend class SatInterface;
847 friend class KnapsackInterface;
848
849 // Debugging: verify that the given MPVariable* belongs to this solver.
850 bool OwnsVariable(const MPVariable* var) const;
851
852 private:
853 // Computes the size of the constraint with the largest number of
854 // coefficients with index in [min_constraint_index,
855 // max_constraint_index)
856 int ComputeMaxConstraintSize(int min_constraint_index,
857 int max_constraint_index) const;
858
859 // Returns true if the model has constraints with lower bound > upper bound.
860 bool HasInfeasibleConstraints() const;
861
862 // Returns true if the model has at least 1 integer variable.
863 bool HasIntegerVariables() const;
864
865 // Generates the map from variable names to their indices.
866 void GenerateVariableNameIndex() const;
867
868 // Generates the map from constraint names to their indices.
869 void GenerateConstraintNameIndex() const;
870
871 // The name of the linear programming problem.
872 const std::string name_;
873
874 // The type of the linear programming problem.
875 const OptimizationProblemType problem_type_;
876
877 // The solver interface.
878 std::unique_ptr<MPSolverInterface> interface_;
879
880 // The vector of variables in the problem.
881 std::vector<MPVariable*> variables_;
882 // A map from a variable's name to its index in variables_.
883 mutable absl::optional<absl::flat_hash_map<std::string, int> >
884 variable_name_to_index_;
885 // Whether variables have been extracted to the underlying interface.
886 std::vector<bool> variable_is_extracted_;
887
888 // The vector of constraints in the problem.
889 std::vector<MPConstraint*> constraints_;
890 // A map from a constraint's name to its index in constraints_.
891 mutable absl::optional<absl::flat_hash_map<std::string, int> >
892 constraint_name_to_index_;
893 // Whether constraints have been extracted to the underlying interface.
894 std::vector<bool> constraint_is_extracted_;
895
896 // The linear objective function.
897 std::unique_ptr<MPObjective> objective_;
898
899 // Initial values for all or some of the problem variables that can be
900 // exploited as a starting hint by a solver.
901 //
902 // Note(user): as of 05/05/2015, we can't use >> because of some SWIG errors.
903 //
904 // TODO(user): replace by two vectors, a std::vector<bool> to indicate if a
905 // hint is provided and a std::vector<double> for the hint value.
906 std::vector<std::pair<const MPVariable*, double> > solution_hint_;
907
908 absl::Duration time_limit_ = absl::InfiniteDuration(); // Default = No limit.
909
910 const absl::Time construction_time_;
911
912 // Permanent storage for the number of threads.
913 int num_threads_ = 1;
914
915 // Permanent storage for SetSolverSpecificParametersAsString().
916 std::string solver_specific_parameter_string_;
917
918 static absl::Mutex global_count_mutex_;
919 #ifndef SWIG
920 static int64_t global_num_variables_ ABSL_GUARDED_BY(global_count_mutex_);
921 static int64_t global_num_constraints_ ABSL_GUARDED_BY(global_count_mutex_);
922 #endif
923
924 MPSolverResponseStatus LoadModelFromProtoInternal(
925 const MPModelProto& input_model, bool clear_names,
926 bool check_model_validity, std::string* error_message);
927
928 DISALLOW_COPY_AND_ASSIGN(MPSolver);
929 };
930
SolverTypeIsMip(MPSolver::OptimizationProblemType solver_type)931 inline bool SolverTypeIsMip(MPSolver::OptimizationProblemType solver_type) {
932 return SolverTypeIsMip(static_cast<MPModelRequest::SolverType>(solver_type));
933 }
934
935 const absl::string_view ToString(
936 MPSolver::OptimizationProblemType optimization_problem_type);
937
938 inline std::ostream& operator<<(
939 std::ostream& os,
940 MPSolver::OptimizationProblemType optimization_problem_type) {
941 return os << ToString(optimization_problem_type);
942 }
943
944 inline std::ostream& operator<<(std::ostream& os,
945 MPSolver::ResultStatus status) {
946 return os << ProtoEnumToString<MPSolverResponseStatus>(
947 static_cast<MPSolverResponseStatus>(status));
948 }
949
950 bool AbslParseFlag(absl::string_view text,
951 MPSolver::OptimizationProblemType* solver_type,
952 std::string* error);
953
AbslUnparseFlag(MPSolver::OptimizationProblemType solver_type)954 inline std::string AbslUnparseFlag(
955 MPSolver::OptimizationProblemType solver_type) {
956 return std::string(ToString(solver_type));
957 }
958
959 /// A class to express a linear objective.
960 class MPObjective {
961 public:
962 /**
963 * Clears the offset, all variables and coefficients, and the optimization
964 * direction.
965 */
966 void Clear();
967
968 /**
969 * Sets the coefficient of the variable in the objective.
970 *
971 * If the variable does not belong to the solver, the function just returns,
972 * or crashes in non-opt mode.
973 */
974 void SetCoefficient(const MPVariable* const var, double coeff);
975
976 /**
977 * Gets the coefficient of a given variable in the objective
978 *
979 * It returns 0 if the variable does not appear in the objective).
980 */
981 double GetCoefficient(const MPVariable* const var) const;
982
983 /**
984 * Returns a map from variables to their coefficients in the objective.
985 *
986 * If a variable is not present in the map, then its coefficient is zero.
987 */
terms()988 const absl::flat_hash_map<const MPVariable*, double>& terms() const {
989 return coefficients_;
990 }
991
992 /// Sets the constant term in the objective.
993 void SetOffset(double value);
994
995 /// Gets the constant term in the objective.
offset()996 double offset() const { return offset_; }
997
998 /**
999 * Resets the current objective to take the value of linear_expr, and sets the
1000 * objective direction to maximize if "is_maximize", otherwise minimizes.
1001 */
1002 void OptimizeLinearExpr(const LinearExpr& linear_expr, bool is_maximization);
1003
1004 /// Resets the current objective to maximize linear_expr.
MaximizeLinearExpr(const LinearExpr & linear_expr)1005 void MaximizeLinearExpr(const LinearExpr& linear_expr) {
1006 OptimizeLinearExpr(linear_expr, true);
1007 }
1008 /// Resets the current objective to minimize linear_expr.
MinimizeLinearExpr(const LinearExpr & linear_expr)1009 void MinimizeLinearExpr(const LinearExpr& linear_expr) {
1010 OptimizeLinearExpr(linear_expr, false);
1011 }
1012
1013 /// Adds linear_expr to the current objective, does not change the direction.
1014 void AddLinearExpr(const LinearExpr& linear_expr);
1015
1016 /// Sets the optimization direction (maximize: true or minimize: false).
1017 void SetOptimizationDirection(bool maximize);
1018
1019 /// Sets the optimization direction to minimize.
SetMinimization()1020 void SetMinimization() { SetOptimizationDirection(false); }
1021
1022 /// Sets the optimization direction to maximize.
SetMaximization()1023 void SetMaximization() { SetOptimizationDirection(true); }
1024
1025 /// Is the optimization direction set to maximize?
1026 bool maximization() const;
1027
1028 /// Is the optimization direction set to minimize?
1029 bool minimization() const;
1030
1031 /**
1032 * Returns the objective value of the best solution found so far.
1033 *
1034 * It is the optimal objective value if the problem has been solved to
1035 * optimality.
1036 *
1037 * Note: the objective value may be slightly different than what you could
1038 * compute yourself using \c MPVariable::solution_value(); please use the
1039 * --verify_solution flag to gain confidence about the numerical stability of
1040 * your solution.
1041 */
1042 double Value() const;
1043
1044 /**
1045 * Returns the best objective bound.
1046 *
1047 * In case of minimization, it is a lower bound on the objective value of the
1048 * optimal integer solution. Only available for discrete problems.
1049 */
1050 double BestBound() const;
1051
1052 private:
1053 friend class MPSolver;
1054 friend class MPSolverInterface;
1055 friend class CBCInterface;
1056 friend class CLPInterface;
1057 friend class GLPKInterface;
1058 friend class SCIPInterface;
1059 friend class SLMInterface;
1060 friend class GurobiInterface;
1061 friend class CplexInterface;
1062 friend class XpressInterface;
1063 friend class GLOPInterface;
1064 friend class BopInterface;
1065 friend class SatInterface;
1066 friend class KnapsackInterface;
1067
1068 // Constructor. An objective points to a single MPSolverInterface
1069 // that is specified in the constructor. An objective cannot belong
1070 // to several models.
1071 // At construction, an MPObjective has no terms (which is equivalent
1072 // on having a coefficient of 0 for all variables), and an offset of 0.
MPObjective(MPSolverInterface * const interface_in)1073 explicit MPObjective(MPSolverInterface* const interface_in)
1074 : interface_(interface_in), coefficients_(1), offset_(0.0) {}
1075
1076 MPSolverInterface* const interface_;
1077
1078 // Mapping var -> coefficient.
1079 absl::flat_hash_map<const MPVariable*, double> coefficients_;
1080 // Constant term.
1081 double offset_;
1082
1083 DISALLOW_COPY_AND_ASSIGN(MPObjective);
1084 };
1085
1086 /// The class for variables of a Mathematical Programming (MP) model.
1087 class MPVariable {
1088 public:
1089 /// Returns the name of the variable.
name()1090 const std::string& name() const { return name_; }
1091
1092 /// Sets the integrality requirement of the variable.
1093 void SetInteger(bool integer);
1094
1095 /// Returns the integrality requirement of the variable.
integer()1096 bool integer() const { return integer_; }
1097
1098 /**
1099 * Returns the value of the variable in the current solution.
1100 *
1101 * If the variable is integer, then the value will always be an integer (the
1102 * underlying solver handles floating-point values only, but this function
1103 * automatically rounds it to the nearest integer; see: man 3 round).
1104 */
1105 double solution_value() const;
1106
1107 /// Returns the index of the variable in the MPSolver::variables_.
index()1108 int index() const { return index_; }
1109
1110 /// Returns the lower bound.
lb()1111 double lb() const { return lb_; }
1112
1113 /// Returns the upper bound.
ub()1114 double ub() const { return ub_; }
1115
1116 /// Sets the lower bound.
SetLB(double lb)1117 void SetLB(double lb) { SetBounds(lb, ub_); }
1118
1119 /// Sets the upper bound.
SetUB(double ub)1120 void SetUB(double ub) { SetBounds(lb_, ub); }
1121
1122 /// Sets both the lower and upper bounds.
1123 void SetBounds(double lb, double ub);
1124
1125 /**
1126 * Advanced usage: unrounded solution value.
1127 *
1128 * The returned value won't be rounded to the nearest integer even if the
1129 * variable is integer.
1130 */
1131 double unrounded_solution_value() const;
1132
1133 /**
1134 * Advanced usage: returns the reduced cost of the variable in the current
1135 * solution (only available for continuous problems).
1136 */
1137 double reduced_cost() const;
1138
1139 /**
1140 * Advanced usage: returns the basis status of the variable in the current
1141 * solution (only available for continuous problems).
1142 *
1143 * @see MPSolver::BasisStatus.
1144 */
1145 MPSolver::BasisStatus basis_status() const;
1146
1147 /**
1148 * Advanced usage: Certain MIP solvers (e.g. Gurobi or SCIP) allow you to set
1149 * a per-variable priority for determining which variable to branch on.
1150 *
1151 * A value of 0 is treated as default, and is equivalent to not setting the
1152 * branching priority. The solver looks first to branch on fractional
1153 * variables in higher priority levels. As of 2019-05, only Gurobi and SCIP
1154 * support setting branching priority; all other solvers will simply ignore
1155 * this annotation.
1156 */
branching_priority()1157 int branching_priority() const { return branching_priority_; }
1158 void SetBranchingPriority(int priority);
1159
1160 protected:
1161 friend class MPSolver;
1162 friend class MPSolverInterface;
1163 friend class CBCInterface;
1164 friend class CLPInterface;
1165 friend class GLPKInterface;
1166 friend class SCIPInterface;
1167 friend class SLMInterface;
1168 friend class GurobiInterface;
1169 friend class CplexInterface;
1170 friend class XpressInterface;
1171 friend class GLOPInterface;
1172 friend class MPVariableSolutionValueTest;
1173 friend class BopInterface;
1174 friend class SatInterface;
1175 friend class KnapsackInterface;
1176
1177 // Constructor. A variable points to a single MPSolverInterface that
1178 // is specified in the constructor. A variable cannot belong to
1179 // several models.
MPVariable(int index,double lb,double ub,bool integer,const std::string & name,MPSolverInterface * const interface_in)1180 MPVariable(int index, double lb, double ub, bool integer,
1181 const std::string& name, MPSolverInterface* const interface_in)
1182 : index_(index),
1183 lb_(lb),
1184 ub_(ub),
1185 integer_(integer),
1186 name_(name.empty() ? absl::StrFormat("auto_v_%09d", index) : name),
1187 solution_value_(0.0),
1188 reduced_cost_(0.0),
1189 interface_(interface_in) {}
1190
set_solution_value(double value)1191 void set_solution_value(double value) { solution_value_ = value; }
set_reduced_cost(double reduced_cost)1192 void set_reduced_cost(double reduced_cost) { reduced_cost_ = reduced_cost; }
1193
1194 private:
1195 const int index_;
1196 double lb_;
1197 double ub_;
1198 bool integer_;
1199 const std::string name_;
1200 double solution_value_;
1201 double reduced_cost_;
1202 int branching_priority_ = 0;
1203 MPSolverInterface* const interface_;
1204 DISALLOW_COPY_AND_ASSIGN(MPVariable);
1205 };
1206
1207 /**
1208 * The class for constraints of a Mathematical Programming (MP) model.
1209 *
1210 * A constraint is represented as a linear equation or inequality.
1211 */
1212 class MPConstraint {
1213 public:
1214 /// Returns the name of the constraint.
name()1215 const std::string& name() const { return name_; }
1216
1217 /// Clears all variables and coefficients. Does not clear the bounds.
1218 void Clear();
1219
1220 /**
1221 * Sets the coefficient of the variable on the constraint.
1222 *
1223 * If the variable does not belong to the solver, the function just returns,
1224 * or crashes in non-opt mode.
1225 */
1226 void SetCoefficient(const MPVariable* const var, double coeff);
1227
1228 /**
1229 * Gets the coefficient of a given variable on the constraint (which is 0 if
1230 * the variable does not appear in the constraint).
1231 */
1232 double GetCoefficient(const MPVariable* const var) const;
1233
1234 /**
1235 * Returns a map from variables to their coefficients in the constraint.
1236 *
1237 * If a variable is not present in the map, then its coefficient is zero.
1238 */
terms()1239 const absl::flat_hash_map<const MPVariable*, double>& terms() const {
1240 return coefficients_;
1241 }
1242
1243 /// Returns the lower bound.
lb()1244 double lb() const { return lb_; }
1245
1246 /// Returns the upper bound.
ub()1247 double ub() const { return ub_; }
1248
1249 /// Sets the lower bound.
SetLB(double lb)1250 void SetLB(double lb) { SetBounds(lb, ub_); }
1251
1252 /// Sets the upper bound.
SetUB(double ub)1253 void SetUB(double ub) { SetBounds(lb_, ub); }
1254
1255 /// Sets both the lower and upper bounds.
1256 void SetBounds(double lb, double ub);
1257
1258 /// Advanced usage: returns true if the constraint is "lazy" (see below).
is_lazy()1259 bool is_lazy() const { return is_lazy_; }
1260
1261 /**
1262 * Advanced usage: sets the constraint "laziness".
1263 *
1264 * <em>This is only supported for SCIP and has no effect on other
1265 * solvers.</em>
1266 *
1267 * When \b laziness is true, the constraint is only considered by the Linear
1268 * Programming solver if its current solution violates the constraint. In this
1269 * case, the constraint is definitively added to the problem. This may be
1270 * useful in some MIP problems, and may have a dramatic impact on performance.
1271 *
1272 * For more info see: http://tinyurl.com/lazy-constraints.
1273 */
set_is_lazy(bool laziness)1274 void set_is_lazy(bool laziness) { is_lazy_ = laziness; }
1275
indicator_variable()1276 const MPVariable* indicator_variable() const { return indicator_variable_; }
indicator_value()1277 bool indicator_value() const { return indicator_value_; }
1278
1279 /// Returns the index of the constraint in the MPSolver::constraints_.
index()1280 int index() const { return index_; }
1281
1282 /**
1283 * Advanced usage: returns the dual value of the constraint in the current
1284 * solution (only available for continuous problems).
1285 */
1286 double dual_value() const;
1287
1288 /**
1289 * Advanced usage: returns the basis status of the constraint.
1290 *
1291 * It is only available for continuous problems).
1292 *
1293 * Note that if a constraint "linear_expression in [lb, ub]" is transformed
1294 * into "linear_expression + slack = 0" with slack in [-ub, -lb], then this
1295 * status is the same as the status of the slack variable with AT_UPPER_BOUND
1296 * and AT_LOWER_BOUND swapped.
1297 *
1298 * @see MPSolver::BasisStatus.
1299 */
1300 MPSolver::BasisStatus basis_status() const;
1301
1302 protected:
1303 friend class MPSolver;
1304 friend class MPSolverInterface;
1305 friend class CBCInterface;
1306 friend class CLPInterface;
1307 friend class GLPKInterface;
1308 friend class SCIPInterface;
1309 friend class SLMInterface;
1310 friend class GurobiInterface;
1311 friend class CplexInterface;
1312 friend class XpressInterface;
1313 friend class GLOPInterface;
1314 friend class BopInterface;
1315 friend class SatInterface;
1316 friend class KnapsackInterface;
1317
1318 // Constructor. A constraint points to a single MPSolverInterface
1319 // that is specified in the constructor. A constraint cannot belong
1320 // to several models.
MPConstraint(int index,double lb,double ub,const std::string & name,MPSolverInterface * const interface_in)1321 MPConstraint(int index, double lb, double ub, const std::string& name,
1322 MPSolverInterface* const interface_in)
1323 : coefficients_(1),
1324 index_(index),
1325 lb_(lb),
1326 ub_(ub),
1327 name_(name.empty() ? absl::StrFormat("auto_c_%09d", index) : name),
1328 is_lazy_(false),
1329 indicator_variable_(nullptr),
1330 dual_value_(0.0),
1331 interface_(interface_in) {}
1332
set_dual_value(double dual_value)1333 void set_dual_value(double dual_value) { dual_value_ = dual_value; }
1334
1335 private:
1336 // Returns true if the constraint contains variables that have not
1337 // been extracted yet.
1338 bool ContainsNewVariables();
1339
1340 // Mapping var -> coefficient.
1341 absl::flat_hash_map<const MPVariable*, double> coefficients_;
1342
1343 const int index_; // See index().
1344
1345 // The lower bound for the linear constraint.
1346 double lb_;
1347
1348 // The upper bound for the linear constraint.
1349 double ub_;
1350
1351 // Name.
1352 const std::string name_;
1353
1354 // True if the constraint is "lazy", i.e. the constraint is added to the
1355 // underlying Linear Programming solver only if it is violated.
1356 // By default this parameter is 'false'.
1357 bool is_lazy_;
1358
1359 // If given, this constraint is only active if `indicator_variable_`'s value
1360 // is equal to `indicator_value_`.
1361 const MPVariable* indicator_variable_;
1362 bool indicator_value_;
1363
1364 double dual_value_;
1365 MPSolverInterface* const interface_;
1366 DISALLOW_COPY_AND_ASSIGN(MPConstraint);
1367 };
1368
1369 /**
1370 * This class stores parameter settings for LP and MIP solvers. Some parameters
1371 * are marked as advanced: do not change their values unless you know what you
1372 * are doing!
1373 *
1374 * For developers: how to add a new parameter:
1375 * - Add the new Foo parameter in the DoubleParam or IntegerParam enum.
1376 * - If it is a categorical param, add a FooValues enum.
1377 * - Decide if the wrapper should define a default value for it: yes
1378 * if it controls the properties of the solution (example:
1379 * tolerances) or if it consistently improves performance, no
1380 * otherwise. If yes, define kDefaultFoo.
1381 * - Add a foo_value_ member and, if no default value is defined, a
1382 * foo_is_default_ member.
1383 * - Add code to handle Foo in Set...Param, Reset...Param,
1384 * Get...Param, Reset and the constructor.
1385 * - In class MPSolverInterface, add a virtual method SetFoo, add it
1386 * to SetCommonParameters or SetMIPParameters, and implement it for
1387 * each solver. Sometimes, parameters need to be implemented
1388 * differently, see for example the INCREMENTALITY implementation.
1389 * - Add a test in linear_solver_test.cc.
1390 *
1391 * TODO(user): store the parameter values in a protocol buffer
1392 * instead. We need to figure out how to deal with the subtleties of
1393 * the default values.
1394 */
1395 class MPSolverParameters {
1396 public:
1397 /// Enumeration of parameters that take continuous values.
1398 enum DoubleParam {
1399 /// Limit for relative MIP gap.
1400 RELATIVE_MIP_GAP = 0,
1401
1402 /**
1403 * Advanced usage: tolerance for primal feasibility of basic solutions.
1404 *
1405 * This does not control the integer feasibility tolerance of integer
1406 * solutions for MIP or the tolerance used during presolve.
1407 */
1408 PRIMAL_TOLERANCE = 1,
1409 /// Advanced usage: tolerance for dual feasibility of basic solutions.
1410 DUAL_TOLERANCE = 2
1411 };
1412
1413 /// Enumeration of parameters that take integer or categorical values.
1414 enum IntegerParam {
1415 /// Advanced usage: presolve mode.
1416 PRESOLVE = 1000,
1417 /// Algorithm to solve linear programs.
1418 LP_ALGORITHM = 1001,
1419 /// Advanced usage: incrementality from one solve to the next.
1420 INCREMENTALITY = 1002,
1421 /// Advanced usage: enable or disable matrix scaling.
1422 SCALING = 1003
1423 };
1424
1425 /// For each categorical parameter, enumeration of possible values.
1426 enum PresolveValues {
1427 /// Presolve is off.
1428 PRESOLVE_OFF = 0,
1429 /// Presolve is on.
1430 PRESOLVE_ON = 1
1431 };
1432
1433 /// LP algorithm to use.
1434 enum LpAlgorithmValues {
1435 /// Dual simplex.
1436 DUAL = 10,
1437 /// Primal simplex.
1438 PRIMAL = 11,
1439 /// Barrier algorithm.
1440 BARRIER = 12
1441 };
1442
1443 /// Advanced usage: Incrementality options.
1444 enum IncrementalityValues {
1445 /// Start solve from scratch.
1446 INCREMENTALITY_OFF = 0,
1447
1448 /**
1449 * Reuse results from previous solve as much as the underlying solver
1450 * allows.
1451 */
1452 INCREMENTALITY_ON = 1
1453 };
1454
1455 /// Advanced usage: Scaling options.
1456 enum ScalingValues {
1457 /// Scaling is off.
1458 SCALING_OFF = 0,
1459 /// Scaling is on.
1460 SCALING_ON = 1
1461 };
1462
1463 // Placeholder value to indicate that a parameter is set to
1464 // the default value defined in the wrapper.
1465 static const double kDefaultDoubleParamValue;
1466 static const int kDefaultIntegerParamValue;
1467
1468 // Placeholder value to indicate that a parameter is unknown.
1469 static const double kUnknownDoubleParamValue;
1470 static const int kUnknownIntegerParamValue;
1471
1472 // Default values for parameters. Only parameters that define the
1473 // properties of the solution returned need to have a default value
1474 // (that is the same for all solvers). You can also define a default
1475 // value for performance parameters when you are confident it is a
1476 // good choice (example: always turn presolve on).
1477 static const double kDefaultRelativeMipGap;
1478 static const double kDefaultPrimalTolerance;
1479 static const double kDefaultDualTolerance;
1480 static const PresolveValues kDefaultPresolve;
1481 static const IncrementalityValues kDefaultIncrementality;
1482
1483 /// The constructor sets all parameters to their default value.
1484 MPSolverParameters();
1485
1486 /// Sets a double parameter to a specific value.
1487 void SetDoubleParam(MPSolverParameters::DoubleParam param, double value);
1488
1489 /// Sets a integer parameter to a specific value.
1490 void SetIntegerParam(MPSolverParameters::IntegerParam param, int value);
1491
1492 /**
1493 * Sets a double parameter to its default value (default value defined in
1494 * MPSolverParameters if it exists, otherwise the default value defined in
1495 * the underlying solver).
1496 */
1497 void ResetDoubleParam(MPSolverParameters::DoubleParam param);
1498
1499 /**
1500 * Sets an integer parameter to its default value (default value defined in
1501 * MPSolverParameters if it exists, otherwise the default value defined in
1502 * the underlying solver).
1503 */
1504 void ResetIntegerParam(MPSolverParameters::IntegerParam param);
1505
1506 /// Sets all parameters to their default value.
1507 void Reset();
1508
1509 /// Returns the value of a double parameter.
1510 double GetDoubleParam(MPSolverParameters::DoubleParam param) const;
1511
1512 /// Returns the value of an integer parameter.
1513 int GetIntegerParam(MPSolverParameters::IntegerParam param) const;
1514
1515 private:
1516 // Parameter value for each parameter.
1517 // @see DoubleParam
1518 // @see IntegerParam
1519 double relative_mip_gap_value_;
1520 double primal_tolerance_value_;
1521 double dual_tolerance_value_;
1522 int presolve_value_;
1523 int scaling_value_;
1524 int lp_algorithm_value_;
1525 int incrementality_value_;
1526
1527 // Boolean value indicating whether each parameter is set to the
1528 // solver's default value. Only parameters for which the wrapper
1529 // does not define a default value need such an indicator.
1530 bool lp_algorithm_is_default_;
1531
1532 DISALLOW_COPY_AND_ASSIGN(MPSolverParameters);
1533 };
1534
1535 // Whether the given MPSolverResponseStatus (of a solve) would yield an RPC
1536 // error when happening on the linear solver stubby server, see
1537 // ./linear_solver_service.proto.
1538 // Note that RPC errors forbid to carry a response to the client, who can only
1539 // see the RPC error itself (error code + error message).
1540 bool MPSolverResponseStatusIsRpcError(MPSolverResponseStatus status);
1541
1542 // This class wraps the actual mathematical programming solvers. Each
1543 // solver (GLOP, CLP, CBC, GLPK, SCIP) has its own interface class that
1544 // derives from this abstract class. This class is never directly
1545 // accessed by the user.
1546 // @see glop_interface.cc
1547 // @see cbc_interface.cc
1548 // @see clp_interface.cc
1549 // @see glpk_interface.cc
1550 // @see scip_interface.cc
1551 class MPSolverInterface {
1552 public:
1553 enum SynchronizationStatus {
1554 // The underlying solver (CLP, GLPK, ...) and MPSolver are not in
1555 // sync for the model nor for the solution.
1556 MUST_RELOAD,
1557 // The underlying solver and MPSolver are in sync for the model
1558 // but not for the solution: the model has changed since the
1559 // solution was computed last.
1560 MODEL_SYNCHRONIZED,
1561 // The underlying solver and MPSolver are in sync for the model and
1562 // the solution.
1563 SOLUTION_SYNCHRONIZED
1564 };
1565
1566 // When the underlying solver does not provide the number of simplex
1567 // iterations.
1568 static constexpr int64_t kUnknownNumberOfIterations = -1;
1569 // When the underlying solver does not provide the number of
1570 // branch-and-bound nodes.
1571 static constexpr int64_t kUnknownNumberOfNodes = -1;
1572
1573 // Constructor. The user will access the MPSolverInterface through the
1574 // MPSolver passed as argument.
1575 explicit MPSolverInterface(MPSolver* const solver);
1576 virtual ~MPSolverInterface();
1577
1578 // ----- Solve -----
1579 // Solves problem with specified parameter values. Returns true if the
1580 // solution is optimal.
1581 virtual MPSolver::ResultStatus Solve(const MPSolverParameters& param) = 0;
1582
1583 // Attempts to directly solve a MPModelRequest, bypassing the MPSolver data
1584 // structures entirely. Like MPSolver::SolveWithProto(), optionally takes in
1585 // an 'interrupt' boolean.
1586 // Returns {} (eg. absl::nullopt) if direct-solve is not supported by the
1587 // underlying solver (possibly because interrupt != nullptr), in which case
1588 // the user should fall back to using MPSolver.
DirectlySolveProto(const MPModelRequest & request,std::atomic<bool> * interrupt)1589 virtual absl::optional<MPSolutionResponse> DirectlySolveProto(
1590 const MPModelRequest& request,
1591 // `interrupt` is non-const because the internal
1592 // solver may set it to true itself, in some cases.
1593 std::atomic<bool>* interrupt) {
1594 return absl::nullopt;
1595 }
1596
1597 // Writes the model using the solver internal write function. Currently only
1598 // available for GurobiInterface.
1599 virtual void Write(const std::string& filename);
1600
1601 // ----- Model modifications and extraction -----
1602 // Resets extracted model.
1603 virtual void Reset() = 0;
1604
1605 // Sets the optimization direction (min/max).
1606 virtual void SetOptimizationDirection(bool maximize) = 0;
1607
1608 // Modifies bounds of an extracted variable.
1609 virtual void SetVariableBounds(int index, double lb, double ub) = 0;
1610
1611 // Modifies integrality of an extracted variable.
1612 virtual void SetVariableInteger(int index, bool integer) = 0;
1613
1614 // Modify bounds of an extracted variable.
1615 virtual void SetConstraintBounds(int index, double lb, double ub) = 0;
1616
1617 // Adds a linear constraint.
1618 virtual void AddRowConstraint(MPConstraint* const ct) = 0;
1619
1620 // Adds an indicator constraint. Returns true if the feature is supported by
1621 // the underlying solver.
AddIndicatorConstraint(MPConstraint * const ct)1622 virtual bool AddIndicatorConstraint(MPConstraint* const ct) {
1623 LOG(ERROR) << "Solver doesn't support indicator constraints.";
1624 return false;
1625 }
1626
1627 // Add a variable.
1628 virtual void AddVariable(MPVariable* const var) = 0;
1629
1630 // Changes a coefficient in a constraint.
1631 virtual void SetCoefficient(MPConstraint* const constraint,
1632 const MPVariable* const variable,
1633 double new_value, double old_value) = 0;
1634
1635 // Clears a constraint from all its terms.
1636 virtual void ClearConstraint(MPConstraint* const constraint) = 0;
1637
1638 // Changes a coefficient in the linear objective.
1639 virtual void SetObjectiveCoefficient(const MPVariable* const variable,
1640 double coefficient) = 0;
1641
1642 // Changes the constant term in the linear objective.
1643 virtual void SetObjectiveOffset(double value) = 0;
1644
1645 // Clears the objective from all its terms.
1646 virtual void ClearObjective() = 0;
1647
BranchingPriorityChangedForVariable(int var_index)1648 virtual void BranchingPriorityChangedForVariable(int var_index) {}
1649 // ------ Query statistics on the solution and the solve ------
1650 // Returns the number of simplex iterations. The problem must be discrete,
1651 // otherwise it crashes, or returns kUnknownNumberOfIterations in NDEBUG mode.
1652 virtual int64_t iterations() const = 0;
1653 // Returns the number of branch-and-bound nodes. The problem must be discrete,
1654 // otherwise it crashes, or returns kUnknownNumberOfNodes in NDEBUG mode.
1655 virtual int64_t nodes() const = 0;
1656 // Returns the best objective bound. The problem must be discrete, otherwise
1657 // it crashes, or returns trivial bound (+/- inf) in NDEBUG mode.
1658 double best_objective_bound() const;
1659 // Returns the objective value of the best solution found so far.
1660 double objective_value() const;
1661
1662 // Returns the basis status of a row.
1663 virtual MPSolver::BasisStatus row_status(int constraint_index) const = 0;
1664 // Returns the basis status of a constraint.
1665 virtual MPSolver::BasisStatus column_status(int variable_index) const = 0;
1666
1667 // Checks whether the solution is synchronized with the model, i.e. whether
1668 // the model has changed since the solution was computed last.
1669 // If it isn't, it crashes in NDEBUG, and returns false othwerwise.
1670 bool CheckSolutionIsSynchronized() const;
1671 // Checks whether a feasible solution exists. The behavior is similar to
1672 // CheckSolutionIsSynchronized() above.
1673 virtual bool CheckSolutionExists() const;
1674 // Handy shortcut to do both checks above (it is often used).
CheckSolutionIsSynchronizedAndExists()1675 bool CheckSolutionIsSynchronizedAndExists() const {
1676 return CheckSolutionIsSynchronized() && CheckSolutionExists();
1677 }
1678
1679 // ----- Misc -----
1680 // Queries problem type. For simplicity, the distinction between
1681 // continuous and discrete is based on the declaration of the user
1682 // when the solver is created (example: GLPK_LINEAR_PROGRAMMING
1683 // vs. GLPK_MIXED_INTEGER_PROGRAMMING), not on the actual content of
1684 // the model.
1685 // Returns true if the problem is continuous.
1686 virtual bool IsContinuous() const = 0;
1687 // Returns true if the problem is continuous and linear.
1688 virtual bool IsLP() const = 0;
1689 // Returns true if the problem is discrete and linear.
1690 virtual bool IsMIP() const = 0;
1691
1692 // Returns the index of the last variable extracted.
last_variable_index()1693 int last_variable_index() const { return last_variable_index_; }
1694
variable_is_extracted(int var_index)1695 bool variable_is_extracted(int var_index) const {
1696 return solver_->variable_is_extracted_[var_index];
1697 }
set_variable_as_extracted(int var_index,bool extracted)1698 void set_variable_as_extracted(int var_index, bool extracted) {
1699 solver_->variable_is_extracted_[var_index] = extracted;
1700 }
constraint_is_extracted(int ct_index)1701 bool constraint_is_extracted(int ct_index) const {
1702 return solver_->constraint_is_extracted_[ct_index];
1703 }
set_constraint_as_extracted(int ct_index,bool extracted)1704 void set_constraint_as_extracted(int ct_index, bool extracted) {
1705 solver_->constraint_is_extracted_[ct_index] = extracted;
1706 }
1707
1708 // Returns the boolean indicating the verbosity of the solver output.
quiet()1709 bool quiet() const { return quiet_; }
1710 // Sets the boolean indicating the verbosity of the solver output.
set_quiet(bool quiet_value)1711 void set_quiet(bool quiet_value) { quiet_ = quiet_value; }
1712
1713 // Returns the result status of the last solve.
result_status()1714 MPSolver::ResultStatus result_status() const {
1715 CheckSolutionIsSynchronized();
1716 return result_status_;
1717 }
1718
1719 // Returns a string describing the underlying solver and its version.
1720 virtual std::string SolverVersion() const = 0;
1721
1722 // Returns the underlying solver.
1723 virtual void* underlying_solver() = 0;
1724
1725 // Computes exact condition number. Only available for continuous
1726 // problems and only implemented in GLPK.
1727 virtual double ComputeExactConditionNumber() const;
1728
1729 // See MPSolver::SetStartingLpBasis().
SetStartingLpBasis(const std::vector<MPSolver::BasisStatus> & variable_statuses,const std::vector<MPSolver::BasisStatus> & constraint_statuses)1730 virtual void SetStartingLpBasis(
1731 const std::vector<MPSolver::BasisStatus>& variable_statuses,
1732 const std::vector<MPSolver::BasisStatus>& constraint_statuses) {
1733 LOG(FATAL) << "Not supported by this solver.";
1734 }
1735
InterruptSolve()1736 virtual bool InterruptSolve() { return false; }
1737
1738 // See MPSolver::NextSolution() for contract.
NextSolution()1739 virtual bool NextSolution() { return false; }
1740
1741 // See MPSolver::SetCallback() for details.
SetCallback(MPCallback * mp_callback)1742 virtual void SetCallback(MPCallback* mp_callback) {
1743 LOG(FATAL) << "Callbacks not supported for this solver.";
1744 }
1745
SupportsCallbacks()1746 virtual bool SupportsCallbacks() const { return false; }
1747
1748 friend class MPSolver;
1749
1750 // To access the maximize_ bool and the MPSolver.
1751 friend class MPConstraint;
1752 friend class MPObjective;
1753
1754 protected:
1755 MPSolver* const solver_;
1756 // Indicates whether the model and the solution are synchronized.
1757 SynchronizationStatus sync_status_;
1758 // Indicates whether the solve has reached optimality,
1759 // infeasibility, a limit, etc.
1760 MPSolver::ResultStatus result_status_;
1761 // Optimization direction.
1762 bool maximize_;
1763
1764 // Index in MPSolver::variables_ of last constraint extracted.
1765 int last_constraint_index_;
1766 // Index in MPSolver::constraints_ of last variable extracted.
1767 int last_variable_index_;
1768
1769 // The value of the objective function.
1770 double objective_value_;
1771
1772 // The value of the best objective bound. Used only for MIP solvers.
1773 double best_objective_bound_;
1774
1775 // Boolean indicator for the verbosity of the solver output.
1776 bool quiet_;
1777
1778 // Index of dummy variable created for empty constraints or the
1779 // objective offset.
1780 static const int kDummyVariableIndex;
1781
1782 // Extracts model stored in MPSolver.
1783 void ExtractModel();
1784 // Extracts the variables that have not been extracted yet.
1785 virtual void ExtractNewVariables() = 0;
1786 // Extracts the constraints that have not been extracted yet.
1787 virtual void ExtractNewConstraints() = 0;
1788 // Extracts the objective.
1789 virtual void ExtractObjective() = 0;
1790 // Resets the extraction information.
1791 void ResetExtractionInformation();
1792 // Change synchronization status from SOLUTION_SYNCHRONIZED to
1793 // MODEL_SYNCHRONIZED. To be used for model changes.
1794 void InvalidateSolutionSynchronization();
1795
1796 // Sets parameters common to LP and MIP in the underlying solver.
1797 void SetCommonParameters(const MPSolverParameters& param);
1798 // Sets MIP specific parameters in the underlying solver.
1799 void SetMIPParameters(const MPSolverParameters& param);
1800 // Sets all parameters in the underlying solver.
1801 virtual void SetParameters(const MPSolverParameters& param) = 0;
1802 // Sets an unsupported double parameter.
1803 void SetUnsupportedDoubleParam(MPSolverParameters::DoubleParam param);
1804 // Sets an unsupported integer parameter.
1805 virtual void SetUnsupportedIntegerParam(
1806 MPSolverParameters::IntegerParam param);
1807 // Sets a supported double parameter to an unsupported value.
1808 void SetDoubleParamToUnsupportedValue(MPSolverParameters::DoubleParam param,
1809 double value);
1810 // Sets a supported integer parameter to an unsupported value.
1811 virtual void SetIntegerParamToUnsupportedValue(
1812 MPSolverParameters::IntegerParam param, int value);
1813 // Sets each parameter in the underlying solver.
1814 virtual void SetRelativeMipGap(double value) = 0;
1815 virtual void SetPrimalTolerance(double value) = 0;
1816 virtual void SetDualTolerance(double value) = 0;
1817 virtual void SetPresolveMode(int value) = 0;
1818
1819 // Sets the number of threads to be used by the solver.
1820 virtual absl::Status SetNumThreads(int num_threads);
1821
1822 // Pass solver specific parameters in text format. The format is
1823 // solver-specific and is the same as the corresponding solver configuration
1824 // file format. Returns true if the operation was successful.
1825 //
1826 // Default implementation returns true if the input is empty. It returns false
1827 // and logs a WARNING if the input is not empty.
1828 virtual bool SetSolverSpecificParametersAsString(
1829 const std::string& parameters);
1830
1831 // Sets the scaling mode.
1832 virtual void SetScalingMode(int value) = 0;
1833 virtual void SetLpAlgorithm(int value) = 0;
1834 };
1835
1836 } // namespace operations_research
1837
1838 #endif // OR_TOOLS_LINEAR_SOLVER_LINEAR_SOLVER_H_
1839