C#的线性规划代码

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

namespace Dispatch.linprog
{
    /***************
    单纯形法是解线性规划问题的一个重要方法。 
    其原理的基本框架为： 
    第一步：将LP线性规划变标准型，确定一个初始可行解（顶点）。 
    第二步：对初始基可行解最优性判别，若最优，停止；否则转下一步。 
    第三步：从初始基可行解向相邻的基可行解（顶点）转换，且使目标值有所改善—目标函数值增加，重复第二和第三步直到找到最优解。 
    Max z=CX;
    AX <= b
    X=(x1,x2,...xn)>=(0,0,...,0)
    b=(b1,b2,...,bn)
    C = (c1, c2, ... ,cn)
    A= a11, a12, ..., a1n
       a21, a22, ..., a2n
       a31, a32, ..., a3n
          ......
       am1, am2, ..., amn
    **********************/

     
    class SimplexLP
    {
        public double INT_MAX = 1000000000000;
        public double[][] Matrix;
        public int CNum, BNum;
        //public double Z;
        List<int> PList;

        public bool OutputFlag = false;

        //标准系数矩阵，设 变量 （CNum维向量） 为X，前BNum维是基变量
        //第一行为求Max的函数系数序列。
        //最后一列为线性规划问题的值向量b
        //除第一行最后一列的系数矩阵为 约束矩阵A
        public bool InitData(double[][] matrix, int bn = 0,  int cn = 0)
        {
            if (matrix != null && matrix.Length > 0 && matrix[0].Length > 0)
            {
                Matrix = matrix;
                this.BNum = bn;
                this.CNum = cn;

                if (bn == 0)
                    BNum = matrix.Length;
                if (cn == 0)
                    CNum = matrix[0].Length;

                //不等式矩阵的第一行为Max函数的系数，PList包含 列数-行数-2 个非基变量。
                PList = new List<int>();
                for (int i = 0; i < BNum - 1; i++)
                    PList.Add(CNum - i - 2);

                return true;
            }
            return false;
        }

        bool Pivot(ref int x, ref int y)//返回0表示所有的非轴元素都小于0
        {
            double cmax = -INT_MAX;
            List<double> clist = Matrix[0].ToList();
            List<double> blist = new List<double>();

            x = 0;
            y = 0;
            for (int i = 0; i < BNum; i++)//值序列
            {
                blist.Add(Matrix[i][CNum - 1]);
            }

            for (int i = 0; i < clist.Count; i++)//在第一行的非轴元素系数中找最大的c
            {
                if (cmax < clist[i] && !PList.Contains(i))
                {
                    cmax = clist[i];
                    y = i;
                }
            }
            if (cmax < 0)//未找到正系数
                return false;

            double bmin = INT_MAX;
            for (int i = 1; i < BNum; i++)
            {
                double tmp = blist[i] / Matrix[i][y];
                if (Matrix[i][y] != 0 && bmin > tmp)
                {
                    bmin = tmp;
                    x = i;
                }
            }

            for(int i =0; i <  PList.Count; i ++)
            {
                if (Matrix[x][PList[i]] != 0)
                {
                    PList.Remove(PList[i]);
                    break;
                }
            }
            //把已经处理过的列标加入
            PList.Add(y);
            return true;
        }

        void showMaxtrix()
        {
            for (int i = 0; i < Matrix.Length; i++)
            {
                for (int j = 0; j < Matrix[0].Length; j++)
                {
                    Console.Write($"{Matrix[i][j]}\t");
                }
                Console.Write("\n");
            }
            Console.Write($"result z: {-Matrix[0][CNum - 1]}\n");
        }

        void Gaussian(int x, int y)//行变换
        {
            double norm = Matrix[x][y];
            for (int i = 0; i < CNum; i++)//主行归一化
            {
                Matrix[x][i] /= norm;
            }
            for (int i = 0; i < BNum && i != x; i++)
            {
                if (Matrix[i][y] != 0)
                {
                    double tmpnorm = Matrix[i][y];
                    for (int j = 0; j < CNum; j++)
                    {
                        Matrix[i][j] = Matrix[i][j] - (tmpnorm * Matrix[x][j]);
                    }
                }
            }
        }

        public void Solve()
        {
            int x = 0;
            int y = 0;
            while (true)
            {
                if(OutputFlag)
                    showMaxtrix();
                if (Pivot(ref x, ref y) == false)
                {
                    return;
                }

                if (OutputFlag)
                {
                    Console.Write($"{x} {y}\n");
                    foreach (int i in PList)
                        Console.Write($"{i} ");
                    Console.Write("\n");
                }

                Gaussian(x, y);
            }
        }
    }
}