source: stable/1.17/Clp/src/ClpInterior.hpp @ 2439

Last change on this file since 2439 was 2385, checked in by unxusr, 10 months ago

formatting

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 17.2 KB
Line 
1/* $Id: ClpInterior.hpp 2385 2019-01-06 19:43:06Z stefan $ */
2// Copyright (C) 2003, International Business Machines
3// Corporation and others.  All Rights Reserved.
4// This code is licensed under the terms of the Eclipse Public License (EPL).
5/*
6   Authors
7
8   John Tomlin (pdco)
9   John Forrest (standard predictor-corrector)
10
11   Note JJF has added arrays - this takes more memory but makes
12   flow easier to understand and hopefully easier to extend
13
14 */
15#ifndef ClpInterior_H
16#define ClpInterior_H
17
18#include <iostream>
19#include <cfloat>
20#include "ClpModel.hpp"
21#include "ClpMatrixBase.hpp"
22#include "ClpSolve.hpp"
23#include "CoinDenseVector.hpp"
24class ClpLsqr;
25class ClpPdcoBase;
26/// ******** DATA to be moved into protected section of ClpInterior
27typedef struct {
28  double atolmin;
29  double r3norm;
30  double LSdamp;
31  double *deltay;
32} Info;
33/// ******** DATA to be moved into protected section of ClpInterior
34
35typedef struct {
36  double atolold;
37  double atolnew;
38  double r3ratio;
39  int istop;
40  int itncg;
41} Outfo;
42/// ******** DATA to be moved into protected section of ClpInterior
43
44typedef struct {
45  double gamma;
46  double delta;
47  int MaxIter;
48  double FeaTol;
49  double OptTol;
50  double StepTol;
51  double x0min;
52  double z0min;
53  double mu0;
54  int LSmethod; // 1=Cholesky    2=QR    3=LSQR
55  int LSproblem; // See below
56  int LSQRMaxIter;
57  double LSQRatol1; // Initial  atol
58  double LSQRatol2; // Smallest atol (unless atol1 is smaller)
59  double LSQRconlim;
60  int wait;
61} Options;
62class Lsqr;
63class ClpCholeskyBase;
64// ***** END
65/** This solves LPs using interior point methods
66
67    It inherits from ClpModel and all its arrays are created at
68    algorithm time.
69
70*/
71
72class ClpInterior : public ClpModel {
73  friend void ClpInteriorUnitTest(const std::string &mpsDir,
74    const std::string &netlibDir);
75
76public:
77  /**@name Constructors and destructor and copy */
78  //@{
79  /// Default constructor
80  ClpInterior();
81
82  /// Copy constructor.
83  ClpInterior(const ClpInterior &);
84  /// Copy constructor from model.
85  ClpInterior(const ClpModel &);
86  /** Subproblem constructor.  A subset of whole model is created from the
87         row and column lists given.  The new order is given by list order and
88         duplicates are allowed.  Name and integer information can be dropped
89     */
90  ClpInterior(const ClpModel *wholeModel,
91    int numberRows, const int *whichRows,
92    int numberColumns, const int *whichColumns,
93    bool dropNames = true, bool dropIntegers = true);
94  /// Assignment operator. This copies the data
95  ClpInterior &operator=(const ClpInterior &rhs);
96  /// Destructor
97  ~ClpInterior();
98  // Ones below are just ClpModel with some changes
99  /** Loads a problem (the constraints on the
100           rows are given by lower and upper bounds). If a pointer is 0 then the
101           following values are the default:
102           <ul>
103             <li> <code>colub</code>: all columns have upper bound infinity
104             <li> <code>collb</code>: all columns have lower bound 0
105             <li> <code>rowub</code>: all rows have upper bound infinity
106             <li> <code>rowlb</code>: all rows have lower bound -infinity
107         <li> <code>obj</code>: all variables have 0 objective coefficient
108           </ul>
109       */
110  void loadProblem(const ClpMatrixBase &matrix,
111    const double *collb, const double *colub,
112    const double *obj,
113    const double *rowlb, const double *rowub,
114    const double *rowObjective = NULL);
115  void loadProblem(const CoinPackedMatrix &matrix,
116    const double *collb, const double *colub,
117    const double *obj,
118    const double *rowlb, const double *rowub,
119    const double *rowObjective = NULL);
120
121  /** Just like the other loadProblem() method except that the matrix is
122       given in a standard column major ordered format (without gaps). */
123  void loadProblem(const int numcols, const int numrows,
124    const CoinBigIndex *start, const int *index,
125    const double *value,
126    const double *collb, const double *colub,
127    const double *obj,
128    const double *rowlb, const double *rowub,
129    const double *rowObjective = NULL);
130  /// This one is for after presolve to save memory
131  void loadProblem(const int numcols, const int numrows,
132    const CoinBigIndex *start, const int *index,
133    const double *value, const int *length,
134    const double *collb, const double *colub,
135    const double *obj,
136    const double *rowlb, const double *rowub,
137    const double *rowObjective = NULL);
138  /// Read an mps file from the given filename
139  int readMps(const char *filename,
140    bool keepNames = false,
141    bool ignoreErrors = false);
142  /** Borrow model.  This is so we dont have to copy large amounts
143         of data around.  It assumes a derived class wants to overwrite
144         an empty model with a real one - while it does an algorithm.
145         This is same as ClpModel one. */
146  void borrowModel(ClpModel &otherModel);
147  /** Return model - updates any scalars */
148  void returnModel(ClpModel &otherModel);
149  //@}
150
151  /**@name Functions most useful to user */
152  //@{
153  /** Pdco algorithm - see ClpPdco.hpp for method */
154  int pdco();
155  // ** Temporary version
156  int pdco(ClpPdcoBase *stuff, Options &options, Info &info, Outfo &outfo);
157  /// Primal-Dual Predictor-Corrector barrier
158  int primalDual();
159  //@}
160
161  /**@name most useful gets and sets */
162  //@{
163  /// If problem is primal feasible
164  inline bool primalFeasible() const
165  {
166    return (sumPrimalInfeasibilities_ <= 1.0e-5);
167  }
168  /// If problem is dual feasible
169  inline bool dualFeasible() const
170  {
171    return (sumDualInfeasibilities_ <= 1.0e-5);
172  }
173  /// Current (or last) algorithm
174  inline int algorithm() const
175  {
176    return algorithm_;
177  }
178  /// Set algorithm
179  inline void setAlgorithm(int value)
180  {
181    algorithm_ = value;
182  }
183  /// Sum of dual infeasibilities
184  inline CoinWorkDouble sumDualInfeasibilities() const
185  {
186    return sumDualInfeasibilities_;
187  }
188  /// Sum of primal infeasibilities
189  inline CoinWorkDouble sumPrimalInfeasibilities() const
190  {
191    return sumPrimalInfeasibilities_;
192  }
193  /// dualObjective.
194  inline CoinWorkDouble dualObjective() const
195  {
196    return dualObjective_;
197  }
198  /// primalObjective.
199  inline CoinWorkDouble primalObjective() const
200  {
201    return primalObjective_;
202  }
203  /// diagonalNorm
204  inline CoinWorkDouble diagonalNorm() const
205  {
206    return diagonalNorm_;
207  }
208  /// linearPerturbation
209  inline CoinWorkDouble linearPerturbation() const
210  {
211    return linearPerturbation_;
212  }
213  inline void setLinearPerturbation(CoinWorkDouble value)
214  {
215    linearPerturbation_ = value;
216  }
217  /// projectionTolerance
218  inline CoinWorkDouble projectionTolerance() const
219  {
220    return projectionTolerance_;
221  }
222  inline void setProjectionTolerance(CoinWorkDouble value)
223  {
224    projectionTolerance_ = value;
225  }
226  /// diagonalPerturbation
227  inline CoinWorkDouble diagonalPerturbation() const
228  {
229    return diagonalPerturbation_;
230  }
231  inline void setDiagonalPerturbation(CoinWorkDouble value)
232  {
233    diagonalPerturbation_ = value;
234  }
235  /// gamma
236  inline CoinWorkDouble gamma() const
237  {
238    return gamma_;
239  }
240  inline void setGamma(CoinWorkDouble value)
241  {
242    gamma_ = value;
243  }
244  /// delta
245  inline CoinWorkDouble delta() const
246  {
247    return delta_;
248  }
249  inline void setDelta(CoinWorkDouble value)
250  {
251    delta_ = value;
252  }
253  /// ComplementarityGap
254  inline CoinWorkDouble complementarityGap() const
255  {
256    return complementarityGap_;
257  }
258  //@}
259
260  /**@name most useful gets and sets */
261  //@{
262  /// Largest error on Ax-b
263  inline CoinWorkDouble largestPrimalError() const
264  {
265    return largestPrimalError_;
266  }
267  /// Largest error on basic duals
268  inline CoinWorkDouble largestDualError() const
269  {
270    return largestDualError_;
271  }
272  /// Maximum iterations
273  inline int maximumBarrierIterations() const
274  {
275    return maximumBarrierIterations_;
276  }
277  inline void setMaximumBarrierIterations(int value)
278  {
279    maximumBarrierIterations_ = value;
280  }
281  /// Set cholesky (and delete present one)
282  void setCholesky(ClpCholeskyBase *cholesky);
283  /// Return number fixed to see if worth presolving
284  int numberFixed() const;
285  /** fix variables interior says should be.  If reallyFix false then just
286         set values to exact bounds */
287  void fixFixed(bool reallyFix = true);
288  /// Primal erturbation vector
289  inline CoinWorkDouble *primalR() const
290  {
291    return primalR_;
292  }
293  /// Dual erturbation vector
294  inline CoinWorkDouble *dualR() const
295  {
296    return dualR_;
297  }
298  //@}
299
300protected:
301  /**@name protected methods */
302  //@{
303  /// Does most of deletion
304  void gutsOfDelete();
305  /// Does most of copying
306  void gutsOfCopy(const ClpInterior &rhs);
307  /// Returns true if data looks okay, false if not
308  bool createWorkingData();
309  void deleteWorkingData();
310  /// Sanity check on input rim data
311  bool sanityCheck();
312  ///  This does housekeeping
313  int housekeeping();
314  //@}
315public:
316  /**@name public methods */
317  //@{
318  /// Raw objective value (so always minimize)
319  inline CoinWorkDouble rawObjectiveValue() const
320  {
321    return objectiveValue_;
322  }
323  /// Returns 1 if sequence indicates column
324  inline int isColumn(int sequence) const
325  {
326    return sequence < numberColumns_ ? 1 : 0;
327  }
328  /// Returns sequence number within section
329  inline int sequenceWithin(int sequence) const
330  {
331    return sequence < numberColumns_ ? sequence : sequence - numberColumns_;
332  }
333  /// Checks solution
334  void checkSolution();
335  /** Modifies djs to allow for quadratic.
336         returns quadratic offset */
337  CoinWorkDouble quadraticDjs(CoinWorkDouble *djRegion, const CoinWorkDouble *solution,
338    CoinWorkDouble scaleFactor);
339
340  /// To say a variable is fixed
341  inline void setFixed(int sequence)
342  {
343    status_[sequence] = static_cast< unsigned char >(status_[sequence] | 1);
344  }
345  inline void clearFixed(int sequence)
346  {
347    status_[sequence] = static_cast< unsigned char >(status_[sequence] & ~1);
348  }
349  inline bool fixed(int sequence) const
350  {
351    return ((status_[sequence] & 1) != 0);
352  }
353
354  /// To flag a variable
355  inline void setFlagged(int sequence)
356  {
357    status_[sequence] = static_cast< unsigned char >(status_[sequence] | 2);
358  }
359  inline void clearFlagged(int sequence)
360  {
361    status_[sequence] = static_cast< unsigned char >(status_[sequence] & ~2);
362  }
363  inline bool flagged(int sequence) const
364  {
365    return ((status_[sequence] & 2) != 0);
366  }
367
368  /// To say a variable is fixed OR free
369  inline void setFixedOrFree(int sequence)
370  {
371    status_[sequence] = static_cast< unsigned char >(status_[sequence] | 4);
372  }
373  inline void clearFixedOrFree(int sequence)
374  {
375    status_[sequence] = static_cast< unsigned char >(status_[sequence] & ~4);
376  }
377  inline bool fixedOrFree(int sequence) const
378  {
379    return ((status_[sequence] & 4) != 0);
380  }
381
382  /// To say a variable has lower bound
383  inline void setLowerBound(int sequence)
384  {
385    status_[sequence] = static_cast< unsigned char >(status_[sequence] | 8);
386  }
387  inline void clearLowerBound(int sequence)
388  {
389    status_[sequence] = static_cast< unsigned char >(status_[sequence] & ~8);
390  }
391  inline bool lowerBound(int sequence) const
392  {
393    return ((status_[sequence] & 8) != 0);
394  }
395
396  /// To say a variable has upper bound
397  inline void setUpperBound(int sequence)
398  {
399    status_[sequence] = static_cast< unsigned char >(status_[sequence] | 16);
400  }
401  inline void clearUpperBound(int sequence)
402  {
403    status_[sequence] = static_cast< unsigned char >(status_[sequence] & ~16);
404  }
405  inline bool upperBound(int sequence) const
406  {
407    return ((status_[sequence] & 16) != 0);
408  }
409
410  /// To say a variable has fake lower bound
411  inline void setFakeLower(int sequence)
412  {
413    status_[sequence] = static_cast< unsigned char >(status_[sequence] | 32);
414  }
415  inline void clearFakeLower(int sequence)
416  {
417    status_[sequence] = static_cast< unsigned char >(status_[sequence] & ~32);
418  }
419  inline bool fakeLower(int sequence) const
420  {
421    return ((status_[sequence] & 32) != 0);
422  }
423
424  /// To say a variable has fake upper bound
425  inline void setFakeUpper(int sequence)
426  {
427    status_[sequence] = static_cast< unsigned char >(status_[sequence] | 64);
428  }
429  inline void clearFakeUpper(int sequence)
430  {
431    status_[sequence] = static_cast< unsigned char >(status_[sequence] & ~64);
432  }
433  inline bool fakeUpper(int sequence) const
434  {
435    return ((status_[sequence] & 64) != 0);
436  }
437  //@}
438
439  ////////////////// data //////////////////
440protected:
441  /**@name data.  Many arrays have a row part and a column part.
442      There is a single array with both - columns then rows and
443      then normally two arrays pointing to rows and columns.  The
444      single array is the owner of memory
445     */
446  //@{
447  /// Largest error on Ax-b
448  CoinWorkDouble largestPrimalError_;
449  /// Largest error on basic duals
450  CoinWorkDouble largestDualError_;
451  /// Sum of dual infeasibilities
452  CoinWorkDouble sumDualInfeasibilities_;
453  /// Sum of primal infeasibilities
454  CoinWorkDouble sumPrimalInfeasibilities_;
455  /// Worst complementarity
456  CoinWorkDouble worstComplementarity_;
457  ///
458public:
459  CoinWorkDouble xsize_;
460  CoinWorkDouble zsize_;
461
462protected:
463  /// Working copy of lower bounds (Owner of arrays below)
464  CoinWorkDouble *lower_;
465  /// Row lower bounds - working copy
466  CoinWorkDouble *rowLowerWork_;
467  /// Column lower bounds - working copy
468  CoinWorkDouble *columnLowerWork_;
469  /// Working copy of upper bounds (Owner of arrays below)
470  CoinWorkDouble *upper_;
471  /// Row upper bounds - working copy
472  CoinWorkDouble *rowUpperWork_;
473  /// Column upper bounds - working copy
474  CoinWorkDouble *columnUpperWork_;
475  /// Working copy of objective
476  CoinWorkDouble *cost_;
477
478public:
479  /// Rhs
480  CoinWorkDouble *rhs_;
481  CoinWorkDouble *x_;
482  CoinWorkDouble *y_;
483  CoinWorkDouble *dj_;
484
485protected:
486  /// Pointer to Lsqr object
487  ClpLsqr *lsqrObject_;
488  /// Pointer to stuff
489  ClpPdcoBase *pdcoStuff_;
490  /// Below here is standard barrier stuff
491  /// mu.
492  CoinWorkDouble mu_;
493  /// objectiveNorm.
494  CoinWorkDouble objectiveNorm_;
495  /// rhsNorm.
496  CoinWorkDouble rhsNorm_;
497  /// solutionNorm.
498  CoinWorkDouble solutionNorm_;
499  /// dualObjective.
500  CoinWorkDouble dualObjective_;
501  /// primalObjective.
502  CoinWorkDouble primalObjective_;
503  /// diagonalNorm.
504  CoinWorkDouble diagonalNorm_;
505  /// stepLength
506  CoinWorkDouble stepLength_;
507  /// linearPerturbation
508  CoinWorkDouble linearPerturbation_;
509  /// diagonalPerturbation
510  CoinWorkDouble diagonalPerturbation_;
511  // gamma from Saunders and Tomlin regularized
512  CoinWorkDouble gamma_;
513  // delta from Saunders and Tomlin regularized
514  CoinWorkDouble delta_;
515  /// targetGap
516  CoinWorkDouble targetGap_;
517  /// projectionTolerance
518  CoinWorkDouble projectionTolerance_;
519  /// maximumRHSError.  maximum Ax
520  CoinWorkDouble maximumRHSError_;
521  /// maximumBoundInfeasibility.
522  CoinWorkDouble maximumBoundInfeasibility_;
523  /// maximumDualError.
524  CoinWorkDouble maximumDualError_;
525  /// diagonalScaleFactor.
526  CoinWorkDouble diagonalScaleFactor_;
527  /// scaleFactor.  For scaling objective
528  CoinWorkDouble scaleFactor_;
529  /// actualPrimalStep
530  CoinWorkDouble actualPrimalStep_;
531  /// actualDualStep
532  CoinWorkDouble actualDualStep_;
533  /// smallestInfeasibility
534  CoinWorkDouble smallestInfeasibility_;
535  /// historyInfeasibility.
536#define LENGTH_HISTORY 5
537  CoinWorkDouble historyInfeasibility_[LENGTH_HISTORY];
538  /// complementarityGap.
539  CoinWorkDouble complementarityGap_;
540  /// baseObjectiveNorm
541  CoinWorkDouble baseObjectiveNorm_;
542  /// worstDirectionAccuracy
543  CoinWorkDouble worstDirectionAccuracy_;
544  /// maximumRHSChange
545  CoinWorkDouble maximumRHSChange_;
546  /// errorRegion. i.e. Ax
547  CoinWorkDouble *errorRegion_;
548  /// rhsFixRegion.
549  CoinWorkDouble *rhsFixRegion_;
550  /// upperSlack
551  CoinWorkDouble *upperSlack_;
552  /// lowerSlack
553  CoinWorkDouble *lowerSlack_;
554  /// diagonal
555  CoinWorkDouble *diagonal_;
556  /// solution
557  CoinWorkDouble *solution_;
558  /// work array
559  CoinWorkDouble *workArray_;
560  /// delta X
561  CoinWorkDouble *deltaX_;
562  /// delta Y
563  CoinWorkDouble *deltaY_;
564  /// deltaZ.
565  CoinWorkDouble *deltaZ_;
566  /// deltaW.
567  CoinWorkDouble *deltaW_;
568  /// deltaS.
569  CoinWorkDouble *deltaSU_;
570  CoinWorkDouble *deltaSL_;
571  /// Primal regularization array
572  CoinWorkDouble *primalR_;
573  /// Dual regularization array
574  CoinWorkDouble *dualR_;
575  /// rhs B
576  CoinWorkDouble *rhsB_;
577  /// rhsU.
578  CoinWorkDouble *rhsU_;
579  /// rhsL.
580  CoinWorkDouble *rhsL_;
581  /// rhsZ.
582  CoinWorkDouble *rhsZ_;
583  /// rhsW.
584  CoinWorkDouble *rhsW_;
585  /// rhs C
586  CoinWorkDouble *rhsC_;
587  /// zVec
588  CoinWorkDouble *zVec_;
589  /// wVec
590  CoinWorkDouble *wVec_;
591  /// cholesky.
592  ClpCholeskyBase *cholesky_;
593  /// numberComplementarityPairs i.e. ones with lower and/or upper bounds (not fixed)
594  int numberComplementarityPairs_;
595  /// numberComplementarityItems_ i.e. number of active bounds
596  int numberComplementarityItems_;
597  /// Maximum iterations
598  int maximumBarrierIterations_;
599  /// gonePrimalFeasible.
600  bool gonePrimalFeasible_;
601  /// goneDualFeasible.
602  bool goneDualFeasible_;
603  /// Which algorithm being used
604  int algorithm_;
605  //@}
606};
607//#############################################################################
608/** A function that tests the methods in the ClpInterior class. The
609    only reason for it not to be a member method is that this way it doesn't
610    have to be compiled into the library. And that's a gain, because the
611    library should be compiled with optimization on, but this method should be
612    compiled with debugging.
613
614    It also does some testing of ClpFactorization class
615 */
616void ClpInteriorUnitTest(const std::string &mpsDir,
617  const std::string &netlibDir);
618
619#endif
620
621/* vi: softtabstop=2 shiftwidth=2 expandtab tabstop=2
622*/
Note: See TracBrowser for help on using the repository browser.