xzdcx.h

这是一个求解线性规划问题的程序

/*
  这个文件是本程序的算法最核心的东西,
  定义了一个修正单纯型表类TxzDcx(修正单纯型拼音的首字母)
  和一个单纯型表算法通道TxzDcxChannel
*/
//---------------------------------------------------------------------------
#ifndef xzDcxH
#define xzDcxH

#include "Matrix.h"
#include "Vector.h"
#include "TList.h"
#include "global.h"
#include "TDcxQ.h"
//---------------------------------------------------------------------------
#define TxzDE_OUTOFMEMMORY      0
#define TxzDE_INVALIDOPERATION  1
#define TxzDE_INVALIDPARAM      2
#define TxzDE_OUTOFBOUND        3

//修正单纯型表异常类
class TxzDException
{
public:
    unsigned int mErrorCode;
public:
    TxzDException(){
    }
    TxzDException(unsigned int aE){
        mErrorCode=aE;
    }
    operator unsigned int(void){    //重载(int)强制转换:如 if(e) dosomething;
        return mErrorCode;
    }
};



//修正单纯型表节点数据
class TxzDcx
{
private:
   TCVector *mCVetBx;   //当前基解的正分量
   TMatrix *mMatU;      //当前矩阵的逆阵
   TCVector *mCVety;    //平衡方程的解
   TMVALUE mValz0;      //当前最优值~cx

   KTList<unsigned int> *mLB;    //当前基
public:
   //拷贝构造一个节点数据
   TxzDcx(TxzDcx &axzDcx)
   {
      int i,j;
      //当前基
      mLB = NULL;
      mLB = new KTList<unsigned int>;
      if (!mLB)
         throw TxzDException(TxzDE_OUTOFMEMMORY);
      for (axzDcx.mLB->First();!axzDcx.mLB->IsEof();axzDcx.mLB->Next()){
         mLB->Append(axzDcx.mLB->GetData());
      }

      //当前基矩阵的逆阵
      mMatU = NULL;
      mMatU = new TMatrix(axzDcx.mMatU);
      if (!mMatU)
         throw TxzDException(TxzDE_OUTOFMEMMORY);

      //当前基解的正分量
      mCVetBx = NULL;
      mCVetBx = new TCVector(axzDcx.mCVetBx);
      if (!mCVetBx)
         throw TxzDException(TxzDE_OUTOFMEMMORY);

      //平衡方程的解(~y = ~Bc * U)
      mCVety = NULL;
      mCVety = new TCVector(axzDcx.mCVety);
      if (!mCVety)
         throw TxzDException(TxzDE_OUTOFMEMMORY);

      //当前最优值 (z(0) = Bx * ~Bc)
      mValz0 = axzDcx.mValz0;
   }

   //由原线性规划问题构造修正单纯型表
   TxzDcx(TDcxQuestion *aDcxQ)
   {
      unsigned int i,j;

      unsigned int lRowCount = aDcxQ->mMatA->GetRowCount();
      unsigned int lColCount = aDcxQ->mMatA->GetColCount();

      //当前基
      mLB = NULL;
      mLB = new KTList<unsigned int>;
      if (!mLB)
         throw TxzDException(TxzDE_OUTOFMEMMORY);
      for (i=1;i<=lColCount;i++){
         if ((*aDcxQ->mCVetx)(i) > 0)
             mLB->Append(i);
      }

      //构造基矩阵
      TMatrix lMatM(lRowCount,lRowCount);
      unsigned int k;
      for (k=1,mLB->First();!mLB->IsEof();mLB->Next(),k++) {
         j=mLB->GetData();
         for (i=1;i<=lRowCount;i++)
            lMatM(i,k) = (*aDcxQ->mMatA)(i,j);
      }

      //当前基矩阵的逆阵
      mMatU = NULL;
      mMatU = new TMatrix(lRowCount,lRowCount);
      if (!mMatU)
         throw TxzDException(TxzDE_OUTOFMEMMORY);
      *mMatU = lMatM^-1;

      //当前基解的正分量
      mCVetBx = NULL;
      mCVetBx = new TCVector(lRowCount);
      if (!mCVetBx)
         throw TxzDException(TxzDE_OUTOFMEMMORY);
      for (i=1,mLB->First();!mLB->IsEof();mLB->Next(),i++) {
         (*mCVetBx)(i) = (*aDcxQ->mCVetx)(mLB->GetData());
      }

      //当前基费用向量
      TCVector *lCVetBc;   //当前基费用向量
      lCVetBc = NULL;
      lCVetBc = new TCVector(lRowCount);
      if (!lCVetBc)
         throw TxzDException(TxzDE_OUTOFMEMMORY);
      for (i=1,mLB->First();!mLB->IsEof();mLB->Next(),i++) {
         (*lCVetBc)(i) = (*aDcxQ->mCVetc)(mLB->GetData());
      }

      //平衡方程的解(~y = ~Bc * U)
      mCVety = NULL;
      mCVety = new TCVector(lRowCount);
      if (!mCVety)
         throw TxzDException(TxzDE_OUTOFMEMMORY);
      *mCVety = (~*mMatU) * (*(TMatrix*)lCVetBc);

      //当前最优值 (z(0) = ~Bx * Bc)
      mValz0 = (~*mCVetBx * (*((TMatrix*)lCVetBc)))(1,1);
   }


   ~TxzDcx()
   {
      delete mCVetBx;   //当前基解的正分量
      delete mMatU;     //当前矩阵的逆阵
      delete mCVety;    //平衡方程的解

      delete mLB;      //当前基
   }

   //传进当前问题,计算取得该节点的非基判别行矩阵MatJudge及枢列坐标as、枢列向量VetPivot及枢行坐标位置,
   //返回:
   // 1  - 所有判别行z(j)-c(j)都<=0(即当前已是最优解)
   // -1 - 所有枢列VetPivot(j)都<=0(即无最优解，目标函数值趋向于正无穷)
   // 0  - 正常
   int GetPivot(TDcxQuestion* aDcxQ,TMatrix* aMatJudge,unsigned int* as,TCVector* aCVetPivot,unsigned int *ar)
   {
      unsigned int s,i;
      unsigned int lRowCount,lColCount;
      lRowCount = aDcxQ->mMatA->GetRowCount();
      lColCount = aDcxQ->mMatA->GetColCount();

      //计算非基变量的判别行:z(s)-c(s) = ~y*a(s) - c(s)
     // delete aMatJudge;
     // aMatJudge = new TMatrix(lColCount - lRowCount,2);
      TCVector lCVet(lRowCount);
      for (s=1,i=1;s<=lColCount;s++) {
         if (!mLB->Locate(s)) {
            (*aMatJudge)(i,1) = s;   //下标
            lCVet = aDcxQ->mMatA->CVector(s);
            (*aMatJudge)(i,2) = ((~*mCVety)*(*(TMatrix*)&lCVet))(1,1) - (*aDcxQ->mCVetc)(s); //值
       //?  (*aMatJudge)(i,2) = ((~*mCVety)*(aDcxQ->mMatA->CVector(s)))(1,1) - (*aDcxQ->mCVetc)(s); //值
            i++;
         }
      }

      //找非基变量的判别行的z(j)-c(j)的最大元素lMaxJudge及其位置s即枢列
      TMVALUE lMaxJudge = (*aMatJudge)(1,2);
      *as = 1;//StrToInt((AnsiString)((*aMatJudge)(1,1)));
      for (i=2;i<=lColCount-lRowCount;i++)
         if (lMaxJudge < (*aMatJudge)(i,2)) {
            lMaxJudge = (*aMatJudge)(i,2);
            *as = i;//StrToInt((AnsiString)((*aMatJudge)(i,1)));
         }

      //所有判别行z(j)-c(j)都<=0(即当前已是最优解)
      if (lMaxJudge <= 0)
         return 1;

      //计算枢列各元素(t(i,s) = ~U(i) * A(s))   <= A(s) = M * t(s)
      lCVet = aDcxQ->mMatA->CVector(StrToInt((AnsiString)((*aMatJudge)(*as,1))));
      for (i=1;i<=lRowCount;i++) {
         (*aCVetPivot)(i) = (mMatU->RVector(i) * *(TMatrix*)&lCVet)(1,1);
      }

      //确定枢行r
      TMVALUE lMinPivot;
      TMVALUE lTemp;
      bool lIsFirst=true;
      *ar = 0;
      for (i=1;i<=lRowCount;i++) {
         lTemp = (*aCVetPivot)(i);
         if (lTemp > 0) {
            lTemp = (*mCVetBx)(i) / lTemp;
            if (lIsFirst) {
               lMinPivot = lTemp;
               lIsFirst = false;
               *ar = i;
            }
            else if (lTemp < lMinPivot) {
               lMinPivot = lTemp;
               *ar = i;
            }
         }
      }

      //如果所有枢列VetPivot(j)都<=0(即无最优解，目标函数值趋向于正无穷)
      if (*ar == 0)
         return -1;

      //还无法确定能否达到最优
      return 0;
   }

   //根据枢列向量VetPivot和最大的z(s)-c(s)变换基
   void ChangeBase(TCVector* aCVetPivot,unsigned int ar,TMVALUE aMaxJudge,unsigned int as)
   {
      unsigned int i,j;
      unsigned int lRowCount;
      lRowCount = aCVetPivot->GetDimCount();

      //新i行 = (旧i行) -  (t(i)/t(r)) * (旧r行)
      for (i=1;i<=lRowCount;i++) {
         if (i!=ar) {
            //矩阵U
            mMatU->RAddition(i,ar,-(*aCVetPivot)(i)/(*aCVetPivot)(ar));
            //基解的正分量
            (*mCVetBx)(i) = (*mCVetBx)(i) - ((*aCVetPivot)(i)/(*aCVetPivot)(ar)) * (*mCVetBx)(ar);
         }
      }
      //平衡方程的解y (只有1行)
      for (j=1;j<=lRowCount;j++){
        (*mCVety)(j) = (*mCVety)(j) - (aMaxJudge/(*aCVetPivot)(ar)) * (*mMatU)(ar,j);
      }
      //当前最优值z(0)(只有1行,1列)
      mValz0 = mValz0 - (aMaxJudge/(*aCVetPivot)(ar)) * (*mCVetBx)(ar);

      //新r行 = (旧r行) / t(r)
      //矩阵U
      mMatU->RScale(ar,1/(*aCVetPivot)(ar));
      //基解的正分量
      (*mCVetBx)(ar) = (*mCVetBx)(ar) / (*aCVetPivot)(ar);

      //换基
      //B的第r个元素离基, s入基
      for (i=1,mLB->First();!mLB->IsEof() && i<=lRowCount;mLB->Next(),i++) {
         if (i == ar) {
            mLB->SetData(as);
            break;
         }
      }
   }

   //以整一个单存型表返回
   TMatrix GetMatxzDcx()
   {
      unsigned int lDimCount = mMatU->GetRowCount()+1;
      TMatrix lMatR(lDimCount,lDimCount);

      unsigned int i,j;
      //U
      for (i=1;i<=lDimCount-1;i++)
         for (j=2;j<=lDimCount;j++)
            lMatR(i,j) = (*mMatU)(i,j-1);

      //y-最后一行
      for (j=2;j<=lDimCount;j++)
         lMatR(lDimCount,j) = (*mCVety)(j-1);

      //Bx-第一列
      for (i=1;i<=lDimCount-1;i++)
         lMatR(i,1) = (*mCVetBx)(i);

      //z0-左下角
      lMatR(lDimCount,1) = mValz0;

      return lMatR;
   }


   TMatrix* GetMatU() { return mMatU; }
   TCVector* GetCVety() { return mCVety; }
   TCVector* GetCVetBx() { return mCVetBx; }
   TMVALUE GetValz0() { return mValz0; }
   KTList<unsigned int>* GetLB() { return mLB; }
};
//---------------------------------------------------------------------------

//单纯型表算法通道
class TxzDcxChannel
{
private:
   KTList<TxzDcx*> *mLxzDcx;  //解路径
   TxzDcx *mxzDcx;            //当前数据状态

   TDcxQuestion *mDcxQ;       //原单纯型问题

   //当前节点可得到的状态,中间变量
   TMatrix *mMatJudge;     //非基判别行矩阵MatJudge
   TCVector *mCVetPivot;   //枢列向量
   unsigned int *ms,*mr;   //枢列坐标s,枢行坐标r
   int mResult;            //计算结果 1: 当前已是最优解, -1:当前可判断无最优解，0: 正常

private:
   bool mIsValid;            //当前通道是否有效
   bool mIsRecord;         //是否记录中间状态

public:

   TxzDcxChannel()
   {
      mIsValid = false;
      mIsRecord = true;  //默认记录中间状态
   }

   ~TxzDcxChannel()
   {
      if (mIsValid) {  //如果通道有效(即有数据),先删除数据
         for (mLxzDcx->First();!mLxzDcx->IsEof();mLxzDcx->Next())
            delete mLxzDcx->GetData();    //节点数据

         delete mLxzDcx;    //解路径
         delete mDcxQ;      //原单纯型问题
         delete mMatJudge;  //非基判别行矩阵MatJudge
         delete mCVetPivot; //枢列向量
         delete ms;         //枢列坐标s
         delete mr;         //枢行坐标r
      }
   }

   //设置工作状态
   void SetStates(bool aIsRecord)
   {
      mIsRecord = aIsRecord;
      if (mIsValid) {  //如果通道有效(即有数据),先删除数据
         for (mLxzDcx->First();!mLxzDcx->IsEof();mLxzDcx->Next())
            delete mLxzDcx->GetData();    //节点数据

         delete mLxzDcx;    //解路径
         delete mDcxQ;      //原单纯型问题
         delete mMatJudge;  //非基判别行矩阵MatJudge
         delete mCVetPivot; //枢列向量
         delete ms;         //枢列坐标s
         delete mr;         //枢行坐标r
      }
      mIsValid = false;
   }

   //提交数据,立即执行
   void operator<<(TDcxQuestion* aDcxQ)
   {
      if (mIsValid) {  //如果通道有效(即有数据),先删除数据
         for (mLxzDcx->First();!mLxzDcx->IsEof();mLxzDcx->Next())
            delete mLxzDcx->GetData();    //节点数据

         delete mLxzDcx;    //解路径
         delete mDcxQ;      //原单纯型问题
         delete mMatJudge;  //非基判别行矩阵MatJudge
         delete mCVetPivot; //枢列向量
         delete ms;         //枢列坐标s
         delete mr;         //枢行坐标r
      }

      //拷贝原单纯型问题
      mDcxQ = new TDcxQuestion(*aDcxQ);

      //由原单纯型问题构造一个修正单纯型表节点数据
      mxzDcx = new TxzDcx(mDcxQ);

      mMatJudge = new TMatrix(mDcxQ->mMatA->GetColCount()-mDcxQ->mMatA->GetRowCount(),2);   //非基判别行矩阵MatJudge
      mCVetPivot = new TCVector(mDcxQ->mMatA->GetRowCount());     //枢列向量
      ms = new unsigned int;            //枢列坐标s
      mr = new unsigned int;            //枢行坐标r

      //传进原问题,计算取得非基判别行矩阵MatJudge及枢列坐标as、枢列向量VetPivot及枢行坐标位置,
      mResult = mxzDcx->GetPivot(mDcxQ,mMatJudge,ms,mCVetPivot,mr);

      mLxzDcx = new KTList<TxzDcx*>;  //解路径
      mLxzDcx->Append(mxzDcx);

      if (mIsRecord) {  //要记录中间结果
         while (mResult==0) {
            mxzDcx = new TxzDcx(*mxzDcx);   //拷贝一个新节点
            mxzDcx->ChangeBase(mCVetPivot,*mr,(*mMatJudge)(*ms,2),StrToInt((AnsiString)((*mMatJudge)(*ms,1))));  //换基
            mResult = mxzDcx->GetPivot(mDcxQ,mMatJudge,ms,mCVetPivot,mr);//计算枢行枢列
            mLxzDcx->Append(mxzDcx);
         }
      }
      else {     //不要记录中间结果
         while (mResult==0) {
            mxzDcx->ChangeBase(mCVetPivot,*mr,(*mMatJudge)(*ms,2),StrToInt((AnsiString)((*mMatJudge)(*ms,1))));  //换基
            mResult = mxzDcx->GetPivot(mDcxQ,mMatJudge,ms,mCVetPivot,mr);//计算枢行枢列
         }
      }

      mIsValid = true;
   }

   //用链表取得所有节点数据
   void operator>>(KTList<TxzDcx*>* aLxzDcx)
   {
      TxzDcx *lxzDcx;

      //清空原来的数据
      for (aLxzDcx->First();!aLxzDcx->IsEof();aLxzDcx->Next())
         delete aLxzDcx->GetData();
      aLxzDcx->Empty();
         
      for (mLxzDcx->First();!mLxzDcx->IsEof();mLxzDcx->Next()) {
         lxzDcx = new TxzDcx(*mLxzDcx->GetData());
         aLxzDcx->Append(lxzDcx);
      }
   }


};
//---------------------------------------------------------------------------
#endif