1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /* */
3 /* This file is part of the class library */
4 /* SoPlex --- the Sequential object-oriented simPlex. */
5 /* */
6 /* Copyright (C) 1996-2022 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 /**@file spxequilisc.hpp
17 * @brief Equilibrium row/column scaling.
18 */
19 #include <assert.h>
20
21 #include "soplex/spxequilisc.h"
22 #include "soplex/spxout.h"
23 #include "soplex/spxlpbase.h"
24 #include "soplex/spxlp.h"
|
(1) Event include_recursion: |
#include file "/sapmnt/home1/d029903/my_work/SCIPSoPlex_coverity/scipoptsuite-8.0.1/soplex/src/soplex.h" includes itself: soplex.h -> spxequilisc.h -> spxequilisc.hpp -> soplex.h |
|
(2) Event caretline: |
^ |
25 #include "soplex.h"
26
27 namespace soplex
28 {
29 static inline const char* makename(bool doBoth)
30 {
31 return doBoth ? "bi-Equilibrium" : "uni-Equilibrium";
32 }
33
34 /// maximum ratio between absolute biggest and smallest element in any scaled row/column.
35 template <class R>
36 static R maxPrescaledRatio(const SPxLPBase<R>& lp, const std::vector<R>& coScaleval, bool rowRatio)
37 {
38 R pmax = 0.0;
39 const int n = rowRatio ? lp.nRows() : lp.nCols();
40
41 for(int i = 0; i < n; ++i)
42 {
43 const SVectorBase<R>& vec = rowRatio ? lp.rowVector(i) : lp.colVector(i);
44 R mini = R(infinity);
45 R maxi = 0.0;
46
47 for(int j = 0; j < vec.size(); ++j)
48 {
49 assert(vec.index(j) >= 0);
50 const R x = spxAbs(vec.value(j)) * coScaleval[unsigned(vec.index(j))];
51
52 if(isZero(x))
53 continue;
54
55 if(x < mini)
56 mini = x;
57
58 if(x > maxi)
59 maxi = x;
60 }
61
62 if(mini == R(infinity))
63 continue;
64
65 const R p = maxi / mini;
66
67 if(p > pmax)
68 pmax = p;
69 }
70
71 return pmax;
72 }
73
74 template <class R>
75 void SPxEquiliSC<R>::computeEquiExpVec(const SVSetBase<R>* vecset, const DataArray<int>& coScaleExp,
76 DataArray<int>& scaleExp)
77 {
78 assert(vecset != nullptr);
79
80 for(int i = 0; i < vecset->num(); ++i)
81 {
82 const SVectorBase<R>& vec = (*vecset)[i];
83
84 R maxi = 0.0;
85
86 for(int j = 0; j < vec.size(); ++j)
87 {
88 const R x = spxAbs(spxLdexp(vec.value(j), coScaleExp[vec.index(j)]));
89
90 if(GT(x, maxi))
91 maxi = x;
92 }
93
94 // empty rows/cols are possible
95 if(maxi == 0.0)
96 maxi = 1.0;
97
98 assert(maxi > 0.0);
99
100 spxFrexp(Real(1.0 / maxi), &(scaleExp[i]));
101
102 scaleExp[i] -= 1;
103 }
104 }
105
106 template <class R>
107 void SPxEquiliSC<R>::computeEquiExpVec(const SVSetBase<R>* vecset, const std::vector<R>& coScaleVal,
108 DataArray<int>& scaleExp)
109 {
110 assert(vecset != nullptr);
111
112 for(int i = 0; i < vecset->num(); ++i)
113 {
114 const SVectorBase<R>& vec = (*vecset)[i];
115
116 R maxi = 0.0;
117
118 for(int j = 0; j < vec.size(); ++j)
119 {
120 assert(vec.index(j) >= 0);
121 const R x = spxAbs(vec.value(j) * coScaleVal[unsigned(vec.index(j))]);
122
123 if(GT(x, maxi))
124 maxi = x;
125 }
126
127 // empty rows/cols are possible
128 if(maxi == 0.0)
129 maxi = 1.0;
130
131 assert(maxi > 0.0);
132
133 spxFrexp(Real(1.0 / maxi), &(scaleExp[i]));
134
135 scaleExp[i] -= 1;
136 }
137 }
138
139 template <class R>
140 void SPxEquiliSC<R>::computePostequiExpVecs(const SPxLPBase<R>& lp,
141 const std::vector<R>& preRowscale, const std::vector<R>& preColscale,
142 DataArray<int>& rowscaleExp, DataArray<int>& colscaleExp)
143 {
144 const R colratio = maxPrescaledRatio(lp, preRowscale, false);
145 const R rowratio = maxPrescaledRatio(lp, preColscale, true);
146
147 const bool colFirst = colratio < rowratio;
148
149 // see SPxEquiliSC<R>::scale for reason behind this branch
150 if(colFirst)
151 {
152 computeEquiExpVec(lp.colSet(), preRowscale, colscaleExp);
153 computeEquiExpVec(lp.rowSet(), colscaleExp, rowscaleExp);
154 }
155 else
156 {
157 computeEquiExpVec(lp.rowSet(), preColscale, rowscaleExp);
158 computeEquiExpVec(lp.colSet(), rowscaleExp, colscaleExp);
159 }
160 }
161
162 template <class R>
163 SPxEquiliSC<R>::SPxEquiliSC(bool doBoth)
164 : SPxScaler<R>(makename(doBoth), false, doBoth)
165 {}
166
167 template <class R>
168 SPxEquiliSC<R>::SPxEquiliSC(const SPxEquiliSC<R>& old)
169 : SPxScaler<R>(old)
170 {}
171
172 template <class R>
173 SPxEquiliSC<R>& SPxEquiliSC<R>::operator=(const SPxEquiliSC<R>& rhs)
174 {
175 if(this != &rhs)
176 {
177 SPxScaler<R>::operator=(rhs);
178 }
179
180 return *this;
181 }
182
183 template <class R>
184 void SPxEquiliSC<R>::scale(SPxLPBase<R>& lp, bool persistent)
185 {
186
187 MSG_INFO1((*this->spxout), (*this->spxout) << "Equilibrium scaling LP" <<
188 (persistent ? " (persistent)" : "") << std::endl;)
189
190 this->setup(lp);
191
192 /* We want to do the direction first, which has a lower maximal ratio,
193 * since the lowest value in the scaled matrix is bounded from below by
194 * the inverse of the maximum ratio of the direction that is done first
195 * Example:
196 *
197 * Rowratio
198 * 0.1 1 10
199 * 10 1 10
200 *
201 * Colratio 100 1
202 *
203 * Row first => Col next =>
204 * 0.1 1 0.1 1
205 * 1 0.1 1 0.1
206 *
207 * Col first => Row next =>
208 * 0.01 1 0.01 1
209 * 1 1 1 1
210 *
211 */
212 R colratio = this->maxColRatio(lp);
213 R rowratio = this->maxRowRatio(lp);
214
215 bool colFirst = colratio < rowratio;
216
217 MSG_INFO2((*this->spxout), (*this->spxout) << "before scaling:"
218 << " min= " << lp.minAbsNzo()
219 << " max= " << lp.maxAbsNzo()
220 << " col-ratio= " << colratio
221 << " row-ratio= " << rowratio
222 << std::endl;)
223
224 if(colFirst)
225 {
226 computeEquiExpVec(lp.colSet(), *this->m_activeRowscaleExp, *this->m_activeColscaleExp);
227
228 if(this->m_doBoth)
229 computeEquiExpVec(lp.rowSet(), *this->m_activeColscaleExp, *this->m_activeRowscaleExp);
230 }
231 else
232 {
233 computeEquiExpVec(lp.rowSet(), *this->m_activeColscaleExp, *this->m_activeRowscaleExp);
234
235 if(this->m_doBoth)
236 computeEquiExpVec(lp.colSet(), *this->m_activeRowscaleExp, *this->m_activeColscaleExp);
237 }
238
239 /* scale */
240 this->applyScaling(lp);
241
242 MSG_INFO3((*this->spxout), (*this->spxout) << "Row scaling min= " << this->minAbsRowscale()
243 << " max= " << this->maxAbsRowscale()
244 << std::endl
245 << "Col scaling min= " << this->minAbsColscale()
246 << " max= " << this->maxAbsColscale()
247 << std::endl;)
248
249 MSG_INFO2((*this->spxout), (*this->spxout) << "after scaling: "
250 << " min= " << lp.minAbsNzo(false)
251 << " max= " << lp.maxAbsNzo(false)
252 << " col-ratio= " << this->maxColRatio(lp)
253 << " row-ratio= " << this->maxRowRatio(lp)
254 << std::endl;)
255
256 }
257
258 } // namespace soplex
259