source: trunk/Clp/src/AbcSimplexDual.hpp @ 2385

Last change on this file since 2385 was 2385, checked in by unxusr, 3 months ago

formatting

  • Property svn:keywords set to Id
File size: 11.6 KB
Line 
1/* $Id: AbcSimplexDual.hpp 2385 2019-01-06 19:43:06Z unxusr $ */
2// Copyright (C) 2002, International Business Machines
3// Corporation and others, Copyright (C) 2012, FasterCoin.  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 AbcSimplexDual_H
12#define AbcSimplexDual_H
13
14#include "AbcSimplex.hpp"
15#if 0
16#undef ABC_PARALLEL
17#define ABC_PARALLEL 2
18#undef cilk_for
19#undef cilk_spawn
20#undef cilk_sync
21#include <cilk/cilk.h>
22#endif
23typedef struct {
24  double theta;
25  double totalThru;
26  double useThru;
27  double bestEverPivot;
28  double increaseInObjective;
29  double tentativeTheta;
30  double lastPivotValue;
31  double thisPivotValue;
32  double thruThis;
33  double increaseInThis;
34  int lastSequence;
35  int sequence;
36  int block;
37  int numberSwapped;
38  int numberRemaining;
39  int numberLastSwapped;
40  bool modifyCosts;
41} dualColumnResult;
42/** This solves LPs using the dual simplex method
43   
44    It inherits from AbcSimplex.  It has no data of its own and
45    is never created - only cast from a AbcSimplex object at algorithm time.
46   
47*/
48
49class AbcSimplexDual : public AbcSimplex {
50
51public:
52  /**@name Description of algorithm */
53  //@{
54  /** Dual algorithm
55     
56      Method
57     
58      It tries to be a single phase approach with a weight of 1.0 being
59      given to getting optimal and a weight of updatedDualBound_ being
60      given to getting dual feasible.  In this version I have used the
61      idea that this weight can be thought of as a fake bound.  If the
62      distance between the lower and upper bounds on a variable is less
63      than the feasibility weight then we are always better off flipping
64      to other bound to make dual feasible.  If the distance is greater
65      then we make up a fake bound updatedDualBound_ away from one bound.
66      If we end up optimal or primal infeasible, we check to see if
67      bounds okay.  If so we have finished, if not we increase updatedDualBound_
68      and continue (after checking if unbounded). I am undecided about
69      free variables - there is coding but I am not sure about it.  At
70      present I put them in basis anyway.
71     
72      The code is designed to take advantage of sparsity so arrays are
73      seldom zeroed out from scratch or gone over in their entirety.
74      The only exception is a full scan to find outgoing variable for
75      Dantzig row choice.  For steepest edge we keep an updated list
76      of infeasibilities (actually squares).
77      On easy problems we don't need full scan - just
78      pick first reasonable.
79     
80      One problem is how to tackle degeneracy and accuracy.  At present
81      I am using the modification of costs which I put in OSL and some
82      of what I think is the dual analog of Gill et al.
83      I am still not sure of the exact details.
84     
85      The flow of dual is three while loops as follows:
86     
87      while (not finished) {
88     
89      while (not clean solution) {
90     
91      Factorize and/or clean up solution by flipping variables so
92      dual feasible.  If looks finished check fake dual bounds.
93      Repeat until status is iterating (-1) or finished (0,1,2)
94     
95      }
96     
97      while (status==-1) {
98     
99      Iterate until no pivot in or out or time to re-factorize.
100     
101      Flow is:
102     
103      choose pivot row (outgoing variable).  if none then
104      we are primal feasible so looks as if done but we need to
105      break and check bounds etc.
106     
107      Get pivot row in tableau
108     
109      Choose incoming column.  If we don't find one then we look
110      primal infeasible so break and check bounds etc.  (Also the
111      pivot tolerance is larger after any iterations so that may be
112      reason)
113     
114      If we do find incoming column, we may have to adjust costs to
115      keep going forwards (anti-degeneracy).  Check pivot will be stable
116      and if unstable throw away iteration and break to re-factorize.
117      If minor error re-factorize after iteration.
118     
119      Update everything (this may involve flipping variables to stay
120      dual feasible.
121     
122      }
123     
124      }
125     
126      TODO's (or maybe not)
127     
128      At present we never check we are going forwards.  I overdid that in
129      OSL so will try and make a last resort.
130     
131      Needs partial scan pivot out option.
132     
133      May need other anti-degeneracy measures, especially if we try and use
134      loose tolerances as a way to solve in fewer iterations.
135     
136      I like idea of dynamic scaling.  This gives opportunity to decouple
137      different implications of scaling for accuracy, iteration count and
138      feasibility tolerance.
139     
140      for use of exotic parameter startFinishoptions see Abcsimplex.hpp
141  */
142
143  int dual();
144  /** For strong branching.  On input lower and upper are new bounds
145      while on output they are change in objective function values
146      (>1.0e50 infeasible).
147      Return code is 0 if nothing interesting, -1 if infeasible both
148      ways and +1 if infeasible one way (check values to see which one(s))
149      Solutions are filled in as well - even down, odd up - also
150      status and number of iterations
151  */
152  int strongBranching(int numberVariables, const int *variables,
153    double *newLower, double *newUpper,
154    double **outputSolution,
155    int *outputStatus, int *outputIterations,
156    bool stopOnFirstInfeasible = true,
157    bool alwaysFinish = false,
158    int startFinishOptions = 0);
159  /// This does first part of StrongBranching
160  AbcSimplexFactorization *setupForStrongBranching(char *arrays, int numberRows,
161    int numberColumns, bool solveLp = false);
162  /// This cleans up after strong branching
163  void cleanupAfterStrongBranching(AbcSimplexFactorization *factorization);
164  //@}
165
166  /**@name Functions used in dual */
167  //@{
168  /** This has the flow between re-factorizations
169      Broken out for clarity and will be used by strong branching
170     
171      Reasons to come out:
172      -1 iterations etc
173      -2 inaccuracy
174      -3 slight inaccuracy (and done iterations)
175      +0 looks optimal (might be unbounded - but we will investigate)
176      +1 looks infeasible
177      +3 max iterations
178     
179      If givenPi not NULL then in values pass (copy from ClpSimplexDual)
180  */
181  int whileIteratingSerial();
182#if ABC_PARALLEL == 1
183  int whileIteratingThread();
184#endif
185#if ABC_PARALLEL == 2
186  int whileIteratingCilk();
187#endif
188  void whileIterating2();
189  int whileIteratingParallel(int numberIterations);
190  int whileIterating3();
191  void updatePrimalSolution();
192  int noPivotRow();
193  int noPivotColumn();
194  void dualPivotColumn();
195  /// Create dual pricing vector
196  void createDualPricingVectorSerial();
197  int getTableauColumnFlipAndStartReplaceSerial();
198  void getTableauColumnPart1Serial();
199#if ABC_PARALLEL == 1
200  /// Create dual pricing vector
201  void createDualPricingVectorThread();
202  int getTableauColumnFlipAndStartReplaceThread();
203  void getTableauColumnPart1Thread();
204#endif
205#if ABC_PARALLEL == 2
206  /// Create dual pricing vector
207  void createDualPricingVectorCilk();
208  int getTableauColumnFlipAndStartReplaceCilk();
209  void getTableauColumnPart1Cilk();
210#endif
211  void getTableauColumnPart2();
212  int checkReplace();
213  void replaceColumnPart3();
214  void checkReplacePart1();
215  void checkReplacePart1a();
216  void checkReplacePart1b();
217  /// The duals are updated
218  void updateDualsInDual();
219  /** The duals are updated by the given arrays.
220      This is in values pass - so no changes to primal is made
221  */
222  //void updateDualsInValuesPass(CoinIndexedVector * array,
223  //                           double theta);
224  /** While dualColumn gets flips
225      this does actual flipping.
226      returns number flipped
227  */
228  int flipBounds();
229  /** Undo a flip
230  */
231  void flipBack(int number);
232  /**
233     Array has tableau row (row section)
234     Puts candidates for rows in list
235     Returns guess at upper theta (infinite if no pivot) and may set sequenceIn_ if free
236     Can do all (if tableauRow created)
237  */
238  void dualColumn1(bool doAll = false);
239  /**
240     Array has tableau row (row section)
241     Just does slack part
242     Returns guess at upper theta (infinite if no pivot) and may set sequenceIn_ if free
243  */
244  double dualColumn1A();
245  /// Do all given tableau row
246  double dualColumn1B();
247  /**
248     Chooses incoming
249     Puts flipped ones in list
250     If necessary will modify costs
251  */
252  void dualColumn2();
253  void dualColumn2Most(dualColumnResult &result);
254  void dualColumn2First(dualColumnResult &result);
255  /**
256     Chooses part of incoming
257     Puts flipped ones in list
258     If necessary will modify costs
259  */
260  void dualColumn2(dualColumnResult &result);
261  /**
262     This sees what is best thing to do in branch and bound cleanup
263     If sequenceIn_ < 0 then can't do anything
264  */
265  void checkPossibleCleanup(CoinIndexedVector *array);
266  /**
267     Chooses dual pivot row
268     Would be faster with separate region to scan
269     and will have this (with square of infeasibility) when steepest
270     For easy problems we can just choose one of the first rows we look at
271  */
272  void dualPivotRow();
273  /** Checks if any fake bounds active - if so returns number and modifies
274      updatedDualBound_ and everything.
275      Free variables will be left as free
276      Returns number of bounds changed if >=0
277      Returns -1 if not initialize and no effect
278      fills cost of change vector
279  */
280  int changeBounds(int initialize, double &changeCost);
281  /** As changeBounds but just changes new bounds for a single variable.
282      Returns true if change */
283  bool changeBound(int iSequence);
284  /// Restores bound to original bound
285  void originalBound(int iSequence);
286  /** Checks if tentative optimal actually means unbounded in dual
287      Returns -3 if not, 2 if is unbounded */
288  int checkUnbounded(CoinIndexedVector &ray, double changeCost);
289  /**  Refactorizes if necessary
290       Checks if finished.  Updates status.
291       lastCleaned refers to iteration at which some objective/feasibility
292       cleaning too place.
293       
294       type - 0 initial so set up save arrays etc
295       - 1 normal -if good update save
296       - 2 restoring from saved
297  */
298  void statusOfProblemInDual(int type);
299  /** Fast iterations.  Misses out a lot of initialization.
300      Normally stops on maximum iterations, first re-factorization
301      or tentative optimum.  If looks interesting then continues as
302      normal.  Returns 0 if finished properly, 1 otherwise.
303  */
304  /// Gets tableau column - does flips and checks what to do next
305  /// Knows tableau column in 1, flips in 2 and gets an array for flips (as serial here)
306  int whatNext();
307  /// see if cutoff reached
308  bool checkCutoff(bool computeObjective);
309  /// Does something about fake tolerances
310  int bounceTolerances(int type);
311  /// Perturbs problem
312  void perturb(double factor);
313  /// Perturbs problem B
314  void perturbB(double factor, int type);
315  /// Make non free variables dual feasible by moving to a bound
316  int makeNonFreeVariablesDualFeasible(bool changeCosts = false);
317  int fastDual(bool alwaysFinish = false);
318  /** Checks number of variables at fake bounds.  This is used by fastDual
319      so can exit gracefully before end */
320  int numberAtFakeBound();
321
322  /** Pivot in a variable and choose an outgoing one.  Assumes dual
323      feasible - will not go through a reduced cost.  Returns step length in theta
324      Return codes as before but -1 means no acceptable pivot
325  */
326  int pivotResultPart1();
327  /** Get next free , -1 if none */
328  int nextSuperBasic();
329  /// Startup part of dual
330  void startupSolve();
331  /// Ending part of dual
332  void finishSolve();
333  void gutsOfDual();
334  //int dual2(int ifValuesPass,int startFinishOptions=0);
335  int resetFakeBounds(int type);
336
337  //@}
338};
339#if ABC_PARALLEL == 1
340void *abc_parallelManager(void *simplex);
341#endif
342#endif
343
344/* vi: softtabstop=2 shiftwidth=2 expandtab tabstop=2
345*/
Note: See TracBrowser for help on using the repository browser.