★ FEATURED ARTICLE 简介这份资源是面向算法学习者与C开发者的单纯形法实现程序用于求解线性规划问题的最优解。单纯形法由丹齐格于1947年提出通过迭代在基本可行解之间转换逐步逼近目标函数最优值是运筹学与数值计算中的经典算法。资源包共3个文件包含2个cpp源文件与1个h头文件压缩包约2KB代码围绕单纯形表、基变量与非基变量替换、迭代终止条件等核心环节展开结构紧凑便于直接阅读与二次开发。目前已有742人学习下载说明其在算法教学与练习场景中具有一定参考价值。读者可借助这份代码理解线性规划标准形式的转换思路、单纯形表的矩阵存储方式以及用C容器组织系数矩阵与向量的具体做法并在此基础上尝试性能优化与异常情况处理适合作为课程实验、算法练习或自学单纯形法的实践素材。1. 从一份 C 单纯形算法计算程序说起线性规划求解的工程化落地很多做运筹优化或算法工程的朋友第一次接触线性规划求解往往是从 Python 的 scipy.optimize.linprog 或者 MATLAB 的 linprog 起步的。调包确实快但一旦遇到需要嵌入到 C 后端服务、嵌入式设备或者需要自己控制迭代过程、观察单纯形表每一步变化的场景调包就不够用了。这份 C 实现的单纯形算法计算程序解决的正是这个问题它把线性规划的标准型求解过程完整地用 C 写了出来包含单纯形表的构建、入基变量选择、出基变量比值检验、高斯消元迭代以及最优性判定。适合正在学运筹学课程、需要交算法作业的学生也适合想把求解器嵌入自己 C 项目的开发者。它不依赖任何第三方数学库纯标准 C 实现拿到手编译就能跑。2. 项目结构与核心模块拆解一个线性规划问题是怎么被吃进去的2.1 目录结构与文件职责拿到这份源码第一件事是看目录。典型的工程化组织方式大致是这样simplex-solver/ ├── include/ │ ├── SimplexSolver.h // 求解器类声明 │ ├── LinearProgram.h // 线性规划模型数据结构 │ └── Matrix.h // 轻量矩阵封装 ├── src/ │ ├── SimplexSolver.cpp // 单纯形迭代核心实现 │ ├── LinearProgram.cpp // 模型构建与标准化 │ └── main.cpp // 命令行入口 ├── tests/ │ └── test_cases.txt // 测试用例 ├── CMakeLists.txt └── README.md这个结构不算复杂但职责切得很清楚。LinearProgram负责描述问题本身——目标函数系数、约束矩阵、右端项、变量取值范围SimplexSolver负责算法迭代main.cpp负责读输入、调求解器、打印结果。我一般建议二次开发时不要把所有逻辑塞进一个文件否则改到后面自己都找不到入基变量是在哪一行选的。2.2 线性规划模型的数据结构先看模型是怎么表示的。标准型线性规划通常写成max c^T x s.t. Ax b x 0对应到 C 里核心结构体大概长这样// LinearProgram.h struct LinearProgram { int numVars; // 决策变量个数 int numConstraints; // 约束个数 std::vectordouble c; // 目标函数系数 std::vectorstd::vectordouble A; // 约束系数矩阵 std::vectordouble b; // 右端项 std::vectorint constraintType; // 0: , 1: , 2: };这里有几个参数需要留意。constraintType决定了后续要不要引入松弛变量、剩余变量或人工变量。b向量里的元素如果出现负数标准化阶段需要两边同乘 -1 并翻转不等号方向这一步是很多手写实现翻车的地方。c向量在求最小值问题时要整体取反因为单纯形法默认处理最大化问题。2.3 标准化把各种形式统一成单纯形法能吃的形状单纯形法本身只认一种形式所有约束都是等式右端项非负并且存在一组初始基变量。所以从用户输入的原始问题到可迭代的标准型中间要做几步转换// LinearProgram.cpp 片段 void LinearProgram::toStandardForm() { // 1. 若为 min 问题目标系数取反 if (isMinimization) { for (auto ci : c) ci -ci; } // 2. 对每个约束引入松弛/剩余/人工变量 for (int i 0; i numConstraints; i) { if (constraintType[i] 0) { // 加松弛变量 addSlackVariable(i); } else if (constraintType[i] 1) { // 减剩余变量再加人工变量 addSurplusVariable(i); addArtificialVariable(i); } else { // 直接加人工变量 addArtificialVariable(i); } } // 3. 若 b[i] 0整行取反 for (int i 0; i numConstraints; i) { if (b[i] 0) { for (auto aij : A[i]) aij -aij; b[i] -b[i]; } } }这段代码的逻辑是松弛变量用来把变成系数为 1剩余变量处理系数为 -1人工变量则是为了给等式约束凑出一个初始可行基。人工变量在目标函数里的系数要设成一个大负数大 M 法或者单独走两阶段法。这份实现里用的是大 M 法M 的取值一般设成 1e6 到 1e9 之间太小压不住人工变量太大又容易引起数值不稳定。注意如果原始问题里同时存在和约束标准化之后变量总数会明显增加单纯形表的规模也跟着涨。调试时建议先把标准化后的矩阵打印出来核对一遍。3. 单纯形迭代核心从入基到出基的完整链路3.1 单纯形表的构建与初始基标准化完成后所有变量原始变量 松弛 剩余 人工构成一个大的系数矩阵。单纯形表本质上就是把这个矩阵和目标函数行拼在一起// SimplexSolver.cpp void SimplexSolver::buildTableau() { int totalVars lp.getTotalVars(); tableau.resize(lp.numConstraints 1, std::vectordouble(totalVars 1, 0.0)); // 填充约束行 for (int i 0; i lp.numConstraints; i) { for (int j 0; j totalVars; j) { tableau[i][j] lp.A[i][j]; } tableau[i][totalVars] lp.b[i]; // 最后一列是右端项 } // 填充目标函数行检验数行 for (int j 0; j totalVars; j) { tableau[lp.numConstraints][j] -lp.c[j]; } // 记录基变量索引 basis.resize(lp.numConstraints); for (int i 0; i lp.numConstraints; i) { basis[i] lp.getInitialBasisIndex(i); } }这里tableau的最后一行是检验数行初始时直接取目标系数的相反数。basis数组记录每一行当前对应的基变量是哪个。初始基通常来自松弛变量和人工变量具体是哪些列取决于标准化阶段的处理顺序。3.2 入基变量的选择策略入基变量的选择直接决定迭代效率。最常见的是 Dantzig 规则——选检验数最负的那一列int SimplexSolver::selectEnteringVariable() { int enterCol -1; double minVal -1e-9; // 容差避免浮点误差导致的死循环 int lastRow lp.numConstraints; for (int j 0; j lp.getTotalVars(); j) { if (tableau[lastRow][j] minVal) { minVal tableau[lastRow][j]; enterCol j; } } return enterCol; // 返回 -1 表示已最优 }容差1e-9是个血泪经验。如果不加容差浮点运算产生的极小负数会被当成可以继续优化的信号导致算法在最优解附近来回震荡。Dantzig 规则实现简单但在某些问题上会出现循环实际工程里如果遇到退化问题可以换成 Bland 规则——选下标最小的负检验数对应的变量。3.3 出基变量的比值检验选定入基变量后要在该列的正系数里找比值最小的那一行int SimplexSolver::selectLeavingVariable(int enterCol) { int leaveRow -1; double minRatio 1e18; for (int i 0; i lp.numConstraints; i) { if (tableau[i][enterCol] 1e-9) { // 只考虑正系数 double ratio tableau[i][lp.getTotalVars()] / tableau[i][enterCol]; if (ratio minRatio) { minRatio ratio; leaveRow i; } } } return leaveRow; // 返回 -1 表示问题无界 }比值检验里有两个关键点。第一只考虑入基列系数为正的行系数为零或负的行不能作为出基候选否则会导致右端项变成负数破坏可行性。第二如果所有行的系数都不为正说明问题无界此时应该直接返回并报告无界而不是继续迭代。3.4 高斯消元与基变换确定入基和出基之后做一次行变换把入基列变成单位向量void SimplexSolver::pivot(int enterCol, int leaveRow) { double pivotVal tableau[leaveRow][enterCol]; // 1. 主元行归一化 for (int j 0; j lp.getTotalVars(); j) { tableau[leaveRow][j] / pivotVal; } // 2. 其他行消元 for (int i 0; i lp.numConstraints; i) { if (i leaveRow) continue; double factor tableau[i][enterCol]; if (std::abs(factor) 1e-12) continue; for (int j 0; j lp.getTotalVars(); j) { tableau[i][j] - factor * tableau[leaveRow][j]; } } // 3. 更新基变量索引 basis[leaveRow] enterCol; }pivotVal是主元归一化之后主元行对应的基变量值就是右端项。其他行按照factor倍减去主元行把入基列的其他元素消成零。basis数组同步更新保证后续读取解的时候能对上号。3.5 主循环与终止条件把上面几步串起来主循环的逻辑就很清晰了SolveResult SimplexSolver::solve() { buildTableau(); while (true) { int enterCol selectEnteringVariable(); if (enterCol -1) break; // 最优 int leaveRow selectLeavingVariable(enterCol); if (leaveRow -1) { // 无界 return {UNBOUNDED, {}}; } pivot(enterCol, leaveRow); if (iterations MAX_ITER) { // 防止死循环 return {ITER_LIMIT, {}}; } } return extractSolution(); }MAX_ITER一般设成变量数的几倍就够了设太大反而掩盖了潜在的循环问题。extractSolution从basis和右端项里把原始变量的值读出来人工变量如果非零说明原问题不可行。4. 编译运行与测试从源码到可执行文件4.1 编译环境与 CMake 配置这份代码只依赖标准库编译门槛很低。用 CMake 的话CMakeLists.txt大致是cmake_minimum_required(VERSION 3.10) project(SimplexSolver CXX) set(CMAKE_CXX_STANDARD 17) set(CMAKE_CXX_STANDARD_REQUIRED ON) include_directories(include) add_executable(simplex src/main.cpp src/SimplexSolver.cpp src/LinearProgram.cpp )编译命令mkdir build cd build cmake .. make如果不想用 CMake直接 g 一把梭也行g -stdc17 -O2 -Iinclude src/*.cpp -o simplex-O2建议加上单纯形迭代里矩阵运算比较多开优化后大规模问题的求解速度差别很明显。4.2 输入格式与运行示例程序通常从标准输入或文件读取问题描述。一个典型的输入格式是3 2 2 3 1 1 1 0 2 1 1 4 3 0 0第一行是变量数和约束数第二行是目标函数系数接下来是约束矩阵和右端项最后一行是约束类型。运行./simplex tests/test_cases.txt输出会包含最优目标值、各变量取值、迭代次数。如果问题不可行或无界也会给出对应状态。4.3 测试用例设计测试用例至少要覆盖这几种情况标准最大化问题、最小化问题、含等式约束的问题、含约束的问题、无界问题、不可行问题。每类问题准备两到三个手工算一遍最优解跟程序输出对比。我一般会在tests/下放一个简单的对比脚本把预期结果和实际结果逐行 diff这样改代码之后跑一遍就知道有没有改坏。5. 避坑与常见问题排查5.1 迭代不收敛程序卡死现象程序运行后长时间不输出结果CPU 占用高。原因通常是浮点误差导致检验数行出现极小的负数入基变量反复选中同一列形成循环。解决办法是在selectEnteringVariable里加容差判断小于1e-9的检验数直接视为零如果问题本身退化严重切换到 Bland 规则强制选下标最小的候选变量。5.2 最优解数值偏差大现象程序给出的最优值和手工计算结果差了好几个数量级。原因多半是大 M 法里 M 取值不当或者标准化阶段右端项取反时漏掉了某些行。排查方法是把标准化后的tableau完整打印出来逐行核对系数和右端项符号。M 值建议从1e6起步如果目标函数系数本身很大按比例放大。5.3 无界问题被误判为最优现象问题明明无界程序却返回了一个“最优解”。原因通常是比值检验时没有正确识别所有系数非正的情况或者容差设得太大把本该参与比值计算的正系数过滤掉了。检查selectLeavingVariable里的阈值确保只过滤掉真正接近零的系数而不是把所有小系数都当成零。5.4 人工变量残留在最终基里现象求解结束后extractSolution读出的原始变量值看起来正常但人工变量对应的值不为零。这说明原问题不可行但程序没有正确报告。解决办法是在提取解之后检查所有人工变量对应的基变量值如果任何一个大于容差直接返回INFEASIBLE状态而不是硬给一个解。5.5 大规模问题内存暴涨现象变量数和约束数上去之后程序内存占用迅速增加甚至崩溃。原因是tableau用的是稠密矩阵规模是(m1) x (nm1)当 m 和 n 都到几千时内存开销是平方级的。如果确实要处理大规模稀疏问题需要把稠密矩阵换成稀疏存储或者改用修正单纯形法只维护基矩阵的逆。这份实现定位是教学和小规模求解上生产环境前要评估规模。6. 二次开发与改造思路从能跑到好用这份代码跑通之后真正有意思的是怎么按自己的需求改。我梳理几个常见的改造方向。第一个方向是换入基规则。Dantzig 规则实现简单但迭代次数多可以改成最大改进量规则——选使目标函数下降最多的那一列代价是每轮多算一次比值。也可以上 Bland 规则防循环虽然慢一点但稳定。改的时候只需要动selectEnteringVariable一个函数其他逻辑不用碰。第二个方向是加对偶单纯形。有些场景下初始基可行但检验数不满足最优条件比如在灵敏度分析里改了右端项之后用对偶单纯形从当前基继续迭代比重新求解快得多。实现上就是把入基和出基的判断条件对调主循环结构不变。第三个方向是输出迭代过程。教学场景下每一步的单纯形表、入基变量、出基变量、当前目标值都是有用的信息。可以在pivot之后加一个回调或者日志输出把tableau和basis打出来。注意别在性能敏感的路径上做字符串拼接用条件编译或者日志级别控制。第四个方向是接口封装。如果要把求解器嵌到更大的系统里建议把main.cpp里的输入输出逻辑剥掉只保留LinearProgram和SimplexSolver两个类对外暴露一个solve(lp)接口。输入解析交给调用方求解器只负责算。改造方向改动范围难度适用场景换入基规则selectEnteringVariable低迭代次数优化对偶单纯形主循环 选择逻辑中灵敏度分析迭代日志pivot后加回调低教学演示稀疏矩阵tableau存储结构高大规模问题接口封装剥离 main低嵌入系统验证改造是否正确的办法很直接准备一组已知最优解的测试用例改之前跑一遍记录结果改之后再跑一遍对比。如果结果一致说明改造没有破坏正确性如果迭代次数变了但结果一致说明改造生效了。我一般会额外准备几个退化问题和无界问题专门用来测边界情况。从那以后我每次改求解器核心逻辑都强制走一遍完整测试用例哪怕只改了一行容差。单纯形法的数值稳定性问题往往在特定输入下才暴露靠肉眼 review 根本看不出来。希望这份拆解能帮你少走点弯路把这份 C 单纯形算法程序真正用起来。本文还有配套的精品资源点击获取
阅读完成 · 觉得有帮助？
