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

Last change on this file since 1729 was 1729, checked in by forrest, 9 years ago

messages for Cbc fathoming and allow perturbation and some safety?

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