source: trunk/Clp/src/ClpSimplex.hpp @ 1326

Last change on this file since 1326 was 1326, checked in by forrest, 12 years ago

First attempt at Benders decomposition

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