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

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

extend ranging

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