OR-Tools  8.2
xpress_interface.cc
Go to the documentation of this file.
1 // Copyright 2019 RTE
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 // Initial version of this code was provided by RTE
15 
16 #if defined(USE_XPRESS)
17 
18 #include <algorithm>
19 #include <limits>
20 #include <memory>
21 #include <string>
22 
23 #include "absl/strings/str_format.h"
25 #include "ortools/base/logging.h"
26 #include "ortools/base/timer.h"
28 
29 extern "C" {
30 #include "xprs.h"
31 }
32 
33 #define XPRS_INTEGER 'I'
34 #define XPRS_CONTINUOUS 'C'
35 #define STRINGIFY2(X) #X
36 #define STRINGIFY(X) STRINGIFY2(X)
37 
38 void printError(const XPRSprob& mLp, int line) {
39  char errmsg[512];
40  XPRSgetlasterror(mLp, errmsg);
41  VLOG(0) << absl::StrFormat("Function line %d did not execute correctly: %s\n",
42  line, errmsg);
43  exit(0);
44 }
45 
46 int XPRSgetnumcols(const XPRSprob& mLp) {
47  int nCols = 0;
48  XPRSgetintattrib(mLp, XPRS_COLS, &nCols);
49  return nCols;
50 }
51 
52 int XPRSgetnumrows(const XPRSprob& mLp) {
53  int nRows = 0;
54  XPRSgetintattrib(mLp, XPRS_ROWS, &nRows);
55  return nRows;
56 }
57 
58 int XPRSgetitcnt(const XPRSprob& mLp) {
59  int nIters = 0;
60  XPRSgetintattrib(mLp, XPRS_SIMPLEXITER, &nIters);
61  return nIters;
62 }
63 
64 int XPRSgetnodecnt(const XPRSprob& mLp) {
65  int nNodes = 0;
66  XPRSgetintattrib(mLp, XPRS_NODES, &nNodes);
67  return nNodes;
68 }
69 
70 int XPRSsetobjoffset(const XPRSprob& mLp, double value) {
71  XPRSsetdblcontrol(mLp, XPRS_OBJRHS, value);
72  return 0;
73 }
74 
75 enum XPRS_BASIS_STATUS {
76  XPRS_AT_LOWER = 0,
77  XPRS_BASIC = 1,
78  XPRS_AT_UPPER = 2,
79  XPRS_FREE_SUPER = 3
80 };
81 
82 // In case we need to return a double but don't have a value for that
83 // we just return a NaN.
84 #if !defined(CPX_NAN)
85 #define XPRS_NAN std::numeric_limits<double>::quiet_NaN()
86 #endif
87 
88 // The argument to this macro is the invocation of a XPRS function that
89 // returns a status. If the function returns non-zero the macro aborts
90 // the program with an appropriate error message.
91 #define CHECK_STATUS(s) \
92  do { \
93  int const status_ = s; \
94  CHECK_EQ(0, status_); \
95  } while (0)
96 
97 namespace operations_research {
98 
99 using std::unique_ptr;
100 
101 // For a model that is extracted to an instance of this class there is a
102 // 1:1 corresponence between MPVariable instances and XPRESS columns: the
103 // index of an extracted variable is the column index in the XPRESS model.
104 // Similar for instances of MPConstraint: the index of the constraint in
105 // the model is the row index in the XPRESS model.
106 class XpressInterface : public MPSolverInterface {
107  public:
108  // NOTE: 'mip' specifies the type of the problem (either continuous or
109  // mixed integer. This type is fixed for the lifetime of the
110  // instance. There are no dynamic changes to the model type.
111  explicit XpressInterface(MPSolver* const solver, bool mip);
112  ~XpressInterface();
113 
114  // Sets the optimization direction (min/max).
115  virtual void SetOptimizationDirection(bool maximize);
116 
117  // ----- Solve -----
118  // Solve the problem using the parameter values specified.
119  virtual MPSolver::ResultStatus Solve(MPSolverParameters const& param);
120 
121  // ----- Model modifications and extraction -----
122  // Resets extracted model
123  virtual void Reset();
124 
125  virtual void SetVariableBounds(int var_index, double lb, double ub);
126  virtual void SetVariableInteger(int var_index, bool integer);
127  virtual void SetConstraintBounds(int row_index, double lb, double ub);
128 
129  virtual void AddRowConstraint(MPConstraint* const ct);
130  virtual void AddVariable(MPVariable* const var);
131  virtual void SetCoefficient(MPConstraint* const constraint,
132  MPVariable const* const variable,
133  double new_value, double old_value);
134 
135  // Clear a constraint from all its terms.
136  virtual void ClearConstraint(MPConstraint* const constraint);
137  // Change a coefficient in the linear objective
138  virtual void SetObjectiveCoefficient(MPVariable const* const variable,
139  double coefficient);
140  // Change the constant term in the linear objective.
141  virtual void SetObjectiveOffset(double value);
142  // Clear the objective from all its terms.
143  virtual void ClearObjective();
144 
145  // ------ Query statistics on the solution and the solve ------
146  // Number of simplex iterations
147  virtual int64 iterations() const;
148  // Number of branch-and-bound nodes. Only available for discrete problems.
149  virtual int64 nodes() const;
150 
151  // Returns the basis status of a row.
152  virtual MPSolver::BasisStatus row_status(int constraint_index) const;
153  // Returns the basis status of a column.
154  virtual MPSolver::BasisStatus column_status(int variable_index) const;
155 
156  // ----- Misc -----
157 
158  // Query problem type.
159  // Remember that problem type is a static property that is set
160  // in the constructor and never changed.
161  virtual bool IsContinuous() const { return IsLP(); }
162  virtual bool IsLP() const { return !mMip; }
163  virtual bool IsMIP() const { return mMip; }
164 
165  virtual void ExtractNewVariables();
166  virtual void ExtractNewConstraints();
167  virtual void ExtractObjective();
168 
169  virtual std::string SolverVersion() const;
170 
171  virtual void* underlying_solver() { return reinterpret_cast<void*>(mLp); }
172 
173  virtual double ComputeExactConditionNumber() const {
174  if (!IsContinuous()) {
175  LOG(DFATAL) << "ComputeExactConditionNumber not implemented for"
176  << " XPRESS_MIXED_INTEGER_PROGRAMMING";
177  return 0.0;
178  }
179 
180  // TODO(user,user): Not yet working.
181  LOG(DFATAL) << "ComputeExactConditionNumber not implemented for"
182  << " XPRESS_LINEAR_PROGRAMMING";
183  return 0.0;
184  }
185 
186  protected:
187  // Set all parameters in the underlying solver.
188  virtual void SetParameters(MPSolverParameters const& param);
189  // Set each parameter in the underlying solver.
190  virtual void SetRelativeMipGap(double value);
191  virtual void SetPrimalTolerance(double value);
192  virtual void SetDualTolerance(double value);
193  virtual void SetPresolveMode(int value);
194  virtual void SetScalingMode(int value);
195  virtual void SetLpAlgorithm(int value);
196 
197  virtual bool ReadParameterFile(std::string const& filename);
198  virtual std::string ValidFileExtensionForParameterFile() const;
199 
200  private:
201  // Mark modeling object "out of sync". This implicitly invalidates
202  // solution information as well. It is the counterpart of
203  // MPSolverInterface::InvalidateSolutionSynchronization
204  void InvalidateModelSynchronization() {
205  mCstat = 0;
206  mRstat = 0;
207  sync_status_ = MUST_RELOAD;
208  }
209 
210  // Transform XPRESS basis status to MPSolver basis status.
211  static MPSolver::BasisStatus xformBasisStatus(int xpress_basis_status);
212 
213  private:
214  XPRSprob mLp;
215  bool const mMip;
216  // Incremental extraction.
217  // Without incremental extraction we have to re-extract the model every
218  // time we perform a solve. Due to the way the Reset() function is
219  // implemented, this will lose MIP start or basis information from a
220  // previous solve. On the other hand, if there is a significant changes
221  // to the model then just re-extracting everything is usually faster than
222  // keeping the low-level modeling object in sync with the high-level
223  // variables/constraints.
224  // Note that incremental extraction is particularly expensive in function
225  // ExtractNewVariables() since there we must scan _all_ old constraints
226  // and update them with respect to the new variables.
227  bool const supportIncrementalExtraction;
228 
229  // Use slow and immediate updates or try to do bulk updates.
230  // For many updates to the model we have the option to either perform
231  // the update immediately with a potentially slow operation or to
232  // just mark the low-level modeling object out of sync and re-extract
233  // the model later.
234  enum SlowUpdates {
235  SlowSetCoefficient = 0x0001,
236  SlowClearConstraint = 0x0002,
237  SlowSetObjectiveCoefficient = 0x0004,
238  SlowClearObjective = 0x0008,
239  SlowSetConstraintBounds = 0x0010,
240  SlowSetVariableInteger = 0x0020,
241  SlowSetVariableBounds = 0x0040,
242  SlowUpdatesAll = 0xffff
243  } const slowUpdates;
244  // XPRESS has no method to query the basis status of a single variable.
245  // Hence we query the status only once and cache the array. This is
246  // much faster in case the basis status of more than one row/column
247  // is required.
248  unique_ptr<int[]> mutable mCstat;
249  unique_ptr<int[]> mutable mRstat;
250 
251  // Setup the right-hand side of a constraint from its lower and upper bound.
252  static void MakeRhs(double lb, double ub, double& rhs, char& sense,
253  double& range);
254 };
255 
257 int init_xpress_env(int xpress_oem_license_key = 0) {
258  int code;
259 
260  const char* xpress_from_env = getenv("XPRESS");
261  std::string xpresspath;
262 
263  if (xpress_from_env == nullptr) {
264 #if defined(XPRESS_PATH)
265  std::string path(STRINGIFY(XPRESS_PATH));
266  LOG(WARNING)
267  << "Environment variable XPRESS undefined. Trying compile path "
268  << "'" << path << "'";
269 #if defined(_MSC_VER)
270  // need to remove the enclosing '\"' from the string itself.
271  path.erase(std::remove(path.begin(), path.end(), '\"'), path.end());
272  xpresspath = path + "\\bin";
273 #else // _MSC_VER
274  xpresspath = path + "/bin";
275 #endif // _MSC_VER
276 #else
277  LOG(WARNING)
278  << "XpressInterface Error : Environment variable XPRESS undefined.\n";
279  return -1;
280 #endif
281  } else {
282  xpresspath = xpress_from_env;
283  }
284 
286  if (xpress_oem_license_key == 0) {
287  LOG(WARNING) << "XpressInterface : Initialising xpress-MP with parameter "
288  << xpresspath << std::endl;
289 
290  code = XPRSinit(xpresspath.c_str());
291 
292  if (!code) {
294  char banner[1000];
295  XPRSgetbanner(banner);
296 
297  LOG(WARNING) << "XpressInterface : Xpress banner :\n"
298  << banner << std::endl;
299  return 0;
300  } else {
301  char errmsg[256];
302  XPRSgetlicerrmsg(errmsg, 256);
303 
304  VLOG(0) << "XpressInterface : License error : " << errmsg << std::endl;
305  VLOG(0) << "XpressInterface : XPRSinit returned code : " << code << "\n";
306 
307  char banner[1000];
308  XPRSgetbanner(banner);
309 
310  LOG(ERROR) << "XpressInterface : Xpress banner :\n" << banner << "\n";
311  return -1;
312  }
313  } else {
315  LOG(WARNING) << "XpressInterface : Initialising xpress-MP with OEM key "
316  << xpress_oem_license_key << "\n";
317 
318  int nvalue = 0;
319  int ierr;
320  char slicmsg[256] = "";
321  char errmsg[256];
322 
323  XPRSlicense(&nvalue, slicmsg);
324  VLOG(0) << "XpressInterface : First message from XPRSLicense : " << slicmsg
325  << "\n";
326 
327  nvalue = xpress_oem_license_key - ((nvalue * nvalue) / 19);
328  ierr = XPRSlicense(&nvalue, slicmsg);
329 
330  VLOG(0) << "XpressInterface : Second message from XPRSLicense : " << slicmsg
331  << "\n";
332  if (ierr == 16) {
333  VLOG(0) << "XpressInterface : Optimizer development software detected\n";
334  } else if (ierr != 0) {
336  XPRSgetlicerrmsg(errmsg, 256);
337 
338  LOG(ERROR) << "XpressInterface : " << errmsg << "\n";
339  return -1;
340  }
341 
342  code = XPRSinit(NULL);
343 
344  if (!code) {
345  return 0;
346  } else {
347  LOG(ERROR) << "XPRSinit returned code : " << code << "\n";
348  return -1;
349  }
350  }
351 }
352 
353 // Creates a LP/MIP instance.
354 XpressInterface::XpressInterface(MPSolver* const solver, bool mip)
355  : MPSolverInterface(solver),
356  mLp(0),
357  mMip(mip),
358  supportIncrementalExtraction(false),
359  slowUpdates(static_cast<SlowUpdates>(SlowSetObjectiveCoefficient |
360  SlowClearObjective)),
361  mCstat(),
362  mRstat() {
363  int status = init_xpress_env();
364  CHECK_STATUS(status);
365  status = XPRScreateprob(&mLp);
366  CHECK_STATUS(status);
367  DCHECK(mLp != nullptr); // should not be NULL if status=0
368  CHECK_STATUS(XPRSloadlp(mLp, "newProb", 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0));
369 
370  CHECK_STATUS(
371  XPRSchgobjsense(mLp, maximize_ ? XPRS_OBJ_MAXIMIZE : XPRS_OBJ_MINIMIZE));
372 }
373 
374 XpressInterface::~XpressInterface() {
375  CHECK_STATUS(XPRSdestroyprob(mLp));
377 }
378 
379 std::string XpressInterface::SolverVersion() const {
380  // We prefer XPRSversionnumber() over XPRSversion() since the
381  // former will never pose any encoding issues.
382  int version = 0;
383  CHECK_STATUS(XPRSgetintcontrol(mLp, XPRS_VERSION, &version));
384 
385  int const major = version / 1000000;
386  version -= major * 1000000;
387  int const release = version / 10000;
388  version -= release * 10000;
389  int const mod = version / 100;
390  version -= mod * 100;
391  int const fix = version;
392 
393  return absl::StrFormat("XPRESS library version %d.%02d.%02d.%02d", major,
394  release, mod, fix);
395 }
396 
397 // ------ Model modifications and extraction -----
398 
399 void XpressInterface::Reset() {
400  // Instead of explicitly clearing all modeling objects we
401  // just delete the problem object and allocate a new one.
402  CHECK_STATUS(XPRSdestroyprob(mLp));
403 
404  int status;
405  status = XPRScreateprob(&mLp);
406  CHECK_STATUS(status);
407  DCHECK(mLp != nullptr); // should not be NULL if status=0
408  CHECK_STATUS(XPRSloadlp(mLp, "newProb", 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0));
409 
410  CHECK_STATUS(
411  XPRSchgobjsense(mLp, maximize_ ? XPRS_OBJ_MAXIMIZE : XPRS_OBJ_MINIMIZE));
412 
413  ResetExtractionInformation();
414  mCstat = 0;
415  mRstat = 0;
416 }
417 
418 void XpressInterface::SetOptimizationDirection(bool maximize) {
419  InvalidateSolutionSynchronization();
420  XPRSchgobjsense(mLp, maximize ? XPRS_OBJ_MAXIMIZE : XPRS_OBJ_MINIMIZE);
421 }
422 
423 void XpressInterface::SetVariableBounds(int var_index, double lb, double ub) {
424  InvalidateSolutionSynchronization();
425 
426  // Changing the bounds of a variable is fast. However, doing this for
427  // many variables may still be slow. So we don't perform the update by
428  // default. However, if we support incremental extraction
429  // (supportIncrementalExtraction is true) then we MUST perform the
430  // update here or we will lose it.
431 
432  if (!supportIncrementalExtraction && !(slowUpdates & SlowSetVariableBounds)) {
433  InvalidateModelSynchronization();
434  } else {
435  if (variable_is_extracted(var_index)) {
436  // Variable has already been extracted, so we must modify the
437  // modeling object.
438  DCHECK_LT(var_index, last_variable_index_);
439  char const lu[2] = {'L', 'U'};
440  double const bd[2] = {lb, ub};
441  int const idx[2] = {var_index, var_index};
442  CHECK_STATUS(XPRSchgbounds(mLp, 2, idx, lu, bd));
443  } else {
444  // Variable is not yet extracted. It is sufficient to just mark
445  // the modeling object "out of sync"
446  InvalidateModelSynchronization();
447  }
448  }
449 }
450 
451 // Modifies integrality of an extracted variable.
452 void XpressInterface::SetVariableInteger(int var_index, bool integer) {
453  InvalidateSolutionSynchronization();
454 
455  // NOTE: The type of the model (continuous or mixed integer) is
456  // defined once and for all in the constructor. There are no
457  // dynamic changes to the model type.
458 
459  // Changing the type of a variable should be fast. Still, doing all
460  // updates in one big chunk right before solve() is usually faster.
461  // However, if we support incremental extraction
462  // (supportIncrementalExtraction is true) then we MUST change the
463  // type of extracted variables here.
464 
465  if (!supportIncrementalExtraction && !slowUpdates &&
466  !SlowSetVariableInteger) {
467  InvalidateModelSynchronization();
468  } else {
469  if (mMip) {
470  if (variable_is_extracted(var_index)) {
471  // Variable is extracted. Change the type immediately.
472  // TODO: Should we check the current type and don't do anything
473  // in case the type does not change?
474  DCHECK_LE(var_index, XPRSgetnumcols(mLp));
475  char const type = integer ? XPRS_INTEGER : XPRS_CONTINUOUS;
476  CHECK_STATUS(XPRSchgcoltype(mLp, 1, &var_index, &type));
477  } else {
478  InvalidateModelSynchronization();
479  }
480  } else {
481  LOG(DFATAL)
482  << "Attempt to change variable to integer in non-MIP problem!";
483  }
484  }
485 }
486 
487 // Setup the right-hand side of a constraint.
488 void XpressInterface::MakeRhs(double lb, double ub, double& rhs, char& sense,
489  double& range) {
490  if (lb == ub) {
491  // Both bounds are equal -> this is an equality constraint
492  rhs = lb;
493  range = 0.0;
494  sense = 'E';
495  } else if (lb > XPRS_MINUSINFINITY && ub < XPRS_PLUSINFINITY) {
496  // Both bounds are finite -> this is a ranged constraint
497  // The value of a ranged constraint is allowed to be in
498  // [ rhs[i], rhs[i]+rngval[i] ]
499  // see also the reference documentation for XPRSnewrows()
500  if (ub < lb) {
501  // The bounds for the constraint are contradictory. XPRESS models
502  // a range constraint l <= ax <= u as
503  // ax = l + v
504  // where v is an auxiliary variable the range of which is controlled
505  // by l and u: if l < u then v in [0, u-l]
506  // else v in [u-l, 0]
507  // (the range is specified as the rngval[] argument to XPRSnewrows).
508  // Thus XPRESS cannot represent range constraints with contradictory
509  // bounds and we must error out here.
510  CHECK_STATUS(-1);
511  }
512  rhs = ub;
513  range = ub - lb;
514  sense = 'R';
515  } else if (ub < XPRS_PLUSINFINITY || (std::abs(ub) == XPRS_PLUSINFINITY &&
516  std::abs(lb) > XPRS_PLUSINFINITY)) {
517  // Finite upper, infinite lower bound -> this is a <= constraint
518  rhs = ub;
519  range = 0.0;
520  sense = 'L';
521  } else if (lb > XPRS_MINUSINFINITY || (std::abs(lb) == XPRS_PLUSINFINITY &&
522  std::abs(ub) > XPRS_PLUSINFINITY)) {
523  // Finite lower, infinite upper bound -> this is a >= constraint
524  rhs = lb;
525  range = 0.0;
526  sense = 'G';
527  } else {
528  // Lower and upper bound are both infinite.
529  // This is used for example in .mps files to specify alternate
530  // objective functions.
531  // Note that the case lb==ub was already handled above, so we just
532  // pick the bound with larger magnitude and create a constraint for it.
533  // Note that we replace the infinite bound by XPRS_PLUSINFINITY since
534  // bounds with larger magnitude may cause other XPRESS functions to
535  // fail (for example the export to LP files).
536  DCHECK_GT(std::abs(lb), XPRS_PLUSINFINITY);
537  DCHECK_GT(std::abs(ub), XPRS_PLUSINFINITY);
538  if (std::abs(lb) > std::abs(ub)) {
539  rhs = (lb < 0) ? -XPRS_PLUSINFINITY : XPRS_PLUSINFINITY;
540  sense = 'G';
541  } else {
542  rhs = (ub < 0) ? -XPRS_PLUSINFINITY : XPRS_PLUSINFINITY;
543  sense = 'L';
544  }
545  range = 0.0;
546  }
547 }
548 
549 void XpressInterface::SetConstraintBounds(int index, double lb, double ub) {
550  InvalidateSolutionSynchronization();
551 
552  // Changing rhs, sense, or range of a constraint is not too slow.
553  // Still, doing all the updates in one large operation is faster.
554  // Note however that if we do not want to re-extract the full model
555  // for each solve (supportIncrementalExtraction is true) then we MUST
556  // update the constraint here, otherwise we lose this update information.
557 
558  if (!supportIncrementalExtraction &&
559  !(slowUpdates & SlowSetConstraintBounds)) {
560  InvalidateModelSynchronization();
561  } else {
562  if (constraint_is_extracted(index)) {
563  // Constraint is already extracted, so we must update its bounds
564  // and its type.
565  DCHECK(mLp != NULL);
566  char sense;
567  double range, rhs;
568  MakeRhs(lb, ub, rhs, sense, range);
569  CHECK_STATUS(XPRSchgrhs(mLp, 1, &index, &lb));
570  CHECK_STATUS(XPRSchgrowtype(mLp, 1, &index, &sense));
571  CHECK_STATUS(XPRSchgrhsrange(mLp, 1, &index, &range));
572  } else {
573  // Constraint is not yet extracted. It is sufficient to mark the
574  // modeling object as "out of sync"
575  InvalidateModelSynchronization();
576  }
577  }
578 }
579 
580 void XpressInterface::AddRowConstraint(MPConstraint* const ct) {
581  // This is currently only invoked when a new constraint is created,
582  // see MPSolver::MakeRowConstraint().
583  // At this point we only have the lower and upper bounds of the
584  // constraint. We could immediately call XPRSaddrows() here but it is
585  // usually much faster to handle the fully populated constraint in
586  // ExtractNewConstraints() right before the solve.
587  InvalidateModelSynchronization();
588 }
589 
590 void XpressInterface::AddVariable(MPVariable* const ct) {
591  // This is currently only invoked when a new variable is created,
592  // see MPSolver::MakeVar().
593  // At this point the variable does not appear in any constraints or
594  // the objective function. We could invoke XPRSaddcols() to immediately
595  // create the variable here but it is usually much faster to handle the
596  // fully setup variable in ExtractNewVariables() right before the solve.
597  InvalidateModelSynchronization();
598 }
599 
600 void XpressInterface::SetCoefficient(MPConstraint* const constraint,
601  MPVariable const* const variable,
602  double new_value, double) {
603  InvalidateSolutionSynchronization();
604 
605  // Changing a single coefficient in the matrix is potentially pretty
606  // slow since that coefficient has to be found in the sparse matrix
607  // representation. So by default we don't perform this update immediately
608  // but instead mark the low-level modeling object "out of sync".
609  // If we want to support incremental extraction then we MUST perform
610  // the modification immediately or we will lose it.
611 
612  if (!supportIncrementalExtraction && !(slowUpdates & SlowSetCoefficient)) {
613  InvalidateModelSynchronization();
614  } else {
615  int const row = constraint->index();
616  int const col = variable->index();
617  if (constraint_is_extracted(row) && variable_is_extracted(col)) {
618  // If row and column are both extracted then we can directly
619  // update the modeling object
620  DCHECK_LE(row, last_constraint_index_);
621  DCHECK_LE(col, last_variable_index_);
622  CHECK_STATUS(XPRSchgcoef(mLp, row, col, new_value));
623  } else {
624  // If either row or column is not yet extracted then we can
625  // defer the update to ExtractModel()
626  InvalidateModelSynchronization();
627  }
628  }
629 }
630 
631 void XpressInterface::ClearConstraint(MPConstraint* const constraint) {
632  int const row = constraint->index();
633  if (!constraint_is_extracted(row))
634  // There is nothing to do if the constraint was not even extracted.
635  return;
636 
637  // Clearing a constraint means setting all coefficients in the corresponding
638  // row to 0 (we cannot just delete the row since that would renumber all
639  // the constraints/rows after it).
640  // Modifying coefficients in the matrix is potentially pretty expensive
641  // since they must be found in the sparse matrix representation. That is
642  // why by default we do not modify the coefficients here but only mark
643  // the low-level modeling object "out of sync".
644 
645  if (!(slowUpdates & SlowClearConstraint)) {
646  InvalidateModelSynchronization();
647  } else {
648  InvalidateSolutionSynchronization();
649 
650  int const len = constraint->coefficients_.size();
651  unique_ptr<int[]> rowind(new int[len]);
652  unique_ptr<int[]> colind(new int[len]);
653  unique_ptr<double[]> val(new double[len]);
654  int j = 0;
655  const auto& coeffs = constraint->coefficients_;
656  for (auto it(coeffs.begin()); it != coeffs.end(); ++it) {
657  int const col = it->first->index();
658  if (variable_is_extracted(col)) {
659  rowind[j] = row;
660  colind[j] = col;
661  val[j] = 0.0;
662  ++j;
663  }
664  }
665  if (j)
666  CHECK_STATUS(XPRSchgmcoef(mLp, j, rowind.get(), colind.get(), val.get()));
667  }
668 }
669 
670 void XpressInterface::SetObjectiveCoefficient(MPVariable const* const variable,
671  double coefficient) {
672  int const col = variable->index();
673  if (!variable_is_extracted(col))
674  // Nothing to do if variable was not even extracted
675  return;
676 
677  InvalidateSolutionSynchronization();
678 
679  // The objective function is stored as a dense vector, so updating a
680  // single coefficient is O(1). So by default we update the low-level
681  // modeling object here.
682  // If we support incremental extraction then we have no choice but to
683  // perform the update immediately.
684 
685  if (supportIncrementalExtraction ||
686  (slowUpdates & SlowSetObjectiveCoefficient)) {
687  CHECK_STATUS(XPRSchgobj(mLp, 1, &col, &coefficient));
688  } else {
689  InvalidateModelSynchronization();
690  }
691 }
692 
693 void XpressInterface::SetObjectiveOffset(double value) {
694  // Changing the objective offset is O(1), so we always do it immediately.
695  InvalidateSolutionSynchronization();
696  CHECK_STATUS(XPRSsetobjoffset(mLp, value));
697 }
698 
699 void XpressInterface::ClearObjective() {
700  InvalidateSolutionSynchronization();
701 
702  // Since the objective function is stored as a dense vector updating
703  // it is O(n), so we usually perform the update immediately.
704  // If we want to support incremental extraction then we have no choice
705  // but to perform the update immediately.
706 
707  if (supportIncrementalExtraction || (slowUpdates & SlowClearObjective)) {
708  int const cols = XPRSgetnumcols(mLp);
709  unique_ptr<int[]> ind(new int[cols]);
710  unique_ptr<double[]> zero(new double[cols]);
711  int j = 0;
712  const auto& coeffs = solver_->objective_->coefficients_;
713  for (auto it(coeffs.begin()); it != coeffs.end(); ++it) {
714  int const idx = it->first->index();
715  // We only need to reset variables that have been extracted.
716  if (variable_is_extracted(idx)) {
717  DCHECK_LT(idx, cols);
718  ind[j] = idx;
719  zero[j] = 0.0;
720  ++j;
721  }
722  }
723  if (j > 0) CHECK_STATUS(XPRSchgobj(mLp, j, ind.get(), zero.get()));
724  CHECK_STATUS(XPRSsetobjoffset(mLp, 0.0));
725  } else {
726  InvalidateModelSynchronization();
727  }
728 }
729 
730 // ------ Query statistics on the solution and the solve ------
731 
732 int64 XpressInterface::iterations() const {
733  if (!CheckSolutionIsSynchronized()) return kUnknownNumberOfIterations;
734  return static_cast<int64>(XPRSgetitcnt(mLp));
735 }
736 
737 int64 XpressInterface::nodes() const {
738  if (mMip) {
739  if (!CheckSolutionIsSynchronized()) return kUnknownNumberOfNodes;
740  return static_cast<int64>(XPRSgetnodecnt(mLp));
741  } else {
742  LOG(DFATAL) << "Number of nodes only available for discrete problems";
743  return kUnknownNumberOfNodes;
744  }
745 }
746 
747 // Transform a XPRESS basis status to an MPSolver basis status.
748 MPSolver::BasisStatus XpressInterface::xformBasisStatus(
749  int xpress_basis_status) {
750  switch (xpress_basis_status) {
751  case XPRS_AT_LOWER:
752  return MPSolver::AT_LOWER_BOUND;
753  case XPRS_BASIC:
754  return MPSolver::BASIC;
755  case XPRS_AT_UPPER:
756  return MPSolver::AT_UPPER_BOUND;
757  case XPRS_FREE_SUPER:
758  return MPSolver::FREE;
759  default:
760  LOG(DFATAL) << "Unknown XPRESS basis status";
761  return MPSolver::FREE;
762  }
763 }
764 
765 // Returns the basis status of a row.
766 MPSolver::BasisStatus XpressInterface::row_status(int constraint_index) const {
767  if (mMip) {
768  LOG(FATAL) << "Basis status only available for continuous problems";
769  return MPSolver::FREE;
770  }
771 
772  if (CheckSolutionIsSynchronized()) {
773  if (!mRstat) {
774  int const rows = XPRSgetnumrows(mLp);
775  unique_ptr<int[]> data(new int[rows]);
776  mRstat.swap(data);
777  CHECK_STATUS(XPRSgetbasis(mLp, 0, mRstat.get()));
778  }
779  } else {
780  mRstat = 0;
781  }
782 
783  if (mRstat) {
784  return xformBasisStatus(mRstat[constraint_index]);
785  } else {
786  LOG(FATAL) << "Row basis status not available";
787  return MPSolver::FREE;
788  }
789 }
790 
791 // Returns the basis status of a column.
792 MPSolver::BasisStatus XpressInterface::column_status(int variable_index) const {
793  if (mMip) {
794  LOG(FATAL) << "Basis status only available for continuous problems";
795  return MPSolver::FREE;
796  }
797 
798  if (CheckSolutionIsSynchronized()) {
799  if (!mCstat) {
800  int const cols = XPRSgetnumcols(mLp);
801  unique_ptr<int[]> data(new int[cols]);
802  mCstat.swap(data);
803  CHECK_STATUS(XPRSgetbasis(mLp, mCstat.get(), 0));
804  }
805  } else {
806  mCstat = 0;
807  }
808 
809  if (mCstat) {
810  return xformBasisStatus(mCstat[variable_index]);
811  } else {
812  LOG(FATAL) << "Column basis status not available";
813  return MPSolver::FREE;
814  }
815 }
816 
817 // Extract all variables that have not yet been extracted.
818 void XpressInterface::ExtractNewVariables() {
819  // NOTE: The code assumes that a linear expression can never contain
820  // non-zero duplicates.
821 
822  InvalidateSolutionSynchronization();
823 
824  if (!supportIncrementalExtraction) {
825  // Without incremental extraction ExtractModel() is always called
826  // to extract the full model.
827  CHECK(last_variable_index_ == 0 ||
828  last_variable_index_ == solver_->variables_.size());
829  CHECK(last_constraint_index_ == 0 ||
830  last_constraint_index_ == solver_->constraints_.size());
831  }
832 
833  int const last_extracted = last_variable_index_;
834  int const var_count = solver_->variables_.size();
835  int newcols = var_count - last_extracted;
836  if (newcols > 0) {
837  // There are non-extracted variables. Extract them now.
838 
839  unique_ptr<double[]> obj(new double[newcols]);
840  unique_ptr<double[]> lb(new double[newcols]);
841  unique_ptr<double[]> ub(new double[newcols]);
842  unique_ptr<char[]> ctype(new char[newcols]);
843  unique_ptr<const char*[]> colname(new const char*[newcols]);
844 
845  bool have_names = false;
846  for (int j = 0, varidx = last_extracted; j < newcols; ++j, ++varidx) {
847  MPVariable const* const var = solver_->variables_[varidx];
848  lb[j] = var->lb();
849  ub[j] = var->ub();
850  ctype[j] = var->integer() ? XPRS_INTEGER : XPRS_CONTINUOUS;
851  colname[j] = var->name().empty() ? 0 : var->name().c_str();
852  have_names = have_names || var->name().empty();
853  obj[j] = solver_->objective_->GetCoefficient(var);
854  }
855 
856  // Arrays for modifying the problem are setup. Update the index
857  // of variables that will get extracted now. Updating indices
858  // _before_ the actual extraction makes things much simpler in
859  // case we support incremental extraction.
860  // In case of error we just reset the indices.
861  std::vector<MPVariable*> const& variables = solver_->variables();
862  for (int j = last_extracted; j < var_count; ++j) {
863  CHECK(!variable_is_extracted(variables[j]->index()));
864  set_variable_as_extracted(variables[j]->index(), true);
865  }
866 
867  try {
868  bool use_newcols = true;
869 
870  if (supportIncrementalExtraction) {
871  // If we support incremental extraction then we must
872  // update existing constraints with the new variables.
873  // To do that we use XPRSaddcols() to actually create the
874  // variables. This is supposed to be faster than combining
875  // XPRSnewcols() and XPRSchgcoeflist().
876 
877  // For each column count the size of the intersection with
878  // existing constraints.
879  unique_ptr<int[]> collen(new int[newcols]);
880  for (int j = 0; j < newcols; ++j) collen[j] = 0;
881  int nonzeros = 0;
882  // TODO: Use a bitarray to flag the constraints that actually
883  // intersect new variables?
884  for (int i = 0; i < last_constraint_index_; ++i) {
885  MPConstraint const* const ct = solver_->constraints_[i];
886  CHECK(constraint_is_extracted(ct->index()));
887  const auto& coeffs = ct->coefficients_;
888  for (auto it(coeffs.begin()); it != coeffs.end(); ++it) {
889  int const idx = it->first->index();
890  if (variable_is_extracted(idx) && idx > last_variable_index_) {
891  collen[idx - last_variable_index_]++;
892  ++nonzeros;
893  }
894  }
895  }
896 
897  if (nonzeros > 0) {
898  // At least one of the new variables did intersect with an
899  // old constraint. We have to create the new columns via
900  // XPRSaddcols().
901  use_newcols = false;
902  unique_ptr<int[]> begin(new int[newcols + 2]);
903  unique_ptr<int[]> cmatind(new int[nonzeros]);
904  unique_ptr<double[]> cmatval(new double[nonzeros]);
905 
906  // Here is how cmatbeg[] is setup:
907  // - it is initialized as
908  // [ 0, 0, collen[0], collen[0]+collen[1], ... ]
909  // so that cmatbeg[j+1] tells us where in cmatind[] and
910  // cmatval[] we need to put the next nonzero for column
911  // j
912  // - after nonzeros have been setup the array looks like
913  // [ 0, collen[0], collen[0]+collen[1], ... ]
914  // so that it is the correct input argument for XPRSaddcols
915  int* cmatbeg = begin.get();
916  cmatbeg[0] = 0;
917  cmatbeg[1] = 0;
918  ++cmatbeg;
919  for (int j = 0; j < newcols; ++j)
920  cmatbeg[j + 1] = cmatbeg[j] + collen[j];
921 
922  for (int i = 0; i < last_constraint_index_; ++i) {
923  MPConstraint const* const ct = solver_->constraints_[i];
924  int const row = ct->index();
925  const auto& coeffs = ct->coefficients_;
926  for (auto it(coeffs.begin()); it != coeffs.end(); ++it) {
927  int const idx = it->first->index();
928  if (variable_is_extracted(idx) && idx > last_variable_index_) {
929  int const nz = cmatbeg[idx]++;
930  cmatind[nz] = row;
931  cmatval[nz] = it->second;
932  }
933  }
934  }
935  --cmatbeg;
936  CHECK_STATUS(XPRSaddcols(mLp, newcols, nonzeros, obj.get(), cmatbeg,
937  cmatind.get(), cmatval.get(), lb.get(),
938  ub.get()));
939  }
940  }
941 
942  if (use_newcols) {
943  // Either incremental extraction is not supported or none of
944  // the new variables did intersect an existing constraint.
945  // We can just use XPRSnewcols() to create the new variables.
946  std::vector<int> collen(newcols, 0);
947  std::vector<int> cmatbeg(newcols, 0);
948  unique_ptr<int[]> cmatind(new int[1]);
949  unique_ptr<double[]> cmatval(new double[1]);
950  cmatind[0] = 0;
951  cmatval[0] = 1.0;
952 
953  CHECK_STATUS(XPRSaddcols(mLp, newcols, 0, obj.get(), cmatbeg.data(),
954  cmatind.get(), cmatval.get(), lb.get(),
955  ub.get()));
956  int const cols = XPRSgetnumcols(mLp);
957  unique_ptr<int[]> ind(new int[newcols]);
958  for (int j = 0; j < cols; ++j) ind[j] = j;
959  CHECK_STATUS(
960  XPRSchgcoltype(mLp, cols - last_extracted, ind.get(), ctype.get()));
961  } else {
962  // Incremental extraction: we must update the ctype of the
963  // newly created variables (XPRSaddcols() does not allow
964  // specifying the ctype)
965  if (mMip && XPRSgetnumcols(mLp) > 0) {
966  // Query the actual number of columns in case we did not
967  // manage to extract all columns.
968  int const cols = XPRSgetnumcols(mLp);
969  unique_ptr<int[]> ind(new int[newcols]);
970  for (int j = last_extracted; j < cols; ++j)
971  ind[j - last_extracted] = j;
972  CHECK_STATUS(XPRSchgcoltype(mLp, cols - last_extracted, ind.get(),
973  ctype.get()));
974  }
975  }
976  } catch (...) {
977  // Undo all changes in case of error.
978  int const cols = XPRSgetnumcols(mLp);
979  if (cols > last_extracted) {
980  std::vector<int> colsToDelete;
981  for (int i = last_extracted; i < cols; ++i) colsToDelete.push_back(i);
982  (void)XPRSdelcols(mLp, colsToDelete.size(), colsToDelete.data());
983  }
984  std::vector<MPVariable*> const& variables = solver_->variables();
985  int const size = variables.size();
986  for (int j = last_extracted; j < size; ++j)
987  set_variable_as_extracted(j, false);
988  throw;
989  }
990  }
991 }
992 
993 // Extract constraints that have not yet been extracted.
994 void XpressInterface::ExtractNewConstraints() {
995  // NOTE: The code assumes that a linear expression can never contain
996  // non-zero duplicates.
997  if (!supportIncrementalExtraction) {
998  // Without incremental extraction ExtractModel() is always called
999  // to extract the full model.
1000  CHECK(last_variable_index_ == 0 ||
1001  last_variable_index_ == solver_->variables_.size());
1002  CHECK(last_constraint_index_ == 0 ||
1003  last_constraint_index_ == solver_->constraints_.size());
1004  }
1005 
1006  int const offset = last_constraint_index_;
1007  int const total = solver_->constraints_.size();
1008 
1009  if (total > offset) {
1010  // There are constraints that are not yet extracted.
1011 
1012  InvalidateSolutionSynchronization();
1013 
1014  int newCons = total - offset;
1015  int const cols = XPRSgetnumcols(mLp);
1016  DCHECK_EQ(last_variable_index_, cols);
1017  int const chunk = newCons; // 10; // max number of rows to add in one shot
1018 
1019  // Update indices of new constraints _before_ actually extracting
1020  // them. In case of error we will just reset the indices.
1021  for (int c = offset; c < total; ++c) set_constraint_as_extracted(c, true);
1022 
1023  try {
1024  unique_ptr<int[]> rmatind(new int[cols]);
1025  unique_ptr<double[]> rmatval(new double[cols]);
1026  unique_ptr<int[]> rmatbeg(new int[chunk]);
1027  unique_ptr<char[]> sense(new char[chunk]);
1028  unique_ptr<double[]> rhs(new double[chunk]);
1029  unique_ptr<char const*[]> name(new char const*[chunk]);
1030  unique_ptr<double[]> rngval(new double[chunk]);
1031  unique_ptr<int[]> rngind(new int[chunk]);
1032  bool haveRanges = false;
1033 
1034  // Loop over the new constraints, collecting rows for up to
1035  // CHUNK constraints into the arrays so that adding constraints
1036  // is faster.
1037  for (int c = 0; c < newCons; /* nothing */) {
1038  // Collect up to CHUNK constraints into the arrays.
1039  int nextRow = 0;
1040  int nextNz = 0;
1041  for (/* nothing */; c < newCons && nextRow < chunk; ++c, ++nextRow) {
1042  MPConstraint const* const ct = solver_->constraints_[offset + c];
1043 
1044  // Stop if there is not enough room in the arrays
1045  // to add the current constraint.
1046  if (nextNz + ct->coefficients_.size() > cols) {
1047  DCHECK_GT(nextRow, 0);
1048  break;
1049  }
1050 
1051  // Setup right-hand side of constraint.
1052  MakeRhs(ct->lb(), ct->ub(), rhs[nextRow], sense[nextRow],
1053  rngval[nextRow]);
1054  haveRanges = haveRanges || (rngval[nextRow] != 0.0);
1055  rngind[nextRow] = offset + c;
1056 
1057  // Setup left-hand side of constraint.
1058  rmatbeg[nextRow] = nextNz;
1059  const auto& coeffs = ct->coefficients_;
1060  for (auto it(coeffs.begin()); it != coeffs.end(); ++it) {
1061  int const idx = it->first->index();
1062  if (variable_is_extracted(idx)) {
1063  DCHECK_LT(nextNz, cols);
1064  DCHECK_LT(idx, cols);
1065  rmatind[nextNz] = idx;
1066  rmatval[nextNz] = it->second;
1067  ++nextNz;
1068  }
1069  }
1070 
1071  // Finally the name of the constraint.
1072  name[nextRow] = ct->name().empty() ? 0 : ct->name().c_str();
1073  }
1074  if (nextRow > 0) {
1075  CHECK_STATUS(XPRSaddrows(mLp, nextRow, nextNz, sense.get(), rhs.get(),
1076  rngval.get(), rmatbeg.get(), rmatind.get(),
1077  rmatval.get()));
1078  if (haveRanges) {
1079  CHECK_STATUS(
1080  XPRSchgrhsrange(mLp, nextRow, rngind.get(), rngval.get()));
1081  }
1082  }
1083  }
1084  } catch (...) {
1085  // Undo all changes in case of error.
1086  int const rows = XPRSgetnumrows(mLp);
1087  std::vector<int> rowsToDelete;
1088  for (int i = offset; i < rows; ++i) rowsToDelete.push_back(i);
1089  if (rows > offset)
1090  (void)XPRSdelrows(mLp, rowsToDelete.size(), rowsToDelete.data());
1091  std::vector<MPConstraint*> const& constraints = solver_->constraints();
1092  int const size = constraints.size();
1093  for (int i = offset; i < size; ++i) set_constraint_as_extracted(i, false);
1094  throw;
1095  }
1096  }
1097 }
1098 
1099 // Extract the objective function.
1100 void XpressInterface::ExtractObjective() {
1101  // NOTE: The code assumes that the objective expression does not contain
1102  // any non-zero duplicates.
1103 
1104  int const cols = XPRSgetnumcols(mLp);
1105  DCHECK_EQ(last_variable_index_, cols);
1106 
1107  unique_ptr<int[]> ind(new int[cols]);
1108  unique_ptr<double[]> val(new double[cols]);
1109  for (int j = 0; j < cols; ++j) {
1110  ind[j] = j;
1111  val[j] = 0.0;
1112  }
1113 
1114  const auto& coeffs = solver_->objective_->coefficients_;
1115  for (auto it = coeffs.begin(); it != coeffs.end(); ++it) {
1116  int const idx = it->first->index();
1117  if (variable_is_extracted(idx)) {
1118  DCHECK_LT(idx, cols);
1119  val[idx] = it->second;
1120  }
1121  }
1122 
1123  CHECK_STATUS(XPRSchgobj(mLp, cols, ind.get(), val.get()));
1124  CHECK_STATUS(XPRSsetobjoffset(mLp, solver_->Objective().offset()));
1125 }
1126 
1127 // ------ Parameters -----
1128 
1129 void XpressInterface::SetParameters(const MPSolverParameters& param) {
1130  SetCommonParameters(param);
1131  if (mMip) SetMIPParameters(param);
1132 }
1133 
1134 void XpressInterface::SetRelativeMipGap(double value) {
1135  if (mMip) {
1136  CHECK_STATUS(XPRSsetdblcontrol(mLp, XPRS_MIPRELSTOP, value));
1137  } else {
1138  LOG(WARNING) << "The relative MIP gap is only available "
1139  << "for discrete problems.";
1140  }
1141 }
1142 
1143 void XpressInterface::SetPrimalTolerance(double value) {
1144  CHECK_STATUS(XPRSsetdblcontrol(mLp, XPRS_FEASTOL, value));
1145 }
1146 
1147 void XpressInterface::SetDualTolerance(double value) {
1148  CHECK_STATUS(XPRSsetdblcontrol(mLp, XPRS_OPTIMALITYTOL, value));
1149 }
1150 
1151 void XpressInterface::SetPresolveMode(int value) {
1152  MPSolverParameters::PresolveValues const presolve =
1153  static_cast<MPSolverParameters::PresolveValues>(value);
1154 
1155  switch (presolve) {
1156  case MPSolverParameters::PRESOLVE_OFF:
1157  CHECK_STATUS(XPRSsetintcontrol(mLp, XPRS_PRESOLVE, 0));
1158  return;
1159  case MPSolverParameters::PRESOLVE_ON:
1160  CHECK_STATUS(XPRSsetintcontrol(mLp, XPRS_PRESOLVE, 1));
1161  return;
1162  }
1163  SetIntegerParamToUnsupportedValue(MPSolverParameters::PRESOLVE, value);
1164 }
1165 
1166 // Sets the scaling mode.
1167 void XpressInterface::SetScalingMode(int value) {
1168  MPSolverParameters::ScalingValues const scaling =
1169  static_cast<MPSolverParameters::ScalingValues>(value);
1170 
1171  switch (scaling) {
1172  case MPSolverParameters::SCALING_OFF:
1173  CHECK_STATUS(XPRSsetintcontrol(mLp, XPRS_SCALING, 0));
1174  break;
1175  case MPSolverParameters::SCALING_ON:
1176  CHECK_STATUS(XPRSsetdefaultcontrol(mLp, XPRS_SCALING));
1177  // In Xpress, scaling is not a binary on/off control, but a bit vector
1178  // control setting it to 1 would only enable bit 1. Instead we reset it to
1179  // its default (163 for the current version 8.6) Alternatively, we could
1180  // call CHECK_STATUS(XPRSsetintcontrol(mLp, XPRS_SCALING, 163));
1181  break;
1182  }
1183 }
1184 
1185 // Sets the LP algorithm : primal, dual or barrier. Note that XPRESS offers
1186 // other LP algorithm (e.g. network) and automatic selection
1187 void XpressInterface::SetLpAlgorithm(int value) {
1188  MPSolverParameters::LpAlgorithmValues const algorithm =
1189  static_cast<MPSolverParameters::LpAlgorithmValues>(value);
1190 
1191  int alg = 1;
1192 
1193  switch (algorithm) {
1194  case MPSolverParameters::DUAL:
1195  alg = 2;
1196  break;
1197  case MPSolverParameters::PRIMAL:
1198  alg = 3;
1199  break;
1200  case MPSolverParameters::BARRIER:
1201  alg = 4;
1202  break;
1203  }
1204 
1205  if (alg == XPRS_DEFAULTALG) {
1206  SetIntegerParamToUnsupportedValue(MPSolverParameters::LP_ALGORITHM, value);
1207  } else {
1208  CHECK_STATUS(XPRSsetintcontrol(mLp, XPRS_DEFAULTALG, alg));
1209  }
1210 }
1211 
1212 bool XpressInterface::ReadParameterFile(std::string const& filename) {
1213  // Return true on success and false on error.
1214  LOG(DFATAL) << "ReadParameterFile not implemented for XPRESS interface";
1215  return false;
1216 }
1217 
1218 std::string XpressInterface::ValidFileExtensionForParameterFile() const {
1219  return ".prm";
1220 }
1221 
1222 MPSolver::ResultStatus XpressInterface::Solve(MPSolverParameters const& param) {
1223  int status;
1224 
1225  // Delete chached information
1226  mCstat = 0;
1227  mRstat = 0;
1228 
1229  WallTimer timer;
1230  timer.Start();
1231 
1232  // Set incrementality
1233  MPSolverParameters::IncrementalityValues const inc =
1234  static_cast<MPSolverParameters::IncrementalityValues>(
1235  param.GetIntegerParam(MPSolverParameters::INCREMENTALITY));
1236  switch (inc) {
1237  case MPSolverParameters::INCREMENTALITY_OFF: {
1238  Reset(); // This should not be required but re-extracting everything
1239  // may be faster, so we do it.
1240  break;
1241  }
1242  case MPSolverParameters::INCREMENTALITY_ON: {
1243  XPRSsetintcontrol(mLp, XPRS_CRASH, 0);
1244  break;
1245  }
1246  }
1247 
1248  // Extract the model to be solved.
1249  // If we don't support incremental extraction and the low-level modeling
1250  // is out of sync then we have to re-extract everything. Note that this
1251  // will lose MIP starts or advanced basis information from a previous
1252  // solve.
1253  if (!supportIncrementalExtraction && sync_status_ == MUST_RELOAD) Reset();
1254  ExtractModel();
1255  VLOG(1) << absl::StrFormat("Model build in %.3f seconds.", timer.Get());
1256 
1257  // Set log level.
1258  XPRSsetintcontrol(mLp, XPRS_OUTPUTLOG, quiet() ? 0 : 1);
1259  // Set parameters.
1260  // NOTE: We must invoke SetSolverSpecificParametersAsString() _first_.
1261  // Its current implementation invokes ReadParameterFile() which in
1262  // turn invokes XPRSreadcopyparam(). The latter will _overwrite_
1263  // all current parameter settings in the environment.
1264  solver_->SetSolverSpecificParametersAsString(
1265  solver_->solver_specific_parameter_string_);
1266  SetParameters(param);
1267  if (solver_->time_limit()) {
1268  VLOG(1) << "Setting time limit = " << solver_->time_limit() << " ms.";
1269  // In Xpress, a time limit should usually have a negative sign. With a
1270  // positive sign, the solver will only stop when a solution has been found.
1271  CHECK_STATUS(XPRSsetdblcontrol(mLp, XPRS_MAXTIME,
1272  -1.0 * solver_->time_limit_in_secs()));
1273  }
1274 
1275  // Solve.
1276  // Do not CHECK_STATUS here since some errors (for example CPXERR_NO_MEMORY)
1277  // still allow us to query useful information.
1278  timer.Restart();
1279 
1280  int xpressstat = 0;
1281  if (mMip) {
1282  if (this->maximize_)
1283  status = XPRSmaxim(mLp, "g");
1284  else
1285  status = XPRSminim(mLp, "g");
1286  XPRSgetintattrib(mLp, XPRS_MIPSTATUS, &xpressstat);
1287  } else {
1288  if (this->maximize_)
1289  status = XPRSmaxim(mLp, "");
1290  else
1291  status = XPRSminim(mLp, "");
1292  XPRSgetintattrib(mLp, XPRS_LPSTATUS, &xpressstat);
1293  }
1294 
1295  // Disable screen output right after solve
1296  XPRSsetintcontrol(mLp, XPRS_OUTPUTLOG, 0);
1297 
1298  if (status) {
1299  VLOG(1) << absl::StrFormat("Failed to optimize MIP. Error %d", status);
1300  // NOTE: We do not return immediately since there may be information
1301  // to grab (for example an incumbent)
1302  } else {
1303  VLOG(1) << absl::StrFormat("Solved in %.3f seconds.", timer.Get());
1304  }
1305 
1306  VLOG(1) << absl::StrFormat("XPRESS solution status %d.", xpressstat);
1307 
1308  // Figure out what solution we have.
1309  bool const feasible = (mMip && (xpressstat == XPRS_MIP_OPTIMAL ||
1310  xpressstat == XPRS_MIP_SOLUTION)) ||
1311  (!mMip && xpressstat == XPRS_LP_OPTIMAL);
1312 
1313  // Get problem dimensions for solution queries below.
1314  int const rows = XPRSgetnumrows(mLp);
1315  int const cols = XPRSgetnumcols(mLp);
1316  DCHECK_EQ(rows, solver_->constraints_.size());
1317  DCHECK_EQ(cols, solver_->variables_.size());
1318 
1319  // Capture objective function value.
1320  objective_value_ = XPRS_NAN;
1321  best_objective_bound_ = XPRS_NAN;
1322  if (feasible) {
1323  if (mMip) {
1324  CHECK_STATUS(XPRSgetdblattrib(mLp, XPRS_MIPOBJVAL, &objective_value_));
1325  CHECK_STATUS(XPRSgetdblattrib(mLp, XPRS_BESTBOUND, &best_objective_bound_));
1326  } else {
1327  CHECK_STATUS(XPRSgetdblattrib(mLp, XPRS_LPOBJVAL, &objective_value_));
1328  }
1329  }
1330  VLOG(1) << "objective=" << objective_value_
1331  << ", bound=" << best_objective_bound_;
1332 
1333  // Capture primal and dual solutions
1334  if (mMip) {
1335  // If there is a primal feasible solution then capture it.
1336  if (feasible) {
1337  if (cols > 0) {
1338  unique_ptr<double[]> x(new double[cols]);
1339  CHECK_STATUS(XPRSgetmipsol(mLp, x.get(), 0));
1340  for (int i = 0; i < solver_->variables_.size(); ++i) {
1341  MPVariable* const var = solver_->variables_[i];
1342  var->set_solution_value(x[i]);
1343  VLOG(3) << var->name() << ": value =" << x[i];
1344  }
1345  }
1346  } else {
1347  for (int i = 0; i < solver_->variables_.size(); ++i)
1348  solver_->variables_[i]->set_solution_value(XPRS_NAN);
1349  }
1350 
1351  // MIP does not have duals
1352  for (int i = 0; i < solver_->variables_.size(); ++i)
1353  solver_->variables_[i]->set_reduced_cost(XPRS_NAN);
1354  for (int i = 0; i < solver_->constraints_.size(); ++i)
1355  solver_->constraints_[i]->set_dual_value(XPRS_NAN);
1356  } else {
1357  // Continuous problem.
1358  if (cols > 0) {
1359  unique_ptr<double[]> x(new double[cols]);
1360  unique_ptr<double[]> dj(new double[cols]);
1361  if (feasible) CHECK_STATUS(XPRSgetlpsol(mLp, x.get(), 0, 0, dj.get()));
1362  for (int i = 0; i < solver_->variables_.size(); ++i) {
1363  MPVariable* const var = solver_->variables_[i];
1364  var->set_solution_value(x[i]);
1365  bool value = false, dual = false;
1366 
1367  if (feasible) {
1368  var->set_solution_value(x[i]);
1369  value = true;
1370  } else {
1371  var->set_solution_value(XPRS_NAN);
1372  }
1373  if (feasible) {
1374  var->set_reduced_cost(dj[i]);
1375  dual = true;
1376  } else {
1377  var->set_reduced_cost(XPRS_NAN);
1378  }
1379  VLOG(3) << var->name() << ":"
1380  << (value ? absl::StrFormat(" value = %f", x[i]) : "")
1381  << (dual ? absl::StrFormat(" reduced cost = %f", dj[i]) : "");
1382  }
1383  }
1384 
1385  if (rows > 0) {
1386  unique_ptr<double[]> pi(new double[rows]);
1387  if (feasible) {
1388  CHECK_STATUS(XPRSgetlpsol(mLp, 0, 0, pi.get(), 0));
1389  }
1390  for (int i = 0; i < solver_->constraints_.size(); ++i) {
1391  MPConstraint* const ct = solver_->constraints_[i];
1392  bool dual = false;
1393  if (feasible) {
1394  ct->set_dual_value(pi[i]);
1395  dual = true;
1396  } else {
1397  ct->set_dual_value(XPRS_NAN);
1398  }
1399  VLOG(4) << "row " << ct->index() << ":"
1400  << (dual ? absl::StrFormat(" dual = %f", pi[i]) : "");
1401  }
1402  }
1403  }
1404 
1405  // Map XPRESS status to more generic solution status in MPSolver
1406  if (mMip) {
1407  switch (xpressstat) {
1408  case XPRS_MIP_OPTIMAL:
1409  result_status_ = MPSolver::OPTIMAL;
1410  break;
1411  case XPRS_MIP_INFEAS:
1412  result_status_ = MPSolver::INFEASIBLE;
1413  break;
1414  case XPRS_MIP_UNBOUNDED:
1415  result_status_ = MPSolver::UNBOUNDED;
1416  break;
1417  default:
1418  result_status_ = feasible ? MPSolver::FEASIBLE : MPSolver::ABNORMAL;
1419  break;
1420  }
1421  } else {
1422  switch (xpressstat) {
1423  case XPRS_LP_OPTIMAL:
1424  result_status_ = MPSolver::OPTIMAL;
1425  break;
1426  case XPRS_LP_INFEAS:
1427  result_status_ = MPSolver::INFEASIBLE;
1428  break;
1429  case XPRS_LP_UNBOUNDED:
1430  result_status_ = MPSolver::UNBOUNDED;
1431  break;
1432  default:
1433  result_status_ = feasible ? MPSolver::FEASIBLE : MPSolver::ABNORMAL;
1434  break;
1435  }
1436  }
1437 
1438  sync_status_ = SOLUTION_SYNCHRONIZED;
1439  return result_status_;
1440 }
1441 
1442 MPSolverInterface* BuildXpressInterface(bool mip, MPSolver* const solver) {
1443  return new XpressInterface(solver, mip);
1444 }
1445 
1446 } // namespace operations_research
1447 #endif // #if defined(USE_XPRESS)
#define CHECK(condition)
Definition: base/logging.h:495
#define DCHECK_LE(val1, val2)
Definition: base/logging.h:887
#define DCHECK_GT(val1, val2)
Definition: base/logging.h:890
#define DCHECK_LT(val1, val2)
Definition: base/logging.h:888
#define LOG(severity)
Definition: base/logging.h:420
#define DCHECK(condition)
Definition: base/logging.h:884
#define DCHECK_EQ(val1, val2)
Definition: base/logging.h:885
#define VLOG(verboselevel)
Definition: base/logging.h:978
void Start()
Definition: timer.h:31
void Restart()
Definition: timer.h:35
double Get() const
Definition: timer.h:45
ResultStatus
The status of solving the problem.
BasisStatus
Advanced usage: possible basis status values for a variable and the slack variable of a linear constr...
const std::string name
const Constraint * ct
int64 value
IntVar * var
Definition: expr_array.cc:1858
int64_t int64
A C++ wrapper that provides a simple and unified interface to several linear programming and mixed in...
const int WARNING
Definition: log_severity.h:31
const int ERROR
Definition: log_severity.h:32
const int FATAL
Definition: log_severity.h:32
ColIndex col
Definition: markowitz.cc:176
RowIndex row
Definition: markowitz.cc:175
Definition: cleanup.h:22
void ShutdownGoogleLogging()
CpSolverResponse Solve(const CpModelProto &model_proto)
Solves the given CpModelProto and returns an instance of CpSolverResponse.
The vehicle routing library lets one model and solve generic vehicle routing problems ranging from th...
int index
Definition: pack.cc:508
int64 coefficient
const bool maximize_
Definition: search.cc:2499