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

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

changes for advanced use of Clp

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