source: stable/1.11/Clp/src/ClpSimplex.hpp @ 1458

Last change on this file since 1458 was 1458, checked in by forrest, 11 years ago

Creating new stable branch 1.11 from trunk (rev 1457)

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 55.2 KB
Line 
1/* $Id: ClpSimplex.hpp 1458 2009-11-05 12:34:07Z forrest $ */
2// Copyright (C) 2002, International Business Machines
3// Corporation and others.  All Rights Reserved.
4
5/*
6   Authors
7 
8   John Forrest
9
10 */
11#ifndef ClpSimplex_H
12#define ClpSimplex_H
13
14#include <iostream>
15#include <cfloat>
16#include "ClpModel.hpp"
17#include "ClpMatrixBase.hpp"
18#include "ClpSolve.hpp"
19class ClpDualRowPivot;
20class ClpPrimalColumnPivot;
21class ClpFactorization;
22#include "ClpFactorization.hpp"
23class CoinIndexedVector;
24class ClpNonLinearCost;
25class ClpNodeStuff;
26class CoinStructuredModel;
27class OsiClpSolverInterface;
28class CoinWarmStartBasis;
29class ClpDisasterHandler;
30class ClpConstraint;
31
32/** This solves LPs using the simplex method
33
34    It inherits from ClpModel and all its arrays are created at
35    algorithm time. Originally I tried to work with model arrays
36    but for simplicity of coding I changed to single arrays with
37    structural variables then row variables.  Some coding is still
38    based on old style and needs cleaning up.
39
40    For a description of algorithms:
41
42    for dual see ClpSimplexDual.hpp and at top of ClpSimplexDual.cpp
43    for primal see ClpSimplexPrimal.hpp and at top of ClpSimplexPrimal.cpp
44
45    There is an algorithm data member.  + for primal variations
46    and - for dual variations
47
48*/
49
50class ClpSimplex : public ClpModel {
51  friend void ClpSimplexUnitTest(const std::string & mpsDir);
52
53public:
54  /** enums for status of various sorts.
55      First 4 match CoinWarmStartBasis,
56      isFixed means fixed at lower bound and out of basis
57  */
58  enum Status {
59    isFree = 0x00,
60    basic = 0x01,
61    atUpperBound = 0x02,
62    atLowerBound = 0x03,
63    superBasic = 0x04,
64    isFixed = 0x05
65  };
66  // For Dual
67  enum FakeBound {
68    noFake = 0x00,
69    lowerFake = 0x01,
70    upperFake = 0x02,
71    bothFake = 0x03
72  };
73
74  /**@name Constructors and destructor and copy */
75  //@{
76  /// Default constructor
77    ClpSimplex (bool emptyMessages = false  );
78
79  /** Copy constructor. May scale depending on mode
80      -1 leave mode as is
81      0 -off, 1 equilibrium, 2 geometric, 3, auto, 4 dynamic(later)
82  */
83  ClpSimplex(const ClpSimplex & rhs, int scalingMode =-1);
84  /** Copy constructor from model. May scale depending on mode
85      -1 leave mode as is
86      0 -off, 1 equilibrium, 2 geometric, 3, auto, 4 dynamic(later)
87  */
88  ClpSimplex(const ClpModel & rhs, int scalingMode=-1);
89  /** Subproblem constructor.  A subset of whole model is created from the
90      row and column lists given.  The new order is given by list order and
91      duplicates are allowed.  Name and integer information can be dropped
92      Can optionally modify rhs to take into account variables NOT in list
93      in this case duplicates are not allowed (also see getbackSolution)
94  */
95  ClpSimplex (const ClpModel * wholeModel,
96              int numberRows, const int * whichRows,
97              int numberColumns, const int * whichColumns,
98              bool dropNames=true, bool dropIntegers=true,
99              bool fixOthers=false);
100  /** Subproblem constructor.  A subset of whole model is created from the
101      row and column lists given.  The new order is given by list order and
102      duplicates are allowed.  Name and integer information can be dropped
103      Can optionally modify rhs to take into account variables NOT in list
104      in this case duplicates are not allowed (also see getbackSolution)
105  */
106  ClpSimplex (const ClpSimplex * wholeModel,
107              int numberRows, const int * whichRows,
108              int numberColumns, const int * whichColumns,
109              bool dropNames=true, bool dropIntegers=true,
110              bool fixOthers=false);
111  /** This constructor modifies original ClpSimplex and stores
112      original stuff in created ClpSimplex.  It is only to be used in
113      conjunction with originalModel */
114  ClpSimplex (ClpSimplex * wholeModel,
115              int numberColumns, const int * whichColumns);
116  /** This copies back stuff from miniModel and then deletes miniModel.
117      Only to be used with mini constructor */
118  void originalModel(ClpSimplex * miniModel);
119  /** Array persistence flag
120      If 0 then as now (delete/new)
121      1 then only do arrays if bigger needed
122      2 as 1 but give a bit extra if bigger needed
123  */
124  void setPersistenceFlag(int value);
125  /// Save a copy of model with certain state - normally without cuts
126  void makeBaseModel();
127  /// Switch off base model
128  void deleteBaseModel();
129  /// See if we have base model
130  inline ClpSimplex *  baseModel() const
131  { return baseModel_;}
132  /** Reset to base model (just size and arrays needed)
133      If model NULL use internal copy
134  */
135  void setToBaseModel(ClpSimplex * model=NULL);
136  /// Assignment operator. This copies the data
137    ClpSimplex & operator=(const ClpSimplex & rhs);
138  /// Destructor
139   ~ClpSimplex (  );
140  // Ones below are just ClpModel with some changes
141  /** Loads a problem (the constraints on the
142        rows are given by lower and upper bounds). If a pointer is 0 then the
143        following values are the default:
144        <ul>
145          <li> <code>colub</code>: all columns have upper bound infinity
146          <li> <code>collb</code>: all columns have lower bound 0
147          <li> <code>rowub</code>: all rows have upper bound infinity
148          <li> <code>rowlb</code>: all rows have lower bound -infinity
149          <li> <code>obj</code>: all variables have 0 objective coefficient
150        </ul>
151    */
152  void loadProblem (  const ClpMatrixBase& matrix,
153                     const double* collb, const double* colub,   
154                     const double* obj,
155                     const double* rowlb, const double* rowub,
156                      const double * rowObjective=NULL);
157  void loadProblem (  const CoinPackedMatrix& matrix,
158                     const double* collb, const double* colub,   
159                     const double* obj,
160                     const double* rowlb, const double* rowub,
161                      const double * rowObjective=NULL);
162
163  /** Just like the other loadProblem() method except that the matrix is
164        given in a standard column major ordered format (without gaps). */
165  void loadProblem (  const int numcols, const int numrows,
166                     const CoinBigIndex* start, const int* index,
167                     const double* value,
168                     const double* collb, const double* colub,   
169                     const double* obj,
170                      const double* rowlb, const double* rowub,
171                      const double * rowObjective=NULL);
172  /// This one is for after presolve to save memory
173  void loadProblem (  const int numcols, const int numrows,
174                     const CoinBigIndex* start, const int* index,
175                      const double* value,const int * length,
176                     const double* collb, const double* colub,   
177                     const double* obj,
178                      const double* rowlb, const double* rowub,
179                      const double * rowObjective=NULL);
180  /** This loads a model from a coinModel object - returns number of errors.
181      If keepSolution true and size is same as current then
182      keeps current status and solution
183  */
184  int loadProblem (  CoinModel & modelObject,bool keepSolution=false);
185  /// Read an mps file from the given filename
186  int readMps(const char *filename,
187              bool keepNames=false,
188              bool ignoreErrors = false);
189  /// Read GMPL files from the given filenames
190  int readGMPL(const char *filename,const char * dataName,
191               bool keepNames=false);
192  /// Read file in LP format from file with name filename.
193  /// See class CoinLpIO for description of this format.
194  int readLp(const char *filename, const double epsilon = 1e-5);
195  /** Borrow model.  This is so we dont have to copy large amounts
196      of data around.  It assumes a derived class wants to overwrite
197      an empty model with a real one - while it does an algorithm.
198      This is same as ClpModel one, but sets scaling on etc. */
199  void borrowModel(ClpModel & otherModel);
200  void borrowModel(ClpSimplex & otherModel);
201   /// Pass in Event handler (cloned and deleted at end)
202   void passInEventHandler(const ClpEventHandler * eventHandler);
203  /// Puts solution back into small model
204  void getbackSolution(const ClpSimplex & smallModel,const int * whichRow, const int * whichColumn);
205  /** Load nonlinear part of problem from AMPL info
206      Returns 0 if linear
207      1 if quadratic objective
208      2 if quadratic constraints
209      3 if nonlinear objective
210      4 if nonlinear constraints
211      -1 on failure
212  */
213  int loadNonLinear(void * info, int & numberConstraints, 
214                    ClpConstraint ** & constraints);
215  //@}
216
217  /**@name Functions most useful to user */
218  //@{
219  /** General solve algorithm which can do presolve.
220      See  ClpSolve.hpp for options
221   */
222  int initialSolve(ClpSolve & options);
223  /// Default initial solve
224  int initialSolve();
225  /// Dual initial solve
226  int initialDualSolve();
227  /// Primal initial solve
228  int initialPrimalSolve();
229 /// Barrier initial solve
230  int initialBarrierSolve();
231  /// Barrier initial solve, not to be followed by crossover
232  int initialBarrierNoCrossSolve();
233  /** Dual algorithm - see ClpSimplexDual.hpp for method.
234      ifValuesPass==2 just does values pass and then stops.
235
236      startFinishOptions - bits
237      1 - do not delete work areas and factorization at end
238      2 - use old factorization if same number of rows
239      4 - skip as much initialization of work areas as possible
240          (based on whatsChanged in clpmodel.hpp) ** work in progress
241      maybe other bits later
242  */
243  int dual(int ifValuesPass=0, int startFinishOptions=0);
244  // If using Debug
245  int dualDebug(int ifValuesPass=0, int startFinishOptions=0);
246  /** Primal algorithm - see ClpSimplexPrimal.hpp for method.
247      ifValuesPass==2 just does values pass and then stops.
248
249      startFinishOptions - bits
250      1 - do not delete work areas and factorization at end
251      2 - use old factorization if same number of rows
252      4 - skip as much initialization of work areas as possible
253          (based on whatsChanged in clpmodel.hpp) ** work in progress
254      maybe other bits later
255  */
256  int primal(int ifValuesPass=0, int startFinishOptions=0);
257  /** Solves nonlinear problem using SLP - may be used as crash
258      for other algorithms when number of iterations small.
259      Also exits if all problematical variables are changing
260      less than deltaTolerance
261  */
262  int nonlinearSLP(int numberPasses,double deltaTolerance);
263  /** Solves problem with nonlinear constraints using SLP - may be used as crash
264      for other algorithms when number of iterations small.
265      Also exits if all problematical variables are changing
266      less than deltaTolerance
267  */
268  int nonlinearSLP(int numberConstraints, ClpConstraint ** constraints,
269                   int numberPasses,double deltaTolerance);
270  /** Solves using barrier (assumes you have good cholesky factor code).
271      Does crossover to simplex if asked*/
272  int barrier(bool crossover=true);
273  /** Solves non-linear using reduced gradient.  Phase = 0 get feasible,
274      =1 use solution */
275  int reducedGradient(int phase=0);
276  /// Solve using structure of model and maybe in parallel
277  int solve(CoinStructuredModel * model);
278  /** This loads a model from a CoinStructuredModel object - returns number of errors.
279      If originalOrder then keep to order stored in blocks,
280      otherwise first column/rows correspond to first block - etc.
281      If keepSolution true and size is same as current then
282      keeps current status and solution
283  */
284  int loadProblem (  CoinStructuredModel & modelObject,
285                     bool originalOrder=true,bool keepSolution=false);
286  /**
287     When scaling is on it is possible that the scaled problem
288     is feasible but the unscaled is not.  Clp returns a secondary
289     status code to that effect.  This option allows for a cleanup.
290     If you use it I would suggest 1.
291     This only affects actions when scaled optimal
292     0 - no action
293     1 - clean up using dual if primal infeasibility
294     2 - clean up using dual if dual infeasibility
295     3 - clean up using dual if primal or dual infeasibility
296     11,12,13 - as 1,2,3 but use primal
297
298     return code as dual/primal
299  */
300  int cleanup(int cleanupScaling);
301  /** Dual ranging.
302      This computes increase/decrease in cost for each given variable and corresponding
303      sequence numbers which would change basis.  Sequence numbers are 0..numberColumns
304      and numberColumns.. for artificials/slacks.
305      For non-basic variables the information is trivial to compute and the change in cost is just minus the
306      reduced cost and the sequence number will be that of the non-basic variables.
307      For basic variables a ratio test is between the reduced costs for non-basic variables
308      and the row of the tableau corresponding to the basic variable.
309      The increase/decrease value is always >= 0.0
310
311      Up to user to provide correct length arrays where each array is of length numberCheck.
312      which contains list of variables for which information is desired.  All other
313      arrays will be filled in by function.  If fifth entry in which is variable 7 then fifth entry in output arrays
314      will be information for variable 7.
315
316      If valueIncrease/Decrease not NULL (both must be NULL or both non NULL) then these are filled with
317      the value of variable if such a change in cost were made (the existing bounds are ignored)
318
319      Returns non-zero if infeasible unbounded etc
320  */
321  int dualRanging(int numberCheck,const int * which,
322                  double * costIncrease, int * sequenceIncrease,
323                  double * costDecrease, int * sequenceDecrease,
324                  double * valueIncrease=NULL, double * valueDecrease=NULL);
325  /** Primal ranging.
326      This computes increase/decrease in value for each given variable and corresponding
327      sequence numbers which would change basis.  Sequence numbers are 0..numberColumns
328      and numberColumns.. for artificials/slacks.
329      This should only be used for non-basic variabls as otherwise information is pretty useless
330      For basic variables the sequence number will be that of the basic variables.
331
332      Up to user to provide correct length arrays where each array is of length numberCheck.
333      which contains list of variables for which information is desired.  All other
334      arrays will be filled in by function.  If fifth entry in which is variable 7 then fifth entry in output arrays
335      will be information for variable 7.
336
337      Returns non-zero if infeasible unbounded etc
338  */
339  int primalRanging(int numberCheck,const int * which,
340                  double * valueIncrease, int * sequenceIncrease,
341                  double * valueDecrease, int * sequenceDecrease);
342  /** Write the basis in MPS format to the specified file.
343      If writeValues true writes values of structurals
344      (and adds VALUES to end of NAME card)
345     
346      Row and column names may be null.
347      formatType is
348      <ul>
349      <li> 0 - normal
350      <li> 1 - extra accuracy
351      <li> 2 - IEEE hex (later)
352      </ul>
353     
354      Returns non-zero on I/O error
355  */
356  int writeBasis(const char *filename,
357                 bool writeValues=false,
358                 int formatType=0) const;
359  /** Read a basis from the given filename,
360      returns -1 on file error, 0 if no values, 1 if values */
361  int readBasis(const char *filename);
362  /// Returns a basis (to be deleted by user)
363  CoinWarmStartBasis * getBasis() const;
364  /// Passes in factorization
365  void setFactorization( ClpFactorization & factorization);
366  // Swaps factorization
367  ClpFactorization * swapFactorization( ClpFactorization * factorization);
368  /// Copies in factorization to existing one
369  void copyFactorization( ClpFactorization & factorization);
370  /** Tightens primal bounds to make dual faster.  Unless
371      fixed or doTight>10, bounds are slightly looser than they could be.
372      This is to make dual go faster and is probably not needed
373      with a presolve.  Returns non-zero if problem infeasible.
374
375      Fudge for branch and bound - put bounds on columns of factor *
376      largest value (at continuous) - should improve stability
377      in branch and bound on infeasible branches (0.0 is off)
378  */
379  int tightenPrimalBounds(double factor=0.0,int doTight=0,bool tightIntegers=false);
380  /** Crash - at present just aimed at dual, returns
381      -2 if dual preferred and crash basis created
382      -1 if dual preferred and all slack basis preferred
383       0 if basis going in was not all slack
384       1 if primal preferred and all slack basis preferred
385       2 if primal preferred and crash basis created.
386       
387       if gap between bounds <="gap" variables can be flipped
388       ( If pivot -1 then can be made super basic!)
389
390       If "pivot" is
391       -1 No pivoting - always primal
392       0 No pivoting (so will just be choice of algorithm)
393       1 Simple pivoting e.g. gub
394       2 Mini iterations
395  */
396  int crash(double gap,int pivot);
397  /// Sets row pivot choice algorithm in dual
398  void setDualRowPivotAlgorithm(ClpDualRowPivot & choice);
399  /// Sets column pivot choice algorithm in primal
400  void setPrimalColumnPivotAlgorithm(ClpPrimalColumnPivot & choice);
401  /** For strong branching.  On input lower and upper are new bounds
402      while on output they are change in objective function values
403      (>1.0e50 infeasible).
404      Return code is 0 if nothing interesting, -1 if infeasible both
405      ways and +1 if infeasible one way (check values to see which one(s))
406      Solutions are filled in as well - even down, odd up - also
407      status and number of iterations
408  */
409  int strongBranching(int numberVariables,const int * variables,
410                      double * newLower, double * newUpper,
411                      double ** outputSolution,
412                      int * outputStatus, int * outputIterations,
413                      bool stopOnFirstInfeasible=true,
414                      bool alwaysFinish=false,
415                      int startFinishOptions=0);
416  /// Fathom - 1 if solution
417  int fathom(void * stuff);
418  /** Do up to N deep - returns
419      -1 - no solution nNodes_ valid nodes
420      >= if solution and that node gives solution
421      ClpNode array is 2**N long.  Values for N and
422      array are in stuff (nNodes_ also in stuff) */
423  int fathomMany(void * stuff);
424  /// Double checks OK
425  double doubleCheck();
426  /// Starts Fast dual2
427  int startFastDual2(ClpNodeStuff * stuff);
428  /// Like Fast dual
429  int fastDual2(ClpNodeStuff * stuff);
430  /// Stops Fast dual2
431  void stopFastDual2(ClpNodeStuff * stuff);
432  /** Deals with crunch aspects
433      mode 0 - in
434           1 - out with solution
435           2 - out without solution
436      returns small model or NULL
437  */
438  ClpSimplex * fastCrunch(ClpNodeStuff * stuff, int mode);
439  //@}
440
441  /**@name Needed for functionality of OsiSimplexInterface */
442  //@{
443  /** Pivot in a variable and out a variable.  Returns 0 if okay,
444      1 if inaccuracy forced re-factorization, -1 if would be singular.
445      Also updates primal/dual infeasibilities.
446      Assumes sequenceIn_ and pivotRow_ set and also directionIn and Out.
447  */
448  int pivot();
449
450  /** Pivot in a variable and choose an outgoing one.  Assumes primal
451      feasible - will not go through a bound.  Returns step length in theta
452      Returns ray in ray_ (or NULL if no pivot)
453      Return codes as before but -1 means no acceptable pivot
454  */
455  int primalPivotResult();
456 
457  /** Pivot out a variable and choose an incoing one.  Assumes dual
458      feasible - will not go through a reduced cost. 
459      Returns step length in theta
460      Returns ray in ray_ (or NULL if no pivot)
461      Return codes as before but -1 means no acceptable pivot
462  */
463  int dualPivotResult();
464
465  /** Common bits of coding for dual and primal.  Return 0 if okay,
466      1 if bad matrix, 2 if very bad factorization
467
468      startFinishOptions - bits
469      1 - do not delete work areas and factorization at end
470      2 - use old factorization if same number of rows
471      4 - skip as much initialization of work areas as possible
472          (based on whatsChanged in clpmodel.hpp) ** work in progress
473      maybe other bits later
474     
475  */
476  int startup(int ifValuesPass,int startFinishOptions=0);
477  void finish(int startFinishOptions=0);
478 
479  /** Factorizes and returns true if optimal.  Used by user */
480  bool statusOfProblem(bool initial=false);
481  /// If user left factorization frequency then compute
482  void defaultFactorizationFrequency();
483  //@}
484
485  /**@name most useful gets and sets */
486  //@{
487  /// If problem is primal feasible
488  inline bool primalFeasible() const
489         { return (numberPrimalInfeasibilities_==0);}
490  /// If problem is dual feasible
491  inline bool dualFeasible() const
492         { return (numberDualInfeasibilities_==0);}
493  /// factorization
494  inline ClpFactorization * factorization() const 
495          { return factorization_;}
496  /// Sparsity on or off
497  bool sparseFactorization() const;
498  void setSparseFactorization(bool value);
499  /// Factorization frequency
500  int factorizationFrequency() const;
501  void setFactorizationFrequency(int value);
502  /// Dual bound
503  inline double dualBound() const
504          { return dualBound_;}
505  void setDualBound(double value);
506  /// Infeasibility cost
507  inline double infeasibilityCost() const
508          { return infeasibilityCost_;}
509  void setInfeasibilityCost(double value);
510  /** Amount of print out:
511      0 - none
512      1 - just final
513      2 - just factorizations
514      3 - as 2 plus a bit more
515      4 - verbose
516      above that 8,16,32 etc just for selective debug
517  */
518  /** Perturbation:
519      50  - switch on perturbation
520      100 - auto perturb if takes too long (1.0e-6 largest nonzero)
521      101 - we are perturbed
522      102 - don't try perturbing again
523      default is 100
524      others are for playing
525  */
526  inline int perturbation() const
527    { return perturbation_;}
528  void setPerturbation(int value);
529  /// Current (or last) algorithm
530  inline int algorithm() const 
531  {return algorithm_; } 
532  /// Set algorithm
533  inline void setAlgorithm(int value)
534  {algorithm_=value; } 
535  /// Return true if the objective limit test can be relied upon
536  bool isObjectiveLimitTestValid() const ;
537  /// Sum of dual infeasibilities
538  inline double sumDualInfeasibilities() const 
539          { return sumDualInfeasibilities_;} 
540  inline void setSumDualInfeasibilities(double value)
541          { sumDualInfeasibilities_=value;} 
542  /// Sum of relaxed dual infeasibilities
543  inline double sumOfRelaxedDualInfeasibilities() const 
544          { return sumOfRelaxedDualInfeasibilities_;} 
545  inline void setSumOfRelaxedDualInfeasibilities(double value)
546          { sumOfRelaxedDualInfeasibilities_=value;} 
547  /// Number of dual infeasibilities
548  inline int numberDualInfeasibilities() const 
549          { return numberDualInfeasibilities_;} 
550  inline void setNumberDualInfeasibilities(int value)
551          { numberDualInfeasibilities_=value;} 
552  /// Number of dual infeasibilities (without free)
553  inline int numberDualInfeasibilitiesWithoutFree() const 
554          { return numberDualInfeasibilitiesWithoutFree_;} 
555  /// Sum of primal infeasibilities
556  inline double sumPrimalInfeasibilities() const 
557          { return sumPrimalInfeasibilities_;} 
558  inline void setSumPrimalInfeasibilities(double value)
559          { sumPrimalInfeasibilities_=value;} 
560  /// Sum of relaxed primal infeasibilities
561  inline double sumOfRelaxedPrimalInfeasibilities() const 
562          { return sumOfRelaxedPrimalInfeasibilities_;} 
563  inline void setSumOfRelaxedPrimalInfeasibilities(double value)
564          { sumOfRelaxedPrimalInfeasibilities_=value;} 
565  /// Number of primal infeasibilities
566  inline int numberPrimalInfeasibilities() const 
567          { return numberPrimalInfeasibilities_;} 
568  inline void setNumberPrimalInfeasibilities(int value)
569          { numberPrimalInfeasibilities_=value;} 
570  /** Save model to file, returns 0 if success.  This is designed for
571      use outside algorithms so does not save iterating arrays etc.
572  It does not save any messaging information.
573  Does not save scaling values.
574  It does not know about all types of virtual functions.
575  */
576  int saveModel(const char * fileName);
577  /** Restore model from file, returns 0 if success,
578      deletes current model */
579  int restoreModel(const char * fileName);
580 
581  /** Just check solution (for external use) - sets sum of
582      infeasibilities etc.
583      If setToBounds 0 then primal column values not changed
584      and used to compute primal row activity values.  If 1 or 2
585      then status used - so all nonbasic variables set to
586      indicated bound and if any values changed (or ==2)  basic values re-computed.
587  */
588  void checkSolution(int setToBounds=0);
589  /** Just check solution (for internal use) - sets sum of
590      infeasibilities etc. */
591  void checkSolutionInternal();
592  /// Useful row length arrays (0,1,2,3,4,5)
593  inline CoinIndexedVector * rowArray(int index) const
594  { return rowArray_[index];}
595  /// Useful column length arrays (0,1,2,3,4,5)
596  inline CoinIndexedVector * columnArray(int index) const
597  { return columnArray_[index];}
598  //@}
599
600  /******************** End of most useful part **************/
601  /**@name Functions less likely to be useful to casual user */
602  //@{
603  /** Given an existing factorization computes and checks
604      primal and dual solutions.  Uses input arrays for variables at
605      bounds.  Returns feasibility states */
606  int getSolution (  const double * rowActivities,
607                     const double * columnActivities);
608  /** Given an existing factorization computes and checks
609      primal and dual solutions.  Uses current problem arrays for
610      bounds.  Returns feasibility states */
611  int getSolution ();
612  /** Constructs a non linear cost from list of non-linearities (columns only)
613      First lower of each column is taken as real lower
614      Last lower is taken as real upper and cost ignored
615
616      Returns nonzero if bad data e.g. lowers not monotonic
617  */
618  int createPiecewiseLinearCosts(const int * starts,
619                   const double * lower, const double * gradient);
620  /// dual row pivot choice
621  ClpDualRowPivot * dualRowPivot() const
622  { return dualRowPivot_;}
623  /// Returns true if model looks OK
624  inline bool goodAccuracy() const
625  { return (largestPrimalError_<1.0e-7&&largestDualError_<1.0e-7);}
626  /** Return model - updates any scalars */
627  void returnModel(ClpSimplex & otherModel);
628  /** Factorizes using current basis. 
629      solveType - 1 iterating, 0 initial, -1 external
630      If 10 added then in primal values pass
631      Return codes are as from ClpFactorization unless initial factorization
632      when total number of singularities is returned.
633      Special case is numberRows_+1 -> all slack basis.
634  */
635  int internalFactorize(int solveType);
636  /// Save data
637  ClpDataSave saveData() ;
638  /// Restore data
639  void restoreData(ClpDataSave saved);
640  /// Clean up status
641  void cleanStatus();
642  /// Factorizes using current basis. For external use
643  int factorize();
644  /** Computes duals from scratch. If givenDjs then
645      allows for nonzero basic djs */
646  void computeDuals(double * givenDjs);
647  /// Computes primals from scratch
648  void computePrimals (  const double * rowActivities,
649                     const double * columnActivities);
650  /** Adds multiple of a column into an array */
651  void add(double * array,
652                   int column, double multiplier) const;
653  /**
654     Unpacks one column of the matrix into indexed array
655     Uses sequenceIn_
656     Also applies scaling if needed
657  */
658  void unpack(CoinIndexedVector * rowArray) const ;
659  /**
660     Unpacks one column of the matrix into indexed array
661     Slack if sequence>= numberColumns
662     Also applies scaling if needed
663  */
664  void unpack(CoinIndexedVector * rowArray,int sequence) const;
665  /**
666     Unpacks one column of the matrix into indexed array
667     ** as packed vector
668     Uses sequenceIn_
669     Also applies scaling if needed
670  */
671  void unpackPacked(CoinIndexedVector * rowArray) ;
672  /**
673     Unpacks one column of the matrix into indexed array
674     ** as packed vector
675     Slack if sequence>= numberColumns
676     Also applies scaling if needed
677  */
678  void unpackPacked(CoinIndexedVector * rowArray,int sequence);
679protected: 
680  /**
681      This does basis housekeeping and does values for in/out variables.
682      Can also decide to re-factorize
683  */
684  int housekeeping(double objectiveChange);
685  /** This sets largest infeasibility and most infeasible and sum
686      and number of infeasibilities (Primal) */
687  void checkPrimalSolution(const double * rowActivities=NULL,
688                           const double * columnActivies=NULL);
689  /** This sets largest infeasibility and most infeasible and sum
690      and number of infeasibilities (Dual) */
691  void checkDualSolution();
692  /** This sets sum and number of infeasibilities (Dual and Primal) */
693  void checkBothSolutions();
694  /**  If input negative scales objective so maximum <= -value
695       and returns scale factor used.  If positive unscales and also
696       redoes dual stuff
697  */
698  double scaleObjective(double value);
699  /// Solve using Dantzig-Wolfe decomposition and maybe in parallel
700  int solveDW(CoinStructuredModel * model);
701  /// Solve using Benders decomposition and maybe in parallel
702  int solveBenders(CoinStructuredModel * model);
703public:
704  /** For advanced use.  When doing iterative solves things can get
705      nasty so on values pass if incoming solution has largest
706      infeasibility < incomingInfeasibility throw out variables
707      from basis until largest infeasibility < allowedInfeasibility
708      or incoming largest infeasibility.
709      If allowedInfeasibility>= incomingInfeasibility this is
710      always possible altough you may end up with an all slack basis.
711
712      Defaults are 1.0,10.0
713  */
714  void setValuesPassAction(double incomingInfeasibility,
715                           double allowedInfeasibility);
716  //@}
717  /**@name most useful gets and sets */
718  //@{
719public: 
720  /// Initial value for alpha accuracy calculation (-1.0 off)
721  inline double alphaAccuracy() const
722          { return alphaAccuracy_;} 
723  inline void setAlphaAccuracy(double value)
724          { alphaAccuracy_ = value;} 
725public:
726   /// Objective value
727   //inline double objectiveValue() const {
728  //return (objectiveValue_-bestPossibleImprovement_)*optimizationDirection_ - dblParam_[ClpObjOffset];
729  //}
730  /// Set disaster handler
731  inline void setDisasterHandler(ClpDisasterHandler * handler)
732  { disasterArea_= handler;}
733  /// Get disaster handler
734  inline ClpDisasterHandler * disasterHandler() const
735  { return disasterArea_;}
736  /// Large bound value (for complementarity etc)
737  inline double largeValue() const 
738          { return largeValue_;} 
739  void setLargeValue( double value) ;
740  /// Largest error on Ax-b
741  inline double largestPrimalError() const
742          { return largestPrimalError_;} 
743  /// Largest error on basic duals
744  inline double largestDualError() const
745          { return largestDualError_;} 
746  /// Largest error on Ax-b
747  inline void setLargestPrimalError(double value)
748          { largestPrimalError_=value;} 
749  /// Largest error on basic duals
750  inline void setLargestDualError(double value)
751          { largestDualError_=value;} 
752  /// Get zero tolerance
753  inline double zeroTolerance() const 
754  { return zeroTolerance_;/*factorization_->zeroTolerance();*/} 
755  /// Set zero tolerance
756  inline void setZeroTolerance( double value)
757  { zeroTolerance_ = value;}
758  /// Basic variables pivoting on which rows
759  inline int * pivotVariable() const
760          { return pivotVariable_;}
761  /// If automatic scaling on
762  inline bool automaticScaling() const
763  { return automaticScale_!=0;}
764  inline void setAutomaticScaling(bool onOff)
765  { automaticScale_ = onOff ? 1: 0;} 
766  /// Current dual tolerance
767  inline double currentDualTolerance() const 
768          { return dualTolerance_;} 
769  inline void setCurrentDualTolerance(double value)
770          { dualTolerance_ = value;} 
771  /// Current primal tolerance
772  inline double currentPrimalTolerance() const 
773          { return primalTolerance_;} 
774  inline void setCurrentPrimalTolerance(double value)
775          { primalTolerance_ = value;} 
776  /// How many iterative refinements to do
777  inline int numberRefinements() const 
778          { return numberRefinements_;} 
779  void setNumberRefinements( int value) ;
780  /// Alpha (pivot element) for use by classes e.g. steepestedge
781  inline double alpha() const { return alpha_;}
782  inline void setAlpha(double value) { alpha_ = value;}
783  /// Reduced cost of last incoming for use by classes e.g. steepestedge
784  inline double dualIn() const { return dualIn_;}
785  /// Pivot Row for use by classes e.g. steepestedge
786  inline int pivotRow() const{ return pivotRow_;}
787  inline void setPivotRow(int value) { pivotRow_=value;}
788  /// value of incoming variable (in Dual)
789  double valueIncomingDual() const;
790  //@}
791
792  protected:
793  /**@name protected methods */
794  //@{
795  /** May change basis and then returns number changed.
796      Computation of solutions may be overriden by given pi and solution
797  */
798  int gutsOfSolution ( double * givenDuals,
799                       const double * givenPrimals,
800                       bool valuesPass=false);
801  /// Does most of deletion (0 = all, 1 = most, 2 most + factorization)
802  void gutsOfDelete(int type);
803  /// Does most of copying
804  void gutsOfCopy(const ClpSimplex & rhs);
805  /** puts in format I like (rowLower,rowUpper) also see StandardMatrix
806      1 bit does rows (now and columns), (2 bit does column bounds), 4 bit does objective(s).
807      8 bit does solution scaling in
808      16 bit does rowArray and columnArray indexed vectors
809      and makes row copy if wanted, also sets columnStart_ etc
810      Also creates scaling arrays if needed.  It does scaling if needed.
811      16 also moves solutions etc in to work arrays
812      On 16 returns false if problem "bad" i.e. matrix or bounds bad
813      If startFinishOptions is -1 then called by user in getSolution
814      so do arrays but keep pivotVariable_
815  */
816  bool createRim(int what,bool makeRowCopy=false,int startFinishOptions=0);
817  /// Does rows and columns
818  void createRim1(bool initial);
819  /// Does objective
820  void createRim4(bool initial);
821  /// Does rows and columns and objective
822  void createRim5(bool initial);
823  /** releases above arrays and does solution scaling out.  May also
824      get rid of factorization data -
825      0 get rid of nothing, 1 get rid of arrays, 2 also factorization
826  */
827  void deleteRim(int getRidOfFactorizationData=2);
828  /// Sanity check on input rim data (after scaling) - returns true if okay
829  bool sanityCheck();
830  //@}
831  public:
832  /**@name public methods */
833  //@{
834  /** Return row or column sections - not as much needed as it
835      once was.  These just map into single arrays */
836  inline double * solutionRegion(int section) const
837  { if (!section) return rowActivityWork_; else return columnActivityWork_;}
838  inline double * djRegion(int section) const
839  { if (!section) return rowReducedCost_; else return reducedCostWork_;}
840  inline double * lowerRegion(int section) const
841  { if (!section) return rowLowerWork_; else return columnLowerWork_;}
842  inline double * upperRegion(int section) const
843  { if (!section) return rowUpperWork_; else return columnUpperWork_;}
844  inline double * costRegion(int section) const
845  { if (!section) return rowObjectiveWork_; else return objectiveWork_;}
846  /// Return region as single array
847  inline double * solutionRegion() const
848  { return solution_;}
849  inline double * djRegion() const
850  { return dj_;}
851  inline double * lowerRegion() const
852  { return lower_;}
853  inline double * upperRegion() const
854  { return upper_;}
855  inline double * costRegion() const
856  { return cost_;}
857  inline Status getStatus(int sequence) const
858  {return static_cast<Status> (status_[sequence]&7);}
859  inline void setStatus(int sequence, Status status)
860  {
861    unsigned char & st_byte = status_[sequence];
862    st_byte = static_cast<unsigned char>(st_byte & ~7);
863    st_byte = static_cast<unsigned char>(st_byte | status);
864  }
865  /// Start or reset using maximumRows_ and Columns_ - true if change
866  bool startPermanentArrays();
867  /** Normally the first factorization does sparse coding because
868      the factorization could be singular.  This allows initial dense
869      factorization when it is known to be safe
870  */
871  void setInitialDenseFactorization(bool onOff);
872  bool  initialDenseFactorization() const;
873  /** Return sequence In or Out */
874  inline int sequenceIn() const
875  {return sequenceIn_;}
876  inline int sequenceOut() const
877  {return sequenceOut_;}
878  /** Set sequenceIn or Out */
879  inline void  setSequenceIn(int sequence)
880  { sequenceIn_=sequence;}
881  inline void  setSequenceOut(int sequence)
882  { sequenceOut_=sequence;}
883  /** Return direction In or Out */
884  inline int directionIn() const
885  {return directionIn_;}
886  inline int directionOut() const
887  {return directionOut_;}
888  /** Set directionIn or Out */
889  inline void  setDirectionIn(int direction)
890  { directionIn_=direction;}
891  inline void  setDirectionOut(int direction)
892  { directionOut_=direction;}
893  /// Value of Out variable
894  inline double valueOut() const
895  { return valueOut_;}
896  /// Returns 1 if sequence indicates column
897  inline int isColumn(int sequence) const
898  { return sequence<numberColumns_ ? 1 : 0;}
899  /// Returns sequence number within section
900  inline int sequenceWithin(int sequence) const
901  { return sequence<numberColumns_ ? sequence : sequence-numberColumns_;}
902  /// Return row or column values
903  inline double solution(int sequence)
904  { return solution_[sequence];}
905  /// Return address of row or column values
906  inline double & solutionAddress(int sequence)
907  { return solution_[sequence];}
908  inline double reducedCost(int sequence)
909   { return dj_[sequence];}
910  inline double & reducedCostAddress(int sequence)
911   { return dj_[sequence];}
912  inline double lower(int sequence)
913  { return lower_[sequence];}
914  /// Return address of row or column lower bound
915  inline double & lowerAddress(int sequence)
916  { return lower_[sequence];}
917  inline double upper(int sequence)
918  { return upper_[sequence];}
919  /// Return address of row or column upper bound
920  inline double & upperAddress(int sequence)
921  { return upper_[sequence];}
922  inline double cost(int sequence)
923  { return cost_[sequence];}
924  /// Return address of row or column cost
925  inline double & costAddress(int sequence)
926  { return cost_[sequence];}
927  /// Return original lower bound
928  inline double originalLower(int iSequence) const
929  { if (iSequence<numberColumns_) return columnLower_[iSequence]; else
930    return rowLower_[iSequence-numberColumns_];}
931  /// Return original lower bound
932  inline double originalUpper(int iSequence) const
933  { if (iSequence<numberColumns_) return columnUpper_[iSequence]; else
934    return rowUpper_[iSequence-numberColumns_];}
935  /// Theta (pivot change)
936  inline double theta() const
937  { return theta_;}
938  /** Best possible improvement using djs (primal) or
939      obj change by flipping bounds to make dual feasible (dual) */
940  inline double bestPossibleImprovement() const
941  { return bestPossibleImprovement_;}
942  /// Return pointer to details of costs
943  inline ClpNonLinearCost * nonLinearCost() const
944  { return nonLinearCost_;}
945  /** Return more special options
946      1 bit - if presolve says infeasible in ClpSolve return
947      2 bit - if presolved problem infeasible return
948      4 bit - keep arrays like upper_ around
949      8 bit - if factorization kept can still declare optimal at once
950      16 bit - if checking replaceColumn accuracy before updating
951      32 bit - say optimal if primal feasible!
952      64 bit - give up easily in dual (and say infeasible)
953      128 bit - no objective, 0-1 and in B&B
954      256 bit - in primal from dual or vice versa
955  */
956  inline int moreSpecialOptions() const
957  { return moreSpecialOptions_;}
958  /** Set more special options
959      1 bit - if presolve says infeasible in ClpSolve return
960      2 bit - if presolved problem infeasible return
961      4 bit - keep arrays like upper_ around
962      8 bit - no free or superBasic variables
963      16 bit - if checking replaceColumn accuracy before updating
964      32 bit - say optimal if primal feasible!
965      64 bit - give up easily in dual (and say infeasible)
966      128 bit - no objective, 0-1 and in B&B
967      256 bit - in primal from dual or vice versa
968  */
969  inline void setMoreSpecialOptions(int value)
970  { moreSpecialOptions_ = value;}
971  //@}
972  /**@name status methods */
973  //@{
974  inline void setFakeBound(int sequence, FakeBound fakeBound)
975  {
976    unsigned char & st_byte = status_[sequence];
977    st_byte = static_cast<unsigned char>(st_byte & ~24);
978    st_byte = static_cast<unsigned char>(st_byte | (fakeBound<<3));
979  }
980  inline FakeBound getFakeBound(int sequence) const
981  {return static_cast<FakeBound> ((status_[sequence]>>3)&3);}
982  inline void setRowStatus(int sequence, Status status)
983  {
984    unsigned char & st_byte = status_[sequence+numberColumns_];
985    st_byte = static_cast<unsigned char>(st_byte & ~7);
986    st_byte = static_cast<unsigned char>(st_byte | status);
987  }
988  inline Status getRowStatus(int sequence) const
989  {return static_cast<Status> (status_[sequence+numberColumns_]&7);}
990  inline void setColumnStatus(int sequence, Status status)
991  {
992    unsigned char & st_byte = status_[sequence];
993    st_byte = static_cast<unsigned char>(st_byte & ~7);
994    st_byte = static_cast<unsigned char>(st_byte | status);
995  }
996  inline Status getColumnStatus(int sequence) const
997  {return static_cast<Status> (status_[sequence]&7);}
998  inline void setPivoted( int sequence)
999  { status_[sequence] = static_cast<unsigned char>(status_[sequence] | 32);}
1000  inline void clearPivoted( int sequence)
1001  { status_[sequence] = static_cast<unsigned char>(status_[sequence] & ~32);}
1002  inline bool pivoted(int sequence) const
1003  {return (((status_[sequence]>>5)&1)!=0);}
1004  /// To flag a variable (not inline to allow for column generation)
1005  void setFlagged( int sequence);
1006  inline void clearFlagged( int sequence)
1007  {
1008    status_[sequence] = static_cast<unsigned char>(status_[sequence] & ~64);
1009  }
1010  inline bool flagged(int sequence) const
1011  {return ((status_[sequence]&64)!=0);}
1012  /// To say row active in primal pivot row choice
1013  inline void setActive( int iRow)
1014  {
1015    status_[iRow] = static_cast<unsigned char>(status_[iRow] | 128);
1016  }
1017  inline void clearActive( int iRow)
1018  {
1019    status_[iRow] = static_cast<unsigned char>(status_[iRow] & ~128);
1020  }
1021  inline bool active(int iRow) const
1022  {return ((status_[iRow]&128)!=0);}
1023  /** Set up status array (can be used by OsiClp).
1024      Also can be used to set up all slack basis */
1025  void createStatus() ;
1026  /** Sets up all slack basis and resets solution to
1027      as it was after initial load or readMps */
1028  void allSlackBasis(bool resetSolution=false);
1029   
1030  /// So we know when to be cautious
1031  inline int lastBadIteration() const
1032  {return lastBadIteration_;}
1033  /// Progress flag - at present 0 bit says artificials out
1034  inline int progressFlag() const
1035  {return (progressFlag_&3);}
1036  /// Force re-factorization early
1037  inline void forceFactorization(int value)
1038  { forceFactorization_ = value;}
1039  /// Raw objective value (so always minimize in primal)
1040  inline double rawObjectiveValue() const
1041  { return objectiveValue_;}
1042   /// Compute objective value from solution and put in objectiveValue_
1043  void computeObjectiveValue(bool useWorkingSolution=false);
1044  /// Compute minimization objective value from internal solution without perturbation
1045  double computeInternalObjectiveValue();
1046  /** Number of extra rows.  These are ones which will be dynamically created
1047      each iteration.  This is for GUB but may have other uses.
1048  */
1049  inline int numberExtraRows() const
1050  { return numberExtraRows_;}
1051  /** Maximum number of basic variables - can be more than number of rows if GUB
1052  */
1053  inline int maximumBasic() const
1054  { return maximumBasic_;}
1055  /// Iteration when we entered dual or primal
1056  inline int baseIteration() const
1057  { return baseIteration_;}
1058  /// Create C++ lines to get to current state
1059  void generateCpp( FILE * fp,bool defaultFactor=false);
1060  /// Gets clean and emptyish factorization
1061  ClpFactorization * getEmptyFactorization();
1062  /// May delete or may make clean and emptyish factorization
1063  void setEmptyFactorization();
1064  /// Move status and solution across
1065  void moveInfo(const ClpSimplex & rhs, bool justStatus=false);
1066  //@}
1067
1068  ///@name Basis handling
1069  // These are only to be used using startFinishOptions (ClpSimplexDual, ClpSimplexPrimal)
1070  // *** At present only without scaling
1071  // *** Slacks havve -1.0 element (so == row activity) - take care
1072  ///Get a row of the tableau (slack part in slack if not NULL)
1073  void getBInvARow(int row, double* z, double * slack=NULL);
1074 
1075  ///Get a row of the basis inverse
1076  void getBInvRow(int row, double* z);
1077 
1078  ///Get a column of the tableau
1079  void getBInvACol(int col, double* vec);
1080 
1081  ///Get a column of the basis inverse
1082  void getBInvCol(int col, double* vec);
1083 
1084  /** Get basic indices (order of indices corresponds to the
1085      order of elements in a vector retured by getBInvACol() and
1086      getBInvCol()).
1087  */
1088  void getBasics(int* index);
1089 
1090  //@}
1091    //-------------------------------------------------------------------------
1092    /**@name Changing bounds on variables and constraints */
1093    //@{
1094       /** Set an objective function coefficient */
1095       void setObjectiveCoefficient( int elementIndex, double elementValue );
1096       /** Set an objective function coefficient */
1097       inline void setObjCoeff( int elementIndex, double elementValue )
1098       { setObjectiveCoefficient( elementIndex, elementValue);}
1099
1100      /** Set a single column lower bound<br>
1101          Use -DBL_MAX for -infinity. */
1102       void setColumnLower( int elementIndex, double elementValue );
1103     
1104      /** Set a single column upper bound<br>
1105          Use DBL_MAX for infinity. */
1106       void setColumnUpper( int elementIndex, double elementValue );
1107
1108      /** Set a single column lower and upper bound */
1109      void setColumnBounds( int elementIndex,
1110        double lower, double upper );
1111
1112      /** Set the bounds on a number of columns simultaneously<br>
1113          The default implementation just invokes setColLower() and
1114          setColUpper() over and over again.
1115          @param indexFirst,indexLast pointers to the beginning and after the
1116                 end of the array of the indices of the variables whose
1117                 <em>either</em> bound changes
1118          @param boundList the new lower/upper bound pairs for the variables
1119      */
1120      void setColumnSetBounds(const int* indexFirst,
1121                                   const int* indexLast,
1122                                   const double* boundList);
1123     
1124      /** Set a single column lower bound<br>
1125          Use -DBL_MAX for -infinity. */
1126       inline void setColLower( int elementIndex, double elementValue )
1127       { setColumnLower(elementIndex, elementValue);}
1128      /** Set a single column upper bound<br>
1129          Use DBL_MAX for infinity. */
1130       inline void setColUpper( int elementIndex, double elementValue )
1131       { setColumnUpper(elementIndex, elementValue);}
1132
1133      /** Set a single column lower and upper bound */
1134      inline void setColBounds( int elementIndex,
1135        double lower, double upper )
1136       { setColumnBounds(elementIndex, lower, upper);}
1137
1138      /** Set the bounds on a number of columns simultaneously<br>
1139          @param indexFirst,indexLast pointers to the beginning and after the
1140                 end of the array of the indices of the variables whose
1141                 <em>either</em> bound changes
1142          @param boundList the new lower/upper bound pairs for the variables
1143      */
1144      inline void setColSetBounds(const int* indexFirst,
1145                                   const int* indexLast,
1146                                   const double* boundList)
1147      { setColumnSetBounds(indexFirst, indexLast, boundList);}
1148     
1149      /** Set a single row lower bound<br>
1150          Use -DBL_MAX for -infinity. */
1151      void setRowLower( int elementIndex, double elementValue );
1152     
1153      /** Set a single row upper bound<br>
1154          Use DBL_MAX for infinity. */
1155      void setRowUpper( int elementIndex, double elementValue ) ;
1156   
1157      /** Set a single row lower and upper bound */
1158      void setRowBounds( int elementIndex,
1159                                 double lower, double upper ) ;
1160   
1161      /** Set the bounds on a number of rows simultaneously<br>
1162          @param indexFirst,indexLast pointers to the beginning and after the
1163                 end of the array of the indices of the constraints whose
1164                 <em>either</em> bound changes
1165          @param boundList the new lower/upper bound pairs for the constraints
1166      */
1167      void setRowSetBounds(const int* indexFirst,
1168                                   const int* indexLast,
1169                                   const double* boundList);
1170  /// Resizes rim part of model
1171  void resize (int newNumberRows, int newNumberColumns);
1172   
1173    //@}
1174
1175////////////////// data //////////////////
1176protected:
1177
1178  /**@name data.  Many arrays have a row part and a column part.
1179   There is a single array with both - columns then rows and
1180   then normally two arrays pointing to rows and columns.  The
1181   single array is the owner of memory
1182  */
1183  //@{
1184  /** Best possible improvement using djs (primal) or
1185      obj change by flipping bounds to make dual feasible (dual) */
1186  double bestPossibleImprovement_;
1187  /// Zero tolerance
1188  double zeroTolerance_;
1189  /// Sequence of worst (-1 if feasible)
1190  int columnPrimalSequence_;
1191  /// Sequence of worst (-1 if feasible)
1192  int rowPrimalSequence_;
1193  /// "Best" objective value
1194  double bestObjectiveValue_;
1195  /// More special options - see set for details
1196  int moreSpecialOptions_;
1197  /// Iteration when we entered dual or primal
1198  int baseIteration_;
1199  /// Primal tolerance needed to make dual feasible (<largeTolerance)
1200  double primalToleranceToGetOptimal_;
1201  /// Large bound value (for complementarity etc)
1202  double largeValue_;
1203  /// Largest error on Ax-b
1204  double largestPrimalError_;
1205  /// Largest error on basic duals
1206  double largestDualError_;
1207  /// For computing whether to re-factorize
1208  double alphaAccuracy_;
1209  /// Dual bound
1210  double dualBound_;
1211  /// Alpha (pivot element)
1212  double alpha_;
1213  /// Theta (pivot change)
1214  double theta_;
1215  /// Lower Bound on In variable
1216  double lowerIn_;
1217  /// Value of In variable
1218  double valueIn_;
1219  /// Upper Bound on In variable
1220  double upperIn_;
1221  /// Reduced cost of In variable
1222  double dualIn_;
1223  /// Lower Bound on Out variable
1224  double lowerOut_;
1225  /// Value of Out variable
1226  double valueOut_;
1227  /// Upper Bound on Out variable
1228  double upperOut_;
1229  /// Infeasibility (dual) or ? (primal) of Out variable
1230  double dualOut_;
1231  /// Current dual tolerance for algorithm
1232  double dualTolerance_;
1233  /// Current primal tolerance for algorithm
1234  double primalTolerance_;
1235  /// Sum of dual infeasibilities
1236  double sumDualInfeasibilities_;
1237  /// Sum of primal infeasibilities
1238  double sumPrimalInfeasibilities_;
1239  /// Weight assigned to being infeasible in primal
1240  double infeasibilityCost_;
1241  /// Sum of Dual infeasibilities using tolerance based on error in duals
1242  double sumOfRelaxedDualInfeasibilities_;
1243  /// Sum of Primal infeasibilities using tolerance based on error in primals
1244  double sumOfRelaxedPrimalInfeasibilities_;
1245  /// Acceptable pivot value just after factorization
1246  double acceptablePivot_;
1247  /// Working copy of lower bounds (Owner of arrays below)
1248  double * lower_;
1249  /// Row lower bounds - working copy
1250  double * rowLowerWork_;
1251  /// Column lower bounds - working copy
1252  double * columnLowerWork_;
1253  /// Working copy of upper bounds (Owner of arrays below)
1254  double * upper_;
1255  /// Row upper bounds - working copy
1256  double * rowUpperWork_;
1257  /// Column upper bounds - working copy
1258  double * columnUpperWork_;
1259  /// Working copy of objective (Owner of arrays below)
1260  double * cost_;
1261  /// Row objective - working copy
1262  double * rowObjectiveWork_;
1263  /// Column objective - working copy
1264  double * objectiveWork_;
1265  /// Useful row length arrays
1266  CoinIndexedVector * rowArray_[6];
1267  /// Useful column length arrays
1268  CoinIndexedVector * columnArray_[6];
1269  /// Sequence of In variable
1270  int sequenceIn_;
1271  /// Direction of In, 1 going up, -1 going down, 0 not a clude
1272  int directionIn_;
1273  /// Sequence of Out variable
1274  int sequenceOut_;
1275  /// Direction of Out, 1 to upper bound, -1 to lower bound, 0 - superbasic
1276  int directionOut_;
1277  /// Pivot Row
1278  int pivotRow_;
1279  /// Last good iteration (immediately after a re-factorization)
1280  int lastGoodIteration_;
1281  /// Working copy of reduced costs (Owner of arrays below)
1282  double * dj_;
1283  /// Reduced costs of slacks not same as duals (or - duals)
1284  double * rowReducedCost_;
1285  /// Possible scaled reduced costs
1286  double * reducedCostWork_;
1287  /// Working copy of primal solution (Owner of arrays below)
1288  double * solution_;
1289  /// Row activities - working copy
1290  double * rowActivityWork_;
1291  /// Column activities - working copy
1292  double * columnActivityWork_;
1293  /// Number of dual infeasibilities
1294  int numberDualInfeasibilities_;
1295  /// Number of dual infeasibilities (without free)
1296  int numberDualInfeasibilitiesWithoutFree_;
1297  /// Number of primal infeasibilities
1298  int numberPrimalInfeasibilities_;
1299  /// How many iterative refinements to do
1300  int numberRefinements_;
1301  /// dual row pivot choice
1302  ClpDualRowPivot * dualRowPivot_;
1303  /// primal column pivot choice
1304  ClpPrimalColumnPivot * primalColumnPivot_;
1305  /// Basic variables pivoting on which rows
1306  int * pivotVariable_;
1307  /// factorization
1308  ClpFactorization * factorization_;
1309  /// Saved version of solution
1310  double * savedSolution_;
1311  /// Number of times code has tentatively thought optimal
1312  int numberTimesOptimal_;
1313  /// Disaster handler
1314  ClpDisasterHandler * disasterArea_;
1315  /// If change has been made (first attempt at stopping looping)
1316  int changeMade_;
1317  /// Algorithm >0 == Primal, <0 == Dual
1318  int algorithm_;
1319  /** Now for some reliability aids
1320      This forces re-factorization early */
1321  int forceFactorization_;
1322  /** Perturbation:
1323      -50 to +50 - perturb by this power of ten (-6 sounds good)
1324      100 - auto perturb if takes too long (1.0e-6 largest nonzero)
1325      101 - we are perturbed
1326      102 - don't try perturbing again
1327      default is 100
1328  */
1329  int perturbation_;
1330  /// Saved status regions
1331  unsigned char * saveStatus_;
1332  /** Very wasteful way of dealing with infeasibilities in primal.
1333      However it will allow non-linearities and use of dual
1334      analysis.  If it doesn't work it can easily be replaced.
1335  */
1336  ClpNonLinearCost * nonLinearCost_;
1337  /// So we know when to be cautious
1338  int lastBadIteration_;
1339  /// So we know when to open up again
1340  int lastFlaggedIteration_;
1341  /// Can be used for count of fake bounds (dual) or fake costs (primal)
1342  int numberFake_;
1343  /// Can be used for count of changed costs (dual) or changed bounds (primal)
1344  int numberChanged_;
1345  /// Progress flag - at present 0 bit says artificials out, 1 free in
1346  int progressFlag_;
1347  /// First free/super-basic variable (-1 if none)
1348  int firstFree_;
1349  /** Number of extra rows.  These are ones which will be dynamically created
1350      each iteration.  This is for GUB but may have other uses.
1351  */
1352  int numberExtraRows_;
1353  /** Maximum number of basic variables - can be more than number of rows if GUB
1354  */
1355  int maximumBasic_;
1356  /// If may skip final factorize then allow up to this pivots (default 20)
1357  int dontFactorizePivots_;
1358  /** For advanced use.  When doing iterative solves things can get
1359      nasty so on values pass if incoming solution has largest
1360      infeasibility < incomingInfeasibility throw out variables
1361      from basis until largest infeasibility < allowedInfeasibility.
1362      if allowedInfeasibility>= incomingInfeasibility this is
1363      always possible altough you may end up with an all slack basis.
1364
1365      Defaults are 1.0,10.0
1366  */
1367  double incomingInfeasibility_;
1368  double allowedInfeasibility_;
1369  /// Automatic scaling of objective and rhs and bounds
1370  int automaticScale_;
1371  /// Maximum perturbation array size (take out when code rewritten)
1372  int maximumPerturbationSize_;
1373  /// Perturbation array (maximumPerturbationSize_)
1374  double * perturbationArray_;
1375  /// A copy of model with certain state - normally without cuts
1376  ClpSimplex * baseModel_;
1377  /// For dealing with all issues of cycling etc
1378  ClpSimplexProgress progress_;
1379public:
1380  /// Spare int array for passing information [0]!=0 switches on
1381  mutable int spareIntArray_[4];
1382  /// Spare double array for passing information [0]!=0 switches on
1383  mutable double spareDoubleArray_[4];
1384protected:
1385  /// Allow OsiClp certain perks
1386  friend class OsiClpSolverInterface;
1387  //@}
1388};
1389//#############################################################################
1390/** A function that tests the methods in the ClpSimplex class. The
1391    only reason for it not to be a member method is that this way it doesn't
1392    have to be compiled into the library. And that's a gain, because the
1393    library should be compiled with optimization on, but this method should be
1394    compiled with debugging.
1395
1396    It also does some testing of ClpFactorization class
1397 */
1398void
1399ClpSimplexUnitTest(const std::string & mpsDir);
1400
1401// For Devex stuff
1402#define DEVEX_TRY_NORM 1.0e-4
1403#define DEVEX_ADD_ONE 1.0
1404#endif
Note: See TracBrowser for help on using the repository browser.