source: branches/devel/Clp/src/ClpSimplex.hpp @ 911

Last change on this file since 911 was 911, checked in by forrest, 13 years ago

ranging

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