Scippy

SoPlex

Sequential object-oriented simPlex

spxsolve.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-2015 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 #include <iostream>
18 
19 #include "spxdefines.h"
20 #include "rational.h"
21 #include "spxsolver.h"
22 #include "spxpricer.h"
23 #include "spxratiotester.h"
24 #include "spxdefaultrt.h"
25 #include "spxstarter.h"
26 #include "spxout.h"
27 
28 #define MAXCYCLES 400
29 #define MAXSTALLS 10000
30 #define MAXSTALLRECOVERS 10
31 #define MAXREFACPIVOTS 10
32 
33 namespace soplex
34 {
35 
36 /**@todo check separately for ENTER and LEAVE algorithm */
37 bool SPxSolver::precisionReached(Real& newpricertol) const
38 {
39  Real maxViolRedCost;
40  Real sumViolRedCost;
41  Real maxViolBounds;
42  Real sumViolBounds;
43  Real maxViolConst;
44  Real sumViolConst;
45 
46  qualRedCostViolation(maxViolRedCost, sumViolRedCost);
47  qualBoundViolation(maxViolBounds, sumViolBounds);
48  qualConstraintViolation(maxViolConst, sumViolConst);
49 
50  // is the solution good enough ?
51  bool reached = maxViolRedCost < opttol() && maxViolBounds < feastol() && maxViolConst < feastol();
52 
53  if (!reached)
54  {
55  newpricertol = thepricer->epsilon() / 10.0;
56 
57  MSG_INFO3( (*spxout), (*spxout) << "Precision not reached: Pricer tolerance = "
58  << thepricer->epsilon()
59  << " new tolerance = " << newpricertol
60  << std::endl
61  << " maxViolRedCost= " << maxViolRedCost
62  << " maxViolBounds= " << maxViolBounds
63  << " maxViolConst= " << maxViolConst
64  << std::endl
65  << " sumViolRedCost= " << sumViolRedCost
66  << " sumViolBounds= " << sumViolBounds
67  << " sumViolConst= " << sumViolConst
68  << std::endl; );
69  }
70  return reached;
71 }
72 
73 SPxSolver::Status SPxSolver::solve()
74 {
75 
76  SPxId enterId;
77  int leaveNum;
78  int loopCount = 0;
79  Real minShift = infinity;
80  int cycleCount = 0;
81  bool priced = false;
82  Real lastDelta = 1;
83 
84  /* store the last (primal or dual) feasible objective value to recover/abort in case of stalling */
85  Real stallRefValue;
86  Real stallRefShift;
87  int stallRefIter;
88  int stallNumRecovers;
89 
90  if (dim() <= 0 && coDim() <= 0) // no problem loaded
91  {
92  m_status = NO_PROBLEM;
93  throw SPxStatusException("XSOLVE01 No Problem loaded");
94  }
95 
96  if (slinSolver() == 0) // linear system solver is required.
97  {
98  m_status = NO_SOLVER;
99  throw SPxStatusException("XSOLVE02 No Solver loaded");
100  }
101  if (thepricer == 0) // pricer is required.
102  {
103  m_status = NO_PRICER;
104  throw SPxStatusException("XSOLVE03 No Pricer loaded");
105  }
106  if (theratiotester == 0) // ratiotester is required.
107  {
108  m_status = NO_RATIOTESTER;
109  throw SPxStatusException("XSOLVE04 No RatioTester loaded");
110  }
111  theTime->reset();
112  theTime->start();
113 
114  m_numCycle = 0;
115  iterCount = 0;
116  if (!isInitialized())
117  {
118  /*
119  if(SPxBasis::status() <= NO_PROBLEM)
120  SPxBasis::load(this);
121  */
122  /**@todo != REGULAR is not enough. Also OPTIMAL/DUAL/PRIMAL should
123  * be tested and acted accordingly.
124  */
125  if (thestarter != 0 && status() != REGULAR) // no basis and no starter.
126  thestarter->generate(*this); // generate start basis.
127 
128  init();
129 
130  // Inna/Tobi: init might fail, if the basis is singular
131  if( !isInitialized() )
132  {
133  assert(SPxBasis::status() == SPxBasis::SINGULAR);
134  m_status = UNKNOWN;
135  return status();
136  }
137  }
138 
139  //setType(type());
140 
141  if (!matrixIsSetup)
142  SPxBasis::load(this);
143 
144  //factorized = false;
145 
146  assert(thepricer->solver() == this);
147  assert(theratiotester->solver() == this);
148 
149  // maybe this should be done in init() ?
150  thepricer->setType(type());
151  theratiotester->setType(type());
152 
153  MSG_INFO3( (*spxout),
154  (*spxout) << "starting value = " << value() << std::endl
155  << "starting shift = " << shift() << std::endl;
156  )
157 
158  if (SPxBasis::status() == SPxBasis::OPTIMAL)
159  setBasisStatus(SPxBasis::REGULAR);
160 
161  m_status = RUNNING;
162  bool stop = terminate();
163  leaveCount = 0;
164  enterCount = 0;
165  primalCount = 0;
166  boundflips = 0;
167  totalboundflips = 0;
168 
169  stallNumRecovers = 0;
170 
171  /* if we run into a singular basis, we will retry from regulardesc with tighter tolerance in the ratio test */
172  SPxSolver::Type tightenedtype = type();
173  bool tightened = false;
174 
175  while (!stop)
176  {
177  const SPxBasis::Desc regulardesc = desc();
178 
179  // we need to reset these pointers to avoid unnecessary/wrong solves in leave() or enter()
180  solveVector2 = 0;
181  solveVector3 = 0;
182  coSolveVector2 = 0;
183  coSolveVector3 = 0;
184 
185  try
186  {
187 
188  if (type() == ENTER)
189  {
190  forceRecompNonbasicValue();
191 
192  int enterCycleCount = 0;
193  int enterFacPivotCount = 0;
194 
195  instableEnterVal = 0;
196  instableEnterId = SPxId();
197  instableEnter = false;
198 
199  stallRefIter = iteration()-1;
200  stallRefShift = shift();
201  stallRefValue = value();
202 
203  /* in the entering algorithm, entertol() should be maintained by the ratio test and leavetol() should be
204  * reached by the pricer
205  */
206  Real maxpricertol = leavetol();
207  Real minpricertol = 0.01 * maxpricertol;
208 
209  thepricer->setEpsilon(maxpricertol);
210  priced = false;
211 
212  // to avoid shifts we restrict tolerances in the ratio test
213  if( loopCount > 0 )
214  {
215  lastDelta = (lastDelta < entertol()) ? lastDelta : entertol();
216  lastDelta *= 0.01;
217  theratiotester->setDelta(lastDelta);
218  assert(theratiotester->getDelta() > 0);
219  MSG_DEBUG( std::cout << "decreased delta for ratiotest to: " << theratiotester->getDelta() << std::endl; )
220  }
221  else
222  {
223  lastDelta = 1;
224  theratiotester->setDelta(entertol());
225  }
226 
227  printDisplayLine(true);
228  do
229  {
230  printDisplayLine();
231 
232  enterId = thepricer->selectEnter();
233 
234  if (!enterId.isValid() && instableEnterId.isValid() && lastUpdate() == 0)
235  {
236  /* no entering variable was found, but because of valid instableEnterId we know
237  that this is due to the scaling of the test values. Thus, we use
238  instableEnterId and SPxFastRT::selectEnter shall accept even an instable
239  leaving variable. */
240  MSG_INFO3( (*spxout),
241  (*spxout) << " --- trying instable enter iteration" << std::endl;
242  )
243 
244  enterId = instableEnterId;
245  instableEnter = true;
246  // we also need to reset the test() or coTest() value for getEnterVals()
247  assert(instableEnterVal < 0);
248  if( enterId.isSPxColId() )
249  {
250  int idx = number(SPxColId(enterId));
251  if( rep() == COLUMN )
252  {
253  theTest[idx] = instableEnterVal;
254  if( sparsePricingEnterCo && isInfeasibleCo[idx] == SPxPricer::NOT_VIOLATED )
255  {
256  infeasibilitiesCo.addIdx(idx);
257  isInfeasibleCo[idx] = SPxPricer::VIOLATED;
258  }
259  if( hyperPricingEnter )
260  updateViolsCo.addIdx(idx);
261  }
262  else
263  {
264  theCoTest[idx] = instableEnterVal;
265  if( sparsePricingEnter && isInfeasible[idx] == SPxPricer::NOT_VIOLATED )
266  {
267  infeasibilities.addIdx(idx);
268  isInfeasible[idx] = SPxPricer::VIOLATED;
269  }
270  if( hyperPricingEnter )
271  updateViols.addIdx(idx);
272  }
273  }
274  else
275  {
276  int idx = number(SPxRowId(enterId));
277  if( rep() == COLUMN )
278  {
279  theCoTest[idx] = instableEnterVal;
280  if( sparsePricingEnter && isInfeasible[idx] == SPxPricer::NOT_VIOLATED )
281  {
282  infeasibilities.addIdx(idx);
283  isInfeasible[idx] = SPxPricer::VIOLATED;
284  }
285  if( hyperPricingEnter )
286  updateViols.addIdx(idx);
287  }
288  else
289  {
290  theTest[idx] = instableEnterVal;
291  if( sparsePricingEnterCo && isInfeasibleCo[idx] == SPxPricer::NOT_VIOLATED )
292  {
293  infeasibilitiesCo.addIdx(idx);
294  isInfeasibleCo[idx] = SPxPricer::VIOLATED;
295  }
296  if( hyperPricingEnter )
297  updateViolsCo.addIdx(idx);
298  }
299  }
300  }
301  else
302  {
303  instableEnter = false;
304  }
305 
306  if (!enterId.isValid())
307  {
308  // we are not infeasible and have no shift
309  if ( shift() <= epsilon()
310  && ( SPxBasis::status() == SPxBasis::REGULAR
311  || SPxBasis::status() == SPxBasis::DUAL
312  || SPxBasis::status() == SPxBasis::PRIMAL))
313  {
314  Real newpricertol = minpricertol;
315 
316  MSG_INFO2( (*spxout), (*spxout) << " --- checking feasibility and optimality\n")
317  computeTest();
318  computeCoTest();
319 
320  // is the solution good enough ?
321  // max three times reduced
322  if ((thepricer->epsilon() > minpricertol) && !precisionReached(newpricertol))
323  { // no!
324  // we reduce the pricer tolerance. Note that if the pricer does not find a candiate
325  // with the reduced tolerance, we quit, regardless of the violations.
326  if (newpricertol < minpricertol)
327  newpricertol = minpricertol;
328 
329  thepricer->setEpsilon(newpricertol);
330 
331  MSG_INFO2( (*spxout), (*spxout) << " --- setting pricer tolerance to "
332  << thepricer->epsilon()
333  << std::endl; )
334  }
335  // solution seems good, no check whether we are precise enough
336  else if (lastUpdate() == 0)
337  {
338  priced = true;
339  break;
340  }
341  // We have an iterationlimit and everything looks good? Then stop!
342  // 6 is just a number picked.
343  else if (maxIters > 0 && lastUpdate() < 6)
344  {
345  priced = true;
346  break;
347  }
348  }
349  MSG_INFO3( (*spxout), (*spxout) << " --- solve(enter) triggers refactorization" << std::endl; )
350 
351  // if the factorization is not fresh, we better refactorize and call the pricer again; however, this can
352  // create cycling, so it is performed only a limited number of times per ENTER round
353  if( lastUpdate() > 0 && enterFacPivotCount < MAXREFACPIVOTS )
354  {
355  factorize();
356 
357  // if the factorization was found out to be singular, we have to quit
358  if( SPxBasis::status() < SPxBasis::REGULAR )
359  {
360  MSG_ERROR( std::cerr << "Something wrong with factorization, Basis status: " << SPxBasis::status() << std::endl; )
361  stop = true;
362  break;
363  }
364 
365  // call pricer again
366  enterId = thepricer->selectEnter();
367 
368  // count how often the pricer has found something only after refactorizing
369  if( enterId.isValid() )
370  enterFacPivotCount++;
371  }
372 
373  if( !enterId.isValid() )
374  {
375  priced = true;
376  break;
377  }
378  }
379 
380  /* check if we have iterations left */
381  if (maxIters >= 0 && iterations() >= maxIters)
382  {
383  MSG_INFO2( (*spxout), (*spxout) << " --- maximum number of iterations (" << maxIters
384  << ") reached" << std::endl; )
385  m_status = ABORT_ITER;
386  stop = true;
387  break;
388  }
389 
390  enter(enterId);
391  assert((testBounds(), 1));
392  thepricer->entered4(lastEntered(), lastIndex());
393  stop = terminate();
394  clearUpdateVecs();
395 
396  /* if a successful pivot was performed or a nonbasic variable was flipped to its other bound, we reset the
397  * cycle counter
398  */
399  if( lastEntered().isValid() )
400  enterCycleCount = 0;
401  else
402  {
403  enterCycleCount++;
404  if( enterCycleCount > MAXCYCLES )
405  {
406  MSG_INFO2( (*spxout), (*spxout) << " --- abort solving due to cycling in "
407  << "entering algorithm" << std::endl; );
408  m_status = ABORT_CYCLING;
409  stop = true;
410  }
411  }
412 
413  /* only if the basis has really changed, we increase the iterations counter; this is not the case when only
414  * a nonbasic variable was flipped to its other bound
415  */
416  if( lastIndex() >= 0 )
417  {
418  enterCount++;
419  assert(lastEntered().isValid());
420  }
421 
422  /* check every MAXSTALLS iterations whether shift and objective value have not changed */
423  if( (iteration() - stallRefIter) % MAXSTALLS == 0 )
424  {
425  if( spxAbs(value() - stallRefValue) <= epsilon() && spxAbs(shift() - stallRefShift) <= epsilon() )
426  {
427  if( stallNumRecovers < MAXSTALLRECOVERS )
428  {
429  /* try to recover by unshifting/switching algorithm up to MAXSTALLRECOVERS times (just a number picked) */
430  MSG_INFO3( (*spxout), (*spxout) << " --- stalling detected - trying to recover by switching to LEAVING algorithm." << std::endl; )
431 
432  ++stallNumRecovers;
433  break;
434  }
435  else
436  {
437  /* giving up */
438  MSG_INFO2( (*spxout), (*spxout) << " --- abort solving due to stalling in entering algorithm." << std::endl; );
439 
440  m_status = ABORT_CYCLING;
441  stop = true;
442  }
443  }
444  else
445  {
446  /* merely update reference values */
447  stallRefIter = iteration()-1;
448  stallRefShift = shift();
449  stallRefValue = value();
450  }
451  }
452 
453  //@ assert(isConsistent());
454  }
455  while (!stop);
456 
457  MSG_INFO3( (*spxout),
458  (*spxout) << " --- enter finished. iteration: " << iteration()
459  << ", value: " << value()
460  << ", shift: " << shift()
461  << ", epsilon: " << epsilon()
462  << ", feastol: " << feastol()
463  << ", opttol: " << opttol()
464  << std::endl
465  << "ISOLVE56 stop: " << stop
466  << ", basis status: " << SPxBasis::status() << " (" << int(SPxBasis::status()) << ")"
467  << ", solver status: " << m_status << " (" << int(m_status) << ")" << std::endl;
468  )
469 
470  if (!stop)
471  {
472  /**@todo technically it would be ok to finish already when (priced && maxinfeas + shift() <= entertol()) is
473  * satisfied; maybe at least in the case when SoPlex keeps jumping back between ENTER and LEAVE always
474  * shifting (looping), we may relax this condition here;
475  * note also that unShift may increase shift() slightly due to roundoff errors
476  */
477  if (shift() <= epsilon())
478  {
479  // factorize();
480  unShift();
481 
482  Real maxinfeas = maxInfeas();
483 
484  MSG_INFO3( (*spxout),
485  (*spxout) << " --- maxInfeas: " << maxinfeas
486  << ", shift: " << shift()
487  << ", entertol: " << entertol() << std::endl;
488  )
489 
490  if (priced && maxinfeas + shift() <= entertol())
491  {
492  setBasisStatus(SPxBasis::OPTIMAL);
493  m_status = OPTIMAL;
494  break;
495  }
496  }
497  setType(LEAVE);
498  init();
499  thepricer->setType(type());
500  theratiotester->setType(type());
501  }
502  }
503  else
504  {
505  assert(type() == LEAVE);
506 
507  forceRecompNonbasicValue();
508 
509  int leaveCycleCount = 0;
510  int leaveFacPivotCount = 0;
511 
512  instableLeaveNum = -1;
513  instableLeave = false;
514  instableLeaveVal = 0;
515 
516  stallRefIter = iteration()-1;
517  stallRefShift = shift();
518  stallRefValue = value();
519 
520  /* in the leaving algorithm, leavetol() should be maintained by the ratio test and entertol() should be reached
521  * by the pricer
522  */
523  Real maxpricertol = entertol();
524  Real minpricertol = 0.01 * maxpricertol;
525 
526  thepricer->setEpsilon(maxpricertol);
527  priced = false;
528 
529  // to avoid shifts we restrict tolerances in the ratio test
530  if( loopCount > 0 )
531  {
532  lastDelta = (lastDelta < leavetol()) ? lastDelta : leavetol();
533  lastDelta *= 0.01;
534  theratiotester->setDelta(lastDelta);
535  assert(theratiotester->getDelta() > 0);
536  MSG_DEBUG( std::cout << "decreased delta for ratiotest to: " << theratiotester->getDelta() << std::endl; )
537  }
538  else
539  {
540  lastDelta = 1;
541  theratiotester->setDelta(leavetol());
542  }
543 
544  printDisplayLine(true);
545  do
546  {
547  printDisplayLine();
548 
549  leaveNum = thepricer->selectLeave();
550 
551  if (leaveNum < 0 && instableLeaveNum >= 0 && lastUpdate() == 0)
552  {
553  /* no leaving variable was found, but because of instableLeaveNum >= 0 we know
554  that this is due to the scaling of theCoTest[...]. Thus, we use
555  instableLeaveNum and SPxFastRT::selectEnter shall accept even an instable
556  entering variable. */
557  MSG_INFO3( (*spxout),
558  (*spxout) << " --- trying instable leave iteration" << std::endl;
559  )
560 
561  leaveNum = instableLeaveNum;
562  instableLeave = true;
563  // we also need to reset the fTest() value for getLeaveVals()
564  assert(instableLeaveVal < 0);
565  theCoTest[instableLeaveNum] = instableLeaveVal;
566 
567  if ( sparsePricingLeave )
568  {
569  if ( isInfeasible[instableLeaveNum] == SPxPricer::NOT_VIOLATED )
570  {
571  infeasibilities.addIdx(instableLeaveNum);
572  isInfeasible[instableLeaveNum] = SPxPricer::VIOLATED;
573  }
574  if( hyperPricingLeave )
575  updateViols.addIdx(instableLeaveNum);
576  }
577  }
578  else
579  {
580  instableLeave = false;
581  }
582 
583  if (leaveNum < 0)
584  {
585  // we are not infeasible and have no shift
586  if ( shift() <= epsilon()
587  && ( SPxBasis::status() == SPxBasis::REGULAR
588  || SPxBasis::status() == SPxBasis::DUAL
589  || SPxBasis::status() == SPxBasis::PRIMAL))
590  {
591  Real newpricertol = minpricertol;
592 
593  MSG_INFO2( (*spxout), (*spxout) << " --- checking feasibility and optimality\n")
594  computeFtest();
595 
596  // is the solution good enough ?
597  // max three times reduced
598  if ((thepricer->epsilon() > minpricertol) && !precisionReached(newpricertol))
599  { // no
600  // we reduce the pricer tolerance. Note that if the pricer does not find a candiate
601  // with the reduced pricer tolerance, we quit, regardless of the violations.
602  if (newpricertol < minpricertol)
603  newpricertol = minpricertol;
604 
605  thepricer->setEpsilon(newpricertol);
606 
607  MSG_INFO2( (*spxout), (*spxout) << " --- setting pricer tolerance to "
608  << thepricer->epsilon()
609  << std::endl; );
610  }
611  // solution seems good, no check whether we are precise enough
612  else if (lastUpdate() == 0)
613  {
614  priced = true;
615  break;
616  }
617  // We have an iteration limit and everything looks good? Then stop!
618  // 6 is just a number picked.
619  else if (maxIters > 0 && lastUpdate() < 6)
620  {
621  priced = true;
622  break;
623  }
624  }
625  MSG_INFO3( (*spxout), (*spxout) << " --- solve(leave) triggers refactorization" << std::endl; )
626 
627  // if the factorization is not fresh, we better refactorize and call the pricer again; however, this can
628  // create cycling, so it is performed only a limited number of times per LEAVE round
629  if( lastUpdate() > 0 && leaveFacPivotCount < MAXREFACPIVOTS )
630  {
631  factorize();
632 
633  // Inna/Tobi: if the factorization was found out to be singular, we have to quit
634  if (SPxBasis::status() < SPxBasis::REGULAR)
635  {
636  MSG_ERROR( std::cerr << "Something wrong with factorization, Basis status: " << SPxBasis::status() << std::endl; )
637  stop = true;
638  break;
639  }
640 
641  // call pricer again
642  leaveNum = thepricer->selectLeave();
643 
644  // count how often the pricer has found something only after refactorizing
645  if( leaveNum >= 0 )
646  leaveFacPivotCount++;
647  }
648 
649  if (leaveNum < 0)
650  {
651  priced = true;
652  break;
653  }
654  }
655 
656  /* check if we have iterations left */
657  if (maxIters >= 0 && iterations() >= maxIters)
658  {
659  MSG_INFO2( (*spxout), (*spxout) << " --- maximum number of iterations (" << maxIters
660  << ") reached" << std::endl; )
661  m_status = ABORT_ITER;
662  stop = true;
663  break;
664  }
665 
666  leave(leaveNum);
667  assert((testBounds(), 1));
668  thepricer->left4(lastIndex(), lastLeft());
669  stop = terminate();
670  clearUpdateVecs();
671 
672  /* if a successful pivot was performed or a nonbasic variable was flipped to its other bound, we reset the
673  * cycle counter
674  */
675  if( lastIndex() >= 0 )
676  leaveCycleCount = 0;
677  else
678  {
679  leaveCycleCount++;
680  if( leaveCycleCount > MAXCYCLES )
681  {
682  MSG_INFO2( (*spxout), (*spxout) << " --- abort solving due to cycling in leaving algorithm" << std::endl; );
683  m_status = ABORT_CYCLING;
684  stop = true;
685  }
686  }
687 
688  /* only if the basis has really changed, we increase the iterations counter; this is not the case when only
689  * a nonbasic variable was flipped to its other bound
690  */
691  if( lastEntered().isValid() )
692  {
693  leaveCount++;
694  assert(lastIndex() >= 0);
695  }
696 
697  /* check every MAXSTALLS iterations whether shift and objective value have not changed */
698  if( (iteration() - stallRefIter) % MAXSTALLS == 0 )
699  {
700  if( spxAbs(value() - stallRefValue) <= epsilon() && spxAbs(shift() - stallRefShift) <= epsilon() )
701  {
702  if( stallNumRecovers < MAXSTALLRECOVERS )
703  {
704  /* try to recover by switching algorithm up to MAXSTALLRECOVERS times */
705  MSG_INFO3( (*spxout), (*spxout) << " --- stalling detected - trying to recover by switching to ENTERING algorithm." << std::endl; )
706 
707  ++stallNumRecovers;
708  break;
709  }
710  else
711  {
712  /* giving up */
713  MSG_INFO2( (*spxout), (*spxout) << " --- abort solving due to stalling in leaving algorithm" << std::endl; );
714 
715  m_status = ABORT_CYCLING;
716  stop = true;
717  }
718  }
719  else
720  {
721  /* merely update reference values */
722  stallRefIter = iteration()-1;
723  stallRefShift = shift();
724  stallRefValue = value();
725  }
726  }
727 
728  //@ assert(isConsistent());
729  }
730  while (!stop);
731 
732  MSG_INFO3( (*spxout),
733  (*spxout) << " --- leave finished. iteration: " << iteration()
734  << ", value: " << value()
735  << ", shift: " << shift()
736  << ", epsilon: " << epsilon()
737  << ", feastol: " << feastol()
738  << ", opttol: " << opttol()
739  << std::endl
740  << "ISOLVE57 stop: " << stop
741  << ", basis status: " << SPxBasis::status() << " (" << int(SPxBasis::status()) << ")"
742  << ", solver status: " << m_status << " (" << int(m_status) << ")" << std::endl;
743  )
744 
745  if (!stop)
746  {
747  if( shift() < minShift )
748  {
749  minShift = shift();
750  cycleCount = 0;
751  }
752  else
753  {
754  cycleCount++;
755  if( cycleCount > MAXCYCLES )
756  {
757  m_status = ABORT_CYCLING;
758  throw SPxStatusException("XSOLVE13 Abort solving due to cycling");
759  }
760  MSG_INFO3( (*spxout),
761  (*spxout) << " --- maxInfeas: " << maxInfeas()
762  << ", shift: " << shift()
763  << ", leavetol: " << leavetol()
764  << ", cycle count: " << cycleCount << std::endl;
765  )
766  }
767 
768  /**@todo technically it would be ok to finish already when (priced && maxinfeas + shift() <= entertol()) is
769  * satisfied; maybe at least in the case when SoPlex keeps jumping back between ENTER and LEAVE always
770  * shifting (looping), we may relax this condition here;
771  * note also that unShift may increase shift() slightly due to roundoff errors
772  */
773  if (shift() <= epsilon())
774  {
775  cycleCount = 0;
776  // factorize();
777  unShift();
778 
779  Real maxinfeas = maxInfeas();
780 
781  MSG_INFO3( (*spxout),
782  (*spxout) << " --- maxInfeas: " << maxinfeas
783  << ", shift: " << shift()
784  << ", leavetol: " << leavetol() << std::endl;
785  )
786 
787  // We stop if we are indeed optimal, or if we have already been
788  // two times at this place. In this case it seems futile to
789  // continue.
790  if (loopCount > 2)
791  {
792  m_status = ABORT_CYCLING;
793  throw SPxStatusException("XSOLVE14 Abort solving due to looping");
794  }
795  else if (priced && maxinfeas + shift() <= leavetol())
796  {
797  setBasisStatus(SPxBasis::OPTIMAL);
798  m_status = OPTIMAL;
799  break;
800  }
801  loopCount++;
802  }
803  setType(ENTER);
804  init();
805  thepricer->setType(type());
806  theratiotester->setType(type());
807  }
808  }
809  assert(m_status != SINGULAR);
810 
811  }
812  catch( const SPxException& E )
813  {
814  // if we stopped due to a singular basis, we reload the original basis and try again with tighter
815  // tolerance (only once)
816  if (m_status == SINGULAR && !tightened)
817  {
818  tightenedtype = type();
819 
820  if( tightenedtype == ENTER )
821  {
822  m_entertol = 0.01 * m_entertol;
823 
824  MSG_INFO2( (*spxout), (*spxout) << " --- basis singular: reloading basis and solving with tighter ratio test tolerance " << m_entertol << std::endl; )
825  }
826  else
827  {
828  m_leavetol = 0.01 * m_leavetol;
829 
830  MSG_INFO2( (*spxout), (*spxout) << " --- basis singular: reloading basis and solving with tighter ratio test tolerance " << m_leavetol << std::endl; )
831  }
832 
833  // load original basis
834  int niters = iterations();
835  loadBasis(regulardesc);
836 
837  // remember iteration count
838  iterCount = niters;
839 
840  // try initializing basis (might fail if starting basis was already singular)
841  try
842  {
843  init();
844  theratiotester->setType(type());
845  }
846  catch( const SPxException& Ex )
847  {
848  MSG_INFO2( (*spxout), (*spxout) << " --- reloaded basis singular, resetting original tolerances" << std::endl; )
849 
850  if( tightenedtype == ENTER )
851  m_entertol = 100.0 * m_entertol;
852  else
853  m_leavetol = 100.0 * m_leavetol;
854 
855  theratiotester->setType(type());
856 
857  throw Ex;
858  }
859 
860  // reset status and counters
861  m_status = RUNNING;
862  m_numCycle = 0;
863  leaveCount = 0;
864  enterCount = 0;
865  stallNumRecovers = 0;
866 
867  // continue
868  stop = false;
869  tightened = true;
870  }
871  // reset tolerance to its original value and pass on the exception
872  else if (tightened)
873  {
874  if( tightenedtype == ENTER )
875  m_entertol = 100.0 * m_entertol;
876  else
877  m_leavetol = 100.0 * m_leavetol;
878 
879  theratiotester->setType(type());
880 
881  throw E;
882  }
883  // pass on the exception
884  else
885  throw E;
886  }
887  }
888 
889  // reset tolerance to its original value
890  if (tightened)
891  {
892  if( tightenedtype == ENTER )
893  m_entertol = 100.0 * m_entertol;
894  else
895  m_leavetol = 100.0 * m_leavetol;
896 
897  theratiotester->setType(type());
898  }
899 
900  theTime->stop();
901  theCumulativeTime += time();
902 
903  if (m_status == RUNNING)
904  {
905  m_status = ERROR;
906  throw SPxStatusException("XSOLVE05 Status is still RUNNING when it shouldn't be");
907  }
908 
909  MSG_INFO3( (*spxout),
910  (*spxout) << "Finished solving (status=" << status()
911  << ", iters=" << iterCount
912  << ", leave=" << leaveCount
913  << ", enter=" << enterCount
914  << ", flips=" << totalboundflips;
915  if( status() == OPTIMAL )
916  (*spxout) << ", objValue=" << value();
917  (*spxout) << ")" << std::endl;
918  )
919 
920 #ifdef ENABLE_ADDITIONAL_CHECKS
921  /* check if solution is really feasible */
922  if( status() == OPTIMAL )
923  {
924  int c;
925  Real val;
926  DVector sol( nCols() );
927 
928  getPrimal( sol );
929 
930  for(int row = 0; row < nRows(); ++row )
931  {
932  const SVector& rowvec = rowVector( row );
933  val = 0.0;
934  for( c = 0; c < rowvec.size(); ++c )
935  val += rowvec.value( c ) * sol[rowvec.index( c )];
936 
937  if( LT( val, lhs( row ), feastol() ) ||
938  GT( val, rhs( row ), feastol() ) )
939  {
940  // Minor rhs violations happen frequently, so print these
941  // warnings only with verbose level INFO2 and higher.
942  MSG_INFO2( (*spxout), (*spxout) << "WSOLVE88 Warning! Constraint " << row
943  << " is violated by solution" << std::endl
944  << " lhs:" << lhs( row )
945  << " <= val:" << val
946  << " <= rhs:" << rhs( row ) << std::endl; )
947 
948  if( type() == LEAVE && isRowBasic( row ) )
949  {
950  // find basis variable
951  for( c = 0; c < nRows(); ++c )
952  if (basis().baseId(c).isSPxRowId()
953  && (number(basis().baseId(c)) == row))
954  break;
955 
956  assert( c < nRows() );
957 
958  MSG_WARNING( (*spxout), (*spxout) << "WSOLVE90 basis idx:" << c
959  << " fVec:" << fVec()[c]
960  << " fRhs:" << fRhs()[c]
961  << " fTest:" << fTest()[c] << std::endl; )
962  }
963  }
964  }
965  for(int col = 0; col < nCols(); ++col )
966  {
967  if( LT( sol[col], lower( col ), feastol() ) ||
968  GT( sol[col], upper( col ), feastol() ) )
969  {
970  // Minor bound violations happen frequently, so print these
971  // warnings only with verbose level INFO2 and higher.
972  MSG_INFO2( (*spxout), (*spxout) << "WSOLVE91 Warning! Bound for column " << col
973  << " is violated by solution" << std::endl
974  << " lower:" << lower( col )
975  << " <= val:" << sol[col]
976  << " <= upper:" << upper( col ) << std::endl; )
977 
978  if( type() == LEAVE && isColBasic( col ) )
979  {
980  for( c = 0; c < nRows() ; ++c)
981  if ( basis().baseId( c ).isSPxColId()
982  && ( number( basis().baseId( c ) ) == col ))
983  break;
984 
985  assert( c < nRows() );
986  MSG_WARNING( (*spxout), (*spxout) << "WSOLVE92 basis idx:" << c
987  << " fVec:" << fVec()[c]
988  << " fRhs:" << fRhs()[c]
989  << " fTest:" << fTest()[c] << std::endl; )
990  }
991  }
992  }
993  }
994 #endif // ENABLE_ADDITIONAL_CHECKS
995 
996  primalCount = ( rep() == SPxSolver::COLUMN )
997  ? enterCount
998  : leaveCount;
999 
1000  printDisplayLine(true);
1001  return status();
1002 }
1003 
1004 void SPxSolver::testVecs()
1005 {
1006 
1007  assert(SPxBasis::status() > SPxBasis::SINGULAR);
1008 
1009  DVector tmp(dim());
1010 
1011  tmp = *theCoPvec;
1012  multWithBase(tmp);
1013  tmp -= *theCoPrhs;
1014  if (tmp.length() > leavetol())
1015  {
1016  MSG_INFO3( (*spxout), (*spxout) << "ISOLVE93 " << iteration() << ":\tcoP error = \t"
1017  << tmp.length() << std::endl; )
1018 
1019  tmp.clear();
1020  SPxBasis::coSolve(tmp, *theCoPrhs);
1021  multWithBase(tmp);
1022  tmp -= *theCoPrhs;
1023  MSG_INFO3( (*spxout), (*spxout) << "ISOLVE94\t\t" << tmp.length() << std::endl; )
1024 
1025  tmp.clear();
1026  SPxBasis::coSolve(tmp, *theCoPrhs);
1027  tmp -= *theCoPvec;
1028  MSG_INFO3( (*spxout), (*spxout) << "ISOLVE95\t\t" << tmp.length() << std::endl; )
1029  }
1030 
1031  tmp = *theFvec;
1032  multBaseWith(tmp);
1033  tmp -= *theFrhs;
1034  if (tmp.length() > entertol())
1035  {
1036  MSG_INFO3( (*spxout), (*spxout) << "ISOLVE96 " << iteration() << ":\t F error = \t"
1037  << tmp.length() << std::endl; )
1038 
1039  tmp.clear();
1040  SPxBasis::solve(tmp, *theFrhs);
1041  tmp -= *theFvec;
1042  MSG_INFO3( (*spxout), (*spxout) << "ISOLVE97\t\t" << tmp.length() << std::endl; )
1043  }
1044 
1045  if (type() == ENTER)
1046  {
1047  for (int i = 0; i < dim(); ++i)
1048  {
1049  if (theCoTest[i] < -leavetol() && isCoBasic(i))
1050  {
1051  /// @todo Error message "this shalt not be": shalt this be an assert (also below)?
1052  MSG_ERROR( std::cerr << "ESOLVE98 testVecs: theCoTest: this shalt not be!"
1053  << std::endl
1054  << " i=" << i
1055  << ", theCoTest[i]=" << theCoTest[i]
1056  << ", leavetol()=" << leavetol() << std::endl; )
1057  }
1058  }
1059 
1060  for (int i = 0; i < coDim(); ++i)
1061  {
1062  if (theTest[i] < -leavetol() && isBasic(i))
1063  {
1064  MSG_ERROR( std::cerr << "ESOLVE99 testVecs: theTest: this shalt not be!"
1065  << std::endl
1066  << " i=" << i
1067  << ", theTest[i]=" << theTest[i]
1068  << ", leavetol()=" << leavetol() << std::endl; )
1069  }
1070  }
1071  }
1072 }
1073 
1074 
1075 /// print display line of flying table
1076 void SPxSolver::printDisplayLine(const bool force)
1077 {
1078  MSG_INFO1( (*spxout),
1079  if( displayLine % (displayFreq*30) == 0 )
1080  {
1081  (*spxout) << "type | time | iters | facts | shift | value\n";
1082  }
1083  if( force || (displayLine % displayFreq == 0) )
1084  {
1085  (type() == LEAVE) ? (*spxout) << " L |" : (*spxout) << " E |";
1086  (*spxout) << std::fixed << std::setw(7) << std::setprecision(1) << time() << " |";
1087  (*spxout) << std::scientific << std::setprecision(2);
1088  (*spxout) << std::setw(8) << iteration() << " | "
1089  << std::setw(5) << slinSolver()->getFactorCount() << " | "
1090  << shift() << " | "
1091  << std::setprecision(8) << value() + objOffset()
1092  << std::endl;
1093  }
1094  displayLine++;
1095  );
1096 }
1097 
1098 
1099 bool SPxSolver::terminate()
1100 {
1101 #ifdef ENABLE_ADDITIONAL_CHECKS
1102  if (SPxBasis::status() > SPxBasis::SINGULAR)
1103  testVecs();
1104 #endif
1105 
1106  int redo = dim();
1107  if (redo < 1000)
1108  redo = 1000;
1109 
1110  if (iteration() > 10 && iteration() % redo == 0)
1111  {
1112 #ifdef ENABLE_ADDITIONAL_CHECKS
1113  DVector cr(*theCoPrhs);
1114  DVector fr(*theFrhs);
1115 #endif
1116 
1117  if (type() == ENTER)
1118  computeEnterCoPrhs();
1119  else
1120  computeLeaveCoPrhs();
1121 
1122  computeFrhs();
1123 
1124 #ifdef ENABLE_ADDITIONAL_CHECKS
1125  cr -= *theCoPrhs;
1126  fr -= *theFrhs;
1127  if (cr.length() > leavetol())
1128  MSG_WARNING( (*spxout), (*spxout) << "WSOLVE50 unexpected change of coPrhs "
1129  << cr.length() << std::endl; )
1130  if (fr.length() > entertol())
1131  MSG_WARNING( (*spxout), (*spxout) << "WSOLVE51 unexpected change of Frhs "
1132  << fr.length() << std::endl; )
1133 #endif
1134 
1135  if (updateCount > 1)
1136  {
1137  MSG_INFO3( (*spxout), (*spxout) << " --- terminate triggers refactorization"
1138  << std::endl; )
1139  factorize();
1140  }
1141  SPxBasis::coSolve(*theCoPvec, *theCoPrhs);
1142  SPxBasis::solve (*theFvec, *theFrhs);
1143 
1144  if (pricing() == FULL)
1145  {
1146  computePvec();
1147  if (type() == ENTER)
1148  computeTest();
1149  }
1150 
1151  if (shift() > 0.0)
1152  unShift();
1153  }
1154 
1155  // check time limit and objective limit only for non-terminal bases
1156  if( SPxBasis::status() >= SPxBasis::OPTIMAL ||
1157  SPxBasis::status() <= SPxBasis::SINGULAR )
1158  {
1159  m_status = UNKNOWN;
1160  return true;
1161  }
1162 
1163  if ( isTimeLimitReached() )
1164  {
1165  MSG_INFO2( (*spxout), (*spxout) << " --- timelimit (" << maxTime
1166  << ") reached" << std::endl; )
1167  m_status = ABORT_TIME;
1168  return true;
1169  }
1170 
1171  // objLimit is set and we are running DUAL:
1172  // - objLimit is set if objLimit < infinity
1173  // - DUAL is running if rep() * type() > 0 == DUAL (-1 == PRIMAL)
1174  //
1175  // In this case we have given a objective value limit, e.g, through a
1176  // MIP solver, and we want stop solving the LP if we figure out that the
1177  // optimal value of the current LP can not be better then this objective
1178  // limit. More precisely:
1179  // - MINIMIZATION Problem
1180  // We want stop the solving process if
1181  // objLimit <= current objective value of the DUAL LP
1182  // - MAXIMIZATION Problem
1183  // We want stop the solving process if
1184  // objLimit >= current objective value of the DUAL LP
1185  if (objLimit < infinity && type() * rep() > 0)
1186  {
1187  // We have no bound shifts; therefore, we can trust the current
1188  // objective value.
1189  // It might be even possible to use this termination value in case of
1190  // bound violations (shifting) but in this case it is quite difficult
1191  // to determine if we already reached the limit.
1192  if( shift() < epsilon() && maxInfeas() + shift() <= opttol() )
1193  {
1194  // SPxSense::MINIMIZE == -1, so we have sign = 1 on minimizing
1195  if( spxSense() * value() <= spxSense() * objLimit )
1196  {
1197  MSG_INFO2( (*spxout), (*spxout) << " --- objective value limit (" << objLimit
1198  << ") reached" << std::endl; )
1199  MSG_DEBUG(
1200  (*spxout) << " --- objective value limit reached" << std::endl
1201  << " (value: " << value()
1202  << ", limit: " << objLimit << ")" << std::endl
1203  << " (spxSense: " << int(spxSense())
1204  << ", rep: " << int(rep())
1205  << ", type: " << int(type()) << ")" << std::endl;
1206  )
1207 
1208  m_status = ABORT_VALUE;
1209  return true;
1210  }
1211  }
1212  }
1213 
1214  return false;
1215 }
1216 
1217 SPxSolver::Status SPxSolver::getPrimal (Vector& p_vector) const
1218 {
1219 
1220  if (!isInitialized())
1221  {
1222  /* exit if presolving/simplifier cleared the problem */
1223  if (status() == NO_PROBLEM)
1224  return status();
1225  throw SPxStatusException("XSOLVE06 Not Initialized");
1226  }
1227  if (rep() == ROW)
1228  p_vector = coPvec();
1229  else
1230  {
1231  const SPxBasis::Desc& ds = desc();
1232 
1233  for (int i = 0; i < nCols(); ++i)
1234  {
1235  switch (ds.colStatus(i))
1236  {
1237  case SPxBasis::Desc::P_ON_LOWER :
1238  p_vector[i] = SPxLP::lower(i);
1239  break;
1240  case SPxBasis::Desc::P_ON_UPPER :
1241  case SPxBasis::Desc::P_FIXED :
1242  p_vector[i] = SPxLP::upper(i);
1243  break;
1244  case SPxBasis::Desc::P_FREE :
1245  p_vector[i] = 0;
1246  break;
1247  case SPxBasis::Desc::D_FREE :
1248  case SPxBasis::Desc::D_ON_UPPER :
1249  case SPxBasis::Desc::D_ON_LOWER :
1250  case SPxBasis::Desc::D_ON_BOTH :
1251  case SPxBasis::Desc::D_UNDEFINED :
1252  break;
1253  default:
1254  throw SPxInternalCodeException("XSOLVE07 This should never happen.");
1255  }
1256  }
1257  for (int j = 0; j < dim(); ++j)
1258  {
1259  if (baseId(j).isSPxColId())
1260  p_vector[ number(SPxColId(baseId(j))) ] = fVec()[j];
1261  }
1262  }
1263  return status();
1264 }
1265 
1266 SPxSolver::Status SPxSolver::getDual (Vector& p_vector) const
1267 {
1268 
1269  assert(isInitialized());
1270 
1271  if (!isInitialized())
1272  {
1273  /* exit if presolving/simplifier cleared the problem */
1274  if (status() == NO_PROBLEM)
1275  return status();
1276  throw SPxStatusException("XSOLVE08 No Problem loaded");
1277  }
1278 
1279  if (rep() == ROW)
1280  {
1281  int i;
1282  p_vector = maxRowObj();
1283  for (i = nCols() - 1; i >= 0; --i)
1284  {
1285  if (baseId(i).isSPxRowId())
1286  p_vector[ number(SPxRowId(baseId(i))) ] = fVec()[i];
1287  }
1288  }
1289  else
1290  p_vector = coPvec();
1291 
1292  p_vector *= Real(spxSense());
1293 
1294  return status();
1295 }
1296 
1297 SPxSolver::Status SPxSolver::getRedCost (Vector& p_vector) const
1298 {
1299 
1300  assert(isInitialized());
1301 
1302  if (!isInitialized())
1303  {
1304  throw SPxStatusException("XSOLVE09 No Problem loaded");
1305  // return NOT_INIT;
1306  }
1307 
1308  if (rep() == ROW)
1309  {
1310  int i;
1311  p_vector.clear();
1312  if (spxSense() == SPxLP::MINIMIZE)
1313  {
1314  for (i = dim() - 1; i >= 0; --i)
1315  {
1316  if (baseId(i).isSPxColId())
1317  p_vector[ number(SPxColId(baseId(i))) ] = -fVec()[i];
1318  }
1319  }
1320  else
1321  {
1322  for (i = dim() - 1; i >= 0; --i)
1323  {
1324  if (baseId(i).isSPxColId())
1325  p_vector[ number(SPxColId(baseId(i))) ] = fVec()[i];
1326  }
1327  }
1328  }
1329  else
1330  {
1331  p_vector = maxObj();
1332  p_vector -= pVec();
1333  if (spxSense() == SPxLP::MINIMIZE)
1334  p_vector *= -1.0;
1335  }
1336 
1337  return status();
1338 }
1339 
1340 SPxSolver::Status SPxSolver::getPrimalray (Vector& p_vector) const
1341 {
1342 
1343  assert(isInitialized());
1344 
1345  if (!isInitialized())
1346  {
1347  throw SPxStatusException("XSOLVE10 No Problem loaded");
1348  // return NOT_INIT;
1349  }
1350 
1351  assert(SPxBasis::status() == SPxBasis::UNBOUNDED);
1352  p_vector.clear();
1353  p_vector = primalRay;
1354 
1355  return status();
1356 }
1357 
1358 SPxSolver::Status SPxSolver::getDualfarkas (Vector& p_vector) const
1359 {
1360 
1361  assert(isInitialized());
1362 
1363  if (!isInitialized())
1364  {
1365  throw SPxStatusException("XSOLVE10 No Problem loaded");
1366  // return NOT_INIT;
1367  }
1368 
1369  assert(SPxBasis::status() == SPxBasis::INFEASIBLE);
1370  p_vector.clear();
1371  p_vector = dualFarkas;
1372 
1373  return status();
1374 }
1375 
1376 SPxSolver::Status SPxSolver::getSlacks (Vector& p_vector) const
1377 {
1378 
1379  assert(isInitialized());
1380 
1381  if (!isInitialized())
1382  {
1383  throw SPxStatusException("XSOLVE11 No Problem loaded");
1384  // return NOT_INIT;
1385  }
1386 
1387  if (rep() == COLUMN)
1388  {
1389  int i;
1390  const SPxBasis::Desc& ds = desc();
1391  for (i = nRows() - 1; i >= 0; --i)
1392  {
1393  switch (ds.rowStatus(i))
1394  {
1395  case SPxBasis::Desc::P_ON_LOWER :
1396  p_vector[i] = lhs(i);
1397  break;
1398  case SPxBasis::Desc::P_ON_UPPER :
1399  case SPxBasis::Desc::P_FIXED :
1400  p_vector[i] = rhs(i);
1401  break;
1402  case SPxBasis::Desc::P_FREE :
1403  p_vector[i] = 0;
1404  break;
1405  case SPxBasis::Desc::D_FREE :
1406  case SPxBasis::Desc::D_ON_UPPER :
1407  case SPxBasis::Desc::D_ON_LOWER :
1408  case SPxBasis::Desc::D_ON_BOTH :
1409  case SPxBasis::Desc::D_UNDEFINED :
1410  break;
1411  default:
1412  throw SPxInternalCodeException("XSOLVE12 This should never happen.");
1413  }
1414  }
1415  for (i = dim() - 1; i >= 0; --i)
1416  {
1417  if (baseId(i).isSPxRowId())
1418  p_vector[ number(SPxRowId(baseId(i))) ] = -(*theFvec)[i];
1419  }
1420  }
1421  else
1422  p_vector = pVec();
1423 
1424  return status();
1425 }
1426 
1427 void SPxSolver::setPrimal(Vector& p_vector)
1428 {
1429 
1430  if (!isInitialized())
1431  {
1432  throw SPxStatusException("XSOLVE20 Not Initialized");
1433  }
1434 
1435  if (rep() == ROW)
1436  coPvec() = p_vector;
1437  else
1438  {
1439  for (int j = 0; j < dim(); ++j)
1440  {
1441  if (baseId(j).isSPxColId())
1442  fVec()[j] = p_vector[ number(SPxColId(baseId(j))) ];
1443  }
1444  }
1445 }
1446 
1447 void SPxSolver::setDual(Vector& p_vector)
1448 {
1449 
1450  assert(isInitialized());
1451 
1452  if (!isInitialized())
1453  {
1454  throw SPxStatusException("XSOLVE21 Not Initialized");
1455  }
1456 
1457  if (rep() == ROW)
1458  {
1459  for (int i = nCols() - 1; i >= 0; --i)
1460  {
1461  if (baseId(i).isSPxRowId())
1462  {
1463  if (spxSense() == SPxLP::MAXIMIZE)
1464  fVec()[i] = p_vector[ number(SPxRowId(baseId(i))) ];
1465  else
1466  fVec()[i] = -p_vector[ number(SPxRowId(baseId(i))) ];
1467  }
1468  }
1469  }
1470  else
1471  {
1472  coPvec() = p_vector;
1473  if (spxSense() == SPxLP::MINIMIZE)
1474  coPvec() *= -1.0;
1475  }
1476 }
1477 
1478 void SPxSolver::setRedCost(Vector& p_vector)
1479 {
1480 
1481  assert(isInitialized());
1482 
1483  if (!isInitialized())
1484  {
1485  throw SPxStatusException("XSOLVE22 Not Initialized");
1486  }
1487 
1488  if (rep() == ROW)
1489  {
1490  for (int i = dim() - 1; i >= 0; --i)
1491  {
1492  if (baseId(i).isSPxColId())
1493  {
1494  if (spxSense() == SPxLP::MINIMIZE)
1495  fVec()[i] = -p_vector[ number(SPxColId(baseId(i))) ];
1496  else
1497  fVec()[i] = p_vector[ number(SPxColId(baseId(i))) ];
1498  }
1499  }
1500  }
1501  else
1502  {
1503  pVec() = maxObj();
1504 
1505  if (spxSense() == SPxLP::MINIMIZE)
1506  pVec() += p_vector;
1507  else
1508  pVec() -= p_vector;
1509  }
1510 }
1511 
1512 void SPxSolver::setSlacks(Vector& p_vector)
1513 {
1514 
1515  assert(isInitialized());
1516 
1517  if (!isInitialized())
1518  {
1519  throw SPxStatusException("XSOLVE23 Not Initialized");
1520  }
1521 
1522  if (rep() == COLUMN)
1523  {
1524  for (int i = dim() - 1; i >= 0; --i)
1525  {
1526  if (baseId(i).isSPxRowId())
1527  (*theFvec)[i] = -p_vector[ number(SPxRowId(baseId(i))) ];
1528  }
1529  }
1530  else
1531  pVec() = p_vector;
1532 }
1533 
1534 SPxSolver::Status SPxSolver::status() const
1535 {
1536  switch( m_status )
1537  {
1538  case UNKNOWN :
1539  switch (SPxBasis::status())
1540  {
1541  case SPxBasis::NO_PROBLEM :
1542  return NO_PROBLEM;
1543  case SPxBasis::SINGULAR :
1544  return SINGULAR;
1545  case SPxBasis::REGULAR :
1546  case SPxBasis::DUAL :
1547  case SPxBasis::PRIMAL :
1548  return UNKNOWN;
1549  case SPxBasis::OPTIMAL :
1550  return OPTIMAL;
1551  case SPxBasis::UNBOUNDED :
1552  return UNBOUNDED;
1553  case SPxBasis::INFEASIBLE :
1554  return INFEASIBLE;
1555  default:
1556  return ERROR;
1557  }
1558  case SINGULAR :
1559  return m_status;
1560  case OPTIMAL :
1561  assert( SPxBasis::status() == SPxBasis::OPTIMAL );
1562  /*lint -fallthrough*/
1563  case ABORT_CYCLING :
1564  case ABORT_TIME :
1565  case ABORT_ITER :
1566  case ABORT_VALUE :
1567  case RUNNING :
1568  case REGULAR :
1569  case NOT_INIT :
1570  case NO_SOLVER :
1571  case NO_PRICER :
1572  case NO_RATIOTESTER :
1573  case ERROR:
1574  return m_status;
1575  default:
1576  return ERROR;
1577  }
1578 }
1579 
1580 SPxSolver::Status SPxSolver::getResult(
1581  Real* p_value,
1582  Vector* p_primal,
1583  Vector* p_slacks,
1584  Vector* p_dual,
1585  Vector* reduCosts)
1586 {
1587  if (p_value)
1588  *p_value = this->value();
1589  if (p_primal)
1590  getPrimal(*p_primal);
1591  if (p_slacks)
1592  getSlacks(*p_slacks);
1593  if (p_dual)
1594  getDual(*p_dual);
1595  if (reduCosts)
1596  getRedCost(*reduCosts);
1597  return status();
1598 }
1599 } // namespace soplex