1. 项目概述与核心问题拆解
做电力市场优化方向的同学们,看到这个标题的第一反应应该和我一样——这又是一个典型的“双层决策 + 风险度量 + 强对偶转化”的组合问题。先说结论:这个项目本质上解决的是省间交易商在“省间现货市场 + 省内市场”两级环境下,如何制定最优购电策略,同时把市场价格波动带来的风险纳入决策框架的问题。
先说清楚这个模型的应用场景。我国电力市场正在从省级市场向省间市场延伸,省间交易商(可以理解为跨省购电的售电公司或大用户)面临的不再是单一市场,而是“先在省间市场买电,再在省内市场卖出或部分转售”的双层结构。这里有两个核心难点:一是两级市场的价格形成机制不同、时序不同、约束不同;二是交易商的购电决策会影响省间市场的出清结果,这种“决策影响出清、出清反过来约束决策”的循环关系,用数学语言描述就是一个双层优化问题。
项目的应用价值很实在:给交易商提供一个科学的购电决策工具,不是拍脑袋买电,而是基于市场机制、风险偏好和价格预测,求解一个最优的购电组合。适合的读者群体包括电力经济方向的研究生、电网公司交易中心的策略研究人员,以及做电力市场软件开发的技术人员。
我拿到这个题目后首先梳理了几个关键词:MATLAB意味着整个模型用MATLAB编程实现;Cplex是求解器,说明问题最终要转化成可求解的数学规划模型;强对偶是转化的核心工具,用来把双层问题压平为单层问题。这个技术路线是电力市场优化建模里的经典套路,下文我会把每一环都拆开讲透。
2. 两级市场框架与风险建模:先把场景立起来
2.1 “两级市场”到底指什么
省间和省内不是简单的时间先后关系,而是空间维度和时间维度的双重嵌套。省间现货市场一般是日前或日内滚动出清,交易标的是跨省跨区的电能量;省内市场则是各省自己组织的现货或中长期市场。交易商在省间市场中标后,拿到的电量和电价,成为省内市场中投标或售电的成本基础。
用生活化的类比来理解——你把北京的水果运到上海卖,你在北京批发市场的采购价和采购量,会影响北京市场当天的批发行情;而你到了上海,面对的是另一个供需决定的价格体系。你作为中间商,赚的是两个市场的价差,但同时也承担着两边价格波动的双重风险。
建模时,两级市场的耦合点在于:省间购电成本进入省内市场的成本项,而省内市场的需求响应约束会反过来限制省间购电量。常见的处理方式是建立两个出清模型,其中省内出清作为下层问题,省间购电决策作为上层变量。
2.2 风险从哪来,怎么度量
这个项目标题里特别强调“计及风险”。这里的风险主要是价格不确定性——省间市场价格受新能源出力波动、负荷预测偏差、通道阻塞等因素影响,可能在一天内出现剧烈波动。如果不把风险纳入决策,最优解可能是一个“期望收益最大但方差极大”的方案,实际执行中可能大幅亏损。
风险建模的主流选择是条件风险价值(CVaR,Conditional Value at Risk),也叫平均超额损失。相比传统的方差或VaR(风险价值,Value at Risk),CVaR具有次可加性、凸性,适合嵌入优化模型。直观解释:CVaR回答的是“在最差的5%情景下,平均损失是多少”,它对尾部风险的刻画更准确。
数学形式上,假设有S个价格场景,每个场景的概率为π_s,决策对应的总成本为f_s(x),置信水平为β(比如95%),CVaR的线性化形式是:
CVaR_β = ζ + (1/(1-β)) * Σ_s π_s * max(f_s(x) - ζ, 0)
引入辅助变量η_s ≥ f_s(x) - ζ,η_s ≥ 0,就能转化为线性约束。这一步也是Cplex最擅长的处理范围——线性规划加整数变量基本都能高效求解。
2.3 目标函数怎么搭
交易商的优化目标可以拆成两层逻辑。上层目标函数是总购电成本最小化(或综合收益最大化),下层目标函数是省内市场的社会福利最大化或购电成本最小化。合并之后,单层化的目标函数一般写成:
min 期望购电成本 + 风险惩罚项
展开来看,期望购电成本 = Σ_s π_s * (省间购电费用 + 省内交易费用),风险项就是上面提到的CVaR乘以权重系数λ。λ的大小代表交易商的风险厌恶程度,λ越大,方案越保守,越偏向于在价格波动大的时段减少购电或增加价格锁定。
这里有一个容易踩坑的地方:很多初学者直接把期望成本和CVaR相加,但没有统一量纲。价差超过一定范围时,风险项在数值上可能压过期望项,导致模型过度保守。实际处理中建议对风险项做归一化或敏感性分析,先跑几组λ值看解的变化趋势,再确定合理区间。
3. 强对偶理论:双层问题如何压平成单层
3.1 为什么要用强对偶
双层优化在计算上非常棘手,因为下层问题的最优解必须嵌入上层的约束条件。如果下层问题是线性规划且满足强对偶条件,就可以用“原始最优值 = 对偶最优值”这一等式,把下层问题的求解条件转化为上层问题的约束,从而把双层问题等价转化成一个单层的、带互补约束的优化问题(MPEC,Mathematical Program with Equilibrium Constraints)。
可以这样理解强对偶:你把下层市场出清问题看作一个“黑盒”,这个黑盒的输入是上层决策变量(比如购电量),输出是市场价格和出清量。强对偶告诉你,这个黑盒的行为可以用一组“价格信号”和对偶变量来描述,而不需要真的迭代求解下层问题。
使用强对偶的前提是下层问题必须满足强对偶定理的条件,即线性规划可行且有界,或者满足Slater条件的凸优化问题。电力市场出清模型通常是线性规划,天然满足这个要求,这也是为什么这种方法在电力市场文献中占主导地位。
3.2 强对偶转化的标准步骤
完整转化分四步走。
第一步:写下层问题的标准形式。对省内市场出清问题,决策变量是省内的机组出力、购电量、联络线潮流等,目标函数是购电成本或社会福利最大化,约束包括功率平衡、线路潮流约束、机组上下限等。
第二步:写出对偶问题。按拉格朗日函数写出对偶变量,包括功率平衡约束的对偶变量(这就是节点电价)、线路约束的对偶变量(阻塞影子价格)等。
第三步:应用强对偶等式。原始目标函数在最优解处的值等于对偶目标函数值,把两个表达式相等作为一个约束加入到上层模型中。
第四步:加入互补松弛条件或直接依赖强对偶等式。如果加了互补条件,模型会变成MPEC,需要MCP(混合互补问题)求解器或特殊处理;如果只依赖强对偶等式,则有可能丢失解的可行性信息,需要额外验证。
我在实际项目中的经验是:直接加强对偶等式是最稳定的做法,但不加互补条件时可能出现“虽然强对偶成立但原始解不可行”的情况。一个补救措施是额外加入原问题的全约束集,确保最终解在原始可行域内,这样才严格。
3.3 强对偶在代码里长什么样
以省内出清为线性规划为例,假设原问题形式是:
min c'T * y
s.t. Ay = b(功率平衡)
Dy ≤ h(线路约束)
y ≥ 0
对偶变量记为λ(等式约束)、μ(不等式约束,μ≥0),则对偶问题为:
max b'T * λ + h'T * μ
s.t. A'T * λ + D'T * μ ≤ c
μ ≥ 0
强对偶等式是 c'T * y = b'T * λ + h'T * μ。将这个等式作为上层问题的附加约束,同时保留下层问题的全部原始约束和变量,完成压平。
在MATLAB代码中实现时,主要工作是用Cplex的矩阵接口把这些约束填进去。注意对偶变量在Cplex中要显式声明为自由变量和受限变量:λ是自由变量(因为等式约束),μ是非负变量(因为不等式约束小于等于)。
4. MATLAB + Cplex实战:从建模到代码落地
4.1 环境配置:Cplex在MATLAB中的三种调用方式
很多同学卡在第一步——Cplex装好了,但不知道如何在MATLAB里调用。我梳理了三种常用方式,按推荐程度排序。
第一种是用Cplex的MATLAB接口类,这是最底层也最灵活的方式。安装IBM ILOG CPLEX Optimization Studio后,在MATLAB里添加路径即可调用Cplex类,直接构建模型对象。这个方式的缺点是代码量大,每个约束都要写一行添加语句,适合对Cplex本身比较熟悉的用户。
第二种是使用YALMIP工具箱,这是我现在最推荐的方式。YALMIP是一个建模层的工具,用符号化的方式描述变量和约束,最后调用Cplex求解。优点是代码可读性强,修改模型方便,初学者半天就能上手。缺点是引入了一层抽象,求解效率稍有损耗,但对中小规模问题完全够用。
第三种是使用MATLAB自带的优化工具箱配合linprog或intlinprog,求解器选为Cplex——这种方法局限性较大,MATLAB的优化工具箱无法完全发挥Cplex的高级功能,比如冲突分析、剪枝策略调节等。
我个人的习惯是:中期研究阶段用YALMIP做原型验证,确认模型正确后,再考虑是否用Cplex接口重写以提高求解效率。如果项目直接要求“Cplex代码实现”,你就直接用Cplex接口类,这样在报告和论文里更有说服力。
4.2 数据准备:场景生成与参数设置
购电模型需要输入的数据包括:省间现货市场的预测价格曲线、省内市场的负荷预测、机组参数、线路输电容量的关键参数等。价格场景的生成可以用历史数据分析加蒙特卡洛抽样,也可以用正态分布拟合。
在MATLAB中生成场景的一个示例代码框架如下:
% 生成S个价格场景,假设价格服从正态分布 % mu_price是各时段预测价格向量,sigma_price是标准差向量 S = 100; % 场景数 T = 24; % 时段数 price_scenarios = zeros(S, T); for t = 1:T price_scenarios(:, t) = normrnd(mu_price(t), sigma_price(t), S, 1); end % 等概率场景 prob = ones(S, 1) / S;这个场景集之后会直接用于构建期望成本和CVaR项。需要注意:场景数越多,模型规模越大,求解时间越长。建议先用S=50或100做调试,确认能解之后再增加到500甚至1000。如果求解时间不可接受,考虑使用场景缩减技术(比如同步回代消除法),保留有代表性的少量场景。
4.3 YALMIP版本的完整模型框架
下面给出一个基于YALMIP的模型框架,完整展示如何搭建两级市场购电模型。虽然Cplex接口版代码会更长,但逻辑完全一致。
% 决策变量 x_s = sdpvar(T, 1); % 省间购电量 y_g = sdpvar(G, 1); % 省内机组出力(简化示例,实际要按场景定义) % 对偶变量 lambda = sdpvar(T, 1); % 功率平衡对偶变量(电价) mu_lin = sdpvar(L, 1); % 线路约束对偶变量 % 目标函数:期望购电成本 + CVaR风险项 % 先构建每个场景下的购电成本表达式 cost_scn = zeros(S, 1); for s = 1:S cost_scn(s) = price_scenarios(s, :) * x_s ... % 省间购电成本(简化为线性) + c_g' * y_g; % 省内购电成本 end VaR = sdpvar(1, 1); % CVaR辅助变量ζ eta = sdpvar(S, 1); % 超额损失变量 % 期望成本 exp_cost = prob' * cost_scn; % CVaR项 cvar_term = VaR + (1/(1-beta)) * prob' * eta; % 约束条件 Constraints = []; Constraints = [Constraints, A_eq * [x_s; y_g] == b_eq]; % 功率平衡 Constraints = [Constraints, D_ineq * [x_s; y_g] <= h_ineq]; % 线路约束 % CVaR约束 for s = 1:S Constraints = [Constraints, eta(s) >= cost_scn(s) - VaR]; end Constraints = [Constraints, eta >= 0]; % 强对偶等式 Constraints = [Constraints, c' * [x_s; y_g] == b_eq' * lambda + h_ineq' * mu_lin]; % 对偶变量约束 Constraints = [Constraints, mu_lin >= 0]; % 求解 options = sdpsettings('solver', 'cplex', 'verbose', 2); optimize(Constraints, exp_cost + lambda_risk * cvar_term, options);注意这个框架是一个极简示意,实际代码里省内市场的约束和变量必须按场景展开——因为每个场景下省内出清量不一样,如果只用一个y_g表示,模型就会失真。
4.4 使用Cplex接口类的几个关键细节
如果用Cplex接口类直接写,有几个细节特别重要。
第一,变量索引要提前规划好。Cplex使用连续的变量索引,所以要把所有变量(购电量、机组出力、对偶变量、CVaR辅助变量)编排好索引范围,用一个结构体记录每个变量的起止索引。别图省事用动态增长的方式,模型稍微大一点就会卡死。
第二,约束的稀疏结构提前构建。Cplex添加约束时,用稀疏矩阵方式传入性能最佳。在MATLAB中可以用sparse函数构建A矩阵和b向量,然后一次性addRows,不要一行一行add,效率差距在几百倍以上。
第三,求解器参数要调。核心问题是有无整数变量。如果是纯线性规划,默认设置直接求解即可;如果加入了0-1整数变量(比如机组启停变量),要主动设置相对间隙和MIP容忍度。比如设定相对对偶间隙为1e-4,以解决其保守问题。MIP时间限制也要设置,防止意外情况把求解进程挂住。
第四,注意Cplex对不可行问题的报告。模型构建错误时,Cplex会报“infeasible”,这时可以用refineConflict方法做冲突分析,找出是哪组约束导致无解。这个功能在调试阶段是救命神器。
5. 常见问题与排查技巧实录
5.1 模型总是无解?先查双层转化是否丢条件
这是项目里最常遇到的问题,几乎所有人第一次跑通之前都会卡在这里。排查顺序是这样的:先不看强对偶等式,只用上层约束和价格场景跑一遍主问题,看求解器是否能找到可行解。如果主问题无解,问题出在上层约束本身,比如购电量范围、线路容量设置不合理。
如果主问题有解,加上强对偶等式后无解,那问题几乎一定出在对偶变量域或强对偶等式的正负号处理上。我接手过好几个类似的项目,发现一半以上的错误来自对偶变量的符号设置,特别是“≤”约束对应的对偶变量应该是非负,但代码里写成自由变量。
还有一个隐蔽问题:下层问题是最大化问题时,强对偶等式右侧的目标函数方向容易搞反。建议把所有对偶推导在纸上完整做一遍再写代码,不要凭记忆写。
5.2 Cplex求解器报数值问题:精度设置与Big-M处理
当模型中出现互补约束或大系数时,Cplex会警告“数值问题”或“ill-conditioned”,求解结果可能不可靠。这种情况多出现在MPEC模型中,因为互补条件引入了双线性项。
处理手段有三个:第一,如果可以避免互补条件,就只用强对偶等式,模型就变成线性的,整体难度下降;第二,如果必须加互补条件,用Big-M方法线性化,M的选取要尽量小但又够覆盖变量范围。经验值是M取变量上界的1.5到2倍,不要一味求大;第三,调整Cplex的数值容忍度参数,比如将NumericalEmphasis设置为1,让求解器走更稳健的数值路径。
5.3 求解时间过长:场景缩减与分解策略
当模型规模大了,直接求解时间会指数上升,500个场景可能就要跑上几个小时。除了场景缩减,还可以用的是Benders分解或L型算法:把场景相关的子问题从主问题中分解出去,用割平面迭代逼近。
在MATLAB里实现Benders分解不算太复杂:主问题只包含第一阶段变量(购电策略)和强对偶相关约束;子问题是每个场景下的省内出清问题,返回最优值和灵敏度信息,生成可行性割或最优性割。这个方案代码量会增加一倍,但对于大规模实际系统是必须走的路。
另一个实用技巧是修改Cplex的并行线程数和MIP策略。Cplex默认会探测环境决定线程数,但有时它选得并不合理,比如在集群上默认线程数过高。手动设置CPX_PARAM_THREADS为8或者16,往往能在保证稳定性的前提下显著提速。
5.4 MATLAB调用Cplex时的坑:路径与许可证问题
最常见的一个坑:本机安装了多个版本的MATLAB,Cplex的接口只注册到了其中一个版本,导致另外的版本调用失败。解决办法是重新运行Cplex安装目录下的MATLAB接口配置脚本,或者手动addpath指向对应的bin文件夹。
许可证问题也高频出现。Cplex需要许可证文件,如果用的是学生版或学术版,要确认许可证类型支持MATLAB接口。有些扩展功能(比如MPEC求解器)需要额外模块,许可证不包含时,Cplex会报“unsupported feature”错误。做模型调试前先跑一个简单的Cplex测试实例,确认环境完全正常,再开始大数据计算。
6. 我的实操心得与后续扩展建议
这个项目做完一遍,我对“强对偶转化”的理解深入了一个层次。纸上推导是一回事,把对偶变量放进Cplex求解器里是另一回事——你会发现变量数量翻倍、约束数量翻倍、问题的数值特性也变了,一个小小的正负号错误就能让模型从可解变成不可解。我的建议是:建模阶段用YALMIP快速验证正确性,确定模型表达无误后,再用Cplex接口精调性能,两条腿走路最稳。
CVaR部分也有一些经验想分享。风险厌恶系数λ的选取是决策者和模型之间的“翻译器”,不要一上来就拍脑袋定值。我习惯的做法是跑一个λ从0到1的灵敏度扫描,画出“期望成本-CVaR”前沿曲线,让决策者直观看到风险与成本的权衡。这比给一个单一最优解更有说服力,在实际项目汇报中效果也更好。
后续如果想扩展,可以考虑三个方向:一是把省间市场价格的不确定性从场景抽样升级为随机变量的解析分布,用机会约束规划替代CVaR,体系会更严谨;二是把单一时段的静态模型扩展到多时段的动态决策,引入储能或需求响应资源;三是用深度强化学习替换部分求解环节,特别是当场景数量爆炸导致Cplex难以收敛时,离线训练一个策略网络做在线决策,已经是这个方向的前沿热点。
最后分享一个实际教训:切勿在未做灵敏度分析的情况下直接采纳模型的“最优购电策略”。市场参数稍有变化,最优解可能剧烈跳变。给决策者的输出,一定要带上是哪组场景假设、哪组价格预测得到的结论,这样即便市场变了,也能快速重新求解而不是推翻整套体系。