1 /* $Id: ClpSimplex.hpp 2385 2019-01-06 19:43:06Z unxusr $ */
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"
19 #include "ClpConfig.h"
20 #include "CoinIndexedVector.hpp"
21 class ClpDualRowPivot;
22 class ClpPrimalColumnPivot;
23 class ClpFactorization;
24 class CoinFactorization;
25 class CoinIndexedVector;
26 class ClpNonLinearCost;
27 class ClpNodeStuff;
28 class CoinStructuredModel;
29 class OsiClpSolverInterface;
30 class CoinWarmStartBasis;
31 class ClpDisasterHandler;
32 class ClpConstraint;
33 /*
34 May want to use Clp defaults so that with ABC defined but not used
35 it behaves as Clp (and ABC used will be different than if not defined)
36 */
37 #ifdef ABC_INHERIT
38 #ifndef CLP_INHERIT_MODE
39 #define CLP_INHERIT_MODE 1
40 #endif
41 #ifndef ABC_CLP_DEFAULTS
42 #define ABC_CLP_DEFAULTS 0
43 #endif
44 #else
45 #undef ABC_CLP_DEFAULTS
46 #define ABC_CLP_DEFAULTS 1
47 #endif
48 #ifdef CLP_HAS_ABC
49 #include "AbcCommon.hpp"
50 class AbcTolerancesEtc;
51 class AbcSimplex;
52 #include "CoinAbcCommon.hpp"
53 #endif
54 #ifndef ABC_INHERIT
55 #if ABOCA_LITE
56 #ifndef FAKE_CILK
57 #include <cilk/cilk.h>
58 #else
59 #undef cilk_for
60 #undef cilk_spawn
61 #undef cilk_sync
62 #define cilk_for for
63 #define cilk_spawn
64 #define cilk_sync
65 #endif
66 #ifndef LONG_REGION_2
67 #define LONG_REGION_2 1
68 #endif
69 #define SHORT_REGION 1
70 #else
71 #define cilk_spawn
72 #define cilk_sync
73 #endif
74 #ifdef LONG_REGION_2
75 #define SHORT_REGION 1
76 #else
77 #define SHORT_REGION 2
78 #endif
79 // for now keep simple
80 #undef LONG_REGION_2
81 #undef SHORT_REGION
82 #define SHORT_REGION 2
83 #else
84 //ABC_INHERIT
85 #define LONG_REGION_2 1
86 #define SHORT_REGION 1
87 #endif
88 /** This solves LPs using the simplex method
89
90 It inherits from ClpModel and all its arrays are created at
91 algorithm time. Originally I tried to work with model arrays
92 but for simplicity of coding I changed to single arrays with
93 structural variables then row variables. Some coding is still
94 based on old style and needs cleaning up.
95
96 For a description of algorithms:
97
98 for dual see ClpSimplexDual.hpp and at top of ClpSimplexDual.cpp
99 for primal see ClpSimplexPrimal.hpp and at top of ClpSimplexPrimal.cpp
100
101 There is an algorithm data member. + for primal variations
102 and - for dual variations
103
104 */
105
106 class ClpSimplex : public ClpModel {
107 friend void ClpSimplexUnitTest(const std::string &mpsDir);
108
109 public:
110 /** enums for status of various sorts.
111 First 4 match CoinWarmStartBasis,
112 isFixed means fixed at lower bound and out of basis
113 */
114 enum Status {
115 isFree = 0x00,
116 basic = 0x01,
117 atUpperBound = 0x02,
118 atLowerBound = 0x03,
119 superBasic = 0x04,
120 isFixed = 0x05
121 };
122 // For Dual
123 enum FakeBound {
124 noFake = 0x00,
125 lowerFake = 0x01,
126 upperFake = 0x02,
127 bothFake = 0x03
128 };
129
130 /**@name Constructors and destructor and copy */
131 //@{
132 /// Default constructor
133 ClpSimplex(bool emptyMessages = false);
134
135 /** Copy constructor. May scale depending on mode
136 -1 leave mode as is
137 0 -off, 1 equilibrium, 2 geometric, 3, auto, 4 dynamic(later)
138 */
139 ClpSimplex(const ClpSimplex &rhs, int scalingMode = -1);
140 /** Copy constructor from model. May scale depending on mode
141 -1 leave mode as is
142 0 -off, 1 equilibrium, 2 geometric, 3, auto, 4 dynamic(later)
143 */
144 ClpSimplex(const ClpModel &rhs, int scalingMode = -1);
145 /** Subproblem constructor. A subset of whole model is created from the
146 row and column lists given. The new order is given by list order and
147 duplicates are allowed. Name and integer information can be dropped
148 Can optionally modify rhs to take into account variables NOT in list
149 in this case duplicates are not allowed (also see getbackSolution)
150 */
151 ClpSimplex(const ClpModel *wholeModel,
152 int numberRows, const int *whichRows,
153 int numberColumns, const int *whichColumns,
154 bool dropNames = true, bool dropIntegers = true,
155 bool fixOthers = false);
156 /** Subproblem constructor. A subset of whole model is created from the
157 row and column lists given. The new order is given by list order and
158 duplicates are allowed. Name and integer information can be dropped
159 Can optionally modify rhs to take into account variables NOT in list
160 in this case duplicates are not allowed (also see getbackSolution)
161 */
162 ClpSimplex(const ClpSimplex *wholeModel,
163 int numberRows, const int *whichRows,
164 int numberColumns, const int *whichColumns,
165 bool dropNames = true, bool dropIntegers = true,
166 bool fixOthers = false);
167 /** This constructor modifies original ClpSimplex and stores
168 original stuff in created ClpSimplex. It is only to be used in
169 conjunction with originalModel */
170 ClpSimplex(ClpSimplex *wholeModel,
171 int numberColumns, const int *whichColumns);
172 /** This copies back stuff from miniModel and then deletes miniModel.
173 Only to be used with mini constructor */
174 void originalModel(ClpSimplex *miniModel);
175 #ifdef ABC_INHERIT
abcState() const176 inline int abcState() const
177 {
178 return abcState_;
179 }
setAbcState(int state)180 inline void setAbcState(int state)
181 {
182 abcState_ = state;
183 }
abcSimplex() const184 inline AbcSimplex *abcSimplex() const
185 {
186 return abcSimplex_;
187 }
setAbcSimplex(AbcSimplex * simplex)188 inline void setAbcSimplex(AbcSimplex *simplex)
189 {
190 abcSimplex_ = simplex;
191 }
192 /// Returns 0 if dual can be skipped
193 int doAbcDual();
194 /// Returns 0 if primal can be skipped
195 int doAbcPrimal(int ifValuesPass);
196 #endif
197 /** Array persistence flag
198 If 0 then as now (delete/new)
199 1 then only do arrays if bigger needed
200 2 as 1 but give a bit extra if bigger needed
201 */
202 void setPersistenceFlag(int value);
203 /// Save a copy of model with certain state - normally without cuts
204 void makeBaseModel();
205 /// Switch off base model
206 void deleteBaseModel();
207 /// See if we have base model
baseModel() const208 inline ClpSimplex *baseModel() const
209 {
210 return baseModel_;
211 }
212 /** Reset to base model (just size and arrays needed)
213 If model NULL use internal copy
214 */
215 void setToBaseModel(ClpSimplex *model = NULL);
216 /// Assignment operator. This copies the data
217 ClpSimplex &operator=(const ClpSimplex &rhs);
218 /// Destructor
219 ~ClpSimplex();
220 // Ones below are just ClpModel with some changes
221 /** Loads a problem (the constraints on the
222 rows are given by lower and upper bounds). If a pointer is 0 then the
223 following values are the default:
224 <ul>
225 <li> <code>colub</code>: all columns have upper bound infinity
226 <li> <code>collb</code>: all columns have lower bound 0
227 <li> <code>rowub</code>: all rows have upper bound infinity
228 <li> <code>rowlb</code>: all rows have lower bound -infinity
229 <li> <code>obj</code>: all variables have 0 objective coefficient
230 </ul>
231 */
232 void loadProblem(const ClpMatrixBase &matrix,
233 const double *collb, const double *colub,
234 const double *obj,
235 const double *rowlb, const double *rowub,
236 const double *rowObjective = NULL);
237 void loadProblem(const CoinPackedMatrix &matrix,
238 const double *collb, const double *colub,
239 const double *obj,
240 const double *rowlb, const double *rowub,
241 const double *rowObjective = NULL);
242
243 /** Just like the other loadProblem() method except that the matrix is
244 given in a standard column major ordered format (without gaps). */
245 void loadProblem(const int numcols, const int numrows,
246 const CoinBigIndex *start, const int *index,
247 const double *value,
248 const double *collb, const double *colub,
249 const double *obj,
250 const double *rowlb, const double *rowub,
251 const double *rowObjective = NULL);
252 /// This one is for after presolve to save memory
253 void loadProblem(const int numcols, const int numrows,
254 const CoinBigIndex *start, const int *index,
255 const double *value, const int *length,
256 const double *collb, const double *colub,
257 const double *obj,
258 const double *rowlb, const double *rowub,
259 const double *rowObjective = NULL);
260 /** This loads a model from a coinModel object - returns number of errors.
261 If keepSolution true and size is same as current then
262 keeps current status and solution
263 */
264 int loadProblem(CoinModel &modelObject, bool keepSolution = false);
265 /// Read an mps file from the given filename
266 int readMps(const char *filename,
267 bool keepNames = false,
268 bool ignoreErrors = false);
269 /// Read GMPL files from the given filenames
270 int readGMPL(const char *filename, const char *dataName,
271 bool keepNames = false);
272 /// Read file in LP format from file with name filename.
273 /// See class CoinLpIO for description of this format.
274 int readLp(const char *filename, const double epsilon = 1e-5);
275 /** Write the problem into an Lp file of the given filename.
276 If objSense is non zero then -1.0 forces the code to write a
277 maximization objective and +1.0 to write a minimization one.
278 If 0.0 then solver can do what it wants.*/
279 void writeLp(const char *filename,
280 const char *extension = "lp",
281 double epsilon = 1e-5,
282 int numberAcross = 10,
283 int decimals = 5,
284 double objSense = 0.0,
285 bool useRowNames = true) const;
286 /** Borrow model. This is so we dont have to copy large amounts
287 of data around. It assumes a derived class wants to overwrite
288 an empty model with a real one - while it does an algorithm.
289 This is same as ClpModel one, but sets scaling on etc. */
290 void borrowModel(ClpModel &otherModel);
291 void borrowModel(ClpSimplex &otherModel);
292 /// Pass in Event handler (cloned and deleted at end)
293 void passInEventHandler(const ClpEventHandler *eventHandler);
294 /// Puts solution back into small model
295 void getbackSolution(const ClpSimplex &smallModel, const int *whichRow, const int *whichColumn);
296 /** Load nonlinear part of problem from AMPL info
297 Returns 0 if linear
298 1 if quadratic objective
299 2 if quadratic constraints
300 3 if nonlinear objective
301 4 if nonlinear constraints
302 -1 on failure
303 */
304 int loadNonLinear(void *info, int &numberConstraints,
305 ClpConstraint **&constraints);
306 #ifdef ABC_INHERIT
307 /// Loads tolerances etc
308 void loadTolerancesEtc(const AbcTolerancesEtc &data);
309 /// Unloads tolerances etc
310 void unloadTolerancesEtc(AbcTolerancesEtc &data);
311 #endif
312 //@}
313
314 /**@name Functions most useful to user */
315 //@{
316 /** General solve algorithm which can do presolve.
317 See ClpSolve.hpp for options
318 */
319 int initialSolve(ClpSolve &options);
320 /// Default initial solve
321 int initialSolve();
322 /// Dual initial solve
323 int initialDualSolve();
324 /// Primal initial solve
325 int initialPrimalSolve();
326 /// Barrier initial solve
327 int initialBarrierSolve();
328 /// Barrier initial solve, not to be followed by crossover
329 int initialBarrierNoCrossSolve();
330 /** Dual algorithm - see ClpSimplexDual.hpp for method.
331 ifValuesPass==2 just does values pass and then stops.
332
333 startFinishOptions - bits
334 1 - do not delete work areas and factorization at end
335 2 - use old factorization if same number of rows
336 4 - skip as much initialization of work areas as possible
337 (based on whatsChanged in clpmodel.hpp) ** work in progress
338 maybe other bits later
339 */
340 int dual(int ifValuesPass = 0, int startFinishOptions = 0);
341 // If using Debug
342 int dualDebug(int ifValuesPass = 0, int startFinishOptions = 0);
343 /** Primal algorithm - see ClpSimplexPrimal.hpp for method.
344 ifValuesPass==2 just does values pass and then stops.
345
346 startFinishOptions - bits
347 1 - do not delete work areas and factorization at end
348 2 - use old factorization if same number of rows
349 4 - skip as much initialization of work areas as possible
350 (based on whatsChanged in clpmodel.hpp) ** work in progress
351 maybe other bits later
352 */
353 int primal(int ifValuesPass = 0, int startFinishOptions = 0);
354 /** Solves nonlinear problem using SLP - may be used as crash
355 for other algorithms when number of iterations small.
356 Also exits if all problematical variables are changing
357 less than deltaTolerance
358 */
359 int nonlinearSLP(int numberPasses, double deltaTolerance);
360 /** Solves problem with nonlinear constraints using SLP - may be used as crash
361 for other algorithms when number of iterations small.
362 Also exits if all problematical variables are changing
363 less than deltaTolerance
364 */
365 int nonlinearSLP(int numberConstraints, ClpConstraint **constraints,
366 int numberPasses, double deltaTolerance);
367 /** Solves using barrier (assumes you have good cholesky factor code).
368 Does crossover to simplex if asked*/
369 int barrier(bool crossover = true);
370 /** Solves non-linear using reduced gradient. Phase = 0 get feasible,
371 =1 use solution */
372 int reducedGradient(int phase = 0);
373 /// Solve using structure of model and maybe in parallel
374 int solve(CoinStructuredModel *model);
375 #ifdef ABC_INHERIT
376 /** solvetype 0 for dual, 1 for primal
377 startup 1 for values pass
378 interrupt whether to pass across interrupt handler
379 add 10 to return AbcSimplex
380 */
381 AbcSimplex *dealWithAbc(int solveType, int startUp, bool interrupt = false);
382 //void dealWithAbc(int solveType,int startUp,bool interrupt=false);
383 #endif
384 /** This loads a model from a CoinStructuredModel object - returns number of errors.
385 If originalOrder then keep to order stored in blocks,
386 otherwise first column/rows correspond to first block - etc.
387 If keepSolution true and size is same as current then
388 keeps current status and solution
389 */
390 int loadProblem(CoinStructuredModel &modelObject,
391 bool originalOrder = true, bool keepSolution = false);
392 /**
393 When scaling is on it is possible that the scaled problem
394 is feasible but the unscaled is not. Clp returns a secondary
395 status code to that effect. This option allows for a cleanup.
396 If you use it I would suggest 1.
397 This only affects actions when scaled optimal
398 0 - no action
399 1 - clean up using dual if primal infeasibility
400 2 - clean up using dual if dual infeasibility
401 3 - clean up using dual if primal or dual infeasibility
402 11,12,13 - as 1,2,3 but use primal
403
404 return code as dual/primal
405 */
406 int cleanup(int cleanupScaling);
407 /** Clean primal solution
408 If you expect solution to only have exact multiples of "exactMultiple" then
409 this tries moving solution values to nearest multiple. If still feasible
410 then the solution is replaced.
411
412 This is designed for the case where values should be integral, but Clp may
413 have values at e.g. 1.0e-13
414 Returns 0 if successful, n if n rhs violated
415 The dual version may be written if this gets used.
416 */
417 int cleanPrimalSolution(double exactMultiple);
418 /** Dual ranging.
419 This computes increase/decrease in cost for each given variable and corresponding
420 sequence numbers which would change basis. Sequence numbers are 0..numberColumns
421 and numberColumns.. for artificials/slacks.
422 For non-basic variables the information is trivial to compute and the change in cost is just minus the
423 reduced cost and the sequence number will be that of the non-basic variables.
424 For basic variables a ratio test is between the reduced costs for non-basic variables
425 and the row of the tableau corresponding to the basic variable.
426 The increase/decrease value is always >= 0.0
427
428 Up to user to provide correct length arrays where each array is of length numberCheck.
429 which contains list of variables for which information is desired. All other
430 arrays will be filled in by function. If fifth entry in which is variable 7 then fifth entry in output arrays
431 will be information for variable 7.
432
433 If valueIncrease/Decrease not NULL (both must be NULL or both non NULL) then these are filled with
434 the value of variable if such a change in cost were made (the existing bounds are ignored)
435
436 Returns non-zero if infeasible unbounded etc
437 */
438 int dualRanging(int numberCheck, const int *which,
439 double *costIncrease, int *sequenceIncrease,
440 double *costDecrease, int *sequenceDecrease,
441 double *valueIncrease = NULL, double *valueDecrease = NULL);
442 /** Primal ranging.
443 This computes increase/decrease in value for each given variable and corresponding
444 sequence numbers which would change basis. Sequence numbers are 0..numberColumns
445 and numberColumns.. for artificials/slacks.
446 This should only be used for non-basic variabls as otherwise information is pretty useless
447 For basic variables the sequence number will be that of the basic variables.
448
449 Up to user to provide correct length arrays where each array is of length numberCheck.
450 which contains list of variables for which information is desired. All other
451 arrays will be filled in by function. If fifth entry in which is variable 7 then fifth entry in output arrays
452 will be information for variable 7.
453
454 Returns non-zero if infeasible unbounded etc
455 */
456 int primalRanging(int numberCheck, const int *which,
457 double *valueIncrease, int *sequenceIncrease,
458 double *valueDecrease, int *sequenceDecrease);
459 /**
460 Modifies coefficients etc and if necessary pivots in and out.
461 All at same status will be done (basis may go singular).
462 User can tell which others have been done (i.e. if status matches).
463 If called from outside will change status and return 0.
464 If called from event handler returns non-zero if user has to take action.
465 indices>=numberColumns are slacks (obviously no coefficients)
466 status array is (char) Status enum
467 */
468 int modifyCoefficientsAndPivot(int number,
469 const int *which,
470 const CoinBigIndex *start,
471 const int *row,
472 const double *newCoefficient,
473 const unsigned char *newStatus = NULL,
474 const double *newLower = NULL,
475 const double *newUpper = NULL,
476 const double *newObjective = NULL);
477 /** Take out duplicate rows (includes scaled rows and intersections).
478 On exit whichRows has rows to delete - return code is number can be deleted
479 or -1 if would be infeasible.
480 If tolerance is -1.0 use primalTolerance for equality rows and infeasibility
481 If cleanUp not zero then spend more time trying to leave more stable row
482 and make row bounds exact multiple of cleanUp if close enough
483 */
484 int outDuplicateRows(int numberLook, int *whichRows, bool noOverlaps = false, double tolerance = -1.0,
485 double cleanUp = 0.0);
486 /** Try simple crash like techniques to get closer to primal feasibility
487 returns final sum of infeasibilities */
488 double moveTowardsPrimalFeasible();
489 /** Try simple crash like techniques to remove super basic slacks
490 but only if > threshold */
491 void removeSuperBasicSlacks(int threshold = 0);
492 /** Mini presolve (faster)
493 Char arrays must be numberRows and numberColumns long
494 on entry second part must be filled in as follows -
495 0 - possible
496 >0 - take out and do something (depending on value - TBD)
497 -1 row/column can't vanish but can have entries removed/changed
498 -2 don't touch at all
499 on exit <=0 ones will be in presolved problem
500 struct will be created and will be long enough
501 (information on length etc in first entry)
502 user must delete struct
503 */
504 ClpSimplex *miniPresolve(char *rowType, char *columnType, void **info);
505 /// After mini presolve
506 void miniPostsolve(const ClpSimplex *presolvedModel, void *info);
507 /// mini presolve and solve
508 void miniSolve(char *rowType, char *columnType, int algorithm, int startUp);
509 /** Write the basis in MPS format to the specified file.
510 If writeValues true writes values of structurals
511 (and adds VALUES to end of NAME card)
512
513 Row and column names may be null.
514 formatType is
515 <ul>
516 <li> 0 - normal
517 <li> 1 - extra accuracy
518 <li> 2 - IEEE hex (later)
519 </ul>
520
521 Returns non-zero on I/O error
522 */
523 int writeBasis(const char *filename,
524 bool writeValues = false,
525 int formatType = 0) const;
526 /** Read a basis from the given filename,
527 returns -1 on file error, 0 if no values, 1 if values */
528 int readBasis(const char *filename);
529 /// Returns a basis (to be deleted by user)
530 CoinWarmStartBasis *getBasis() const;
531 /// Passes in factorization
532 void setFactorization(ClpFactorization &factorization);
533 // Swaps factorization
534 ClpFactorization *swapFactorization(ClpFactorization *factorization);
535 /// Copies in factorization to existing one
536 void copyFactorization(ClpFactorization &factorization);
537 /** Tightens primal bounds to make dual faster. Unless
538 fixed or doTight>10, bounds are slightly looser than they could be.
539 This is to make dual go faster and is probably not needed
540 with a presolve. Returns non-zero if problem infeasible.
541
542 Fudge for branch and bound - put bounds on columns of factor *
543 largest value (at continuous) - should improve stability
544 in branch and bound on infeasible branches (0.0 is off)
545 */
546 int tightenPrimalBounds(double factor = 0.0, int doTight = 0, bool tightIntegers = false);
547 /** Crash - at present just aimed at dual, returns
548 -2 if dual preferred and crash basis created
549 -1 if dual preferred and all slack basis preferred
550 0 if basis going in was not all slack
551 1 if primal preferred and all slack basis preferred
552 2 if primal preferred and crash basis created.
553
554 if gap between bounds <="gap" variables can be flipped
555 ( If pivot -1 then can be made super basic!)
556
557 If "pivot" is
558 -1 No pivoting - always primal
559 0 No pivoting (so will just be choice of algorithm)
560 1 Simple pivoting e.g. gub
561 2 Mini iterations
562 */
563 int crash(double gap, int pivot);
564 /// Sets row pivot choice algorithm in dual
565 void setDualRowPivotAlgorithm(ClpDualRowPivot &choice);
566 /// Sets column pivot choice algorithm in primal
567 void setPrimalColumnPivotAlgorithm(ClpPrimalColumnPivot &choice);
568 /// Create a hotstart point of the optimization process
569 void markHotStart(void *&saveStuff);
570 /// Optimize starting from the hotstart
571 void solveFromHotStart(void *saveStuff);
572 /// Delete the snapshot
573 void unmarkHotStart(void *saveStuff);
574 /** For strong branching. On input lower and upper are new bounds
575 while on output they are change in objective function values
576 (>1.0e50 infeasible).
577 Return code is 0 if nothing interesting, -1 if infeasible both
578 ways and +1 if infeasible one way (check values to see which one(s))
579 Solutions are filled in as well - even down, odd up - also
580 status and number of iterations
581 */
582 int strongBranching(int numberVariables, const int *variables,
583 double *newLower, double *newUpper,
584 double **outputSolution,
585 int *outputStatus, int *outputIterations,
586 bool stopOnFirstInfeasible = true,
587 bool alwaysFinish = false,
588 int startFinishOptions = 0);
589 /// Fathom - 1 if solution
590 int fathom(void *stuff);
591 /** Do up to N deep - returns
592 -1 - no solution nNodes_ valid nodes
593 >= if solution and that node gives solution
594 ClpNode array is 2**N long. Values for N and
595 array are in stuff (nNodes_ also in stuff) */
596 int fathomMany(void *stuff);
597 /// Double checks OK
598 double doubleCheck();
599 /// Starts Fast dual2
600 int startFastDual2(ClpNodeStuff *stuff);
601 /// Like Fast dual
602 int fastDual2(ClpNodeStuff *stuff);
603 /// Stops Fast dual2
604 void stopFastDual2(ClpNodeStuff *stuff);
605 /** Deals with crunch aspects
606 mode 0 - in
607 1 - out with solution
608 2 - out without solution
609 returns small model or NULL
610 */
611 ClpSimplex *fastCrunch(ClpNodeStuff *stuff, int mode);
612 //@}
613
614 /**@name Needed for functionality of OsiSimplexInterface */
615 //@{
616 /** Pivot in a variable and out a variable. Returns 0 if okay,
617 1 if inaccuracy forced re-factorization, -1 if would be singular.
618 Also updates primal/dual infeasibilities.
619 Assumes sequenceIn_ and pivotRow_ set and also directionIn and Out.
620 */
621 int pivot();
622
623 /** Pivot in a variable and choose an outgoing one. Assumes primal
624 feasible - will not go through a bound. Returns step length in theta
625 Returns ray in ray_ (or NULL if no pivot)
626 Return codes as before but -1 means no acceptable pivot
627 */
628 int primalPivotResult();
629
630 /** Pivot out a variable and choose an incoing one. Assumes dual
631 feasible - will not go through a reduced cost.
632 Returns step length in theta
633 Return codes as before but -1 means no acceptable pivot
634 */
635 int dualPivotResultPart1();
636 /** Do actual pivot
637 state is 0 if need tableau column, 1 if in rowArray_[1]
638 */
639 int pivotResultPart2(int algorithm, int state);
640
641 /** Common bits of coding for dual and primal. Return 0 if okay,
642 1 if bad matrix, 2 if very bad factorization
643
644 startFinishOptions - bits
645 1 - do not delete work areas and factorization at end
646 2 - use old factorization if same number of rows
647 4 - skip as much initialization of work areas as possible
648 (based on whatsChanged in clpmodel.hpp) ** work in progress
649 maybe other bits later
650
651 */
652 int startup(int ifValuesPass, int startFinishOptions = 0);
653 void finish(int startFinishOptions = 0);
654
655 /** Factorizes and returns true if optimal. Used by user */
656 bool statusOfProblem(bool initial = false);
657 /// If user left factorization frequency then compute
658 void defaultFactorizationFrequency();
659 /// Copy across enabled stuff from one solver to another
660 void copyEnabledStuff(const ClpSimplex *rhs);
661 //@}
662
663 /**@name most useful gets and sets */
664 //@{
665 /// If problem is primal feasible
primalFeasible() const666 inline bool primalFeasible() const
667 {
668 return (numberPrimalInfeasibilities_ == 0);
669 }
670 /// If problem is dual feasible
dualFeasible() const671 inline bool dualFeasible() const
672 {
673 return (numberDualInfeasibilities_ == 0);
674 }
675 /// factorization
factorization() const676 inline ClpFactorization *factorization() const
677 {
678 return factorization_;
679 }
680 /// Sparsity on or off
681 bool sparseFactorization() const;
682 void setSparseFactorization(bool value);
683 /// Factorization frequency
684 int factorizationFrequency() const;
685 void setFactorizationFrequency(int value);
686 /// Dual bound
dualBound() const687 inline double dualBound() const
688 {
689 return dualBound_;
690 }
691 void setDualBound(double value);
692 /// Infeasibility cost
infeasibilityCost() const693 inline double infeasibilityCost() const
694 {
695 return infeasibilityCost_;
696 }
697 void setInfeasibilityCost(double value);
698 /** Amount of print out:
699 0 - none
700 1 - just final
701 2 - just factorizations
702 3 - as 2 plus a bit more
703 4 - verbose
704 above that 8,16,32 etc just for selective debug
705 */
706 /** Perturbation:
707 50 - switch on perturbation
708 100 - auto perturb if takes too long (1.0e-6 largest nonzero)
709 101 - we are perturbed
710 102 - don't try perturbing again
711 default is 100
712 others are for playing
713 */
perturbation() const714 inline int perturbation() const
715 {
716 return perturbation_;
717 }
718 void setPerturbation(int value);
719 /// Current (or last) algorithm
algorithm() const720 inline int algorithm() const
721 {
722 return algorithm_;
723 }
724 /// Set algorithm
setAlgorithm(int value)725 inline void setAlgorithm(int value)
726 {
727 algorithm_ = value;
728 }
729 /// Return true if the objective limit test can be relied upon
730 bool isObjectiveLimitTestValid() const;
731 /// Sum of dual infeasibilities
sumDualInfeasibilities() const732 inline double sumDualInfeasibilities() const
733 {
734 return sumDualInfeasibilities_;
735 }
setSumDualInfeasibilities(double value)736 inline void setSumDualInfeasibilities(double value)
737 {
738 sumDualInfeasibilities_ = value;
739 }
740 /// Sum of relaxed dual infeasibilities
sumOfRelaxedDualInfeasibilities() const741 inline double sumOfRelaxedDualInfeasibilities() const
742 {
743 return sumOfRelaxedDualInfeasibilities_;
744 }
setSumOfRelaxedDualInfeasibilities(double value)745 inline void setSumOfRelaxedDualInfeasibilities(double value)
746 {
747 sumOfRelaxedDualInfeasibilities_ = value;
748 }
749 /// Number of dual infeasibilities
numberDualInfeasibilities() const750 inline int numberDualInfeasibilities() const
751 {
752 return numberDualInfeasibilities_;
753 }
setNumberDualInfeasibilities(int value)754 inline void setNumberDualInfeasibilities(int value)
755 {
756 numberDualInfeasibilities_ = value;
757 }
758 /// Number of dual infeasibilities (without free)
numberDualInfeasibilitiesWithoutFree() const759 inline int numberDualInfeasibilitiesWithoutFree() const
760 {
761 return numberDualInfeasibilitiesWithoutFree_;
762 }
763 /// Sum of primal infeasibilities
sumPrimalInfeasibilities() const764 inline double sumPrimalInfeasibilities() const
765 {
766 return sumPrimalInfeasibilities_;
767 }
setSumPrimalInfeasibilities(double value)768 inline void setSumPrimalInfeasibilities(double value)
769 {
770 sumPrimalInfeasibilities_ = value;
771 }
772 /// Sum of relaxed primal infeasibilities
sumOfRelaxedPrimalInfeasibilities() const773 inline double sumOfRelaxedPrimalInfeasibilities() const
774 {
775 return sumOfRelaxedPrimalInfeasibilities_;
776 }
setSumOfRelaxedPrimalInfeasibilities(double value)777 inline void setSumOfRelaxedPrimalInfeasibilities(double value)
778 {
779 sumOfRelaxedPrimalInfeasibilities_ = value;
780 }
781 /// Number of primal infeasibilities
numberPrimalInfeasibilities() const782 inline int numberPrimalInfeasibilities() const
783 {
784 return numberPrimalInfeasibilities_;
785 }
setNumberPrimalInfeasibilities(int value)786 inline void setNumberPrimalInfeasibilities(int value)
787 {
788 numberPrimalInfeasibilities_ = value;
789 }
790 /** Save model to file, returns 0 if success. This is designed for
791 use outside algorithms so does not save iterating arrays etc.
792 It does not save any messaging information.
793 Does not save scaling values.
794 It does not know about all types of virtual functions.
795 */
796 int saveModel(const char *fileName);
797 /** Restore model from file, returns 0 if success,
798 deletes current model */
799 int restoreModel(const char *fileName);
800
801 /** Just check solution (for external use) - sets sum of
802 infeasibilities etc.
803 If setToBounds 0 then primal column values not changed
804 and used to compute primal row activity values. If 1 or 2
805 then status used - so all nonbasic variables set to
806 indicated bound and if any values changed (or ==2) basic values re-computed.
807 */
808 void checkSolution(int setToBounds = 0);
809 /** Just check solution (for internal use) - sets sum of
810 infeasibilities etc. */
811 void checkSolutionInternal();
812 /// Check unscaled primal solution but allow for rounding error
813 void checkUnscaledSolution();
814 /// Useful row length arrays (0,1,2,3,4,5)
rowArray(int index) const815 inline CoinIndexedVector *rowArray(int index) const
816 {
817 return rowArray_[index];
818 }
819 /// Useful column length arrays (0,1,2,3,4,5)
columnArray(int index) const820 inline CoinIndexedVector *columnArray(int index) const
821 {
822 return columnArray_[index];
823 }
824 //@}
825
826 /******************** End of most useful part **************/
827 /**@name Functions less likely to be useful to casual user */
828 //@{
829 /** Given an existing factorization computes and checks
830 primal and dual solutions. Uses input arrays for variables at
831 bounds. Returns feasibility states */
832 int getSolution(const double *rowActivities,
833 const double *columnActivities);
834 /** Given an existing factorization computes and checks
835 primal and dual solutions. Uses current problem arrays for
836 bounds. Returns feasibility states */
837 int getSolution();
838 /** Constructs a non linear cost from list of non-linearities (columns only)
839 First lower of each column is taken as real lower
840 Last lower is taken as real upper and cost ignored
841
842 Returns nonzero if bad data e.g. lowers not monotonic
843 */
844 int createPiecewiseLinearCosts(const int *starts,
845 const double *lower, const double *gradient);
846 /// dual row pivot choice
dualRowPivot() const847 inline ClpDualRowPivot *dualRowPivot() const
848 {
849 return dualRowPivot_;
850 }
851 /// primal column pivot choice
primalColumnPivot() const852 inline ClpPrimalColumnPivot *primalColumnPivot() const
853 {
854 return primalColumnPivot_;
855 }
856 /// Returns true if model looks OK
goodAccuracy() const857 inline bool goodAccuracy() const
858 {
859 return (largestPrimalError_ < 1.0e-7 && largestDualError_ < 1.0e-7);
860 }
861 /** Return model - updates any scalars */
862 void returnModel(ClpSimplex &otherModel);
863 /** Factorizes using current basis.
864 solveType - 1 iterating, 0 initial, -1 external
865 If 10 added then in primal values pass
866 Return codes are as from ClpFactorization unless initial factorization
867 when total number of singularities is returned.
868 Special case is numberRows_+1 -> all slack basis.
869 */
870 int internalFactorize(int solveType);
871 /// Save data
872 ClpDataSave saveData();
873 /// Restore data
874 void restoreData(ClpDataSave saved);
875 /// Clean up status
876 void cleanStatus();
877 /// Factorizes using current basis. For external use
878 int factorize();
879 /** Computes duals from scratch. If givenDjs then
880 allows for nonzero basic djs */
881 void computeDuals(double *givenDjs);
882 /// Computes primals from scratch
883 void computePrimals(const double *rowActivities,
884 const double *columnActivities);
885 /** Adds multiple of a column into an array */
886 void add(double *array,
887 int column, double multiplier) const;
888 /**
889 Unpacks one column of the matrix into indexed array
890 Uses sequenceIn_
891 Also applies scaling if needed
892 */
893 void unpack(CoinIndexedVector *rowArray) const;
894 /**
895 Unpacks one column of the matrix into indexed array
896 Slack if sequence>= numberColumns
897 Also applies scaling if needed
898 */
899 void unpack(CoinIndexedVector *rowArray, int sequence) const;
900 /**
901 Unpacks one column of the matrix into indexed array
902 ** as packed vector
903 Uses sequenceIn_
904 Also applies scaling if needed
905 */
906 void unpackPacked(CoinIndexedVector *rowArray);
907 /**
908 Unpacks one column of the matrix into indexed array
909 ** as packed vector
910 Slack if sequence>= numberColumns
911 Also applies scaling if needed
912 */
913 void unpackPacked(CoinIndexedVector *rowArray, int sequence);
914 #ifndef CLP_USER_DRIVEN
915 protected:
916 #endif
917 /**
918 This does basis housekeeping and does values for in/out variables.
919 Can also decide to re-factorize
920 */
921 int housekeeping(double objectiveChange);
922 /** This sets largest infeasibility and most infeasible and sum
923 and number of infeasibilities (Primal) */
924 void checkPrimalSolution(const double *rowActivities = NULL,
925 const double *columnActivies = NULL);
926 /** This sets largest infeasibility and most infeasible and sum
927 and number of infeasibilities (Dual) */
928 void checkDualSolution();
929 /** This sets sum and number of infeasibilities (Dual and Primal) */
930 void checkBothSolutions();
931 /** If input negative scales objective so maximum <= -value
932 and returns scale factor used. If positive unscales and also
933 redoes dual stuff
934 */
935 double scaleObjective(double value);
936 /// Solve using Dantzig-Wolfe decomposition and maybe in parallel
937 int solveDW(CoinStructuredModel *model, ClpSolve &options);
938 /// Solve using Benders decomposition and maybe in parallel
939 int solveBenders(CoinStructuredModel *model, ClpSolve &options);
940
941 public:
942 /** For advanced use. When doing iterative solves things can get
943 nasty so on values pass if incoming solution has largest
944 infeasibility < incomingInfeasibility throw out variables
945 from basis until largest infeasibility < allowedInfeasibility
946 or incoming largest infeasibility.
947 If allowedInfeasibility>= incomingInfeasibility this is
948 always possible altough you may end up with an all slack basis.
949
950 Defaults are 1.0,10.0
951 */
952 void setValuesPassAction(double incomingInfeasibility,
953 double allowedInfeasibility);
954 /** Get a clean factorization - i.e. throw out singularities
955 may do more later */
956 int cleanFactorization(int ifValuesPass);
957 //@}
958 /**@name most useful gets and sets */
959 //@{
960 public:
961 /// Initial value for alpha accuracy calculation (-1.0 off)
alphaAccuracy() const962 inline double alphaAccuracy() const
963 {
964 return alphaAccuracy_;
965 }
setAlphaAccuracy(double value)966 inline void setAlphaAccuracy(double value)
967 {
968 alphaAccuracy_ = value;
969 }
970
971 public:
972 /// Objective value
973 //inline double objectiveValue() const {
974 //return (objectiveValue_-bestPossibleImprovement_)*optimizationDirection_ - dblParam_[ClpObjOffset];
975 //}
976 /// Set disaster handler
setDisasterHandler(ClpDisasterHandler * handler)977 inline void setDisasterHandler(ClpDisasterHandler *handler)
978 {
979 disasterArea_ = handler;
980 }
981 /// Get disaster handler
disasterHandler() const982 inline ClpDisasterHandler *disasterHandler() const
983 {
984 return disasterArea_;
985 }
986 /// Large bound value (for complementarity etc)
largeValue() const987 inline double largeValue() const
988 {
989 return largeValue_;
990 }
991 void setLargeValue(double value);
992 /// Largest error on Ax-b
largestPrimalError() const993 inline double largestPrimalError() const
994 {
995 return largestPrimalError_;
996 }
997 /// Largest error on basic duals
largestDualError() const998 inline double largestDualError() const
999 {
1000 return largestDualError_;
1001 }
1002 /// Largest error on Ax-b
setLargestPrimalError(double value)1003 inline void setLargestPrimalError(double value)
1004 {
1005 largestPrimalError_ = value;
1006 }
1007 /// Largest error on basic duals
setLargestDualError(double value)1008 inline void setLargestDualError(double value)
1009 {
1010 largestDualError_ = value;
1011 }
1012 /// Get zero tolerance
zeroTolerance() const1013 inline double zeroTolerance() const
1014 {
1015 return zeroTolerance_; /*factorization_->zeroTolerance();*/
1016 }
1017 /// Set zero tolerance
setZeroTolerance(double value)1018 inline void setZeroTolerance(double value)
1019 {
1020 zeroTolerance_ = value;
1021 }
1022 /// Basic variables pivoting on which rows
pivotVariable() const1023 inline int *pivotVariable() const
1024 {
1025 return pivotVariable_;
1026 }
1027 /// If automatic scaling on
automaticScaling() const1028 inline bool automaticScaling() const
1029 {
1030 return automaticScale_ != 0;
1031 }
setAutomaticScaling(bool onOff)1032 inline void setAutomaticScaling(bool onOff)
1033 {
1034 automaticScale_ = onOff ? 1 : 0;
1035 }
1036 /// Current dual tolerance
currentDualTolerance() const1037 inline double currentDualTolerance() const
1038 {
1039 return dualTolerance_;
1040 }
setCurrentDualTolerance(double value)1041 inline void setCurrentDualTolerance(double value)
1042 {
1043 dualTolerance_ = value;
1044 }
1045 /// Current primal tolerance
currentPrimalTolerance() const1046 inline double currentPrimalTolerance() const
1047 {
1048 return primalTolerance_;
1049 }
setCurrentPrimalTolerance(double value)1050 inline void setCurrentPrimalTolerance(double value)
1051 {
1052 primalTolerance_ = value;
1053 }
1054 /// How many iterative refinements to do
numberRefinements() const1055 inline int numberRefinements() const
1056 {
1057 return numberRefinements_;
1058 }
1059 void setNumberRefinements(int value);
1060 /// Alpha (pivot element) for use by classes e.g. steepestedge
alpha() const1061 inline double alpha() const
1062 {
1063 return alpha_;
1064 }
setAlpha(double value)1065 inline void setAlpha(double value)
1066 {
1067 alpha_ = value;
1068 }
1069 /// Reduced cost of last incoming for use by classes e.g. steepestedge
dualIn() const1070 inline double dualIn() const
1071 {
1072 return dualIn_;
1073 }
1074 /// Set reduced cost of last incoming to force error
setDualIn(double value)1075 inline void setDualIn(double value)
1076 {
1077 dualIn_ = value;
1078 }
1079 /// Pivot Row for use by classes e.g. steepestedge
pivotRow() const1080 inline int pivotRow() const
1081 {
1082 return pivotRow_;
1083 }
setPivotRow(int value)1084 inline void setPivotRow(int value)
1085 {
1086 pivotRow_ = value;
1087 }
1088 /// value of incoming variable (in Dual)
1089 double valueIncomingDual() const;
1090 //@}
1091
1092 #ifndef CLP_USER_DRIVEN
1093 protected:
1094 #endif
1095 /**@name protected methods */
1096 //@{
1097 /** May change basis and then returns number changed.
1098 Computation of solutions may be overriden by given pi and solution
1099 */
1100 int gutsOfSolution(double *givenDuals,
1101 const double *givenPrimals,
1102 bool valuesPass = false);
1103 /// Does most of deletion (0 = all, 1 = most, 2 most + factorization)
1104 void gutsOfDelete(int type);
1105 /// Does most of copying
1106 void gutsOfCopy(const ClpSimplex &rhs);
1107 /** puts in format I like (rowLower,rowUpper) also see StandardMatrix
1108 1 bit does rows (now and columns), (2 bit does column bounds), 4 bit does objective(s).
1109 8 bit does solution scaling in
1110 16 bit does rowArray and columnArray indexed vectors
1111 and makes row copy if wanted, also sets columnStart_ etc
1112 Also creates scaling arrays if needed. It does scaling if needed.
1113 16 also moves solutions etc in to work arrays
1114 On 16 returns false if problem "bad" i.e. matrix or bounds bad
1115 If startFinishOptions is -1 then called by user in getSolution
1116 so do arrays but keep pivotVariable_
1117 */
1118 bool createRim(int what, bool makeRowCopy = false, int startFinishOptions = 0);
1119 /// Does rows and columns
1120 void createRim1(bool initial);
1121 /// Does objective
1122 void createRim4(bool initial);
1123 /// Does rows and columns and objective
1124 void createRim5(bool initial);
1125 /** releases above arrays and does solution scaling out. May also
1126 get rid of factorization data -
1127 0 get rid of nothing, 1 get rid of arrays, 2 also factorization
1128 */
1129 void deleteRim(int getRidOfFactorizationData = 2);
1130 /// Sanity check on input rim data (after scaling) - returns true if okay
1131 bool sanityCheck();
1132 //@}
1133 public:
1134 /**@name public methods */
1135 //@{
1136 /** Return row or column sections - not as much needed as it
1137 once was. These just map into single arrays */
solutionRegion(int section) const1138 inline double *solutionRegion(int section) const
1139 {
1140 if (!section)
1141 return rowActivityWork_;
1142 else
1143 return columnActivityWork_;
1144 }
djRegion(int section) const1145 inline double *djRegion(int section) const
1146 {
1147 if (!section)
1148 return rowReducedCost_;
1149 else
1150 return reducedCostWork_;
1151 }
lowerRegion(int section) const1152 inline double *lowerRegion(int section) const
1153 {
1154 if (!section)
1155 return rowLowerWork_;
1156 else
1157 return columnLowerWork_;
1158 }
upperRegion(int section) const1159 inline double *upperRegion(int section) const
1160 {
1161 if (!section)
1162 return rowUpperWork_;
1163 else
1164 return columnUpperWork_;
1165 }
costRegion(int section) const1166 inline double *costRegion(int section) const
1167 {
1168 if (!section)
1169 return rowObjectiveWork_;
1170 else
1171 return objectiveWork_;
1172 }
1173 /// Return region as single array
solutionRegion() const1174 inline double *solutionRegion() const
1175 {
1176 return solution_;
1177 }
djRegion() const1178 inline double *djRegion() const
1179 {
1180 return dj_;
1181 }
lowerRegion() const1182 inline double *lowerRegion() const
1183 {
1184 return lower_;
1185 }
upperRegion() const1186 inline double *upperRegion() const
1187 {
1188 return upper_;
1189 }
costRegion() const1190 inline double *costRegion() const
1191 {
1192 return cost_;
1193 }
getStatus(int sequence) const1194 inline Status getStatus(int sequence) const
1195 {
1196 return static_cast< Status >(status_[sequence] & 7);
1197 }
setStatus(int sequence,Status newstatus)1198 inline void setStatus(int sequence, Status newstatus)
1199 {
1200 unsigned char &st_byte = status_[sequence];
1201 st_byte = static_cast< unsigned char >(st_byte & ~7);
1202 st_byte = static_cast< unsigned char >(st_byte | newstatus);
1203 }
1204 /// Start or reset using maximumRows_ and Columns_ - true if change
1205 bool startPermanentArrays();
1206 /** Normally the first factorization does sparse coding because
1207 the factorization could be singular. This allows initial dense
1208 factorization when it is known to be safe
1209 */
1210 void setInitialDenseFactorization(bool onOff);
1211 bool initialDenseFactorization() const;
1212 /** Return sequence In or Out */
sequenceIn() const1213 inline int sequenceIn() const
1214 {
1215 return sequenceIn_;
1216 }
sequenceOut() const1217 inline int sequenceOut() const
1218 {
1219 return sequenceOut_;
1220 }
1221 /** Set sequenceIn or Out */
setSequenceIn(int sequence)1222 inline void setSequenceIn(int sequence)
1223 {
1224 sequenceIn_ = sequence;
1225 }
setSequenceOut(int sequence)1226 inline void setSequenceOut(int sequence)
1227 {
1228 sequenceOut_ = sequence;
1229 }
1230 /** Return direction In or Out */
directionIn() const1231 inline int directionIn() const
1232 {
1233 return directionIn_;
1234 }
directionOut() const1235 inline int directionOut() const
1236 {
1237 return directionOut_;
1238 }
1239 /** Set directionIn or Out */
setDirectionIn(int direction)1240 inline void setDirectionIn(int direction)
1241 {
1242 directionIn_ = direction;
1243 }
setDirectionOut(int direction)1244 inline void setDirectionOut(int direction)
1245 {
1246 directionOut_ = direction;
1247 }
1248 /// Value of Out variable
valueOut() const1249 inline double valueOut() const
1250 {
1251 return valueOut_;
1252 }
1253 /// Lower of out variable
lowerOut() const1254 inline double lowerOut() const
1255 {
1256 return lowerOut_;
1257 }
1258 /// Upper of out variable
upperOut() const1259 inline double upperOut() const
1260 {
1261 return upperOut_;
1262 }
1263 /// Set value of out variable
setValueOut(double value)1264 inline void setValueOut(double value)
1265 {
1266 valueOut_ = value;
1267 }
1268 /// Dual value of Out variable
dualOut() const1269 inline double dualOut() const
1270 {
1271 return dualOut_;
1272 }
1273 /// Set dual value of out variable
setDualOut(double value)1274 inline void setDualOut(double value)
1275 {
1276 dualOut_ = value;
1277 }
1278 /// Set lower of out variable
setLowerOut(double value)1279 inline void setLowerOut(double value)
1280 {
1281 lowerOut_ = value;
1282 }
1283 /// Set upper of out variable
setUpperOut(double value)1284 inline void setUpperOut(double value)
1285 {
1286 upperOut_ = value;
1287 }
1288 /// Set theta of out variable
setTheta(double value)1289 inline void setTheta(double value)
1290 {
1291 theta_ = value;
1292 }
1293 /// Returns 1 if sequence indicates column
isColumn(int sequence) const1294 inline int isColumn(int sequence) const
1295 {
1296 return sequence < numberColumns_ ? 1 : 0;
1297 }
1298 /// Returns sequence number within section
sequenceWithin(int sequence) const1299 inline int sequenceWithin(int sequence) const
1300 {
1301 return sequence < numberColumns_ ? sequence : sequence - numberColumns_;
1302 }
1303 /// Return row or column values
solution(int sequence)1304 inline double solution(int sequence)
1305 {
1306 return solution_[sequence];
1307 }
1308 /// Return address of row or column values
solutionAddress(int sequence)1309 inline double &solutionAddress(int sequence)
1310 {
1311 return solution_[sequence];
1312 }
reducedCost(int sequence)1313 inline double reducedCost(int sequence)
1314 {
1315 return dj_[sequence];
1316 }
reducedCostAddress(int sequence)1317 inline double &reducedCostAddress(int sequence)
1318 {
1319 return dj_[sequence];
1320 }
lower(int sequence)1321 inline double lower(int sequence)
1322 {
1323 return lower_[sequence];
1324 }
1325 /// Return address of row or column lower bound
lowerAddress(int sequence)1326 inline double &lowerAddress(int sequence)
1327 {
1328 return lower_[sequence];
1329 }
upper(int sequence)1330 inline double upper(int sequence)
1331 {
1332 return upper_[sequence];
1333 }
1334 /// Return address of row or column upper bound
upperAddress(int sequence)1335 inline double &upperAddress(int sequence)
1336 {
1337 return upper_[sequence];
1338 }
cost(int sequence)1339 inline double cost(int sequence)
1340 {
1341 return cost_[sequence];
1342 }
1343 /// Return address of row or column cost
costAddress(int sequence)1344 inline double &costAddress(int sequence)
1345 {
1346 return cost_[sequence];
1347 }
1348 /// Return original lower bound
originalLower(int iSequence) const1349 inline double originalLower(int iSequence) const
1350 {
1351 if (iSequence < numberColumns_)
1352 return columnLower_[iSequence];
1353 else
1354 return rowLower_[iSequence - numberColumns_];
1355 }
1356 /// Return original lower bound
originalUpper(int iSequence) const1357 inline double originalUpper(int iSequence) const
1358 {
1359 if (iSequence < numberColumns_)
1360 return columnUpper_[iSequence];
1361 else
1362 return rowUpper_[iSequence - numberColumns_];
1363 }
1364 /// Theta (pivot change)
theta() const1365 inline double theta() const
1366 {
1367 return theta_;
1368 }
1369 /// Lower Bound on In variable
lowerIn() const1370 inline double lowerIn() const
1371 {
1372 return lowerIn_;
1373 }
1374 /// Value of In variable
valueIn() const1375 inline double valueIn() const
1376 {
1377 return valueIn_;
1378 }
1379 /// Upper Bound on In variable
upperIn() const1380 inline double upperIn() const
1381 {
1382 return upperIn_;
1383 }
1384 /** Best possible improvement using djs (primal) or
1385 obj change by flipping bounds to make dual feasible (dual) */
bestPossibleImprovement() const1386 inline double bestPossibleImprovement() const
1387 {
1388 return bestPossibleImprovement_;
1389 }
1390 /// Return pointer to details of costs
nonLinearCost() const1391 inline ClpNonLinearCost *nonLinearCost() const
1392 {
1393 return nonLinearCost_;
1394 }
1395 /// Set pointer to details of costs
1396 void setNonLinearCost(ClpNonLinearCost &nonLinearCost);
1397 /** Return more special options
1398 1 bit - if presolve says infeasible in ClpSolve return
1399 2 bit - if presolved problem infeasible return
1400 4 bit - keep arrays like upper_ around
1401 8 bit - no free or superBasic variables
1402 16 bit - if checking replaceColumn accuracy before updating
1403 32 bit - say optimal if primal feasible!
1404 64 bit - give up easily in dual (and say infeasible)
1405 128 bit - no objective, 0-1 and in B&B
1406 256 bit - in primal from dual or vice versa
1407 512 bit - alternative use of solveType_
1408 1024 bit - don't do row copy of factorization
1409 2048 bit - perturb in complete fathoming
1410 4096 bit - try more for complete fathoming
1411 8192 bit - don't even think of using primal if user asks for dual (and vv)
1412 16384 bit - in initialSolve so be more flexible
1413 32768 bit - don't swap algorithms from dual if small infeasibility
1414 65536 bit - perturb in postsolve cleanup (even if < 10000 rows)
1415 131072 bit (*3) initial stateDualColumn
1416 524288 bit - stop when primal feasible
1417 1048576 bit - stop when primal feasible after n-1000000 iterations
1418 2097152 bit - no primal in fastDual2 if feasible
1419 4194304 bit - tolerances have been changed by code
1420 8388608 bit - tolerances are dynamic (at first)
1421 16777216 bit - if factorization kept can still declare optimal at once
1422 */
moreSpecialOptions() const1423 inline int moreSpecialOptions() const
1424 {
1425 return moreSpecialOptions_;
1426 }
1427 /// Get vector mode
vectorMode() const1428 inline int vectorMode() const
1429 {
1430 return vectorMode_;
1431 }
1432 /** Set more special options
1433 1 bit - if presolve says infeasible in ClpSolve return
1434 2 bit - if presolved problem infeasible return
1435 4 bit - keep arrays like upper_ around
1436 8 bit - no free or superBasic variables
1437 16 bit - if checking replaceColumn accuracy before updating
1438 32 bit - say optimal if primal feasible!
1439 64 bit - give up easily in dual (and say infeasible)
1440 128 bit - no objective, 0-1 and in B&B
1441 256 bit - in primal from dual or vice versa
1442 512 bit - alternative use of solveType_
1443 1024 bit - don't do row copy of factorization
1444 2048 bit - perturb in complete fathoming
1445 4096 bit - try more for complete fathoming
1446 8192 bit - don't even think of using primal if user asks for dual (and vv)
1447 16384 bit - in initialSolve so be more flexible
1448 32768 bit - don't swap algorithms from dual if small infeasibility
1449 65536 bit - perturb in postsolve cleanup (even if < 10000 rows)
1450 131072 bit (*3) initial stateDualColumn
1451 524288 bit - stop when primal feasible
1452 1048576 bit - don't perturb even if long time
1453 2097152 bit - no primal in fastDual2 if feasible
1454 4194304 bit - tolerances have been changed by code
1455 8388608 bit - tolerances are dynamic (at first)
1456 16777216 bit - if factorization kept can still declare optimal at once
1457 */
setMoreSpecialOptions(int value)1458 inline void setMoreSpecialOptions(int value)
1459 {
1460 moreSpecialOptions_ = value;
1461 }
1462 /// Set vector mode
setVectorMode(int value)1463 inline void setVectorMode(int value)
1464 {
1465 vectorMode_ = value;
1466 }
1467 //@}
1468 /**@name status methods */
1469 //@{
setFakeBound(int sequence,FakeBound fakeBound)1470 inline void setFakeBound(int sequence, FakeBound fakeBound)
1471 {
1472 unsigned char &st_byte = status_[sequence];
1473 st_byte = static_cast< unsigned char >(st_byte & ~24);
1474 st_byte = static_cast< unsigned char >(st_byte | (fakeBound << 3));
1475 }
getFakeBound(int sequence) const1476 inline FakeBound getFakeBound(int sequence) const
1477 {
1478 return static_cast< FakeBound >((status_[sequence] >> 3) & 3);
1479 }
setRowStatus(int sequence,Status newstatus)1480 inline void setRowStatus(int sequence, Status newstatus)
1481 {
1482 unsigned char &st_byte = status_[sequence + numberColumns_];
1483 st_byte = static_cast< unsigned char >(st_byte & ~7);
1484 st_byte = static_cast< unsigned char >(st_byte | newstatus);
1485 }
getRowStatus(int sequence) const1486 inline Status getRowStatus(int sequence) const
1487 {
1488 return static_cast< Status >(status_[sequence + numberColumns_] & 7);
1489 }
setColumnStatus(int sequence,Status newstatus)1490 inline void setColumnStatus(int sequence, Status newstatus)
1491 {
1492 unsigned char &st_byte = status_[sequence];
1493 st_byte = static_cast< unsigned char >(st_byte & ~7);
1494 st_byte = static_cast< unsigned char >(st_byte | newstatus);
1495 }
getColumnStatus(int sequence) const1496 inline Status getColumnStatus(int sequence) const
1497 {
1498 return static_cast< Status >(status_[sequence] & 7);
1499 }
setPivoted(int sequence)1500 inline void setPivoted(int sequence)
1501 {
1502 status_[sequence] = static_cast< unsigned char >(status_[sequence] | 32);
1503 }
clearPivoted(int sequence)1504 inline void clearPivoted(int sequence)
1505 {
1506 status_[sequence] = static_cast< unsigned char >(status_[sequence] & ~32);
1507 }
pivoted(int sequence) const1508 inline bool pivoted(int sequence) const
1509 {
1510 return (((status_[sequence] >> 5) & 1) != 0);
1511 }
1512 /// To flag a variable (not inline to allow for column generation)
1513 void setFlagged(int sequence);
clearFlagged(int sequence)1514 inline void clearFlagged(int sequence)
1515 {
1516 status_[sequence] = static_cast< unsigned char >(status_[sequence] & ~64);
1517 }
flagged(int sequence) const1518 inline bool flagged(int sequence) const
1519 {
1520 return ((status_[sequence] & 64) != 0);
1521 }
1522 /// To say row active in primal pivot row choice
setActive(int iRow)1523 inline void setActive(int iRow)
1524 {
1525 status_[iRow] = static_cast< unsigned char >(status_[iRow] | 128);
1526 }
clearActive(int iRow)1527 inline void clearActive(int iRow)
1528 {
1529 status_[iRow] = static_cast< unsigned char >(status_[iRow] & ~128);
1530 }
active(int iRow) const1531 inline bool active(int iRow) const
1532 {
1533 return ((status_[iRow] & 128) != 0);
1534 }
1535 /// To say perturbed
setPerturbed(int iSequence)1536 inline void setPerturbed(int iSequence)
1537 {
1538 status_[iSequence] = static_cast< unsigned char >(status_[iSequence] | 128);
1539 }
clearPerturbed(int iSequence)1540 inline void clearPerturbed(int iSequence)
1541 {
1542 status_[iSequence] = static_cast< unsigned char >(status_[iSequence] & ~128);
1543 }
perturbed(int iSequence) const1544 inline bool perturbed(int iSequence) const
1545 {
1546 return ((status_[iSequence] & 128) != 0);
1547 }
1548 /** Set up status array (can be used by OsiClp).
1549 Also can be used to set up all slack basis */
1550 void createStatus();
1551 /** Sets up all slack basis and resets solution to
1552 as it was after initial load or readMps */
1553 void allSlackBasis(bool resetSolution = false);
1554
1555 /// So we know when to be cautious
lastBadIteration() const1556 inline int lastBadIteration() const
1557 {
1558 return lastBadIteration_;
1559 }
1560 /// Set so we know when to be cautious
setLastBadIteration(int value)1561 inline void setLastBadIteration(int value)
1562 {
1563 lastBadIteration_ = value;
1564 }
1565 /// Progress flag - at present 0 bit says artificials out
progressFlag() const1566 inline int progressFlag() const
1567 {
1568 return (progressFlag_ & 3);
1569 }
1570 /// For dealing with all issues of cycling etc
progress()1571 inline ClpSimplexProgress *progress()
1572 {
1573 return &progress_;
1574 }
1575 /// Force re-factorization early value
forceFactorization() const1576 inline int forceFactorization() const
1577 {
1578 return forceFactorization_;
1579 }
1580 /// Force re-factorization early
forceFactorization(int value)1581 inline void forceFactorization(int value)
1582 {
1583 forceFactorization_ = value;
1584 }
1585 /// Raw objective value (so always minimize in primal)
rawObjectiveValue() const1586 inline double rawObjectiveValue() const
1587 {
1588 return objectiveValue_;
1589 }
1590 /// Compute objective value from solution and put in objectiveValue_
1591 void computeObjectiveValue(bool useWorkingSolution = false);
1592 /// Compute minimization objective value from internal solution without perturbation
1593 double computeInternalObjectiveValue();
1594 /** Infeasibility/unbounded ray (NULL returned if none/wrong)
1595 Up to user to use delete [] on these arrays. */
1596 double *infeasibilityRay(bool fullRay = false) const;
1597 /** Number of extra rows. These are ones which will be dynamically created
1598 each iteration. This is for GUB but may have other uses.
1599 */
numberExtraRows() const1600 inline int numberExtraRows() const
1601 {
1602 return numberExtraRows_;
1603 }
1604 /** Maximum number of basic variables - can be more than number of rows if GUB
1605 */
maximumBasic() const1606 inline int maximumBasic() const
1607 {
1608 return maximumBasic_;
1609 }
1610 /// Iteration when we entered dual or primal
baseIteration() const1611 inline int baseIteration() const
1612 {
1613 return baseIteration_;
1614 }
1615 /// Create C++ lines to get to current state
1616 void generateCpp(FILE *fp, bool defaultFactor = false);
1617 /// Gets clean and emptyish factorization
1618 ClpFactorization *getEmptyFactorization();
1619 /// May delete or may make clean and emptyish factorization
1620 void setEmptyFactorization();
1621 /// Move status and solution across
1622 void moveInfo(const ClpSimplex &rhs, bool justStatus = false);
1623 //@}
1624
1625 ///@name Basis handling
1626 // These are only to be used using startFinishOptions (ClpSimplexDual, ClpSimplexPrimal)
1627 // *** At present only without scaling
1628 // *** Slacks havve -1.0 element (so == row activity) - take care
1629 ///Get a row of the tableau (slack part in slack if not NULL)
1630 void getBInvARow(int row, double *z, double *slack = NULL);
1631
1632 ///Get a row of the basis inverse
1633 void getBInvRow(int row, double *z);
1634
1635 ///Get a column of the tableau
1636 void getBInvACol(int col, double *vec);
1637
1638 ///Get a column of the basis inverse
1639 void getBInvCol(int col, double *vec);
1640
1641 /** Get basic indices (order of indices corresponds to the
1642 order of elements in a vector retured by getBInvACol() and
1643 getBInvCol()).
1644 */
1645 void getBasics(int *index);
1646
1647 //@}
1648 //-------------------------------------------------------------------------
1649 /**@name Changing bounds on variables and constraints */
1650 //@{
1651 /** Set an objective function coefficient */
1652 void setObjectiveCoefficient(int elementIndex, double elementValue);
1653 /** Set an objective function coefficient */
setObjCoeff(int elementIndex,double elementValue)1654 inline void setObjCoeff(int elementIndex, double elementValue)
1655 {
1656 setObjectiveCoefficient(elementIndex, elementValue);
1657 }
1658
1659 /** Set a single column lower bound<br>
1660 Use -DBL_MAX for -infinity. */
1661 void setColumnLower(int elementIndex, double elementValue);
1662
1663 /** Set a single column upper bound<br>
1664 Use DBL_MAX for infinity. */
1665 void setColumnUpper(int elementIndex, double elementValue);
1666
1667 /** Set a single column lower and upper bound */
1668 void setColumnBounds(int elementIndex,
1669 double lower, double upper);
1670
1671 /** Set the bounds on a number of columns simultaneously<br>
1672 The default implementation just invokes setColLower() and
1673 setColUpper() over and over again.
1674 @param indexFirst,indexLast pointers to the beginning and after the
1675 end of the array of the indices of the variables whose
1676 <em>either</em> bound changes
1677 @param boundList the new lower/upper bound pairs for the variables
1678 */
1679 void setColumnSetBounds(const int *indexFirst,
1680 const int *indexLast,
1681 const double *boundList);
1682
1683 /** Set a single column lower bound<br>
1684 Use -DBL_MAX for -infinity. */
setColLower(int elementIndex,double elementValue)1685 inline void setColLower(int elementIndex, double elementValue)
1686 {
1687 setColumnLower(elementIndex, elementValue);
1688 }
1689 /** Set a single column upper bound<br>
1690 Use DBL_MAX for infinity. */
setColUpper(int elementIndex,double elementValue)1691 inline void setColUpper(int elementIndex, double elementValue)
1692 {
1693 setColumnUpper(elementIndex, elementValue);
1694 }
1695
1696 /** Set a single column lower and upper bound */
setColBounds(int elementIndex,double newlower,double newupper)1697 inline void setColBounds(int elementIndex,
1698 double newlower, double newupper)
1699 {
1700 setColumnBounds(elementIndex, newlower, newupper);
1701 }
1702
1703 /** Set the bounds on a number of columns simultaneously<br>
1704 @param indexFirst,indexLast pointers to the beginning and after the
1705 end of the array of the indices of the variables whose
1706 <em>either</em> bound changes
1707 @param boundList the new lower/upper bound pairs for the variables
1708 */
setColSetBounds(const int * indexFirst,const int * indexLast,const double * boundList)1709 inline void setColSetBounds(const int *indexFirst,
1710 const int *indexLast,
1711 const double *boundList)
1712 {
1713 setColumnSetBounds(indexFirst, indexLast, boundList);
1714 }
1715
1716 /** Set a single row lower bound<br>
1717 Use -DBL_MAX for -infinity. */
1718 void setRowLower(int elementIndex, double elementValue);
1719
1720 /** Set a single row upper bound<br>
1721 Use DBL_MAX for infinity. */
1722 void setRowUpper(int elementIndex, double elementValue);
1723
1724 /** Set a single row lower and upper bound */
1725 void setRowBounds(int elementIndex,
1726 double lower, double upper);
1727
1728 /** Set the bounds on a number of rows simultaneously<br>
1729 @param indexFirst,indexLast pointers to the beginning and after the
1730 end of the array of the indices of the constraints whose
1731 <em>either</em> bound changes
1732 @param boundList the new lower/upper bound pairs for the constraints
1733 */
1734 void setRowSetBounds(const int *indexFirst,
1735 const int *indexLast,
1736 const double *boundList);
1737 /// Resizes rim part of model
1738 void resize(int newNumberRows, int newNumberColumns);
1739
1740 //@}
1741
1742 ////////////////// data //////////////////
1743 protected:
1744 /**@name data. Many arrays have a row part and a column part.
1745 There is a single array with both - columns then rows and
1746 then normally two arrays pointing to rows and columns. The
1747 single array is the owner of memory
1748 */
1749 //@{
1750 /** Best possible improvement using djs (primal) or
1751 obj change by flipping bounds to make dual feasible (dual) */
1752 double bestPossibleImprovement_;
1753 /// Zero tolerance
1754 double zeroTolerance_;
1755 /// Sequence of worst (-1 if feasible)
1756 int columnPrimalSequence_;
1757 /// Sequence of worst (-1 if feasible)
1758 int rowPrimalSequence_;
1759 /// "Best" objective value
1760 double bestObjectiveValue_;
1761 /// More special options - see set for details
1762 int moreSpecialOptions_;
1763 /// Iteration when we entered dual or primal
1764 int baseIteration_;
1765 /// Vector mode - try and use vector instructions
1766 int vectorMode_;
1767 /// Primal tolerance needed to make dual feasible (<largeTolerance)
1768 double primalToleranceToGetOptimal_;
1769 /// Large bound value (for complementarity etc)
1770 double largeValue_;
1771 /// Largest error on Ax-b
1772 double largestPrimalError_;
1773 /// Largest error on basic duals
1774 double largestDualError_;
1775 /// For computing whether to re-factorize
1776 double alphaAccuracy_;
1777 /// Dual bound
1778 double dualBound_;
1779 /// Alpha (pivot element)
1780 double alpha_;
1781 /// Theta (pivot change)
1782 double theta_;
1783 /// Lower Bound on In variable
1784 double lowerIn_;
1785 /// Value of In variable
1786 double valueIn_;
1787 /// Upper Bound on In variable
1788 double upperIn_;
1789 /// Reduced cost of In variable
1790 double dualIn_;
1791 /// Lower Bound on Out variable
1792 double lowerOut_;
1793 /// Value of Out variable
1794 double valueOut_;
1795 /// Upper Bound on Out variable
1796 double upperOut_;
1797 /// Infeasibility (dual) or ? (primal) of Out variable
1798 double dualOut_;
1799 /// Current dual tolerance for algorithm
1800 double dualTolerance_;
1801 /// Current primal tolerance for algorithm
1802 double primalTolerance_;
1803 /// Sum of dual infeasibilities
1804 double sumDualInfeasibilities_;
1805 /// Sum of primal infeasibilities
1806 double sumPrimalInfeasibilities_;
1807 /// Weight assigned to being infeasible in primal
1808 double infeasibilityCost_;
1809 /// Sum of Dual infeasibilities using tolerance based on error in duals
1810 double sumOfRelaxedDualInfeasibilities_;
1811 /// Sum of Primal infeasibilities using tolerance based on error in primals
1812 double sumOfRelaxedPrimalInfeasibilities_;
1813 /// Acceptable pivot value just after factorization
1814 double acceptablePivot_;
1815 /// Minimum primal tolerance
1816 double minimumPrimalTolerance_;
1817 /// Last few infeasibilities
1818 #define CLP_INFEAS_SAVE 5
1819 double averageInfeasibility_[CLP_INFEAS_SAVE];
1820 /// Working copy of lower bounds (Owner of arrays below)
1821 double *lower_;
1822 /// Row lower bounds - working copy
1823 double *rowLowerWork_;
1824 /// Column lower bounds - working copy
1825 double *columnLowerWork_;
1826 /// Working copy of upper bounds (Owner of arrays below)
1827 double *upper_;
1828 /// Row upper bounds - working copy
1829 double *rowUpperWork_;
1830 /// Column upper bounds - working copy
1831 double *columnUpperWork_;
1832 /// Working copy of objective (Owner of arrays below)
1833 double *cost_;
1834 /// Row objective - working copy
1835 double *rowObjectiveWork_;
1836 /// Column objective - working copy
1837 double *objectiveWork_;
1838 /// Useful row length arrays
1839 CoinIndexedVector *rowArray_[6];
1840 /// Useful column length arrays
1841 CoinIndexedVector *columnArray_[6];
1842 /// Sequence of In variable
1843 int sequenceIn_;
1844 /// Direction of In, 1 going up, -1 going down, 0 not a clue
1845 int directionIn_;
1846 /// Sequence of Out variable
1847 int sequenceOut_;
1848 /// Direction of Out, 1 to upper bound, -1 to lower bound, 0 - superbasic
1849 int directionOut_;
1850 /// Pivot Row
1851 int pivotRow_;
1852 /// Last good iteration (immediately after a re-factorization)
1853 int lastGoodIteration_;
1854 /// Working copy of reduced costs (Owner of arrays below)
1855 double *dj_;
1856 /// Reduced costs of slacks not same as duals (or - duals)
1857 double *rowReducedCost_;
1858 /// Possible scaled reduced costs
1859 double *reducedCostWork_;
1860 /// Working copy of primal solution (Owner of arrays below)
1861 double *solution_;
1862 /// Row activities - working copy
1863 double *rowActivityWork_;
1864 /// Column activities - working copy
1865 double *columnActivityWork_;
1866 /// Number of dual infeasibilities
1867 int numberDualInfeasibilities_;
1868 /// Number of dual infeasibilities (without free)
1869 int numberDualInfeasibilitiesWithoutFree_;
1870 /// Number of primal infeasibilities
1871 int numberPrimalInfeasibilities_;
1872 /// How many iterative refinements to do
1873 int numberRefinements_;
1874 /// dual row pivot choice
1875 ClpDualRowPivot *dualRowPivot_;
1876 /// primal column pivot choice
1877 ClpPrimalColumnPivot *primalColumnPivot_;
1878 /// Basic variables pivoting on which rows
1879 int *pivotVariable_;
1880 /// factorization
1881 ClpFactorization *factorization_;
1882 /// Saved version of solution
1883 double *savedSolution_;
1884 /// Number of times code has tentatively thought optimal
1885 int numberTimesOptimal_;
1886 /// Disaster handler
1887 ClpDisasterHandler *disasterArea_;
1888 /// If change has been made (first attempt at stopping looping)
1889 int changeMade_;
1890 /// Algorithm >0 == Primal, <0 == Dual
1891 int algorithm_;
1892 /** Now for some reliability aids
1893 This forces re-factorization early */
1894 int forceFactorization_;
1895 /** Perturbation:
1896 -50 to +50 - perturb by this power of ten (-6 sounds good)
1897 100 - auto perturb if takes too long (1.0e-6 largest nonzero)
1898 101 - we are perturbed
1899 102 - don't try perturbing again
1900 default is 100
1901 */
1902 int perturbation_;
1903 /// Saved status regions
1904 unsigned char *saveStatus_;
1905 /** Very wasteful way of dealing with infeasibilities in primal.
1906 However it will allow non-linearities and use of dual
1907 analysis. If it doesn't work it can easily be replaced.
1908 */
1909 ClpNonLinearCost *nonLinearCost_;
1910 /// So we know when to be cautious
1911 int lastBadIteration_;
1912 /// So we know when to open up again
1913 int lastFlaggedIteration_;
1914 /// Can be used for count of fake bounds (dual) or fake costs (primal)
1915 int numberFake_;
1916 /// Can be used for count of changed costs (dual) or changed bounds (primal)
1917 int numberChanged_;
1918 /// Progress flag - at present 0 bit says artificials out, 1 free in
1919 int progressFlag_;
1920 /// First free/super-basic variable (-1 if none)
1921 int firstFree_;
1922 /** Number of extra rows. These are ones which will be dynamically created
1923 each iteration. This is for GUB but may have other uses.
1924 */
1925 int numberExtraRows_;
1926 /** Maximum number of basic variables - can be more than number of rows if GUB
1927 */
1928 int maximumBasic_;
1929 /// If may skip final factorize then allow up to this pivots (default 20)
1930 int dontFactorizePivots_;
1931 /** For advanced use. When doing iterative solves things can get
1932 nasty so on values pass if incoming solution has largest
1933 infeasibility < incomingInfeasibility throw out variables
1934 from basis until largest infeasibility < allowedInfeasibility.
1935 if allowedInfeasibility>= incomingInfeasibility this is
1936 always possible altough you may end up with an all slack basis.
1937
1938 Defaults are 1.0,10.0
1939 */
1940 double incomingInfeasibility_;
1941 double allowedInfeasibility_;
1942 /// Automatic scaling of objective and rhs and bounds
1943 int automaticScale_;
1944 /// Maximum perturbation array size (take out when code rewritten)
1945 int maximumPerturbationSize_;
1946 /// Perturbation array (maximumPerturbationSize_)
1947 double *perturbationArray_;
1948 /// A copy of model with certain state - normally without cuts
1949 ClpSimplex *baseModel_;
1950 /// For dealing with all issues of cycling etc
1951 ClpSimplexProgress progress_;
1952 #ifdef ABC_INHERIT
1953 AbcSimplex *abcSimplex_;
1954 int abcState_;
1955 #define CLP_ABC_WANTED 1
1956 #define CLP_ABC_WANTED_PARALLEL 2
1957 #define CLP_ABC_FULL_DONE 8
1958 // bits 256,512,1024 for crash
1959 #endif
1960 #define CLP_ABC_BEEN_FEASIBLE 65536
1961 /// Number of degenerate pivots since last perturbed
1962 int numberDegeneratePivots_;
1963
1964 public:
1965 /// Spare int array for passing information [0]!=0 switches on
1966 mutable int spareIntArray_[4];
1967 /// Spare double array for passing information [0]!=0 switches on
1968 mutable double spareDoubleArray_[4];
1969
1970 protected:
1971 /// Allow OsiClp certain perks
1972 friend class OsiClpSolverInterface;
1973 /// And OsiCLP
1974 friend class OsiCLPSolverInterface;
1975 //@}
1976 };
1977 //#############################################################################
1978 /** A function that tests the methods in the ClpSimplex class. The
1979 only reason for it not to be a member method is that this way it doesn't
1980 have to be compiled into the library. And that's a gain, because the
1981 library should be compiled with optimization on, but this method should be
1982 compiled with debugging.
1983
1984 It also does some testing of ClpFactorization class
1985 */
1986 void ClpSimplexUnitTest(const std::string &mpsDir);
1987
1988 // For Devex stuff
1989 #define DEVEX_TRY_NORM 1.0e-4
1990 #define DEVEX_ADD_ONE 1.0
1991 #if defined(ABC_INHERIT) || defined(CBC_THREAD) || defined(THREADS_IN_ANALYZE)
1992 // Use pthreads
1993 #include <pthread.h>
1994 typedef struct {
1995 double result;
1996 //const CoinIndexedVector * constVector; // can get rid of
1997 //CoinIndexedVector * vectors[2]; // can get rid of
1998 void *extraInfo;
1999 void *extraInfo2;
2000 int status;
2001 int stuff[4];
2002 } CoinThreadInfo;
2003 class CoinPthreadStuff {
2004 public:
2005 /**@name Constructors and destructor and copy */
2006 //@{
2007 /** Main constructor
2008 */
2009 CoinPthreadStuff(int numberThreads = 0,
2010 void *parallelManager(void *stuff) = NULL);
2011 /// Assignment operator. This copies the data
2012 CoinPthreadStuff &operator=(const CoinPthreadStuff &rhs);
2013 /// Destructor
2014 ~CoinPthreadStuff();
2015 /// set stop start
setStopStart(int value)2016 inline void setStopStart(int value)
2017 {
2018 stopStart_ = value;
2019 }
2020 #ifndef NUMBER_THREADS
2021 #define NUMBER_THREADS 8
2022 #endif
2023 // For waking up thread
mutexPointer(int which,int thread=0)2024 inline pthread_mutex_t *mutexPointer(int which, int thread = 0)
2025 {
2026 return mutex_ + which + 3 * thread;
2027 }
2028 #ifdef PTHREAD_BARRIER_SERIAL_THREAD
barrierPointer()2029 inline pthread_barrier_t *barrierPointer()
2030 {
2031 return &barrier_;
2032 }
2033 #endif
whichLocked(int thread=0) const2034 inline int whichLocked(int thread = 0) const
2035 {
2036 return locked_[thread];
2037 }
threadInfoPointer(int thread=0)2038 inline CoinThreadInfo *threadInfoPointer(int thread = 0)
2039 {
2040 return threadInfo_ + thread;
2041 }
2042 void startParallelTask(int type, int iThread, void *info = NULL);
2043 int waitParallelTask(int type, int &iThread, bool allowIdle);
2044 void waitAllTasks();
2045 /// so thread can find out which one it is
2046 int whichThread() const;
2047 void sayIdle(int iThread);
2048 //void startThreads(int numberThreads);
2049 //void stopThreads();
2050 // For waking up thread
2051 pthread_mutex_t mutex_[3 * (NUMBER_THREADS + 1)];
2052 #ifdef PTHREAD_BARRIER_SERIAL_THREAD
2053 pthread_barrier_t barrier_;
2054 #endif
2055 CoinThreadInfo threadInfo_[NUMBER_THREADS + 1];
2056 pthread_t abcThread_[NUMBER_THREADS + 1];
2057 int locked_[NUMBER_THREADS + 1];
2058 int stopStart_;
2059 int numberThreads_;
2060 };
2061 void *clp_parallelManager(void *stuff);
2062 #endif
2063 typedef struct {
2064 double upperTheta;
2065 double bestPossible;
2066 double acceptablePivot;
2067 double tolerance;
2068 double dualTolerance;
2069 double theta;
2070 double primalRatio;
2071 double changeObj;
2072 const double *COIN_RESTRICT cost;
2073 double *COIN_RESTRICT solution;
2074 double *COIN_RESTRICT reducedCost;
2075 const double *COIN_RESTRICT lower;
2076 const double *COIN_RESTRICT upper;
2077 double *COIN_RESTRICT work;
2078 int *COIN_RESTRICT index;
2079 double *COIN_RESTRICT spare;
2080 const unsigned char *COIN_RESTRICT status;
2081 int *COIN_RESTRICT which;
2082 double *COIN_RESTRICT infeas;
2083 const int *COIN_RESTRICT pivotVariable;
2084 const double *COIN_RESTRICT element;
2085 const CoinBigIndex *COIN_RESTRICT start;
2086 const int *COIN_RESTRICT row;
2087 int numberAdded;
2088 int numberInfeasibilities;
2089 int numberRemaining;
2090 int startColumn;
2091 int numberToDo;
2092 int numberColumns;
2093 } clpTempInfo;
2094 #ifndef ABC_INHERIT
2095 #if ABOCA_LITE
2096 void moveAndZero(clpTempInfo *info, int type, void *extra);
2097 // 2 is owner of abcState_
2098 #ifdef ABCSTATE_LITE
2099 #if ABCSTATE_LITE == 2
2100 int abcState_ = 0;
2101 #else
2102 extern int abcState_;
2103 #endif
abcState()2104 inline int abcState()
2105 {
2106 return abcState_;
2107 }
setAbcState(int state)2108 inline void setAbcState(int state)
2109 {
2110 abcState_ = state;
2111 }
2112 #endif
2113 #else
2114 #define abcState 0
2115 #endif
2116 #endif
2117 #ifdef CLP_USER_DRIVEN
2118 // expand as needed
2119 typedef struct {
2120 double alpha;
2121 double totalThru;
2122 double rhs;
2123 double value;
2124 double lower;
2125 double upper;
2126 double cost;
2127 int type;
2128 int row;
2129 int sequence;
2130 int printing;
2131 int change;
2132 } clpUserStruct;
2133 #endif
2134 #endif
2135
2136 /* vi: softtabstop=2 shiftwidth=2 expandtab tabstop=2
2137 */
2138