Julia 教程 / Julia 优化建模(JuMP.jl)
Julia 优化建模(JuMP.jl)
JuMP.jl 是 Julia 生态中的数学优化建模语言,支持线性规划、整数规划、非线性规划等多种优化问题。它提供简洁的代数建模语法,可以对接数十种优化求解器。
1. JuMP 概述
using JuMP, HiGHS
# 最简单的线性规划
# minimize 2x + 3y
# subject to: x + y >= 1, x >= 0, y >= 0
model = Model(HiGHS.Optimizer)
@variable(model, x >= 0)
@variable(model, y >= 0)
@objective(model, Min, 2x + 3y)
@constraint(model, x + y >= 1)
optimize!(model)
println("目标值: $(objective_value(model))")
println("x = $(value(x))")
println("y = $(value(y))")
println("求解状态: $(termination_status(model))")
| 核心宏 | 说明 |
|---|---|
@variable(model, x) | 定义变量 |
@objective(model, sense, expr) | 定义目标函数 |
@constraint(model, expr) | 定义约束 |
@expression(model, expr) | 定义表达式(不创建约束) |
@NLconstraint(model, expr) | 非线性约束 |
@NLobjective(model, sense, expr) | 非线性目标 |
2. 线性规划(LP)
using JuMP, HiGHS
# 生产计划问题
# 工厂生产两种产品,最大化利润
model = Model(HiGHS.Optimizer)
# 变量:产品产量
@variable(model, x1 >= 0) # 产品 A
@variable(model, x2 >= 0) # 产品 B
# 目标:最大化利润
@objective(model, Max, 40x1 + 30x2)
# 约束
@constraint(model, x1 + x2 <= 100) # 原材料限制
@constraint(model, 2x1 + x2 <= 150) # 工时限制
@constraint(model, x1 <= 60) # 产品 A 需求上限
@constraint(model, x2 <= 80) # 产品 B 需求上限
optimize!(model)
if termination_status(model) == OPTIMAL
println("最优生产方案:")
println(" 产品 A: $(value(x1)) 件")
println(" 产品 B: $(value(x2)) 件")
println(" 最大利润: $(objective_value(model)) 元")
# 灵敏度分析
println("原材料影子价格: $(shadow_price(constraint_by_name(model, \"x1 + x2 <= 100\")))")
end
矩阵形式
using JuMP, HiGHS, LinearAlgebra
# 标准形式 LP: min c'x, s.t. Ax <= b, x >= 0
n, m = 5, 3
c = rand(n)
A = rand(m, n)
b = rand(m)
model = Model(HiGHS.Optimizer)
@variable(model, x[1:n] >= 0)
@objective(model, Min, dot(c, x))
@constraint(model, A * x .<= b)
optimize!(model)
3. 整数规划(IP)与混合整数规划(MIP)
using JuMP, HiGHS
# 背包问题:在容量限制下最大化价值
n = 10 # 物品数量
weights = rand(1:20, n)
values = rand(10:100, n)
capacity = 50
model = Model(HiGHS.Optimizer)
# 二进制变量:0-1 背包
@variable(model, x[1:n], Bin)
@objective(model, Max, dot(values, x))
@constraint(model, dot(weights, x) <= capacity)
optimize!(model)
selected = [i for i in 1:n if value(x[i]) > 0.5]
println("选中物品: $selected")
println("总价值: $(objective_value(model))")
println("总重量: $(sum(weights[selected]))")
MIP 示例:设施选址
using JuMP, HiGHS
# 设施选址问题
n_facilities = 5 # 候选设施数
n_customers = 10 # 客户数
# 成本参数
f = rand(50:200, n_facilities) # 设施固定成本
c = rand(5:20, n_facilities, n_customers) # 运输成本
d = rand(5:15, n_customers) # 客户需求
cap = rand(50:100, n_facilities) # 设施容量
model = Model(HiGHS.Optimizer)
# 决策变量
@variable(model, y[1:n_facilities], Bin) # 是否开设设施
@variable(model, x[1:n_facilities, 1:n_customers] >= 0) # 运输量
# 目标:最小化总成本
@objective(model, Min,
sum(f[i] * y[i] for i in 1:n_facilities) +
sum(c[i,j] * x[i,j] for i in 1:n_facilities, j in 1:n_customers)
)
# 容量约束
@constraint(model, [i in 1:n_facilities],
sum(x[i,j] for j in 1:n_customers) <= cap[i] * y[i])
# 需求约束
@constraint(model, [j in 1:n_customers],
sum(x[i,j] for i in 1:n_facilities) >= d[j])
optimize!(model)
opened = [i for i in 1:n_facilities if value(y[i]) > 0.5]
println("开设设施: $opened")
println("总成本: $(objective_value(model))")
4. 约束定义
using JuMP, HiGHS
model = Model(HiGHS.Optimizer)
@variable(model, x[1:5] >= 0)
# 单个约束
@constraint(model, sum(x) == 10)
# 范围约束
@constraint(model, 1 <= x[1] + x[2] <= 5)
# 索引约束族
@constraint(model, [i in 1:5], x[i] <= 3)
# 条件约束
@constraint(model, [i in 1:3, j in 4:5], x[i] + x[j] <= 5)
# 二次约束
@variable(model, y)
@constraint(model, x[1]^2 + x[2]^2 <= 1) # 二阶锥约束
# 半正定约束
@variable(model, Z[1:3, 1:3], PSD) # 半正定矩阵变量
# 命名约束
@constraint(model, budget, sum(x) <= 100)
shadow_price(budget) # 获取影子价格
5. 目标函数
using JuMP, HiGHS
model = Model(HiGHS.Optimizer)
@variable(model, x[1:3] >= 0)
# 线性目标
@objective(model, Min, 2x[1] + 3x[2] + x[3])
@objective(model, Max, x[1] - x[2])
# 二次目标(二次规划 QP)
@objective(model, Min, x[1]^2 + 2x[2]^2 + x[3]^2 - x[1]*x[2])
# 无目标(可行性问题)
# 不调用 @objective 即可
# 修改目标
set_objective(model, Min, x[1] + x[2] + x[3])
set_objective_coefficient(model, x[1], 5.0) # 修改 x[1] 的系数
6. 求解器选择
using JuMP
# 线性规划
model = Model(HiGHS.Optimizer) # HiGHS(开源,推荐)
model = Model(GLPK.Optimizer) # GLPK(开源)
model = Model(Clp.Optimizer) # Clp(开源)
# model = Model(Gurobi.Optimizer) # Gurobi(商业)
# model = Model(CPLEX.Optimizer) # CPLEX(商业)
# 整数规划
model = Model(HiGHS.Optimizer) # HiGHS 支持 MIP
model = Model(GLPK.Optimizer) # GLPK 支持 MIP
# model = Model(Gurobi.Optimizer) # Gurobi(最强 MIP 求解器)
# 二次规划(QP)
model = Model(Ipopt.Optimizer) # Ipopt(非线性)
model = Model(OSQP.Optimizer) # OSQP(凸 QP)
# 非线性规划
model = Model(Ipopt.Optimizer) # Ipopt
model = Model(NLopt.Optimizer) # NLopt
# 半正定规划
model = Model(COSMO.Optimizer) # COSMO
| 求解器 | LP | IP/MIP | QP | NLP | SDP | 开源 |
|---|---|---|---|---|---|---|
| HiGHS | ✅ | ✅ | ✅ | ❌ | ❌ | ✅ |
| GLPK | ✅ | ✅ | ❌ | ❌ | ❌ | ✅ |
| Ipopt | ❌ | ❌ | ✅ | ✅ | ❌ | ✅ |
| OSQP | ✅ | ❌ | ✅ | ❌ | ❌ | ✅ |
| COSMO | ✅ | ❌ | ✅ | ❌ | ✅ | ✅ |
| Gurobi | ✅ | ✅ | ✅ | ❌ | ❌ | ❌ |
| CPLEX | ✅ | ✅ | ✅ | ❌ | ❌ | ❌ |
⚠️ 注意: 商业求解器(Gurobi、CPLEX)在 MIP 问题上性能远超开源求解器,但需要许可证。学术用户通常可以免费申请。
💡 提示: HiGHS 是目前最强的开源线性规划求解器,推荐作为默认选择。对于非线性问题,Ipopt 是首选。
7. 灵敏度分析
using JuMP, HiGHS
model = Model(HiGHS.Optimizer)
@variable(model, 0 <= x <= 10)
@variable(model, 0 <= y <= 10)
@objective(model, Max, 5x + 4y)
@constraint(model, con1, 6x + 4y <= 24)
@constraint(model, con2, x + 2y <= 6)
optimize!(model)
# 影子价格(对偶变量)
println("con1 影子价格: $(shadow_price(con1))")
println("con2 影子价格: $(shadow_price(con2))")
# 约束的右端项灵敏度范围
println("con1 RHS: $(normalized_rhs(con1))")
println("x 上界: $(upper_bound(x))")
# reduced cost(简约成本)
for (var, rc) in zip([x, y], [reduced_cost(x), reduced_cost(y)])
println("$(name(var)) 的 reduced cost: $rc")
end
8. 多目标优化
using JuMP, HiGHS
# 加权和方法
model = Model(HiGHS.Optimizer)
@variable(model, x[1:3] >= 0)
# 两个冲突目标:最小化成本 vs 最大化质量
w = 0.6 # 权重
@objective(model, Min,
w * (10x[1] + 15x[2] + 20x[3]) + # 成本
(1 - w) * -(5x[1] + 8x[2] + 12x[3]) # 质量(取负)
)
@constraint(model, sum(x) <= 100)
optimize!(model)
9. 建模技巧
9.1 线性化
# 非线性项 xy(其中 y 是二进制变量)的线性化
# xy 可以替换为 z,增加约束:
# z <= x, z <= M*y, z >= x - M*(1-y), z >= 0
using JuMP, HiGHS
M = 100 # big-M 常数
model = Model(HiGHS.Optimizer)
@variable(model, 0 <= x <= 50)
@variable(model, y, Bin)
@variable(model, z >= 0)
@constraint(model, z <= x)
@constraint(model, z <= M * y)
@constraint(model, z >= x - M * (1 - y))
@objective(model, Max, z)
optimize!(model)
9.2 分段线性函数
# 用 SOS2 约束建模分段线性函数
using JuMP, HiGHS
bp = [0.0, 10.0, 25.0, 50.0]
val = [0.0, 8.0, 20.0, 35.0]
model = Model(HiGHS.Optimizer)
@variable(model, x >= 0)
@variable(model, λ[1:length(bp)] >= 0)
@constraint(model, sum(λ) == 1)
@constraint(model, x == dot(λ, bp))
@constraint(model, λ in MOI.SOS2(collect(1:length(bp))))
@objective(model, Max, dot(λ, val))
optimize!(model)
10. 实际案例
10.1 运输问题
using JuMP, HiGHS
# 运输问题:最小化从仓库到客户的运输成本
n_supply = 3 # 仓库数
n_demand = 4 # 客户数
supply = [100, 150, 200] # 各仓库供应量
demand = [80, 120, 100, 150] # 各客户需求量
cost = [2 3 4 5; # 运输成本矩阵
3 2 5 1;
4 3 2 6]
# 平衡供需
@assert sum(supply) >= sum(demand)
model = Model(HiGHS.Optimizer)
@variable(model, x[1:n_supply, 1:n_demand] >= 0)
@objective(model, Min, sum(cost[i,j] * x[i,j] for i in 1:n_supply, j in 1:n_demand))
# 供应约束
@constraint(model, [i in 1:n_supply], sum(x[i,j] for j in 1:n_demand) <= supply[i])
# 需求约束
@constraint(model, [j in 1:n_demand], sum(x[i,j] for i in 1:n_supply) >= demand[j])
optimize!(model)
println("运输方案:")
for i in 1:n_supply
for j in 1:n_demand
if value(x[i,j]) > 0.01
println(" 仓库 $i → 客户 $j: $(value(x[i,j]))")
end
end
end
println("总运输成本: $(objective_value(model))")
10.2 排班问题
using JuMP, HiGHS
# 员工排班:7 天,每天需要不同人数
n_days = 7
n_employees = 15
required = [5, 4, 6, 7, 8, 6, 3] # 每天需要最少人数
max_days = 5 # 每人每周最多工作 5 天
model = Model(HiGHS.Optimizer)
@variable(model, x[1:n_employees, 1:n_days], Bin)
# 目标:最小化总工作天数
@objective(model, Min, sum(x))
# 每天人手约束
@constraint(model, [d in 1:n_days],
sum(x[e, d] for e in 1:n_employees) >= required[d])
# 每人最多工作天数
@constraint(model, [e in 1:n_employees],
sum(x[e, d] for d in 1:n_days) <= max_days)
# 至少连续休息 1 天
@constraint(model, [e in 1:n_employees, d in 1:n_days-1],
x[e, d] + x[e, d+1] <= 1) # 不能连续两天工作
optimize!(model)
println("排班表:")
for e in 1:n_employees
schedule = [Int(value(x[e,d])) for d in 1:n_days]
println(" 员工 $e: $schedule")
end
10.3 投资组合优化
using JuMP, HiGHS, LinearAlgebra
# Markowitz 均值-方差投资组合
μ = [0.08, 0.12, 0.10, 0.06, 0.15] # 预期收益率
Σ = [0.0064 0.0008 0.0011 0.0002 0.0018; # 协方差矩阵
0.0008 0.0144 0.0036 0.0012 0.0054;
0.0011 0.0036 0.0100 0.0018 0.0045;
0.0002 0.0012 0.0018 0.0049 0.0010;
0.0018 0.0054 0.0045 0.0010 0.0225]
model = Model(HiGHS.Optimizer)
@variable(model, w[1:5] >= 0) # 权重
@objective(model, Min, sum(Σ[i,j]*w[i]*w[j] for i in 1:5, j in 1:5))
@constraint(model, sum(w) == 1)
@constraint(model, dot(μ, w) >= 0.10) # 最低收益 10%
optimize!(model)
println("最优权重: $(round.(value.(w)*100, digits=1))%")
println("预期收益: $(round(dot(μ, value.(w))*100, digits=2))%")
扩展阅读
- JuMP.jl 官方文档
- JuMP 教程集
- MathOptInterface — 求解器接口
- Optimization.jl — 另一种优化建模接口
- Bertsimas & Tsitsiklis, Introduction to Linear Optimization
- NEOS Server — 在线优化求解
- JuliaOpt 求解器列表