source: trunk/Clp/src/ClpSimplexDual.hpp @ 1502

Last change on this file since 1502 was 1502, checked in by forrest, 10 years ago

moving sandbox stuff to trunk

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 12.1 KB
Line 
1/* $Id: ClpSimplexDual.hpp 1502 2010-01-29 14:25:07Z forrest $ */
2// Copyright (C) 2002, International Business Machines
3// Corporation and others.  All Rights Reserved.
4
5/*
6   Authors
7
8   John Forrest
9
10 */
11#ifndef ClpSimplexDual_H
12#define ClpSimplexDual_H
13
14#include "ClpSimplex.hpp"
15
16/** This solves LPs using the dual simplex method
17
18    It inherits from ClpSimplex.  It has no data of its own and
19    is never created - only cast from a ClpSimplex object at algorithm time.
20
21*/
22
23class ClpSimplexDual : public ClpSimplex {
24
25public:
26
27    /**@name Description of algorithm */
28    //@{
29    /** Dual algorithm
30
31        Method
32
33       It tries to be a single phase approach with a weight of 1.0 being
34       given to getting optimal and a weight of updatedDualBound_ being
35       given to getting dual feasible.  In this version I have used the
36       idea that this weight can be thought of as a fake bound.  If the
37       distance between the lower and upper bounds on a variable is less
38       than the feasibility weight then we are always better off flipping
39       to other bound to make dual feasible.  If the distance is greater
40       then we make up a fake bound updatedDualBound_ away from one bound.
41       If we end up optimal or primal infeasible, we check to see if
42       bounds okay.  If so we have finished, if not we increase updatedDualBound_
43       and continue (after checking if unbounded). I am undecided about
44       free variables - there is coding but I am not sure about it.  At
45       present I put them in basis anyway.
46
47       The code is designed to take advantage of sparsity so arrays are
48       seldom zeroed out from scratch or gone over in their entirety.
49       The only exception is a full scan to find outgoing variable for
50       Dantzig row choice.  For steepest edge we keep an updated list
51       of infeasibilities (actually squares).
52       On easy problems we don't need full scan - just
53       pick first reasonable.
54
55       One problem is how to tackle degeneracy and accuracy.  At present
56       I am using the modification of costs which I put in OSL and some
57       of what I think is the dual analog of Gill et al.
58       I am still not sure of the exact details.
59
60       The flow of dual is three while loops as follows:
61
62       while (not finished) {
63
64         while (not clean solution) {
65
66            Factorize and/or clean up solution by flipping variables so
67        dual feasible.  If looks finished check fake dual bounds.
68        Repeat until status is iterating (-1) or finished (0,1,2)
69
70         }
71
72         while (status==-1) {
73
74           Iterate until no pivot in or out or time to re-factorize.
75
76           Flow is:
77
78           choose pivot row (outgoing variable).  if none then
79       we are primal feasible so looks as if done but we need to
80       break and check bounds etc.
81
82       Get pivot row in tableau
83
84           Choose incoming column.  If we don't find one then we look
85       primal infeasible so break and check bounds etc.  (Also the
86       pivot tolerance is larger after any iterations so that may be
87       reason)
88
89           If we do find incoming column, we may have to adjust costs to
90       keep going forwards (anti-degeneracy).  Check pivot will be stable
91       and if unstable throw away iteration and break to re-factorize.
92       If minor error re-factorize after iteration.
93
94       Update everything (this may involve flipping variables to stay
95       dual feasible.
96
97         }
98
99       }
100
101       TODO's (or maybe not)
102
103       At present we never check we are going forwards.  I overdid that in
104       OSL so will try and make a last resort.
105
106       Needs partial scan pivot out option.
107
108       May need other anti-degeneracy measures, especially if we try and use
109       loose tolerances as a way to solve in fewer iterations.
110
111       I like idea of dynamic scaling.  This gives opportunity to decouple
112       different implications of scaling for accuracy, iteration count and
113       feasibility tolerance.
114
115       for use of exotic parameter startFinishoptions see Clpsimplex.hpp
116    */
117
118    int dual(int ifValuesPass, int startFinishOptions = 0);
119    /** For strong branching.  On input lower and upper are new bounds
120        while on output they are change in objective function values
121        (>1.0e50 infeasible).
122        Return code is 0 if nothing interesting, -1 if infeasible both
123        ways and +1 if infeasible one way (check values to see which one(s))
124        Solutions are filled in as well - even down, odd up - also
125        status and number of iterations
126    */
127    int strongBranching(int numberVariables, const int * variables,
128                        double * newLower, double * newUpper,
129                        double ** outputSolution,
130                        int * outputStatus, int * outputIterations,
131                        bool stopOnFirstInfeasible = true,
132                        bool alwaysFinish = false,
133                        int startFinishOptions = 0);
134    /// This does first part of StrongBranching
135    ClpFactorization * setupForStrongBranching(char * arrays, int numberRows,
136            int numberColumns, bool solveLp = false);
137    /// This cleans up after strong branching
138    void cleanupAfterStrongBranching(ClpFactorization * factorization);
139    //@}
140
141    /**@name Functions used in dual */
142    //@{
143    /** This has the flow between re-factorizations
144        Broken out for clarity and will be used by strong branching
145
146        Reasons to come out:
147        -1 iterations etc
148        -2 inaccuracy
149        -3 slight inaccuracy (and done iterations)
150        +0 looks optimal (might be unbounded - but we will investigate)
151        +1 looks infeasible
152        +3 max iterations
153
154        If givenPi not NULL then in values pass
155     */
156    int whileIterating(double * & givenPi, int ifValuesPass);
157    /** The duals are updated by the given arrays.
158        Returns number of infeasibilities.
159        After rowArray and columnArray will just have those which
160        have been flipped.
161        Variables may be flipped between bounds to stay dual feasible.
162        The output vector has movement of primal
163        solution (row length array) */
164    int updateDualsInDual(CoinIndexedVector * rowArray,
165                          CoinIndexedVector * columnArray,
166                          CoinIndexedVector * outputArray,
167                          double theta,
168                          double & objectiveChange,
169                          bool fullRecompute);
170    /** The duals are updated by the given arrays.
171        This is in values pass - so no changes to primal is made
172    */
173    void updateDualsInValuesPass(CoinIndexedVector * rowArray,
174                                 CoinIndexedVector * columnArray,
175                                 double theta);
176    /** While updateDualsInDual sees what effect is of flip
177        this does actual flipping.
178    */
179    void flipBounds(CoinIndexedVector * rowArray,
180                    CoinIndexedVector * columnArray);
181    /**
182        Row array has row part of pivot row
183        Column array has column part.
184        This chooses pivot column.
185        Spare arrays are used to save pivots which will go infeasible
186        We will check for basic so spare array will never overflow.
187        If necessary will modify costs
188        For speed, we may need to go to a bucket approach when many
189        variables are being flipped.
190        Returns best possible pivot value
191    */
192    double dualColumn(CoinIndexedVector * rowArray,
193                      CoinIndexedVector * columnArray,
194                      CoinIndexedVector * spareArray,
195                      CoinIndexedVector * spareArray2,
196                      double accpetablePivot,
197                      CoinBigIndex * dubiousWeights);
198    /// Does first bit of dualColumn
199    int dualColumn0(const CoinIndexedVector * rowArray,
200                    const CoinIndexedVector * columnArray,
201                    CoinIndexedVector * spareArray,
202                    double acceptablePivot,
203                    double & upperReturn, double &bestReturn, double & badFree);
204    /**
205        Row array has row part of pivot row
206        Column array has column part.
207        This sees what is best thing to do in dual values pass
208        if sequenceIn==sequenceOut can change dual on chosen row and leave variable in basis
209    */
210    void checkPossibleValuesMove(CoinIndexedVector * rowArray,
211                                 CoinIndexedVector * columnArray,
212                                 double acceptablePivot);
213    /**
214        Row array has row part of pivot row
215        Column array has column part.
216        This sees what is best thing to do in branch and bound cleanup
217        If sequenceIn_ < 0 then can't do anything
218    */
219    void checkPossibleCleanup(CoinIndexedVector * rowArray,
220                              CoinIndexedVector * columnArray,
221                              double acceptablePivot);
222    /**
223        This sees if we can move duals in dual values pass.
224        This is done before any pivoting
225    */
226    void doEasyOnesInValuesPass(double * givenReducedCosts);
227    /**
228        Chooses dual pivot row
229        Would be faster with separate region to scan
230        and will have this (with square of infeasibility) when steepest
231        For easy problems we can just choose one of the first rows we look at
232
233        If alreadyChosen >=0 then in values pass and that row has been
234        selected
235    */
236    void dualRow(int alreadyChosen);
237    /** Checks if any fake bounds active - if so returns number and modifies
238        updatedDualBound_ and everything.
239        Free variables will be left as free
240        Returns number of bounds changed if >=0
241        Returns -1 if not initialize and no effect
242        Fills in changeVector which can be used to see if unbounded
243        and cost of change vector
244        If 2 sets to original (just changed)
245    */
246    int changeBounds(int initialize, CoinIndexedVector * outputArray,
247                     double & changeCost);
248    /** As changeBounds but just changes new bounds for a single variable.
249        Returns true if change */
250    bool changeBound( int iSequence);
251    /// Restores bound to original bound
252    void originalBound(int iSequence);
253    /** Checks if tentative optimal actually means unbounded in dual
254        Returns -3 if not, 2 if is unbounded */
255    int checkUnbounded(CoinIndexedVector * ray, CoinIndexedVector * spare,
256                       double changeCost);
257    /**  Refactorizes if necessary
258         Checks if finished.  Updates status.
259         lastCleaned refers to iteration at which some objective/feasibility
260         cleaning too place.
261
262         type - 0 initial so set up save arrays etc
263              - 1 normal -if good update save
264          - 2 restoring from saved
265    */
266    void statusOfProblemInDual(int & lastCleaned, int type,
267                               double * givenDjs, ClpDataSave & saveData,
268                               int ifValuesPass);
269    /** Perturbs problem (method depends on perturbation())
270        returns nonzero if should go to dual */
271    int perturb();
272    /** Fast iterations.  Misses out a lot of initialization.
273        Normally stops on maximum iterations, first re-factorization
274        or tentative optimum.  If looks interesting then continues as
275        normal.  Returns 0 if finished properly, 1 otherwise.
276    */
277    int fastDual(bool alwaysFinish = false);
278    /** Checks number of variables at fake bounds.  This is used by fastDual
279        so can exit gracefully before end */
280    int numberAtFakeBound();
281
282    /** Pivot in a variable and choose an outgoing one.  Assumes dual
283        feasible - will not go through a reduced cost.  Returns step length in theta
284        Returns ray in ray_ (or NULL if no pivot)
285        Return codes as before but -1 means no acceptable pivot
286    */
287    int pivotResult();
288    /** Get next free , -1 if none */
289    int nextSuperBasic();
290    /** Startup part of dual (may be extended to other algorithms)
291        returns 0 if good, 1 if bad */
292    int startupSolve(int ifValuesPass, double * saveDuals, int startFinishOptions);
293    void finishSolve(int startFinishOptions);
294    void gutsOfDual(int ifValuesPass, double * & saveDuals, int initialStatus,
295                    ClpDataSave & saveData);
296    //int dual2(int ifValuesPass,int startFinishOptions=0);
297    void resetFakeBounds(int type);
298
299    //@}
300};
301#endif
Note: See TracBrowser for help on using the repository browser.