运筹优化入门:手把手教你用YALMIP+CPLEX在MATLAB里解第一个线性规划问题
运筹优化实战:从零构建YALMIP+CPLEX线性规划模型
第一次打开MATLAB准备尝试运筹优化时,面对空白的命令窗口,很多人会陷入"工具箱已安装但不知从何下手"的困境。本文将以一个经典的生产计划问题为例,带你完整走通YALMIP建模、CPLEX求解、结果分析的闭环流程。不同于单纯的安装教程,我们将聚焦于工具链的实战价值——如何将数学公式转化为可执行的优化代码。
1. 环境准备与问题定义
在开始编码前,确保你的MATLAB已正确配置YALMIP和CPLEX。验证方法很简单:
% 验证YALMIP安装 which yalmip % 验证CPLEX可用性 which cplex若这两条命令返回路径信息,说明环境就绪。我们假设要解决以下问题:
某工厂生产两种产品A和B,每吨利润分别为5万元和4万元。生产过程中需要消耗两种原料X和Y,库存分别为24吨和18吨。生产1吨A需要消耗3吨X和1吨Y,生产1吨B需要消耗1吨X和2吨Y。如何安排生产计划使总利润最大?
这个问题可以抽象为:
- 决策变量:产品A产量x₁,产品B产量x₂
- 目标函数:max 5x₁ + 4x₂
- 约束条件:
- 3x₁ + x₂ ≤ 24 (原料X限制)
- x₁ + 2x₂ ≤ 18 (原料Y限制)
- x₁, x₂ ≥ 0 (非负约束)
2. YALMIP建模实战
YALMIP的魅力在于其建模语法与数学表达近乎一致。下面我们分步构建模型:
% 初始化模型 clc; clear; close all; model = yalmip('Model'); % 定义决策变量 x = sdpvar(2,1); % 2x1的决策向量 % 目标函数 profit = [5 4]*x; Objective = -profit; % YALMIP默认求最小,加负号转为求最大 % 约束条件 Constraints = [ 3*x(1) + x(2) <= 24 x(1) + 2*x(2) <= 18 x >= 0 % 非负约束 ]; % 查看模型结构 disp('决策变量:'); disp(x); disp('目标函数:'); disp(Objective); disp('约束条件:'); disp(Constraints);关键点说明:
sdpvar用于声明决策变量,参数(2,1)表示创建2行1列的变量向量- 约束条件直接用MATLAB语法编写,运算符(<=, >=)会被YALMIP重载
- 显示模型结构有助于调试,建议在复杂模型构建时使用
3. 调用CPLEX求解与结果解析
模型构建完成后,只需几行代码即可调用CPLEX求解:
% 配置求解器 options = sdpsettings('solver','cplex','verbose',1); % 求解模型 sol = optimize(Constraints, Objective, options); % 检查求解状态 if sol.problem == 0 disp('最优解找到:'); disp(value(x)); disp(['最大利润: ', num2str(-value(Objective)), '万元']); else disp('求解失败,状态码:'); disp(sol.problem); end求解器输出日志通常包含以下关键信息:
CPLEX 12.10.0.0: optimal solution; objective -52 0 dual simplex iterations (0 in phase I)解读要点:
optimal solution表示找到全局最优解objective -52对应我们加负号的目标函数值,实际最大利润为52万元- 迭代次数显示求解效率
验证结果合理性:
% 验证约束条件 disp('原料X实际消耗:'); disp(3*value(x(1)) + value(x(2))); % 应<=24 disp('原料Y实际消耗:'); disp(value(x(1)) + 2*value(x(2))); % 应<=184. 模型扩展与可视化分析
基础模型运行成功后,我们可以进一步探索:
敏感性分析示例:
% 研究产品A利润变化对结果的影响 profit_A_range = 3:0.1:7; solutions = zeros(length(profit_A_range), 2); for i = 1:length(profit_A_range) Objective = -([profit_A_range(i) 4]*x); optimize(Constraints, Objective, options); solutions(i,:) = value(x)'; end % 可视化结果 figure; plot(profit_A_range, solutions); xlabel('产品A单位利润(万元)'); ylabel('最优生产量(吨)'); legend('产品A','产品B'); grid on;添加整数约束:
% 假设产品必须整吨生产 IntegerConstraints = [Constraints, integer(x)]; sol = optimize(IntegerConstraints, Objective, options); disp('整数解:'); disp(value(x));常见问题排查表:
| 问题现象 | 可能原因 | 解决方案 |
|---|---|---|
Unrecognized function or variable 'sdpvar' | YALMIP路径未添加 | 运行addpath(genpath('yalmip路径')) |
No solver available | CPLEX未正确链接 | 检查which cplex是否返回有效路径 |
Infeasible problem | 约束条件矛盾 | 检查约束条件是否可能同时满足 |
5. 工业级应用技巧
在实际项目中,这些技巧能显著提升效率:
批量建模方法:
% 定义多周期生产计划 n_periods = 7; % 7天计划 x = sdpvar(2, n_periods); % 2产品×7天 % 构建矩阵形式的目标函数 unit_profit = [5; 4]; Objective = -sum(unit_profit' * x); % 添加库存动态约束 initial_stock = [24; 18]; consumption_rate = [3 1; 1 2]; Constraints = [x >= 0]; for t = 1:n_periods if t == 1 stock = initial_stock - consumption_rate * x(:,t); else stock = stock - consumption_rate * x(:,t); end Constraints = [Constraints, stock >= 0]; end求解器参数调优:
% 设置CPLEX高级参数 options = sdpsettings('solver','cplex',... 'cplex.timelimit', 3600,... 'cplex.mip.tolerances.mipgap', 0.001,... 'cplex.parallel', -1); % 使用所有可用线程结果导出与报告生成:
% 将结果保存为结构体 results.Production = value(x); results.Profit = -value(Objective); results.Consumption = consumption_rate * value(x); % 导出到Excel writetable(struct2table(results), 'Production_Plan.xlsx'); % 生成可视化报告 figure; subplot(2,1,1); bar(value(x)'); legend('产品A','产品B'); title('最优生产计划'); subplot(2,1,2); pie([results.Profit, sum(results.Consumption(:))], {'利润','成本'}); title('利润成本分析');遇到性能瓶颈时,可以尝试以下优化策略:
- 使用
sdpsettings('solver','cplex','usex0',1)复用初始解 - 对大规模问题启用
'cplex.preprocessing.presolve'选项 - 考虑将部分约束转化为惰性约束(lazy constraints)
