changesoplex.cpp

Go to the documentation of this file.

/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
/*                                                                           */
/*                  This file is part of the class library                   */
/*       SoPlex --- the Sequential object-oriented simPlex.                  */
/*                                                                           */
/*    Copyright (C) 1996-2019 Konrad-Zuse-Zentrum                            */
/*                            fuer Informationstechnik Berlin                */
/*                                                                           */
/*  SoPlex is distributed under the terms of the ZIB Academic Licence.       */
/*                                                                           */
/*  You should have received a copy of the ZIB Academic License              */
/*  along with SoPlex; see the file COPYING. If not email to soplex@zib.de.  */
/*                                                                           */
/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */

#include <assert.h>
#include <iostream>

#include "soplex/spxdefines.h"
#include "soplex/spxsolver.h"
#include "soplex/spxpricer.h"
#include "soplex/spxratiotester.h"
#include "soplex/exceptions.h"

namespace soplex
{

void SPxSolver::addedRows(int n)
{

   if(n > 0)
   {
      SPxLP::addedRows(n);

      unInit();
      reDim();

      if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
         SPxBasis::addedRows(n);
   }

   /* we must not assert consistency here, since addedCols() might be still necessary to obtain a consistent basis */
}

void SPxSolver::addedCols(int n)
{

   if(n > 0)
   {
      SPxLP::addedCols(n);

      unInit();
      reDim();

      if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
         SPxBasis::addedCols(n);
   }

   /* we must not assert consistency here, since addedRows() might be still necessary to obtain a consistent basis */
}

void SPxSolver::doRemoveRow(int i)
{

   SPxLP::doRemoveRow(i);

   unInit();

   if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
   {
      removedRow(i);

      switch(SPxBasis::status())
      {
      case SPxBasis::DUAL:
      case SPxBasis::INFEASIBLE:
         setBasisStatus(SPxBasis::REGULAR);
         break;

      case SPxBasis::OPTIMAL:
         setBasisStatus(SPxBasis::PRIMAL);
         break;

      default:
         break;
      }
   }
}

void SPxSolver::doRemoveRows(int perm[])
{

   SPxLP::doRemoveRows(perm);

   unInit();

   if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
   {
      removedRows(perm);

      switch(SPxBasis::status())
      {
      case SPxBasis::DUAL:
      case SPxBasis::INFEASIBLE:
         setBasisStatus(SPxBasis::REGULAR);
         break;

      case SPxBasis::OPTIMAL:
         setBasisStatus(SPxBasis::PRIMAL);
         break;

      default:
         break;
      }
   }
}

void SPxSolver::doRemoveCol(int i)
{
   forceRecompNonbasicValue();

   SPxLP::doRemoveCol(i);

   unInit();

   if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
   {
      removedCol(i);

      switch(SPxBasis::status())
      {
      case SPxBasis::PRIMAL:
      case SPxBasis::UNBOUNDED:
         setBasisStatus(SPxBasis::REGULAR);
         break;

      case SPxBasis::OPTIMAL:
         setBasisStatus(SPxBasis::DUAL);
         break;

      default:
         break;
      }
   }
}

void SPxSolver::doRemoveCols(int perm[])
{
   forceRecompNonbasicValue();

   SPxLP::doRemoveCols(perm);

   unInit();

   if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
   {
      removedCols(perm);

      switch(SPxBasis::status())
      {
      case SPxBasis::PRIMAL:
      case SPxBasis::UNBOUNDED:
         setBasisStatus(SPxBasis::REGULAR);
         break;

      case SPxBasis::OPTIMAL:
         setBasisStatus(SPxBasis::DUAL);
         break;

      default:
         break;
      }
   }
}

void SPxSolver::changeObj(const Vector& newObj, bool scale)
{
   forceRecompNonbasicValue();

   SPxLP::changeObj(newObj, scale);

   /**@todo Factorization remains valid, we do not need a reDim()
    *       pricing vectors should be recomputed.
    */
   unInit();
}

void SPxSolver::changeObj(int i, const Real& newVal, bool scale)
{
   forceRecompNonbasicValue();

   SPxLP::changeObj(i, newVal, scale);


   /**@todo Factorization remains valid, we do not need a reDim()
    *       pricing vectors should be recomputed.
    */
   unInit();
}

void SPxSolver::changeMaxObj(const Vector& newObj, bool scale)
{
   forceRecompNonbasicValue();

   SPxLP::changeMaxObj(newObj, scale);

   /**@todo Factorization remains valid, we do not need a reDim()
    *       pricing vectors should be recomputed.
    */
   unInit();
}

void SPxSolver::changeMaxObj(int i, const Real& newVal, bool scale)
{
   forceRecompNonbasicValue();

   SPxLP::changeMaxObj(i, newVal, scale);

   /**@todo Factorization remains valid, we do not need a reDim()
    *       pricing vectors should be recomputed.
    */
   unInit();
}

void SPxSolver::changeRowObj(const Vector& newObj, bool scale)
{
   forceRecompNonbasicValue();

   SPxLP::changeRowObj(newObj, scale);

   /**@todo Factorization remains valid, we do not need a reDim()
    *       pricing vectors should be recomputed.
    */
   unInit();
}

void SPxSolver::changeRowObj(int i, const Real& newVal, bool scale)
{
   forceRecompNonbasicValue();

   SPxLP::changeRowObj(i, newVal, scale);

   /**@todo Factorization remains valid, we do not need a reDim()
    *       pricing vectors should be recomputed.
    */
   unInit();
}

void SPxSolver::changeLowerStatus(int i, Real newLower, Real oldLower)
{
   SPxBasis::Desc::Status& stat      = desc().colStatus(i);
   Real                    currUpper = upper(i);
   Real                    objChange = 0.0;

   MSG_DEBUG(std::cout << "DCHANG01 changeLowerStatus(): col " << i
             << "[" << newLower << ":" << currUpper << "] " << stat;)

   switch(stat)
   {
   case SPxBasis::Desc::P_ON_LOWER:
      if(newLower <= -infinity)
      {
         if(currUpper >= infinity)
         {
            stat = SPxBasis::Desc::P_FREE;

            if(m_nonbasicValueUpToDate && rep() == COLUMN)
               objChange = -theLCbound[i] * oldLower;
         }
         else
         {
            stat = SPxBasis::Desc::P_ON_UPPER;

            if(m_nonbasicValueUpToDate && rep() == COLUMN)
               objChange = (theUCbound[i] * currUpper) - (theLCbound[i] * oldLower);
         }
      }
      else if(EQ(newLower, currUpper))
      {
         stat = SPxBasis::Desc::P_FIXED;

         if(m_nonbasicValueUpToDate && rep() == COLUMN)
            objChange = maxObj(i) * (newLower - oldLower);
      }
      else if(m_nonbasicValueUpToDate && rep() == COLUMN)
         objChange = theLCbound[i] * (newLower - oldLower);

      break;

   case SPxBasis::Desc::P_ON_UPPER:
      if(EQ(newLower, currUpper))
         stat = SPxBasis::Desc::P_FIXED;

      break;

   case SPxBasis::Desc::P_FREE:
      if(newLower > -infinity)
      {
         stat = SPxBasis::Desc::P_ON_LOWER;

         if(m_nonbasicValueUpToDate && rep() == COLUMN)
            objChange = theLCbound[i] * newLower;
      }

      break;

   case SPxBasis::Desc::P_FIXED:
      if(NE(newLower, currUpper))
      {
         stat = SPxBasis::Desc::P_ON_UPPER;

         if(isInitialized())
            theUCbound[i] = maxObj(i);
      }

      break;

   case SPxBasis::Desc::D_FREE:
   case SPxBasis::Desc::D_ON_UPPER:
   case SPxBasis::Desc::D_ON_LOWER:
   case SPxBasis::Desc::D_ON_BOTH:
   case SPxBasis::Desc::D_UNDEFINED:
      if(rep() == ROW && theShift > 0.0)
         forceRecompNonbasicValue();

      stat = dualColStatus(i);
      break;

   default:
      throw SPxInternalCodeException("XCHANG01 This should never happen.");
   }

   MSG_DEBUG(std::cout << " -> " << stat << std::endl;)

   // we only need to update the nonbasic value in column representation (see nonbasicValue() for comparison/explanation)
   if(rep() == COLUMN)
      updateNonbasicValue(objChange);
}

void SPxSolver::changeLower(const Vector& newLower, bool scale)
{
   // we better recompute the nonbasic value when changing all lower bounds
   forceRecompNonbasicValue();

   SPxLP::changeLower(newLower, scale);

   if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
   {
      for(int i = 0; i < newLower.dim(); ++i)
         changeLowerStatus(i, lower(i));

      unInit();
   }
}

void SPxSolver::changeLower(int i, const Real& newLower, bool scale)
{
   if(newLower != (scale ? lowerUnscaled(i) : lower(i)))
   {
      Real oldLower = lower(i);
      // This has to be done before calling changeLowerStatus() because that is calling
      // basis.dualColStatus() which calls lower() and needs the changed value.
      SPxLP::changeLower(i, newLower, scale);

      if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
      {
         changeLowerStatus(i, lower(i), oldLower);
         unInit();
      }
   }
}

void SPxSolver::changeUpperStatus(int i, Real newUpper, Real oldUpper)
{
   SPxBasis::Desc::Status& stat      = desc().colStatus(i);
   Real                    currLower = lower(i);
   Real                    objChange = 0.0;

   MSG_DEBUG(std::cout << "DCHANG02 changeUpperStatus(): col " << i
             << "[" << currLower << ":" << newUpper << "] " << stat;)

   switch(stat)
   {
   case SPxBasis::Desc::P_ON_LOWER:
      if(newUpper == currLower)
         stat = SPxBasis::Desc::P_FIXED;

      break;

   case SPxBasis::Desc::P_ON_UPPER:
      if(newUpper >= infinity)
      {
         if(currLower <= -infinity)
         {
            stat = SPxBasis::Desc::P_FREE;

            if(m_nonbasicValueUpToDate && rep() == COLUMN)
               objChange = -theUCbound[i] * oldUpper;
         }
         else
         {
            stat = SPxBasis::Desc::P_ON_LOWER;

            if(m_nonbasicValueUpToDate && rep() == COLUMN)
               objChange = (theLCbound[i] * currLower) - (theUCbound[i] * oldUpper);
         }
      }
      else if(EQ(newUpper, currLower))
      {
         stat = SPxBasis::Desc::P_FIXED;

         if(m_nonbasicValueUpToDate && rep() == COLUMN)
            objChange = maxObj(i) * (newUpper - oldUpper);
      }
      else if(m_nonbasicValueUpToDate && rep() == COLUMN)
         objChange = theUCbound[i] * (newUpper - oldUpper);

      break;

   case SPxBasis::Desc::P_FREE:
      if(newUpper < infinity)
      {
         stat = SPxBasis::Desc::P_ON_UPPER;

         if(m_nonbasicValueUpToDate && rep() == COLUMN)
            objChange = theUCbound[i] * newUpper;
      }

      break;

   case SPxBasis::Desc::P_FIXED:
      if(NE(newUpper, currLower))
      {
         stat = SPxBasis::Desc::P_ON_LOWER;

         if(isInitialized())
            theLCbound[i] = maxObj(i);
      }

      break;

   case SPxBasis::Desc::D_FREE:
   case SPxBasis::Desc::D_ON_UPPER:
   case SPxBasis::Desc::D_ON_LOWER:
   case SPxBasis::Desc::D_ON_BOTH:
   case SPxBasis::Desc::D_UNDEFINED:
      if(rep() == ROW && theShift > 0.0)
         forceRecompNonbasicValue();

      stat = dualColStatus(i);
      break;

   default:
      throw SPxInternalCodeException("XCHANG02 This should never happen.");
   }

   MSG_DEBUG(std::cout << " -> " << stat << std::endl;);

   // we only need to update the nonbasic value in column representation (see nonbasicValue() for comparison/explanation)
   if(rep() == COLUMN)
      updateNonbasicValue(objChange);
}

void SPxSolver::changeUpper(const Vector& newUpper, bool scale)
{
   // we better recompute the nonbasic value when changing all upper bounds
   forceRecompNonbasicValue();

   SPxLP::changeUpper(newUpper, scale);

   if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
   {
      for(int i = 0; i < newUpper.dim(); ++i)
         changeUpperStatus(i, upper(i));

      unInit();
   }
}

void SPxSolver::changeUpper(int i, const Real& newUpper, bool scale)
{
   if(newUpper != (scale ? upperUnscaled(i) : upper(i)))
   {
      Real oldUpper = upper(i);
      SPxLP::changeUpper(i, newUpper, scale);

      if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
      {
         changeUpperStatus(i, upper(i), oldUpper);
         unInit();
      }
   }
}

void SPxSolver::changeBounds(const Vector& newLower, const Vector& newUpper, bool scale)
{
   changeLower(newLower, scale);
   changeUpper(newUpper, scale);
}

void SPxSolver::changeBounds(int i, const Real& newLower, const Real& newUpper, bool scale)
{
   changeLower(i, newLower, scale);
   changeUpper(i, newUpper, scale);
}

void SPxSolver::changeLhsStatus(int i, Real newLhs, Real oldLhs)
{
   SPxBasis::Desc::Status& stat      = desc().rowStatus(i);
   Real                    currRhs   = rhs(i);
   Real                    objChange = 0.0;

   MSG_DEBUG(std::cout << "DCHANG03 changeLhsStatus()  : row " << i
             << ": " << stat;)

   switch(stat)
   {
   case SPxBasis::Desc::P_ON_LOWER:
      if(newLhs <= -infinity)
      {
         if(currRhs >= infinity)
         {
            stat = SPxBasis::Desc::P_FREE;

            if(m_nonbasicValueUpToDate && rep() == COLUMN)
               objChange = -theURbound[i] * oldLhs;
         }
         else
         {
            stat = SPxBasis::Desc::P_ON_UPPER;

            if(m_nonbasicValueUpToDate && rep() == COLUMN)
               objChange = (theLRbound[i] * currRhs) - (theURbound[i] * oldLhs);
         }
      }
      else if(EQ(newLhs, currRhs))
      {
         stat = SPxBasis::Desc::P_FIXED;

         if(m_nonbasicValueUpToDate && rep() == COLUMN)
            objChange = maxRowObj(i) * (newLhs - oldLhs);
      }
      else if(m_nonbasicValueUpToDate && rep() == COLUMN)
         objChange = theURbound[i] * (newLhs - oldLhs);

      break;

   case SPxBasis::Desc::P_ON_UPPER:
      if(EQ(newLhs, currRhs))
         stat = SPxBasis::Desc::P_FIXED;

      break;

   case SPxBasis::Desc::P_FREE:
      if(newLhs > -infinity)
      {
         stat = SPxBasis::Desc::P_ON_LOWER;

         if(m_nonbasicValueUpToDate && rep() == COLUMN)
            objChange = theURbound[i] * newLhs;
      }

      break;

   case SPxBasis::Desc::P_FIXED:
      if(NE(newLhs, currRhs))
      {
         stat = SPxBasis::Desc::P_ON_UPPER;

         if(isInitialized())
            theLRbound[i] = maxRowObj(i);
      }

      break;

   case SPxBasis::Desc::D_FREE:
   case SPxBasis::Desc::D_ON_UPPER:
   case SPxBasis::Desc::D_ON_LOWER:
   case SPxBasis::Desc::D_ON_BOTH:
   case SPxBasis::Desc::D_UNDEFINED:
      if(rep() == ROW && theShift > 0.0)
         forceRecompNonbasicValue();

      stat = dualRowStatus(i);
      break;

   default:
      throw SPxInternalCodeException("XCHANG03 This should never happen.");
   }

   MSG_DEBUG(std::cout << " -> " << stat << std::endl;)

   // we only need to update the nonbasic value in column representation (see nonbasicValue() for comparison/explanation)
   if(rep() == COLUMN)
      updateNonbasicValue(objChange);
}

void SPxSolver::changeLhs(const Vector& newLhs, bool scale)
{
   // we better recompute the nonbasic value when changing all lhs
   forceRecompNonbasicValue();

   SPxLP::changeLhs(newLhs, scale);

   if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
   {
      for(int i = 0; i < nRows(); ++i)
         changeLhsStatus(i, lhs(i));

      unInit();
   }
}

void SPxSolver::changeLhs(int i, const Real& newLhs, bool scale)
{
   if(newLhs != (scale ? lhsUnscaled(i) : lhs(i)))
   {
      Real oldLhs = lhs(i);
      SPxLP::changeLhs(i, newLhs, scale);

      if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
      {
         changeLhsStatus(i, lhs(i), oldLhs);
         unInit();
      }
   }
}

void SPxSolver::changeRhsStatus(int i, Real newRhs, Real oldRhs)
{
   SPxBasis::Desc::Status& stat      = desc().rowStatus(i);
   Real                    currLhs   = lhs(i);
   Real                    objChange = 0.0;

   MSG_DEBUG(std::cout << "DCHANG04 changeRhsStatus()  : row " << i
             << ": " << stat;)

   switch(stat)
   {
   case SPxBasis::Desc::P_ON_UPPER:
      if(newRhs >= infinity)
      {
         if(currLhs <= -infinity)
         {
            stat = SPxBasis::Desc::P_FREE;

            if(m_nonbasicValueUpToDate && rep() == COLUMN)
               objChange = -theLRbound[i] * oldRhs;
         }
         else
         {
            stat = SPxBasis::Desc::P_ON_LOWER;

            if(m_nonbasicValueUpToDate && rep() == COLUMN)
               objChange = (theURbound[i] * currLhs) - (theLRbound[i] * oldRhs);
         }
      }
      else if(EQ(newRhs, currLhs))
      {
         stat = SPxBasis::Desc::P_FIXED;

         if(m_nonbasicValueUpToDate && rep() == COLUMN)
            objChange = maxRowObj(i) * (newRhs - oldRhs);
      }
      else if(m_nonbasicValueUpToDate && rep() == COLUMN)
         objChange = theLRbound[i] * (newRhs - oldRhs);

      break;

   case SPxBasis::Desc::P_ON_LOWER:
      if(EQ(newRhs, currLhs))
         stat = SPxBasis::Desc::P_FIXED;

      break;

   case SPxBasis::Desc::P_FREE:
      if(newRhs < infinity)
      {
         stat = SPxBasis::Desc::P_ON_UPPER;

         if(m_nonbasicValueUpToDate && rep() == COLUMN)
            objChange = theLRbound[i] * newRhs;
      }

      break;

   case SPxBasis::Desc::P_FIXED:
      if(NE(newRhs, currLhs))
      {
         stat = SPxBasis::Desc::P_ON_LOWER;

         if(isInitialized())
            theURbound[i] = maxRowObj(i);
      }

      break;

   case SPxBasis::Desc::D_FREE:
   case SPxBasis::Desc::D_ON_UPPER:
   case SPxBasis::Desc::D_ON_LOWER:
   case SPxBasis::Desc::D_ON_BOTH:
   case SPxBasis::Desc::D_UNDEFINED:
      if(rep() == ROW && theShift > 0.0)
         forceRecompNonbasicValue();

      stat = dualRowStatus(i);
      break;

   default:
      throw SPxInternalCodeException("XCHANG04 This should never happen.");
   }

   MSG_DEBUG(std::cout << " -> " << stat << std::endl;)

   // we only need to update the nonbasic value in column representation (see nonbasicValue() for comparison/explanation)
   if(rep() == COLUMN)
      updateNonbasicValue(objChange);
}


void SPxSolver::changeRhs(const Vector& newRhs, bool scale)
{
   // we better recompute the nonbasic value when changing all rhs
   forceRecompNonbasicValue();

   SPxLP::changeRhs(newRhs, scale);

   if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
   {
      for(int i = 0; i < nRows(); ++i)
         changeRhsStatus(i, rhs(i));

      unInit();
   }
}

void SPxSolver::changeRhs(int i, const Real& newRhs, bool scale)
{
   if(newRhs != (scale ? rhsUnscaled(i) : rhs(i)))
   {
      Real oldRhs = rhs(i);
      SPxLP::changeRhs(i, newRhs, scale);

      if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
      {
         changeRhsStatus(i, rhs(i), oldRhs);
         unInit();
      }
   }
}

void SPxSolver::changeRange(const Vector& newLhs, const Vector& newRhs, bool scale)
{
   // we better recompute the nonbasic value when changing all ranges
   forceRecompNonbasicValue();

   SPxLP::changeLhs(newLhs, scale);
   SPxLP::changeRhs(newRhs, scale);

   if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
   {
      for(int i = nRows() - 1; i >= 0; --i)
      {
         changeLhsStatus(i, lhs(i));
         changeRhsStatus(i, rhs(i));
      }

      unInit();
   }
}

void SPxSolver::changeRange(int i, const Real& newLhs, const Real& newRhs, bool scale)
{
   Real oldLhs = lhs(i);
   Real oldRhs = rhs(i);

   SPxLP::changeLhs(i, newLhs, scale);
   SPxLP::changeRhs(i, newRhs, scale);

   if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
   {
      changeLhsStatus(i, lhs(i), oldLhs);
      changeRhsStatus(i, rhs(i), oldRhs);
      unInit();
   }
}

void SPxSolver::changeRow(int i, const LPRow& newRow, bool scale)
{
   forceRecompNonbasicValue();

   SPxLP::changeRow(i, newRow, scale);

   if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
      SPxBasis::changedRow(i);

   unInit();
}

void SPxSolver::changeCol(int i, const LPCol& newCol, bool scale)
{
   if(i < 0)
      return;

   forceRecompNonbasicValue();

   SPxLP::changeCol(i, newCol, scale);

   if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
      SPxBasis::changedCol(i);

   unInit();
}

void SPxSolver::changeElement(int i, int j, const Real& val, bool scale)
{
   if(i < 0 || j < 0)
      return;

   forceRecompNonbasicValue();

   SPxLP::changeElement(i, j, val, scale);

   if(SPxBasis::status() > SPxBasis::NO_PROBLEM)
      SPxBasis::changedElement(i, j);

   unInit();
}

void SPxSolver::changeSense(SPxSense sns)
{

   SPxLP::changeSense(sns);
   unInit();
}
} // namespace soplex