Scippy

SoPlex

Sequential object-oriented simPlex

spxweightpr.cpp
Go to the documentation of this file.
1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /* */
3 /* This file is part of the class library */
4 /* SoPlex --- the Sequential object-oriented simPlex. */
5 /* */
6 /* Copyright (C) 1996-2016 Konrad-Zuse-Zentrum */
7 /* fuer Informationstechnik Berlin */
8 /* */
9 /* SoPlex is distributed under the terms of the ZIB Academic Licence. */
10 /* */
11 /* You should have received a copy of the ZIB Academic License */
12 /* along with SoPlex; see the file COPYING. If not email to soplex@zib.de. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15 
16 #include <assert.h>
17 
18 #include "spxdefines.h"
19 #include "spxweightpr.h"
20 #include "exceptions.h"
21 
22 namespace soplex
23 {
24 
26 {
27  if (rep == SPxSolver::ROW)
28  {
31  }
32  else
33  {
36  }
37 }
38 
40 {
41  if (thesolver && tp == SPxSolver::LEAVE)
42  {
45  }
46 }
47 
48 void SPxWeightPR::computeLeavePenalty(int start, int end)
49 {
50  const SPxBasis& basis = solver()->basis();
51 
52  for (int i = start; i < end; ++i)
53  {
54  SPxId id = basis.baseId(i);
55  if (id.type() == SPxId::ROW_ID)
57  else
59  }
60 }
61 
62 void SPxWeightPR::computeRP(int start, int end)
63 {
64  for (int i = start; i < end; ++i)
65  {
66  /**@todo TK04NOV98 here is a bug.
67  * solver()->rowVector(i).length() could be zero, so
68  * solver()->rowVector(i).length2() is also zero and we
69  * get an arithmetic exception.
70  */
71  assert(solver()->rowVector(i).length() > 0);
72 
73  rPenalty[i] = (solver()->rowVector(i) * solver()->maxObj()) * objlength
74  / solver()->rowVector(i).length2();
75  ASSERT_WARN( "WWGTPR01", rPenalty[i] > -1 - solver()->epsilon() );
76  }
77 }
78 
79 void SPxWeightPR::computeCP(int start, int end)
80 {
81  for (int i = start; i < end; ++i)
82  {
83  cPenalty[i] = solver()->maxObj(i) * objlength;
84  ASSERT_WARN( "WWGTPR02", cPenalty[i] > -1 - solver()->epsilon() );
85  }
86 }
87 
89 {
90  thesolver = base;
91 
92  rPenalty.reDim(base->nRows());
93  cPenalty.reDim(base->nCols());
94 
95  objlength = 1 / solver()->maxObj().length();
96  computeCP(0, base->nCols());
97  computeRP(0, base->nRows());
98 }
99 
101 {
102  const Real* test = thesolver->fTest().get_const_ptr();
103  Real type = 1 - 2 * (thesolver->rep() == SPxSolver::COLUMN ? 1 : 0);
104  Real best = type * infinity;
105  int lastIdx = -1;
106  Real x;
107  int i;
108 
109  for (i = solver()->dim() - 1; i >= 0; --i)
110  {
111  x = test[i];
112  if (x < -theeps)
113  {
114  x *= leavePenalty[i];
115  if (type * (x-best) < 0.0)
116  {
117  best = x;
118  lastIdx = i;
119  }
120  }
121  }
122  assert(isConsistent());
123  return lastIdx;
124 }
125 
126 #if 0
127 /**@todo remove this code */
128 // ??? This is the old (buggy) version
130 {
131  const Real* test = thesolver->fTest().get_const_ptr();
132  Real type = 1 - 2 * (thesolver->rep() == SPxSolver::COLUMN);
133  Real best = type * infinity;
134  int lastIdx = -1;
135  Real x;
136  int i;
137 
138  for (i = solver()->dim() - 1; i >= 0; --i)
139  {
140  x = test[i];
141  if (x < -theeps)
142  {
143  x *= leavePenalty[i];
144  if (x < best)
145  {
146  best = x;
147  lastIdx = i;
148  }
149  }
150  }
151  assert(isConsistent());
152  return lastIdx;
153 }
154 #endif
155 
157 {
158  const Vector& rTest = (solver()->rep() == SPxSolver::ROW)
159  ? solver()->test() : solver()->coTest();
160  const Vector& cTest = (solver()->rep() == SPxSolver::ROW)
161  ? solver()->coTest() : solver()->test();
162  const SPxBasis::Desc& ds = solver()->basis().desc();
163  Real best = infinity;
164  SPxId lastId;
165  Real x;
166  int i;
167 
168  for (i = solver()->nRows() - 1; i >= 0; --i)
169  {
170  x = rTest[i];
171  if (x < -theeps)
172  {
173  x *= -x;
174  switch (ds.rowStatus(i))
175  {
178  x *= 1 + rPenalty[i];
179  break;
182  x *= 1 - rPenalty[i];
183  break;
186  return SPxId(solver()->rId(i));
188  if (solver()->pVec()[i] > solver()->upBound()[i])
189  x *= 1 + rPenalty[i];
190  else
191  x *= 1 - rPenalty[i];
192  break;
195  default:
196  throw SPxInternalCodeException("XWGTPR01 This should never happen.");
197  }
198  if (x < best)
199  {
200  best = x;
201  lastId = solver()->rId(i);
202  }
203  }
204  }
205 
206  for (i = solver()->nCols() - 1; i >= 0; --i)
207  {
208  x = cTest[i];
209  if (x < -theeps)
210  {
211  x *= -x;
212  switch (ds.colStatus(i))
213  {
216  x *= 1 + cPenalty[i];
217  break;
220  x *= 1 - cPenalty[i];
221  break;
224  return SPxId(solver()->cId(i));
226  if (solver()->coPvec()[i] > solver()->ucBound()[i])
227  x *= 1 + cPenalty[i];
228  else
229  x *= 1 - cPenalty[i];
230  break;
233  default:
234  throw SPxInternalCodeException("XWGTPR02 This should never happen.");
235  }
236  if (x < best)
237  {
238  best = x;
239  lastId = solver()->cId(i);
240  }
241  }
242  }
243  assert(isConsistent());
244  return lastId;
245 }
246 
248 {
249  if (solver()->rep() == SPxSolver::ROW)
250  {
251  int start = rPenalty.dim();
252  rPenalty.reDim(solver()->nRows());
253  computeRP(start, solver()->nRows());
254  }
255  else
256  {
257  int start = cPenalty.dim();
258  cPenalty.reDim(solver()->nCols());
259  computeCP(start, solver()->nCols());
260  }
261  if (solver()->type() == SPxSolver::LEAVE)
262  {
263  int start = leavePenalty.dim();
264  leavePenalty.reDim( solver()->dim() );
265  computeLeavePenalty( start, solver()->dim() );
266  }
267 }
268 
270 {
271  if (solver()->rep() == SPxSolver::COLUMN)
272  {
273  int start = rPenalty.dim();
274  rPenalty.reDim(solver()->nRows());
275  computeRP(start, solver()->nRows());
276  }
277  else
278  {
279  int start = cPenalty.dim();
280  cPenalty.reDim(solver()->nCols());
281  computeCP(start, solver()->nCols());
282  }
283  if (solver()->type() == SPxSolver::LEAVE)
284  {
285  int start = leavePenalty.dim();
286  leavePenalty.reDim( solver()->dim() );
287  computeLeavePenalty( start, solver()->dim() );
288  }
289 }
290 
292 {
293  assert(solver() != 0);
294 
295  if (solver()->rep() == SPxSolver::ROW)
296  {
297  rPenalty[i] = rPenalty[rPenalty.dim()];
298  rPenalty.reDim(solver()->nRows());
299  }
300  else
301  {
302  cPenalty[i] = cPenalty[cPenalty.dim()];
303  cPenalty.reDim(solver()->nCols());
304  }
305 }
306 
307 void SPxWeightPR::removedVecs(const int perm[])
308 {
309  assert(solver() != 0);
310 
311  if (solver()->rep() == SPxSolver::ROW)
312  {
313  int j = rPenalty.dim();
314  for (int i = 0; i < j; ++i)
315  {
316  if (perm[i] >= 0)
317  rPenalty[perm[i]] = rPenalty[i];
318  }
319  rPenalty.reDim(solver()->nRows());
320  }
321  else
322  {
323  int j = cPenalty.dim();
324  for (int i = 0; i < j; ++i)
325  {
326  if (perm[i] >= 0)
327  cPenalty[perm[i]] = cPenalty[i];
328  }
329  cPenalty.reDim(solver()->nCols());
330  }
331 }
332 
334 {
335  assert(solver() != 0);
336 
337  if (solver()->rep() == SPxSolver::COLUMN)
338  {
339  rPenalty[i] = rPenalty[rPenalty.dim()];
340  rPenalty.reDim(solver()->nRows());
341  }
342  else
343  {
344  cPenalty[i] = cPenalty[cPenalty.dim()];
345  cPenalty.reDim(solver()->nCols());
346  }
347 }
348 
349 void SPxWeightPR::removedCoVecs(const int perm[])
350 {
351  assert(solver() != 0);
352 
353  if (solver()->rep() == SPxSolver::COLUMN)
354  {
355  int j = rPenalty.dim();
356  for (int i = 0; i < j; ++i)
357  {
358  if (perm[i] >= 0)
359  rPenalty[perm[i]] = rPenalty[i];
360  }
361  rPenalty.reDim(solver()->nRows());
362  }
363  else
364  {
365  int j = cPenalty.dim();
366  for (int i = 0; i < j; ++i)
367  {
368  if (perm[i] >= 0)
369  cPenalty[perm[i]] = cPenalty[i];
370  }
371  cPenalty.reDim(solver()->nCols());
372  }
373 }
374 
376 {
377 #ifdef ENABLE_CONSISTENCY_CHECKS
378  if (solver() != 0)
379  {
380  if (rPenalty.dim() != solver()->nRows())
381  return MSGinconsistent("SPxWeightPR");
382  if (cPenalty.dim() != solver()->nCols())
383  return MSGinconsistent("SPxWeightPR");
384  }
385 #endif
386 
387  return true;
388 }
389 } // namespace soplex
int dim() const
dimension of basis matrix.
Definition: spxsolver.h:936
Exception class for things that should NEVER happen.This class is derived from the SoPlex exception b...
Definition: exceptions.h:109
Representation rep() const
return the current basis representation.
Definition: spxsolver.h:406
void reDim(int newdim, const bool setZero=true)
Resets DVectorBase&#39;s dimension to newdim.
Definition: dvectorbase.h:249
primal variable is fixed to both bounds
Definition: spxbasis.h:190
const VectorBase< R > & maxObj() const
Returns objective vector for maximization problem.
Definition: spxlpbase.h:384
primal or dual variable has no status
Definition: spxbasis.h:195
Type
Algorithmic type.
Definition: spxsolver.h:124
const Real * penalty
Definition: spxweightpr.h:55
const Vector & fTest() const
Violations of fVec.
Definition: spxsolver.h:1246
void setType(SPxSolver::Type tp)
set entering/leaving algorithm
Definition: spxweightpr.cpp:39
Exception classes for SoPlex.
Status & rowStatus(int i)
Definition: spxbasis.h:239
Representation
LP basis representation.
Definition: spxsolver.h:105
int number(const SPxRowId &id) const
Returns the row number of the row with identifier id.
Definition: spxlpbase.h:450
const R * get_const_ptr() const
Conversion to C-style pointer.
Definition: vectorbase.h:412
#define ASSERT_WARN(prefix, expr)
Macro to turn some assertions into warnings.
Definition: spxdefines.h:72
void computeCP(int start, int end)
compute weights for columns.
Definition: spxweightpr.cpp:79
rowwise representation.
Definition: spxsolver.h:107
dual variable is left free, but unset
Definition: spxbasis.h:191
const Real * coPenalty
Definition: spxweightpr.h:57
primal variable is set to its upper bound
Definition: spxbasis.h:188
Generic Ids for LP rows or columns.Both SPxColIds and SPxRowIds may be treated uniformly as SPxIds: ...
Definition: spxid.h:85
Leaving Simplex.
Definition: spxsolver.h:143
double Real
SOPLEX_DEBUG.
Definition: spxdefines.h:200
const Vector & test() const
Violations of pVec.
Definition: spxsolver.h:1392
virtual void removedVecs(const int perm[])
n vectors have been removed from the loaded LP.
DVector cPenalty
column penalties
Definition: spxweightpr.h:49
SPxRowId rId(int n) const
Returns the row identifier for row n.
Definition: spxlpbase.h:470
int nCols() const
Returns number of columns in LP.
Definition: spxlpbase.h:133
const SPxBasis & basis() const
Return current basis.
Definition: spxsolver.h:1616
dual variable is set to its upper bound
Definition: spxbasis.h:192
virtual void addedVecs(int n)
n vectors have been added to the loaded LP.
primal variable is left free, but unset
Definition: spxbasis.h:189
SPxSolver * thesolver
the solver
Definition: spxpricer.h:56
virtual void removedCoVecs(const int perm[])
n covectors have been removed from the loaded LP.
virtual SPxSolver * solver() const
returns loaded SPxSolver object.
Definition: spxpricer.h:129
Basis descriptor.
Definition: spxbasis.h:104
virtual void removedCoVec(int i)
the i&#39;th covector has been removed from the loaded LP.
int dim() const
Dimension of vector.
Definition: vectorbase.h:174
virtual SPxId selectEnter()
virtual void removedVec(int i)
the i&#39;th vector has been removed from the loaded LP.
row identifier.
Definition: spxid.h:97
DVector leavePenalty
penalties for leaving alg
Definition: spxweightpr.h:53
Status & colStatus(int i)
Definition: spxbasis.h:254
SPxId & baseId(int i)
Definition: spxbasis.h:497
const Vector & coTest() const
violations of coPvec.
Definition: spxsolver.h:1326
Real length() const
Floating point approximation of euclidian norm (without any approximation guarantee).
Definition: vectorbase.h:356
Debugging, floating point type and parameter definitions.
Simplex basis.Consider the linear program as provided from class SPxLP: where , and ...
Definition: spxbasis.h:82
int nRows() const
Returns number of rows in LP.
Definition: spxlpbase.h:127
Sequential object-oriented SimPlex.SPxSolver is an LP solver class using the revised Simplex algorith...
Definition: spxsolver.h:84
virtual bool isConsistent() const
checks for consistency
Everything should be within this namespace.
virtual void addedCoVecs(int n)
n covectors have been added to the loaded LP.
virtual Real epsilon() const
returns violation bound theeps.
Definition: spxpricer.h:135
R length2() const
Squared norm.
Definition: svectorbase.h:477
void computeLeavePenalty(int start, int end)
compute leave penalties.
Definition: spxweightpr.cpp:48
SPxColId cId(int n) const
Returns the column identifier for column n.
Definition: spxlpbase.h:476
primal variable is set to its lower bound
Definition: spxbasis.h:187
Real theeps
violation bound
Definition: spxpricer.h:58
Real objlength
length of objective vector.
Definition: spxweightpr.h:59
dual variable is set to its lower bound
Definition: spxbasis.h:193
void computeRP(int start, int end)
compute weights for rows.
Definition: spxweightpr.cpp:62
const Real infinity
Definition: spxdefines.cpp:26
dual variable has two bounds
Definition: spxbasis.h:194
virtual int selectLeave()
const SVectorBase< R > & rowVector(int i) const
Gets row vector of row i.
Definition: spxlpbase.h:212
const Desc & desc() const
Definition: spxbasis.h:457
#define MSGinconsistent(name)
Definition: spxdefines.h:121
void setRep(SPxSolver::Representation rep)
set row/column representation
Definition: spxweightpr.cpp:25
virtual void load(SPxSolver *base)
sets the solver
Definition: spxweightpr.cpp:88
DVector rPenalty
row penalties
Definition: spxweightpr.h:51
Weighted pricing.
columnwise representation.
Definition: spxsolver.h:108