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

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