1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */ 2 /* */ 3 /* This file is part of the class library */ 4 /* SoPlex --- the Sequential object-oriented simPlex. */ 5 /* */ 6 /* Copyright (c) 1996-2023 Zuse Institute Berlin (ZIB) */ 7 /* */ 8 /* Licensed under the Apache License, Version 2.0 (the "License"); */ 9 /* you may not use this file except in compliance with the License. */ 10 /* You may obtain a copy of the License at */ 11 /* */ 12 /* http://www.apache.org/licenses/LICENSE-2.0 */ 13 /* */ 14 /* Unless required by applicable law or agreed to in writing, software */ 15 /* distributed under the License is distributed on an "AS IS" BASIS, */ 16 /* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. */ 17 /* See the License for the specific language governing permissions and */ 18 /* limitations under the License. */ 19 /* */ 20 /* You should have received a copy of the Apache-2.0 license */ 21 /* along with SoPlex; see the file LICENSE. If not email to soplex@zib.de. */ 22 /* */ 23 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */ 24 25 #include <assert.h> 26 #include <iostream> 27 28 #include "soplex/spxdefines.h" 29 #include "soplex/spxsolver.h" 30 31 namespace soplex 32 { 33 /** Setting up the feasiblity bound for normal primal variables is 34 straightforward. However, slack variables need some more details 35 on how we treat them. This is slightly different from usual 36 textbook versions. Let \f$l_i \le A_i^T x \le u_i\f$. This can be 37 transformed to \f$A_i^Tx + s_i = 0\f$, with \f$-u_i \le s_i \le 38 -l_i\f$. Hence, with this definition of slack variables \f$s_i\f$, we 39 can directly use vectors \f$l\f$ and \f$u\f$ as feasibility bounds. 40 */ 41 template <class R> 42 void SPxSolverBase<R>::setPrimalBounds() 43 { 44 45 theUCbound = SPxLPBase<R>::upper(); 46 theLCbound = SPxLPBase<R>::lower(); 47 48 if(rep() == ROW) 49 { 50 theURbound = this->rhs(); 51 theLRbound = this->lhs(); 52 } 53 else 54 { 55 theURbound = this->lhs(); 56 theLRbound = this->rhs(); 57 theURbound *= -1.0; 58 theLRbound *= -1.0; 59 } 60 } 61 62 /** Setting up the basis for dual simplex requires to install upper and lower 63 feasibility bounds for dual variables (|Lbound| and |Ubound|). Here is a 64 list of how these must be set for inequalities of type \f$l \le a^Tx \le u\f$: 65 66 \f[ 67 \begin{tabular}{cccc} 68 $l$ & $u$ & |Lbound| & |Ubound| \\ 69 \hline 70 $-\infty=l$ & $u=\infty$ & 0 & 0 \\ 71 $-\infty<l$ & $u=\infty$ & 0 & $\infty$ \\ 72 $-\infty=l$ & $u<\infty$ & $-\infty$ & 0 \\ 73 \multicolumn{2}{c}{ 74 $-\infty<l \ne u<\infty$} & 0 & 0 \\ 75 \multicolumn{2}{c}{ 76 $-\infty<l = u<\infty$} & $-\infty$ & $\infty$ \\ 77 \end{tabular} 78 \f] 79 The case \f$l = -\infty\f$, \f$u = \infty\f$ occurs for unbounded primal variables. 80 Such must be treated differently from the general case. 81 82 Given possible upper and lower bounds to a dual variable with |Status stat|, 83 this function clears the bounds according to |stat| by setting them to 84 \f$\infty\f$ or \f$-\infty\f$, respectively. 85 */ 86 template <class R> 87 void SPxSolverBase<R>::clearDualBounds( 88 typename SPxBasisBase<R>::Desc::Status stat, 89 R& upp, 90 R& lw) const 91 { 92 93 switch(stat) 94 { 95 case SPxBasisBase<R>::Desc::P_ON_UPPER + SPxBasisBase<R>::Desc::P_ON_LOWER : 96 case SPxBasisBase<R>::Desc::D_FREE : 97 upp = R(infinity); 98 lw = R(-infinity); 99 break; 100 101 case SPxBasisBase<R>::Desc::P_ON_UPPER : 102 case SPxBasisBase<R>::Desc::D_ON_LOWER : 103 upp = R(infinity); 104 break; 105 106 case SPxBasisBase<R>::Desc::P_ON_LOWER : 107 case SPxBasisBase<R>::Desc::D_ON_UPPER : 108 lw = R(-infinity); 109 break; 110 111 default: 112 break; 113 } 114 } 115 116 template <class R> 117 void SPxSolverBase<R>::setDualColBounds() 118 { 119 120 assert(rep() == COLUMN); 121 122 const typename SPxBasisBase<R>::Desc& ds = this->desc(); 123 124 int i; 125 126 for(i = 0; i < this->nRows(); ++i) 127 { 128 theURbound[i] = this->maxRowObj(i); 129 theLRbound[i] = this->maxRowObj(i); 130 131 clearDualBounds(ds.rowStatus(i), theURbound[i], theLRbound[i]); 132 } 133 134 for(i = 0; i < this->nCols(); ++i) 135 { 136 theUCbound[i] = -this->maxObj(i); 137 theLCbound[i] = -this->maxObj(i); 138 139 // exchanged ... due to definition of slacks! 140 clearDualBounds(ds.colStatus(i), theLCbound[i], theUCbound[i]); 141 142 theUCbound[i] *= -1.0; 143 theLCbound[i] *= -1.0; 144 } 145 } 146 147 template <class R> 148 void SPxSolverBase<R>::setDualRowBounds() 149 { 150 151 assert(rep() == ROW); 152 153 int i; 154 155 for(i = 0; i < this->nRows(); ++i) 156 { 157 theURbound[i] = 0.0; 158 theLRbound[i] = 0.0; 159 160 clearDualBounds(this->dualRowStatus(i), theURbound[i], theLRbound[i]); 161 } 162 163 for(i = 0; i < this->nCols(); ++i) 164 { 165 theUCbound[i] = 0.0; 166 theLCbound[i] = 0.0; 167 168 clearDualBounds(this->dualColStatus(i), theUCbound[i], theLCbound[i]); 169 } 170 } 171 172 173 /** This sets up the bounds for basic variables for entering simplex algorithms. 174 It requires, that all upper lower feasibility bounds have allready been 175 setup. Method |setEnterBound4Row(i, n)| does so for the |i|-th basis 176 variable being row index |n|. Equivalently, method 177 |setEnterBound4Col(i, n)| does so for the |i|-th basis variable being 178 column index |n|. 179 */ 180 template <class R> 181 void SPxSolverBase<R>::setEnterBound4Row(int i, int n) 182 { 183 assert(this->baseId(i).isSPxRowId()); 184 assert(this->number(SPxRowId(this->baseId(i))) == n); 185 186 switch(this->desc().rowStatus(n)) 187 { 188 case SPxBasisBase<R>::Desc::P_ON_LOWER : 189 theLBbound[i] = R(-infinity); 190 theUBbound[i] = theURbound[n]; 191 break; 192 193 case SPxBasisBase<R>::Desc::P_ON_UPPER : 194 theLBbound[i] = theLRbound[n]; 195 theUBbound[i] = R(infinity); 196 break; 197 198 case SPxBasisBase<R>::Desc::P_FIXED: 199 theLBbound[i] = R(-infinity); 200 theUBbound[i] = R(infinity); 201 break; 202 203 default: 204 theUBbound[i] = theURbound[n]; 205 theLBbound[i] = theLRbound[n]; 206 break; 207 } 208 } 209 210 template <class R> 211 void SPxSolverBase<R>::setEnterBound4Col(int i, int n) 212 { 213 assert(this->baseId(i).isSPxColId()); 214 assert(this->number(SPxColId(this->baseId(i))) == n); 215 216 switch(this->desc().colStatus(n)) 217 { 218 case SPxBasisBase<R>::Desc::P_ON_LOWER : 219 theLBbound[i] = R(-infinity); 220 theUBbound[i] = theUCbound[n]; 221 break; 222 223 case SPxBasisBase<R>::Desc::P_ON_UPPER : 224 theLBbound[i] = theLCbound[n]; 225 theUBbound[i] = R(infinity); 226 break; 227 228 case SPxBasisBase<R>::Desc::P_FIXED: 229 theLBbound[i] = R(-infinity); 230 theUBbound[i] = R(infinity); 231 break; 232 233 default: 234 theUBbound[i] = theUCbound[n]; 235 theLBbound[i] = theLCbound[n]; 236 break; 237 } 238 } 239 template <class R> 240 void SPxSolverBase<R>::setEnterBounds() 241 { 242 243 for(int i = 0; i < dim(); ++i) 244 { 245 SPxId base_id = this->baseId(i); 246 247 if(base_id.isSPxRowId()) 248 setEnterBound4Row(i, this->number(SPxRowId(base_id))); 249 else 250 setEnterBound4Col(i, this->number(SPxColId(base_id))); 251 } 252 } 253 254 255 /** This sets up the bounds for basic variables for leaving simplex algorithms. 256 While method |setLeaveBound4Row(i,n)| does so for the |i|-th basic variable 257 being the |n|-th row, |setLeaveBound4Col(i,n)| does so for the |i|-th basic 258 variable being the |n|-th column. 259 */ 260 template <class R> 261 void SPxSolverBase<R>::setLeaveBound4Row(int i, int n) 262 { 263 assert(this->baseId(i).isSPxRowId()); 264 assert(this->number(SPxRowId(this->baseId(i))) == n); 265 266 switch(this->desc().rowStatus(n)) 267 { 268 case SPxBasisBase<R>::Desc::P_ON_LOWER : 269 theLBbound[i] = R(-infinity); 270 theUBbound[i] = -this->maxRowObj(n); 271 break; 272 273 case SPxBasisBase<R>::Desc::P_ON_UPPER : 274 theLBbound[i] = -this->maxRowObj(n); 275 theUBbound[i] = R(infinity); 276 break; 277 278 case SPxBasisBase<R>::Desc::P_ON_UPPER + SPxBasisBase<R>::Desc::P_ON_LOWER : 279 theLBbound[i] = R(-infinity); 280 theUBbound[i] = R(infinity); 281 break; 282 283 case SPxBasisBase<R>::Desc::P_FREE : 284 theLBbound[i] = -this->maxRowObj(n); 285 theUBbound[i] = -this->maxRowObj(n); 286 break; 287 288 default: 289 assert(rep() == COLUMN); 290 theLBbound[i] = -this->rhs(n); // slacks ! 291 theUBbound[i] = -this->lhs(n); // slacks ! 292 break; 293 } 294 } 295 template <class R> 296 void SPxSolverBase<R>::setLeaveBound4Col(int i, int n) 297 { 298 299 assert(this->baseId(i).isSPxColId()); 300 assert(this->number(SPxColId(this->baseId(i))) == n); 301 302 switch(this->desc().colStatus(n)) 303 { 304 case SPxBasisBase<R>::Desc::P_ON_LOWER : 305 theLBbound[i] = R(-infinity); 306 theUBbound[i] = 0; 307 break; 308 309 case SPxBasisBase<R>::Desc::P_ON_UPPER : 310 theLBbound[i] = 0; 311 theUBbound[i] = R(infinity); 312 break; 313 314 case SPxBasisBase<R>::Desc::P_FIXED : 315 theLBbound[i] = R(-infinity); 316 theUBbound[i] = R(infinity); 317 break; 318 319 case SPxBasisBase<R>::Desc::P_FREE : 320 theLBbound[i] = theUBbound[i] = 0; 321 break; 322 323 default: 324 theUBbound[i] = SPxLPBase<R>::upper(n); 325 theLBbound[i] = SPxLPBase<R>::lower(n); 326 break; 327 } 328 } 329 330 template <class R> 331 void SPxSolverBase<R>::setLeaveBounds() 332 { 333 334 for(int i = 0; i < dim(); ++i) 335 { 336 SPxId base_id = this->baseId(i); 337 338 if(base_id.isSPxRowId()) 339 setLeaveBound4Row(i, this->number(SPxRowId(base_id))); 340 else 341 setLeaveBound4Col(i, this->number(SPxColId(base_id))); 342 } 343 } 344 345 template <class R> 346 void SPxSolverBase<R>::testBounds() const 347 { 348 349 if(type() == ENTER) 350 { 351 R viol_max = (1 + this->iterCount) * entertol(); 352 int nlinesprinted = 0; 353 int m = dim(); 354 355 for(int i = 0; i < m; ++i) 356 { 357 // Minor bound violations happen frequently, so print these 358 // warnings only with verbose level INFO2 and higher. 359 if((*theFvec)[i] > theUBbound[i] + viol_max) //@ && theUBbound[i] != theLBbound[i]) 360 { 361 SPX_MSG_INFO2((*this->spxout), (*this->spxout) << "WBOUND01 Invalid upper enter bound " << i 362 << " Fvec: " << (*theFvec)[i] 363 << " UBbound: " << theUBbound[i] 364 << " tolerance: " << viol_max 365 << " violation: " << (*theFvec)[i] - theUBbound[i] << std::endl;) 366 nlinesprinted++; 367 } 368 369 if((*theFvec)[i] < theLBbound[i] - viol_max) //@ && theUBbound[i] != theLBbound[i]) 370 { 371 SPX_MSG_INFO2((*this->spxout), (*this->spxout) << "WBOUND02 Invalid lower enter bound " << i 372 << " Fvec: " << (*theFvec)[i] 373 << " LBbound: " << theLBbound[i] 374 << " tolerance: " << viol_max 375 << " violation: " << theLBbound[i] - (*theFvec)[i] << std::endl;) 376 nlinesprinted++; 377 } 378 379 if(nlinesprinted >= 3) 380 { 381 SPX_MSG_INFO2((*this->spxout), (*this->spxout) << 382 "WBOUND10 suppressing further warnings of type WBOUND{01,02} in this round" << std::endl); 383 break; 384 } 385 } 386 } 387 else 388 { 389 R viol_max = (1 + this->iterCount) * leavetol(); 390 int nlinesprinted = 0; 391 int m = dim(); 392 int n = coDim(); 393 394 for(int i = 0; i < m; ++i) 395 { 396 if((*theCoPvec)[i] > (*theCoUbound)[i] + viol_max) // && (*theCoUbound)[i] != (*theCoLbound)[i]) 397 { 398 SPX_MSG_INFO2((*this->spxout), (*this->spxout) << "WBOUND03 Invalid upper cobound " << i 399 << " CoPvec: " << (*theCoPvec)[i] 400 << " CoUbound: " << (*theCoUbound)[i] 401 << " tolerance: " << viol_max 402 << " violation: " << (*theCoPvec)[i] - (*theCoUbound)[i] << std::endl;) 403 nlinesprinted++; 404 } 405 406 if((*theCoPvec)[i] < (*theCoLbound)[i] - viol_max) // && (*theCoUbound)[i] != (*theCoLbound)[i]) 407 { 408 SPX_MSG_INFO2((*this->spxout), (*this->spxout) << "WBOUND04 Invalid lower cobound " << i 409 << " CoPvec: " << (*theCoPvec)[i] 410 << " CoLbound: " << (*theCoLbound)[i] 411 << " tolerance: " << viol_max 412 << " violation: " << (*theCoLbound)[i] - (*theCoPvec)[i] << std::endl;) 413 nlinesprinted++; 414 } 415 416 if(nlinesprinted >= 3) 417 { 418 SPX_MSG_INFO2((*this->spxout), (*this->spxout) << 419 "WBOUND11 suppressing further warnings of type WBOUND{03,04} in this round" << std::endl); 420 break; 421 } 422 } 423 424 nlinesprinted = 0; 425 426 for(int i = 0; i < n; ++i) 427 { 428 if((*thePvec)[i] > (*theUbound)[i] + viol_max) // && (*theUbound)[i] != (*theLbound)[i]) 429 { 430 SPX_MSG_INFO2((*this->spxout), (*this->spxout) << "WBOUND05 Invalid upper bound " << i 431 << " Pvec: " << (*thePvec)[i] 432 << " Ubound: " << (*theUbound)[i] 433 << " tolerance: " << viol_max 434 << " violation: " << (*thePvec)[i] - (*theUbound)[i] << std::endl;) 435 nlinesprinted++; 436 } 437 438 if((*thePvec)[i] < (*theLbound)[i] - viol_max) // && (*theUbound)[i] != (*theLbound)[i]) 439 { 440 SPX_MSG_INFO2((*this->spxout), (*this->spxout) << "WBOUND06 Invalid lower bound " << i 441 << " Pvec: " << (*thePvec)[i] 442 << " Lbound: " << (*theLbound)[i] 443 << " tolerance: " << viol_max 444 << " violation: " << (*theLbound)[i] - (*thePvec)[i] << std::endl;) 445 nlinesprinted++; 446 } 447 448 if(nlinesprinted >= 3) 449 { 450 SPX_MSG_INFO2((*this->spxout), (*this->spxout) << 451 "WBOUND12 suppressing further warnings of type WBOUND{05,06} in this round" << std::endl); 452 break; 453 } 454 } 455 } 456 } 457 } // namespace soplex 458