C#实现单纯形法

前言

由于列车运行图问题涉及整数分支定界算法,而该算法需要单纯形法计算每个分支最优解,故此处运用C#实现单纯形法。希望大佬批评指点。

说明

Winform窗体,只有一个Button按钮;其中,第一部分代码为算法实现部分,第二部分代码为调用部分,调用部分c为原问题目标函数系数,a为约束系数,b为资源系数,d较为特殊,本文规定,若为<=约束,d值为-1,若为=约束,d值为0,若为>=约束,d值为1。适用问题为最小值问题,若为最大值,请将目标函数乘以-1即可。可直接粘贴调试。

第一部分

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Threading.Tasks;

namespace SimplexMethod
{
    /// <summary>
    /// 单纯形算法
    /// </summary>
    class SimplexMethodAlgorithm
    {
        double[] CMat;//系数矩阵
        double[][] AMat;//约束矩阵
        double[] bMat;//资源矩阵

        int MConstraintCont;//需要添加人工变量的约束个数

        double M = 5000000;//大M法中的M值
        int[] XBIndex;//定义存储基变量标号的一维数组
        int[] XNIndex;//定义存储非基变量标号的一维数组
        double[] XBValue;//定义基变量值的一维数组
        double[] XNEnter;//定义进基变量的判别数
        double[] XBExit;//定义出基变量的判别数
        double[] CBValue;//定义基变量的目标函数系数的向量
        double[] CNValue;//定义非基变量的目标函数系数的向量
        double solutionF = double.MaxValue;//目标函数值
        //double[] CBValueMultiplicationInverseofB;//定义CB与B-1的乘积值
        //double[] InverseofBMultiplicationP;//定义B-1与P的乘积值
        //double[,] B;//定义存储基变量在约束表达式中的系数矩阵


        int t;//定义当前循环次数
        int enterIndex;//定义当前进基序号
        int exitIndex;//定义当前出基序号
        //int[] exitIndexJudgeSolutionDegeneracy;//定义当前出基序号及判断是否退化

        int[] judgeSymbolDirection;//判断约束的符号方向,"<="为-1,“=”为0,“>="为1

        public List<double> solutionSet = new List<double>();//存储最优解列表集

        bool isContinue = true;

        /// <summary>
        /// 算法执行函数
        /// </summary>
        public void ExecuteAlgorithm()
        {
            Console.WriteLine("\r\n\r\n\r\n");
            Console.WriteLine("算法开始迭代:");
            for (t = 1; t < 10; t++)//第t次迭代
            {
                Console.WriteLine("第  " + t + "  次迭代!");
                GetXNEnter(CBValue, XNIndex);//获取非基变量检验数
                GetEnterIndex(XNEnter);//计算进基变量下标
                if (enterIndex != -1)//当算法需继续迭代时
                {
                    GetXBExit(enterIndex);//获取基变量检验数
                    if (isContinue)//重新计算AMat、bMat
                    {
                        GetExitIndex(XBExit);//计算出基变量下标
                        //判断是否存在退化情况
                        List<int> exitMinIndex = new List<int>(); //存储最小比值下标
                        for (int i = 0; i < XBExit.Length; i++)
                        {
                            if (XBExit[i] == XBExit.Min())
                            {
                                exitMinIndex.Add(XBIndex[i]);
                            }
                        }
                        if (exitMinIndex.Count > 1)//如果出现退化情况--运用Bland法则重新选择出入基变量元素
                        {
                            for (int j = 0; j < XBIndex.Length; j++)//查找索引下标最小元素为出基变量
                            {
                                if (XBIndex[j] == exitMinIndex.Min())
                                {
                                    exitIndex = j;
                                    break;
                                }
                            }
                            for (int j = 0; j < XNEnter.Length; j++)//查找索引下标最小元素为入基变量
                            {
                                if (XNEnter[j] < 0)
                                {
                                    enterIndex = XNIndex[j];
                                    break;
                                }
                            }
                        }
                        UpdateAmatandbMatData();//更新A、b数据
                        XBIndex[exitIndex] = enterIndex;//更新基变量下标数组
                        CBValue[exitIndex] = CMat[enterIndex];//更新基变量对应系数值
                        int xnIndex = 0;
                        for (int i = 0; i < CMat.Length; i++)//更新非基变量下标数组
                        {
                            if (!XBIndex.Contains(i))
                            {
                                XNIndex[xnIndex++] = i;
                            }
                        }
                        for (int i = 0; i < XNIndex.Length; i++)//更新CNValue数组
                        {
                            CNValue[i] = CMat[XNIndex[i]];
                        }
                        for (int i = 0; i < bMat.Length; i++)//更新基变量值数组
                        {
                            XBValue[i] = bMat[i];
                        }
                    }
                    else
                    {
                        break;
                    }
                }
                else
                {
                    break;
                }
                solutionSet.Add(solutionF);
                solutionF = 0;
                for (int i = 0; i < CBValue.Length; i++)//计算迭代解
                {
                    solutionF += CBValue[i] * bMat[i];
                }
            }
        }

        /// <summary>
        /// 更新约束矩阵AMat、资源矩阵bMat
        /// </summary>
        private void UpdateAmatandbMatData()
        {
            double douValue = AMat[exitIndex][enterIndex];//定位exitIndex行XNIndex[enterIndex]列的值
            for (int i = 0; i < CMat.Length; i++)//第exitIndex行AMat值线性变换
            {
                AMat[exitIndex][i] /= douValue;
            }

            double[][] AMatCopy = new double[bMat.Length][];//此处须存储约束矩阵AMat
            for (int i = 0; i < AMatCopy.GetLength(0); i++)
            {
                AMatCopy[i] = new double[CMat.Length];
                for (int j = 0; j < AMatCopy[i].Length; j++)
                {
                    AMatCopy[i][j] = AMat[i][j];
                }
            }

            bMat[exitIndex] /= douValue;//第exitIndex行bMat值线性变换
            for (int i = 0; i < bMat.Length; i++)//除exitIndex行之外的AMat值线性变换
            {
                if (i != exitIndex)
                {
                    for (int j = 0; j < CMat.Length; j++)//更新每个约束系数值
                    {
                        AMat[i][j] += AMatCopy[exitIndex][j] * (-AMatCopy[i][enterIndex]);
                    }
                    bMat[i] += bMat[exitIndex] * (-AMatCopy[i][enterIndex]);
                }
            }
        }

        /// <summary>
        /// 计算出基变量索引
        /// </summary>
        /// <param name="_xbExit">出基变量检验数</param>
        private void GetExitIndex(double[] _xbExit)
        {
            int minIndex = 0;
            double douMinValue = _xbExit[0];
            for (int i = 1; i < bMat.Length; i++)
            {
                if (douMinValue >= _xbExit[i])
                {
                    douMinValue = _xbExit[i];
                    minIndex = i;
                }
            }
            exitIndex = minIndex;
        }

        /// <summary>
        /// 获取出基变量检验数
        /// </summary>
        /// <param name="_enterIndex">进基变量索引值</param>
        private void GetXBExit(int _enterIndex)
        {
            int judgeExistSolution = 0;
            for (int i = 0; i < bMat.Length; i++)//判断进基约束列向量中小于等于0的元素个数
            {
                if (AMat[i][_enterIndex] <= 0)
                {
                    judgeExistSolution++;
                }
            }
            if (judgeExistSolution == bMat.Length)//如果所有元素均小于等于0
            {
                Console.WriteLine("**********该问题无可行解!**********");
                exitIndex = -1;
                isContinue = false;
                return;
            }
            for (int i = 0; i < bMat.Length; i++)
            {
                if (AMat[i][_enterIndex] <= 0)
                {
                    XBExit[i] = double.MaxValue;
                }
                else
                {
                    XBExit[i] = bMat[i] / AMat[i][_enterIndex];
                }
                isContinue = true;
            }
        }

        /// <summary>
        /// 计算进基变量索引
        /// </summary>
        /// <param name="_xnEnter">非基变量检验数</param>
        private void GetEnterIndex(double[] _xnEnter)
        {
            List<int> minute = new List<int>();//负值检验数标号集合----可进基标号集合
            List<int> positive = new List<int>();//正值检验数标号集合
            List<int> zero = new List<int>();//零值检验数标号集合
            for (int i = 0; i < _xnEnter.Length; i++)
            {
                if (_xnEnter[i] < 0)
                {
                    minute.Add(i);
                }
                else if (_xnEnter[i] == 0)
                {
                    zero.Add(i);
                }
                else
                {
                    positive.Add(i);
                }
            }
            if (minute.Count == 0)//如果非基变量检验数均大于等于0时,则表示当前解就是最优解或无解
            {
                for (int j = 0; j < CBValue.Length; j++)
                {
                    if (CBValue[j] == M)
                    {
                        if (XBValue[j] != 0)
                        {
                            enterIndex = -1;
                        }
                    }
                }
                if (positive.Count == _xnEnter.Length)//如果当前解为最优解
                {
                    for (int i = 0; i < bMat.Length; i++)
                    {
                        Console.Write("X" + (XBIndex[i] + 1).ToString());
                        Console.Write("        ");
                    }
                    Console.WriteLine();
                    for (int i = 0; i < bMat.Length; i++)
                    {
                        Console.Write(bMat[i].ToString());
                        Console.Write("        ");
                    }
                    Console.WriteLine();
                    Console.Write(solutionF);
                    Console.WriteLine();
                    Console.WriteLine("**********找到最优解!**********");
                    enterIndex = -1;
                }
                else
                {
                    for (int i = 0; i < bMat.Length; i++)
                    {
                        Console.Write("X" + (XBIndex[i] + 1).ToString());
                        Console.Write("        ");
                    }
                    Console.WriteLine();
                    for (int i = 0; i < bMat.Length; i++)
                    {
                        Console.Write(bMat[i].ToString());
                        Console.Write("        ");
                    }
                    Console.WriteLine();
                    Console.Write(solutionF);
                    Console.WriteLine();
                    Console.WriteLine("**********有无穷最优解!**********");
                    enterIndex = -1;
                }
            }
            else//如果非基变量检验数存在负值,该解还不是最优解,需要继续迭代,选择检验数最小
            {
                double minValue = _xnEnter[0];
                int minValueIndex = 0;
                for (int j = 1; j < _xnEnter.Length; j++)
                {
                    if (minValue > _xnEnter[j])
                    {
                        minValueIndex = j;
                        minValue = _xnEnter[j];
                    }
                }
                enterIndex = XNIndex[minValueIndex];
            }
        }

        /// <summary>
        /// 计算非基变量检验数
        /// </summary>
        /// <param name="_cbvalue">基变量目标函数系数向量</param>
        /// <param name="_xnindex">非基变量标号向量</param>
        private void GetXNEnter(double[] _cbvalue, int[] _xnindex)
        {
            for (int i = 0; i < _xnindex.Length; i++)//针对每个非基变量资源向量值
            {
                XNEnter[i] = 0;
                for (int j = 0; j < _cbvalue.Length; j++)//针对每个基变量系数值
                {
                    XNEnter[i] += _cbvalue[j] * AMat[j][_xnindex[i]];
                }
                XNEnter[i] = CMat[_xnindex[i]] - XNEnter[i];
            }
        }

        /// <summary>
        /// 构造函数
        /// </summary>
        ///<param name="_cMat">目标函数系数向量</param>
        ///<param name="_aMat">约束方程系数矩阵</param>
        ///<param name="_bMat">约束方程资源向量</param>
        ///<param name="_judgeSymbolDirection">约束条件的符号方向表示向量</param>
        public SimplexMethodAlgorithm(double[] _cMat, double[][] _aMat, double[] _bMat, int[] _judgeSymbolDirection)
        {
            CMat = _cMat;
            AMat = _aMat;
            bMat = _bMat;
            judgeSymbolDirection = _judgeSymbolDirection;

            XBIndex = new int[bMat.Length];

            //OutParamater(CMat, AMat, bMat,XBIndex);
            Console.WriteLine("\r\n\r\n\r\n");
            GetStandardEquation(ref CMat, ref AMat, ref judgeSymbolDirection);//标准化方程
            OutParamater(CMat, AMat, bMat, XBIndex);//输出标准化方程

            XNIndex = new int[CMat.Length - XBIndex.Length];
            int xnIndex = 0;
            for (int i = 0; i < CMat.Length; i++)
            {
                if (!XBIndex.Contains(i))
                {
                    XNIndex[xnIndex++] = i;
                }
            }

            XBValue = new double[XBIndex.Length];
            XNEnter = new double[XNIndex.Length];
            XBExit = new double[XBIndex.Length];
            CBValue = new double[XBIndex.Length];
            CNValue = new double[XNIndex.Length];
            solutionF = 0;
            //CBValueMultiplicationInverseofB = new double[CBValue.Length];
            //InverseofBMultiplicationP = new double[CBValue.Length];

            for (int i = 0; i < XBValue.Length; i++)
            {
                XBValue[i] = bMat[i];
                CBValue[i] = CMat[XBIndex[i]];
            }
            for (int i = 0; i < XNIndex.Length; i++)
            {
                CNValue[i] = CMat[XNIndex[i]];
            }

            for (int i = 0; i < CBValue.Length; i++)
            {
                solutionF += CBValue[i] * bMat[i];
            }

            solutionSet.Add(solutionF);

            t = 0;
            enterIndex = -1;
            exitIndex = -1;
            //exitIndexJudgeSolutionDegeneracy = new int[2];
        }

        /// <summary>
        /// 方程标准化
        /// </summary>
        /// <param name="_cMat">目标函数系数向量</param>
        /// <param name="_aMat">约束条件系数矩阵</param>
        /// <param name="_judgeSymbolDirection">约束条件符号表示向量</param>
        private void GetStandardEquation(ref double[] _cMat, ref double[][] _aMat, ref int[] _judgeSymbolDirection)
        {
            int index = 0;//基向量索引下标
            for (int i = 0; i < _judgeSymbolDirection.Length; i++)
            {
                //小于等于0的约束
                if (_judgeSymbolDirection[i] == -1)
                {
                    //更新A[i]值
                    double[][] ARowCopy = new double[_aMat.GetLength(0)][];//定义拷贝并拓展A值多维数组
                    for (int j = 0; j < ARowCopy.GetLength(0); j++)
                    {
                        ARowCopy[j] = new double[_cMat.Length + 1];
                        if (i == j)
                        {
                            for (int k = 0; k < _cMat.Length; k++)
                            {
                                ARowCopy[j][k] = _aMat[j][k];
                            }
                            ARowCopy[j][_cMat.Length] = 1;
                        }
                        else
                        {
                            for (int k = 0; k < _cMat.Length; k++)
                            {
                                ARowCopy[j][k] = _aMat[j][k];
                            }
                        }
                    }
                    _aMat = ARowCopy;
                    //更新C值
                    double[] CRowCopy = new double[_cMat.Length + 1];
                    for (int j = 0; j < CRowCopy.Length; j++)
                    {
                        if (j < _cMat.Length)
                        {
                            CRowCopy[j] = _cMat[j];
                        }
                        else
                        {
                            CRowCopy[j] = 0;
                        }
                    }
                    _cMat = CRowCopy;
                    //存储基变量下标
                    XBIndex[index++] = _cMat.Length - 1;
                }
                //等于0的约束
                else if (_judgeSymbolDirection[i] == 0)
                {
                    //更新A[i]值
                    double[][] ARowCopy = new double[_aMat.GetLength(0)][];//定义拷贝并拓展A值多维数组
                    for (int j = 0; j < ARowCopy.GetLength(0); j++)
                    {
                        ARowCopy[j] = new double[_cMat.Length + 1];
                        if (i == j)
                        {
                            for (int k = 0; k < _cMat.Length; k++)
                            {
                                ARowCopy[j][k] = _aMat[j][k];
                            }
                            ARowCopy[j][_cMat.Length] = 1;
                        }
                        else
                        {
                            for (int k = 0; k < _cMat.Length; k++)
                            {
                                ARowCopy[j][k] = _aMat[j][k];
                            }
                        }
                    }
                    _aMat = ARowCopy;
                    //更新C值
                    double[] CRowCopy = new double[_cMat.Length + 1];
                    for (int j = 0; j < CRowCopy.Length; j++)
                    {
                        if (j < _cMat.Length)
                        {
                            CRowCopy[j] = _cMat[j];
                        }
                        else
                        {
                            CRowCopy[j] = M;
                        }
                    }
                    _cMat = CRowCopy;
                    //存储基变量下标
                    XBIndex[index++] = _cMat.Length - 1;
                }
                //大于等于0的约束
                else
                {
                    //更新A[i]值
                    double[][] ARowCopy = new double[_aMat.GetLength(0)][];//定义拷贝并拓展A值多维数组
                    for (int j = 0; j < ARowCopy.GetLength(0); j++)
                    {
                        ARowCopy[j] = new double[_cMat.Length + 2];
                        if (i == j)
                        {
                            for (int k = 0; k < _cMat.Length; k++)
                            {
                                ARowCopy[j][k] = _aMat[j][k];
                            }
                            ARowCopy[j][_cMat.Length] = -1;
                            ARowCopy[j][_cMat.Length + 1] = 1;
                        }
                        else
                        {
                            for (int k = 0; k < _cMat.Length; k++)
                            {
                                ARowCopy[j][k] = _aMat[j][k];
                            }
                        }
                    }
                    _aMat = ARowCopy;
                    //更新C值
                    double[] CRowCopy = new double[_cMat.Length + 2];
                    for (int j = 0; j < CRowCopy.Length; j++)
                    {
                        if (j < _cMat.Length)
                        {
                            CRowCopy[j] = _cMat[j];
                        }
                        else if (j == _cMat.Length)
                        {
                            CRowCopy[j] = 0;
                        }
                        else
                        {
                            CRowCopy[j] = M;
                        }
                    }
                    _cMat = CRowCopy;
                    //存储基变量下标
                    XBIndex[index++] = _cMat.Length - 1;
                }
            }
        }

        /// <summary>
        /// 输出模型参数系数
        /// </summary>
        /// <param name="_cMat">目标函数系数向量</param>
        /// <param name="_aMat">约束条件系数矩阵</param>
        /// <param name="_bMat">约束条件资源向量</param>
        private void OutParamater(double[] _cMat, double[][] _aMat, double[] _bMat, int[] _XBIndex)
        {
            Console.WriteLine("目标函数系数向量:");
            for (int i = 0; i < _cMat.Length; i++)
            {
                Console.Write(_cMat[i] + "   ");
            }
            Console.WriteLine();

            Console.WriteLine("约束条件系数矩阵:");
            for (int i = 0; i < _aMat.Length; i++)
            {
                for (int j = 0; j < _aMat[i].Length; j++)
                {
                    Console.Write(AMat[i][j] + "   ");
                }
                Console.Write("   " + _bMat[i]);
                Console.WriteLine();
            }
            Console.WriteLine("基向量索引下标:");
            for (int i = 0; i < _XBIndex.Length; i++)
            {
                Console.Write(_XBIndex[i] + "   ");
            }
            Console.WriteLine();
        }
    }
}

第二部分

using System;
using System.Collections.Generic;
using System.ComponentModel;
using System.Data;
using System.Drawing;
using System.Linq;
using System.Text;
using System.Threading.Tasks;
using System.Windows.Forms;

namespace SimplexMethod
{
    public partial class Form1 : Form
    {
        public Form1()
        {
            InitializeComponent();
        }

        private void button1_Click(object sender, EventArgs e)
        {
            Console.WriteLine("窗体调用成功!");

            double[] c = new double[] { -2, -1 };
            double[][] a = new double[][] { new double[] { 0, 5 }, new double[] { 6, 2 }, new double[] { 1, 1 } };
            double[] b = new double[] { 15, 24, 5 };
            int[] d = new int[] { -1, -1, -1 };

            //double[] c = new double[] { -2,-3 };
            //double[][] a = new double[][] { new double[] { -1, 1 }, new double[] { 1, 2 }, new double[] { 3, 1 } };
            //double[] b = new double[] { 2,10, 15 };
            //int[] d = new int[] { -1, -1,-1 };

            //double[] c = new double[] { -4, -1 };
            //double[][] a = new double[][] { new double[] { -1, 2 }, new double[] { 2, 3 }, new double[] { 1, -1 } };
            //double[] b = new double[] { 4, 12, 3 };
            //int[] d = new int[] { -1, -1, -1 };

            //double[] c = new double[] { -6, -4 };
            //double[][] a = new double[][] { new double[] { 1, 0 }, new double[] { 0, 1 }, new double[] { 3, 2 } };
            //double[] b = new double[] { 3, 4, 11 };
            //int[] d = new int[] { -1, -1, -1 };



            SimplexMethodAlgorithm simplexMethodAlgorithm = new SimplexMethodAlgorithm(c, a, b,d);
            simplexMethodAlgorithm.ExecuteAlgorithm();

            List<double> solutionList = simplexMethodAlgorithm.solutionSet;
        }
    }
}
  • 1
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值