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

Last change on this file since 1868 was 1868, checked in by forrest, 7 years ago

allow fix to force to stay in dual

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