source: trunk/Clp/src/ClpPrimalColumnDantzig.cpp @ 1732

Last change on this file since 1732 was 1732, checked in by forrest, 8 years ago

various fixes plus slightly weighted pricing plus lagomory switches

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 8.2 KB
Line 
1/* $Id: ClpPrimalColumnDantzig.cpp 1732 2011-05-31 08:09:41Z forrest $ */
2// Copyright (C) 2002, International Business Machines
3// Corporation and others.  All Rights Reserved.
4// This code is licensed under the terms of the Eclipse Public License (EPL).
5
6#include "CoinPragma.hpp"
7
8#include <cstdio>
9
10#include "CoinIndexedVector.hpp"
11
12#include "ClpSimplex.hpp"
13#include "ClpPrimalColumnDantzig.hpp"
14#include "ClpFactorization.hpp"
15#include "ClpPackedMatrix.hpp"
16
17//#############################################################################
18// Constructors / Destructor / Assignment
19//#############################################################################
20
21//-------------------------------------------------------------------
22// Default Constructor
23//-------------------------------------------------------------------
24ClpPrimalColumnDantzig::ClpPrimalColumnDantzig ()
25     : ClpPrimalColumnPivot()
26{
27     type_ = 1;
28}
29
30//-------------------------------------------------------------------
31// Copy constructor
32//-------------------------------------------------------------------
33ClpPrimalColumnDantzig::ClpPrimalColumnDantzig (const ClpPrimalColumnDantzig & source)
34     : ClpPrimalColumnPivot(source)
35{
36
37}
38
39//-------------------------------------------------------------------
40// Destructor
41//-------------------------------------------------------------------
42ClpPrimalColumnDantzig::~ClpPrimalColumnDantzig ()
43{
44
45}
46
47//----------------------------------------------------------------
48// Assignment operator
49//-------------------------------------------------------------------
50ClpPrimalColumnDantzig &
51ClpPrimalColumnDantzig::operator=(const ClpPrimalColumnDantzig& rhs)
52{
53     if (this != &rhs) {
54          ClpPrimalColumnPivot::operator=(rhs);
55     }
56     return *this;
57}
58
59// Returns pivot column, -1 if none
60int
61ClpPrimalColumnDantzig::pivotColumn(CoinIndexedVector * updates,
62                                    CoinIndexedVector * /*spareRow1*/,
63                                    CoinIndexedVector * spareRow2,
64                                    CoinIndexedVector * spareColumn1,
65                                    CoinIndexedVector * spareColumn2)
66{
67     assert(model_);
68     int iSection, j;
69     int number;
70     int * index;
71     double * updateBy;
72     double * reducedCost;
73
74     bool anyUpdates;
75
76     if (updates->getNumElements()) {
77          anyUpdates = true;
78     } else {
79          // sub flip - nothing to do
80          anyUpdates = false;
81     }
82     if (anyUpdates) {
83          model_->factorization()->updateColumnTranspose(spareRow2, updates);
84          // put row of tableau in rowArray and columnArray
85          model_->clpMatrix()->transposeTimes(model_, -1.0,
86                                              updates, spareColumn2, spareColumn1);
87          for (iSection = 0; iSection < 2; iSection++) {
88
89               reducedCost = model_->djRegion(iSection);
90
91               if (!iSection) {
92                    number = updates->getNumElements();
93                    index = updates->getIndices();
94                    updateBy = updates->denseVector();
95               } else {
96                    number = spareColumn1->getNumElements();
97                    index = spareColumn1->getIndices();
98                    updateBy = spareColumn1->denseVector();
99               }
100
101               for (j = 0; j < number; j++) {
102                    int iSequence = index[j];
103                    double value = reducedCost[iSequence];
104                    value -= updateBy[j];
105                    updateBy[j] = 0.0;
106                    reducedCost[iSequence] = value;
107               }
108
109          }
110          updates->setNumElements(0);
111          spareColumn1->setNumElements(0);
112     }
113
114
115     // update of duals finished - now do pricing
116
117     double largest = model_->currentPrimalTolerance();
118     // we can't really trust infeasibilities if there is primal error
119     if (model_->largestDualError() > 1.0e-8)
120          largest *= model_->largestDualError() / 1.0e-8;
121
122
123
124     double bestDj = model_->dualTolerance();
125     int bestSequence = -1;
126
127     double bestFreeDj = model_->dualTolerance();
128     int bestFreeSequence = -1;
129
130     number = model_->numberRows() + model_->numberColumns();
131     int iSequence;
132     reducedCost = model_->djRegion();
133
134#ifndef CLP_PRIMAL_SLACK_MULTIPLIER
135     for (iSequence = 0; iSequence < number; iSequence++) {
136          // check flagged variable
137          if (!model_->flagged(iSequence)) {
138               double value = reducedCost[iSequence];
139               ClpSimplex::Status status = model_->getStatus(iSequence);
140
141               switch(status) {
142
143               case ClpSimplex::basic:
144               case ClpSimplex::isFixed:
145                    break;
146               case ClpSimplex::isFree:
147               case ClpSimplex::superBasic:
148                    if (fabs(value) > bestFreeDj) {
149                         bestFreeDj = fabs(value);
150                         bestFreeSequence = iSequence;
151                    }
152                    break;
153               case ClpSimplex::atUpperBound:
154                    if (value > bestDj) {
155                         bestDj = value;
156                         bestSequence = iSequence;
157                    }
158                    break;
159               case ClpSimplex::atLowerBound:
160                    if (value < -bestDj) {
161                         bestDj = -value;
162                         bestSequence = iSequence;
163                    }
164               }
165          }
166     }
167#else
168     // Columns
169     int numberColumns = model_->numberColumns();
170     for (iSequence = 0; iSequence < numberColumns; iSequence++) {
171          // check flagged variable
172          if (!model_->flagged(iSequence)) {
173               double value = reducedCost[iSequence];
174               ClpSimplex::Status status = model_->getStatus(iSequence);
175
176               switch(status) {
177
178               case ClpSimplex::basic:
179               case ClpSimplex::isFixed:
180                    break;
181               case ClpSimplex::isFree:
182               case ClpSimplex::superBasic:
183                    if (fabs(value) > bestFreeDj) {
184                         bestFreeDj = fabs(value);
185                         bestFreeSequence = iSequence;
186                    }
187                    break;
188               case ClpSimplex::atUpperBound:
189                    if (value > bestDj) {
190                         bestDj = value;
191                         bestSequence = iSequence;
192                    }
193                    break;
194               case ClpSimplex::atLowerBound:
195                    if (value < -bestDj) {
196                         bestDj = -value;
197                         bestSequence = iSequence;
198                    }
199               }
200          }
201     }
202     // Rows
203     for ( ; iSequence < number; iSequence++) {
204          // check flagged variable
205          if (!model_->flagged(iSequence)) {
206               double value = reducedCost[iSequence] * CLP_PRIMAL_SLACK_MULTIPLIER;
207               ClpSimplex::Status status = model_->getStatus(iSequence);
208
209               switch(status) {
210
211               case ClpSimplex::basic:
212               case ClpSimplex::isFixed:
213                    break;
214               case ClpSimplex::isFree:
215               case ClpSimplex::superBasic:
216                    if (fabs(value) > bestFreeDj) {
217                         bestFreeDj = fabs(value);
218                         bestFreeSequence = iSequence;
219                    }
220                    break;
221               case ClpSimplex::atUpperBound:
222                    if (value > bestDj) {
223                         bestDj = value;
224                         bestSequence = iSequence;
225                    }
226                    break;
227               case ClpSimplex::atLowerBound:
228                    if (value < -bestDj) {
229                         bestDj = -value;
230                         bestSequence = iSequence;
231                    }
232               }
233          }
234     }
235#endif
236     // bias towards free
237     if (bestFreeSequence >= 0 && bestFreeDj > 0.1 * bestDj)
238          bestSequence = bestFreeSequence;
239     return bestSequence;
240}
241
242//-------------------------------------------------------------------
243// Clone
244//-------------------------------------------------------------------
245ClpPrimalColumnPivot * ClpPrimalColumnDantzig::clone(bool CopyData) const
246{
247     if (CopyData) {
248          return new ClpPrimalColumnDantzig(*this);
249     } else {
250          return new ClpPrimalColumnDantzig();
251     }
252}
253
Note: See TracBrowser for help on using the repository browser.