source: branches/pre/include/ClpInterior.hpp @ 214

Last change on this file since 214 was 214, checked in by forrest, 16 years ago

More pdco

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 8.8 KB
Line 
1// Copyright (C) 2003, International Business Machines
2// Corporation and others.  All Rights Reserved.
3
4/*
5   Authors
6   
7   John Tomlin (with some help from John Forrest)
8
9 */
10#ifndef ClpInterior_H
11#define ClpInterior_H
12
13#include <iostream>
14#include <cfloat>
15#include "ClpModel.hpp"
16#include "ClpMatrixBase.hpp"
17#include "ClpSolve.hpp"
18class ClpDualRowPivot;
19class ClpPrimalColumnPivot;
20class ClpFactorization;
21class CoinIndexedVector;
22class ClpNonLinearCost;
23class ClpInteriorProgress;
24// ******** DATA to be moved into protected section of ClpInterior
25typedef struct{
26  double  atolmin;
27  double  r3norm;
28  double  LSdamp;
29  double* deltay;
30} Info;
31
32typedef struct{
33  double  atolold;
34  double  atolnew;
35  double  r3ratio;
36  int   istop;
37  int   itncg;
38} Outfo;
39 
40typedef struct{
41double  gamma;
42double  delta;
43int MaxIter;
44double  FeaTol;
45double  OptTol;
46double  StepTol;
47double  x0min;
48double  z0min;
49double  mu0;
50int   LSmethod;   // 1=Cholesky    2=QR    3=LSQR
51int   LSproblem;  // See below
52int LSQRMaxIter;
53double  LSQRatol1; // Initial  atol
54double  LSQRatol2; // Smallest atol (unless atol1 is smaller)
55double  LSQRconlim;
56int  wait;
57} Options;
58class Lsqr;
59// ***** END
60/** This solves LPs using interior point methods
61
62    It inherits from ClpModel and all its arrays are created at
63    algorithm time.
64
65*/
66
67class ClpInterior : public ClpModel {
68   friend void ClpInteriorUnitTest(const std::string & mpsDir,
69                                  const std::string & netlibDir);
70
71public:
72
73  /**@name Constructors and destructor and copy */
74  //@{
75  /// Default constructor
76    ClpInterior (  );
77
78  /// Copy constructor.
79  ClpInterior(const ClpInterior &);
80  /// Copy constructor from model.
81  ClpInterior(const ClpModel &);
82  /** Subproblem constructor.  A subset of whole model is created from the
83      row and column lists given.  The new order is given by list order and
84      duplicates are allowed.  Name and integer information can be dropped
85  */
86  ClpInterior (const ClpModel * wholeModel,
87              int numberRows, const int * whichRows,
88              int numberColumns, const int * whichColumns,
89              bool dropNames=true, bool dropIntegers=true);
90  /// Assignment operator. This copies the data
91    ClpInterior & operator=(const ClpInterior & rhs);
92  /// Destructor
93   ~ClpInterior (  );
94  // Ones below are just ClpModel with some changes
95  /** Loads a problem (the constraints on the
96        rows are given by lower and upper bounds). If a pointer is 0 then the
97        following values are the default:
98        <ul>
99          <li> <code>colub</code>: all columns have upper bound infinity
100          <li> <code>collb</code>: all columns have lower bound 0
101          <li> <code>rowub</code>: all rows have upper bound infinity
102          <li> <code>rowlb</code>: all rows have lower bound -infinity
103          <li> <code>obj</code>: all variables have 0 objective coefficient
104        </ul>
105    */
106  void loadProblem (  const ClpMatrixBase& matrix,
107                     const double* collb, const double* colub,   
108                     const double* obj,
109                     const double* rowlb, const double* rowub,
110                      const double * rowObjective=NULL);
111  void loadProblem (  const CoinPackedMatrix& matrix,
112                     const double* collb, const double* colub,   
113                     const double* obj,
114                     const double* rowlb, const double* rowub,
115                      const double * rowObjective=NULL);
116
117  /** Just like the other loadProblem() method except that the matrix is
118        given in a standard column major ordered format (without gaps). */
119  void loadProblem (  const int numcols, const int numrows,
120                     const CoinBigIndex* start, const int* index,
121                     const double* value,
122                     const double* collb, const double* colub,   
123                     const double* obj,
124                      const double* rowlb, const double* rowub,
125                      const double * rowObjective=NULL);
126  /// This one is for after presolve to save memory
127  void loadProblem (  const int numcols, const int numrows,
128                     const CoinBigIndex* start, const int* index,
129                      const double* value,const int * length,
130                     const double* collb, const double* colub,   
131                     const double* obj,
132                      const double* rowlb, const double* rowub,
133                      const double * rowObjective=NULL);
134  /// Read an mps file from the given filename
135  int readMps(const char *filename,
136              bool keepNames=false,
137              bool ignoreErrors = false);
138  //@}
139
140  /**@name Functions most useful to user */
141  //@{
142  /** Pdco algorithm - see ClpPdco.hpp for method */
143  int pdco();
144  // ** Temporary version
145  int  pdco( Lsqr *lsqr, Options &options, Info &info, Outfo &outfo);
146  //@}
147
148  /**@name most useful gets and sets */
149  //@{
150  /// If problem is primal feasible
151  inline bool primalFeasible() const
152         { return (sumPrimalInfeasibilities_<=1.0e-5);};
153  /// If problem is dual feasible
154  inline bool dualFeasible() const
155         { return (sumDualInfeasibilities_<=1.0e-5);};
156  /// Current (or last) algorithm
157  inline int algorithm() const 
158  {return algorithm_; } ;
159  /// Set algorithm
160  inline void setAlgorithm(int value)
161  {algorithm_=value; } ;
162  /// Sum of dual infeasibilities
163  inline double sumDualInfeasibilities() const 
164          { return sumDualInfeasibilities_;} ;
165  /// Sum of primal infeasibilities
166  inline double sumPrimalInfeasibilities() const 
167          { return sumPrimalInfeasibilities_;} ;
168  //@}
169
170  /******************** End of most useful part **************/
171  /**@name Functions less likely to be useful to casual user */
172  //@{
173  //@}
174  /**@name Matrix times vector methods
175     They can be faster if scalar is +- 1
176     These are covers so user need not worry about scaling
177     Also for simplex I am not using basic/non-basic split */
178  //@{
179    /** Return <code>y + A * x * scalar</code> in <code>y</code>.
180        @pre <code>x</code> must be of size <code>numColumns()</code>
181        @pre <code>y</code> must be of size <code>numRows()</code> */
182   void times(double scalar,
183                       const double * x, double * y) const;
184    /** Return <code>y + x * scalar * A</code> in <code>y</code>.
185        @pre <code>x</code> must be of size <code>numRows()</code>
186        @pre <code>y</code> must be of size <code>numColumns()</code> */
187    void transposeTimes(double scalar,
188                                const double * x, double * y) const ;
189  //@}
190
191  /**@name most useful gets and sets */
192  //@{
193  /// Largest error on Ax-b
194  inline double largestPrimalError() const
195          { return largestPrimalError_;} ;
196  /// Largest error on basic duals
197  inline double largestDualError() const
198          { return largestDualError_;} ;
199  //@}
200
201  protected:
202  /**@name protected methods */
203  //@{
204  /// Does most of deletion
205  void gutsOfDelete();
206  /// Does most of copying
207  void gutsOfCopy(const ClpInterior & rhs);
208  /// Returns true if data looks okay, false if not
209  bool createWorkingData();
210  void deleteWorkingData();
211  /// Sanity check on input rim data
212  bool sanityCheck();
213  ///  This does housekeeping
214  int housekeeping();
215  //@}
216  public:
217  /**@name public methods */
218  //@{
219  /// Raw objective value (so always minimize)
220  inline double rawObjectiveValue() const
221  { return objectiveValue_;};
222  /// Returns 1 if sequence indicates column
223  inline int isColumn(int sequence) const
224  { return sequence<numberColumns_ ? 1 : 0;};
225  /// Returns sequence number within section
226  inline int sequenceWithin(int sequence) const
227  { return sequence<numberColumns_ ? sequence : sequence-numberColumns_;};
228  //@}
229
230////////////////// data //////////////////
231protected:
232
233  /**@name data.  Many arrays have a row part and a column part.
234   There is a single array with both - columns then rows and
235   then normally two arrays pointing to rows and columns.  The
236   single array is the owner of memory
237  */
238  //@{
239  /// Largest error on Ax-b
240  double largestPrimalError_;
241  /// Largest error on basic duals
242  double largestDualError_;
243  /// Sum of dual infeasibilities
244  double sumDualInfeasibilities_;
245  /// Sum of primal infeasibilities
246  double sumPrimalInfeasibilities_;
247  /// Working copy of lower bounds (Owner of arrays below)
248  double * lower_;
249  /// Row lower bounds - working copy
250  double * rowLowerWork_;
251  /// Column lower bounds - working copy
252  double * columnLowerWork_;
253  /// Working copy of upper bounds (Owner of arrays below)
254  double * upper_;
255  /// Row upper bounds - working copy
256  double * rowUpperWork_;
257  /// Column upper bounds - working copy
258  double * columnUpperWork_;
259  /// Working copy of objective
260  double * cost_;
261  /// Which algorithm being used
262  int algorithm_;
263  //@}
264};
265//#############################################################################
266/** A function that tests the methods in the ClpInterior class. The
267    only reason for it not to be a member method is that this way it doesn't
268    have to be compiled into the library. And that's a gain, because the
269    library should be compiled with optimization on, but this method should be
270    compiled with debugging.
271
272    It also does some testing of ClpFactorization class
273 */
274void
275ClpInteriorUnitTest(const std::string & mpsDir,
276                   const std::string & netlibDir);
277
278
279#endif
Note: See TracBrowser for help on using the repository browser.