source: branches/BSP/trunk/Clp/src/ClpSimplex.hpp @ 1087

Last change on this file since 1087 was 1087, checked in by ladanyi, 12 years ago

committing merging w/ trunk of Clp, Osi, Cbc

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