资讯动态

两阶段鲁棒优化在微网容量配置中的Matlab实现与CCG求解

发布时间:2026/10/16 5:06:58 来源:尧图企业网站定制
两阶段鲁棒优化在微网容量配置里这几年是真的火但也是真的让人头大。我在做微网规划项目的时候一开始用确定性方法算出来的容量方案拿到实际运行数据里一测不是光伏出力波动稍微大一点就出现切负荷就是柴油机组被频繁启停折腾得寿命骤降。后来把两阶段鲁棒优化框架搭起来才算是把容量决策和运行调度这两层问题真正捏合到了一起。这篇就围绕我的Matlab实现过程把从建模、求解到代码落地的一系列关键细节讲清楚。项目主场景是含光伏、风电、柴油发电机和储能电池的并网型微网目标是给出既能应对最恶劣风光出力场景、又不会过度浪费投资的电源容量配置方案。如果你想复现一个可运行的鲁棒容量配置Demo或者正在纠结怎么把CCG算法落实到代码里这篇应该能省下你不少调研时间。1. 微网容量配置为什么要引入两阶段鲁棒优化很多刚接触这个领域的人会有个疑问传统的确定性容量配置算得好好的为什么非要用鲁棒优化我在项目里体会最深的一点是微网容量配置本质上是个决策在前、验证在后的问题。你要在规划阶段把光伏板、风机、柴油机、储能的容量定下来但真正运行的时候光照、风速、负荷都是不确定的。如果只用某几个典型场景的期望值来做决策一旦实际风光出力落在比较差的区间容量配置就可能出现系统性缺电。1.1 风光出力的不确定性与容量配置的本质困境微网里的光伏出力和风电出力从来都不是恒定值。某天早上太阳好光伏能做到额定容量的90%以上到了下午云层一遮瞬时出力能掉到20%。风速更是毫无规律的我给项目做数据统计的时候某地冬季风电出力的标准差能到额定容量的30%以上。这种情况下如果你在规划模型里把风光出力当成确定参数算出来的容量配置一定是顾头不顾尾的。容量配置还有一个麻烦是双层时间尺度耦合。投资决策是年尺度甚至二十年尺度的运行调度却是分钟级到小时级的。你不能只算设备装多少还得验证装完之后能不能在每个典型时刻满足功率平衡和备用要求。这就意味着容量配置本质上是一个带运行约束的优化问题计算复杂度比单纯的设备选型高很多。确定性方法之所以不够用是因为它把不确定参数全部用期望值或者典型场景替代了本质上是拿一个平均情况去代表所有可能出现的运行场景自然扛不住恶劣场景的冲击。1.2 随机优化与鲁棒优化一个硬币的两面处理不确定性的主流方法有两类。随机规划需要假设不确定参数的概率分布已知比如光伏出力的预测误差服从正态分布然后用场景树去做期望值优化。它的好处是结果经济性通常比较好不会为极端场景付出过高的容量成本但坏处是你必须对概率分布有足够置信度而且场景数量一多求解规模会爆炸。鲁棒优化走的是另一条路。我不需要知道光伏出力波动的精确概率分布只需要定义它可能波动的范围然后在一个不确定集合里面搜索最坏情况。最坏情况下系统也能安全运行那么其他场景大概率也能扛住。容量配置这种决策问题设备的安装容量一旦买入就是沉没成本你想反悔都不行所以用保守的鲁棒策略反而有天然的适配性。我在项目里对比过随机规划和鲁棒优化的结果随机规划给出的光伏容量会比鲁棒方案低10%左右总成本也低一些但一旦把实际最差出力场景套进去随机规划方案会出现最大20%的甩负荷量。对微网来说供电可靠性往往比一两年的成本差额重要得多鲁棒优化的价值就在这里。1.3 两阶段决策结构容量先行调度随动两阶段鲁棒优化的两阶段不是指时间上的两个步骤而是指决策变量的两个层次。第一阶段变量是现在就要定下来的决策对应到容量配置里就是光伏装机容量、风机装机容量、储能额定容量和功率、柴油机台数。这些变量一旦确定在后续所有场景里都不能变。第二阶段变量是看到不确定性实现之后才能做的决策对应的是每个运行时刻的机组出力、储能充放电功率、弃风弃光量、切负荷量等运行调度量。用数学语言来说就是min投资成本 最坏情况下的运行成本第一阶段决定容量然后自然就产生一个问题容量定了之后在最坏的风光出力场景下运行成本是多少、会不会出现失负荷。这就是内层的max-min问题。外层求容量最优内层搜最坏场景并给出相应调度成本两者嵌套在一起构成了一个典型的两阶段鲁棒优化模型。我在写代码的时候一开始试图把两个阶段合并成一个单层混合整数规划来求解结果变量数直接爆炸根本跑不动。后来老老实实换成CCG迭代求解框架才算是找到了正路。2. 数学模型搭建从投资层到运行层的双层联动模型是整篇代码的地基建模建歪了后面所有求解器调试都是在浪费时间。我这边以典型的光-风-柴-储微网为例把核心变量、目标函数、约束条件和不确定集合逐个拆开讲。你可以直接拿这套模型作为自己项目的骨架根据实际情况增加或删减设备类型。2.1 目标函数投资成本与运行成本的双时间尺度耦合两阶段鲁棒优化的目标函数分两层。第一阶段目标函数是设备投资成本的年化值包含光伏、风机、柴油机组、储能系统四项。每个设备都要考虑单位容量投资成本、等年值系数以及设备的使用寿命。等年值系数就是把一次性投资按照折现率和使用年限换算成每年摊销金额公式是r × (1r)^n / ((1r)^n - 1)其中r是折现率n是设备寿命。我用的是8%折现率光伏和风机按20年寿命算储能按10年柴油机按15年。这套参数也是文献里常用的算出来结果比较有参考价值。第二阶段目标函数是运行成本在给定容量方案和给定风光出力场景下包含柴油机燃料成本、启停成本、储能老化成本、弃风弃光惩罚成本还有切负荷惩罚成本。切负荷惩罚要设得足够高我一般设成10元/kWh以上否则优化器会倾向于用切负荷代替昂贵的燃料成本或储能耗损。两层目标函数通过一个最坏场景连接起来。第二阶段运行成本不是某个固定值而是在不确定性集合内所有可能场景中的最大值。整体形式如下min投资年化成本 max运行成本这个max内层是模型里最核心的部分也是求解难度最大的来源。2.2 约束条件功率平衡、储能电量和机组出力约束条件同样分两个阶段。第一阶段的约束主要是设备容量变量的非负约束和上限约束比如安装面积决定了光伏容量上限变压器容量决定了总装机上限。这些比较简单不多说。第二阶段的运行约束就细了。首先每个时刻都要满足功率平衡光伏出力 风电出力 柴油机出力 储能放电功率 外网购电 负荷 储能充电功率 弃风弃光量注意这里我特意加上了弃风弃光量和切负荷量这两个松弛变量。没有它们模型容易无解因为当风光出力特别高而负荷很低时系统没法把多余的电量完全吸收必须允许弃掉一部分反过来当风光出力特别低时如果所有机组都开满还是不够就需要切负荷来保证功率平衡。松弛变量进入目标函数时带惩罚系数这样优化器不会滥用它们。储能约束是整个模型里最需要小心的地方。储能SOC荷电状态的时序递推关系是SOC(t1) SOC(t) η_c × P_ch(t) - P_dis(t) / η_d其中η_c和η_d分别是充放电效率。充放电功率不能同时为正这个我会通过引入二进制变量来控制也可以直接把充放电功率设成双向变量然后用它的正负和效率函数来处理。柴油机组需要建最小启停时间约束吗理论上需要但实际建模时如果调度周期是24小时引入最小启停时间会带来大量二进制变量大大增加求解难度。我在项目里折中了一下只加了机组出力上下限约束和爬坡约束没有加最小启停时间约束这样既保留了运行调度的关键物理特性又不至于让模型彻底丧失可求解性。2.3 不确定集合的建模盒式、椭球和预算约束两阶段鲁棒优化的核心在于不确定集合的构建。最常见的集合形式是盒式集合就是每个不确定参数各自在某个区间内独立波动。比如光伏出力P_pv可以在预测值上下偏差±ΔP_pv_max的范围内变动风电也是类似。盒式集合最简单但有一个问题它允许所有不确定参数同时达到最坏取值。现实中不会出现所有光伏组件同时出力最低又所有风机同时停转的情况更合理的假设是最多只有一部分不确定参数同时偏离到极端值。这就需要引入预算约束Σ |ΔP_i| / ΔP_i_max ≤ Γ其中Γ是鲁棒控制参数它限制了所有不确定参数总偏离程度的上限。Γ取0就是完全确定性问题Γ取最大允许值就是盒式集合的最保守情况Γ取中间值则可以在保守性和经济性之间调节。我实际用的就是这个带预算约束的盒式集合具体是在每个时段内限制风光出力偏差的累计程度。也有文献用椭球集合它更适合不确定参数之间存在相关性的场景但椭球集合的求解需要二阶锥约束在某些求解器组合下处理起来比较麻烦。做容量配置这种偏工程的项目带预算约束的盒式集合基本够用了而且对Matlab代码实现也更友好。2.4 为什么预算约束能控制保守程度很多刚上手的人对Γ的理解比较抽象。我用一个生活化的类比说明假设你要出门旅行正常情况下路上会有交通拥堵、天气延误、出租车排队这些小意外。如果假设所有最坏情况同时发生你会准备一条几乎不可能用完的巨额应急资金这就是盒式集合的过度保守。但如果假设最多只允许两个坏情况同时发生你准备的应急资金就会合理得多这就是预算约束的价值。在微网容量配置里Γ控制的正是最坏场景有多坏。Γ越小你准备应对的最坏场景越温和算出来的容量成本越低但实际碰到更恶劣场景时的失负荷概率也越高。我在项目里会做一组Γ扫描实验比如从1到最大偏差量观察容量配置结果和总成本的变化曲线然后选择一个成本增长斜率突变的点作为推荐值。我自己的经验是Γ取预算最大值的40%到60%通常是比较合理的区间这个范围内能在可靠性保证和投资成本之间找到不错的性价比。当然不同地区的光照和风速相关性不一样建议你拿到数据后先做统计看看本地区风光出力同时波动的实际规律再定。3. 求解框架列与约束生成算法CCG的原理与实现模型建完之后真正的技术难点在求解。两阶段鲁棒优化是个min-max结构不能直接用现成的优化求解器求解必须用迭代算法。我在项目里用的是列与约束生成算法也就是CCG它比传统的Benders分解在这个问题上收敛快得多代码实现也相对清晰。这一节是我自己从零搭完整套求解流程之后总结出来的能帮你少走不少弯路。3.1 主问题MP与子问题SP的博弈结构CCG的核心思想很简单把原问题拆成主问题和子问题两者互相博弈交替求解直到收敛。主问题MP是一个有限场景的优化问题。原始的不确定集合包含无穷多个可能的场景不可能全部列出来但主问题里只保留已经被识别出的有限个关键场景求一个在这些场景下都可行且成本最低的容量配置和调度方案。主问题求解完得到容量方案x和总成本下界LB。子问题SP的输入就是主问题算出来的容量方案x。子问题要在整个不确定集合内搜索一个最坏场景让给定x下的运行成本最高同时给出这个最坏场景下的最优调度方案。子问题求解完得到最坏场景和对应的运行成本把它和当前容量方案的投资成本相加得到总成本上界UB。如果UB和LB的差距缩小到预设的容差内说明当前容量方案在最坏场景下的成本和主问题内场景的成本已经基本一致迭代就可以终止了。如果差距还大就把子问题找到的最坏场景作为新场景加到主问题里重新求解如此反复。这个过程有点像一个谈判主问题先出一个容量方案子问题立刻指出你这个方案在我手上这个最坏场景里成本很高主问题收到反馈后修正方案子问题再找更坏的场景直到双方达成一致。3.2 子问题的最坏场景搜索对偶、KKT与Big-M子问题的结构是max-min外层要在不确定集合内找最坏的风光出力场景内层要在最坏场景下求最小运行成本。这不是一个可以直接求解的标准数学规划问题需要处理。最常用的处理手法是把内层min问题通过拉格朗日对偶转成max形式然后和外层max合并变成一个单层max问题。这个单层问题里目标函数和约束条件中会出现不确定变量u和调度变量y的乘积项。按照标准鲁棒对偶理论这类乘积项可以通过引入对偶变量来处理最终得到的是一个可以在求解器里求解的优化问题。但这里有个坎如果内层问题存在整数变量比如储能充放电状态的二进制变量那么对偶理论就不直接适用了。不是不能用KKT条件或Big-M线性化来处理但处理难度会上升。我的做法是尽量让内层问题保持为线性规划可解的形式具体办法是储能充放电状态不用二进制变量区分而是把充放电功率拆成两个非负变量并用惩罚项和效率约束来隐式约束不同时充放。这样内层问题保持线性对偶转换就干净了。如果必须保留整数变量怎么办可以把内层问题用KKT条件替换再通过Big-M把互补松弛条件线性化得到一个混合整数规划子问题。这种方法可行但Big-M参数需要仔细选择否则数值稳定性会很差。我在测试中发现M值设成负荷最大值的三到五倍通常能兼顾收敛性和数值稳定。3.3 迭代收敛判定CCG的收敛判定很简单就是持续比较上下界。初始状态下LB可以设成负无穷UB设成正无穷。每轮迭代主问题求得LB子问题求得UB。当(UB-LB)/UB小于预设阈值比如0.01或者0.001就认为收敛。阈值设太大会导致容量方案偏离真实最坏场景太远设太小则迭代次数暴增求解时间从几分钟变成几小时都是可能的。我在代码里用的是相对间隙1e-3配合一个迭代次数上限30次。实测下来对带预测不确定性的容量配置问题一般15到25轮之间能收敛。如果超过30次还没收敛说明模型某个地方大概率有数值问题而不是单纯收敛慢。这时候最有效的排查方式是检查子问题对偶变量的尺度是否相差过大以及Big-M参数是否设得离谱。另外有个小技巧主问题每一轮迭代都会积累新的场景和变量场景变量越来越多主问题规模越来越大。所以必须在主问题里把每一轮迭代生成的最坏场景都保留下来不能只保留最新一轮的。这不是写代码时容易漏掉的地方而且一开始使用MPCC工具时特别容易误以为保留最新场景就够了。3.4 CCG与Benders分解的取舍Benders分解和CCG都能求解两阶段鲁棒问题但我在实际对比中明显感觉到CCG在容量配置这个场景下更优。原因在于CCG每次迭代向主问题添加的是完整的场景变量和约束主问题的约束越来越贴近真实的可行域收敛速度显著更快Benders分解添加的是割平面近似收敛比较慢尤其当最优解落在可行域内部的时候会震荡。Benders分解的一个优势是主问题规模增长速度慢内存占用较小。但对于微网容量配置这种中小规模优化问题CCG增加的那点规模根本无所谓。我做过一个含24时段、5种设备类型的测试算例CCG大概25轮收敛Benders跑了60轮还在震荡结论非常明显。如果你用的是成熟的鲁棒优化工具包可能不需要自己实现CCG。但我还是强烈建议你至少手动实现一遍CCG流程因为只有拆开理解主问题和子问题的博弈关系才能在处理模型的特殊约束时知道该改哪里。4. Matlab实现的关键代码逻辑与工程细节Matlab生态下的组合是YALMIPGurobi或者YALMIPCplex。把CCG算法用YALMIP表达出来好处是可以用数学表达式直接写约束可读性高调试时容易定位问题。坏处是YALMIP本身是建模层的解释器如果模型规模大它的符号运算开销会拖慢整体效率。好在容量配置问题规模不大YALMIP完全够用。4.1 工具链选型YALMIP求解器的搭配逻辑YALMIP负责把优化模型翻译成求解器能理解的标准形式。线性规划、二次规划、混合整数规划都有对应的求解器支持。在我的项目里主问题是一个混合整数线性规划因为储能充放电功率拆成两个变量后不再需要二进制变量但柴油机启停和容量离散选择会产生0-1变量。子问题经过对偶转换后是一个线性规划或者带二阶锥约束的问题。选择具体求解器时Gurobi和Cplex的线性规划、混合整数规划性能都非常强在中小规模问题上差距不大。如果学生党没有商业求解器授权用开源的Cbc、GLPK也能跑通流程只是速度慢一些。我在机器上装的是Gurobi授权的事就不多提了你自己根据条件选择就行。另外需要提醒的是YALMIP版本和求解器版本的兼容性经常出问题。我遇到过YALMIP调用求解器返回空解的情况最后排查下来是YALMIP版本太旧不识别求解器的新参数。建议用较新的YALMIP版本并且求解器也保持更新。4.2 代码骨架与数据结构整个项目的代码结构我分成三块参数初始化脚本、CCG主循环函数、结果后处理脚本。参数初始化脚本里定义所有系统参数包括各类设备的成本参数、效率参数、预测出力曲线、负荷曲线、不确定集合的偏差范围和预算Γ、折现率、寿命年限等。一些重要的心得是把参数集中放到一个结构体里以sysParam命名后面每个函数都传入这个结构体比全局变量好用得多。数据这块最费时间的是准备预测出力曲线。我用的方法是把历史光伏出力和风电出力按小时取平均得到基准预测曲线然后根据历史数据的标准差确定每个时段的不确定偏差范围。如果你手头没有历史数据可以用合成方法生成比如用正弦曲线拟合日照规律再叠加上下偏差但这样出来的结果会偏理想不一定能反映真实波动特征。CCG主循环函数是核心。函数先初始化UB和LB然后循环调用主问题构造函数和子问题构造函数。主问题构造函数接收已经发现的所有最坏场景集合用YALMIP定义容量变量和对应各场景的运行变量求解后返回容量值、目标函数值和割平面信息。子问题构造函数接收主问题给的容量值构建对偶转换后的最坏场景搜索模型求解后返回新的最坏场景和对应运行成本。结果后处理脚本把每个设备的最优容量、总成本、各场景的调度方案等结果输出成表格和曲线方便比对。4.3 核心代码段演示我先展示主问题的核心代码逻辑省略设备参数定义部分只保留建模核心部分function [x_opt, obj_mp, mp_model] build_MP(scenarios, sysParam) % 定义第一阶段容量变量 x_pv sdpvar(1, 1, full); x_wt sdpvar(1, 1, full); x_ess sdpvar(1, 1, full); x_dg sdpvar(1, 1, full); % 投资成本年化值 inv_cost sysParam.cost_pv * x_pv * sysParam.crf_pv ... sysParam.cost_wt * x_wt * sysParam.crf_wt ... sysParam.cost_ess * x_ess * sysParam.crf_ess ... sysParam.cost_dg * x_dg * sysParam.crf_dg; % 循环每个场景定义运行变量和约束 run_cost 0; constraints []; for k 1:length(scenarios) % 场景k的风光出力和负荷数据 p_pv_k scenarios(k).p_pv; p_wt_k scenarios(k).p_wt; p_load_k scenarios(k).p_load; % 第二阶段运行变量 p_dg sdpvar(24, 1, full); p_ess_ch sdpvar(24, 1, full); p_ess_dis sdpvar(24, 1, full); p_curtail sdpvar(24, 1, full); p_shed sdpvar(24, 1, full); soc sdpvar(25, 1, full); u_dg binvar(24, 1); % 柴油机启停状态 constraints [constraints, ... 0 p_dg sysParam.p_dg_max * x_dg, ... 0 p_ess_ch sysParam.p_ess_rated * x_ess, ... 0 p_ess_dis sysParam.p_ess_rated * x_ess, ... 0 p_curtail p_pv_k p_wt_k, ... 0 p_shed p_load_k, ... soc(1) sysParam.soc_init, ... soc(25) sysParam.soc_init, ... soc(2:25) soc(1:24) sysParam.eta_ch * p_ess_ch - p_ess_dis / sysParam.eta_dis, ... sysParam.soc_min soc(2:25) sysParam.soc_max, ... p_pv_k p_wt_k p_dg p_ess_dis p_shed p_curtail p_load_k p_ess_ch]; % 运行成本 run_cost run_cost sum(sysParam.fuel_cost * p_dg) ... sum(sysParam.start_cost * max(0, diff([0; u_dg]))) ... sysParam.penalty_shed * sum(p_shed) ... sysParam.penalty_curtail * sum(p_curtail) ... sysParam.ess_degrad_cost * sum(p_ess_ch p_ess_dis); end % 求解主问题 optimize(constraints, inv_cost run_cost, sdpsettings(solver, gurobi, verbose, 0)); x_opt value([x_pv, x_wt, x_ess, x_dg]); obj_mp value(inv_cost run_cost); end注意这里max(0, diff(...))这种表达在YALMIP里不是标准线性形式实际项目里我会引入辅助变量来表示柴油机启停成本这里写成简单形式是为了展示逻辑。实际代码里应该用线性化约束z_t u_t - u_{t-1}; z_t 0来表示启停动作。4.4 子问题的对偶转换代码实现子问题实现比主问题麻烦需要对内层调度问题做对偶转换。但YALMIP有一个可以偷懒的路径用YALMIP的robust框架或者uncertain声明来直接建模鲁棒优化。我在测试初期试过这个功能不过在带预算约束的不确定集合下YALMIP自带的鲁棒建模功能生成的对偶形式经常不是我想要的而且耦合了自己写的双层结构时容易报错。所以最终我选择手动写对偶转换。做法是先定义内层调度问题的原始变量和对偶变量然后用符号推导写出对偶问题的约束和目标函数再与外层max合并。这个过程最容易出错我建议你在测试小规模问题时先用一个3时段的简化模型验证对偶是否正确再扩展到24时段。对偶转换完成后子问题的目标函数会呈现类似下面的结构function [u_worst, obj_sp, sp_model] build_SP(x_opt, sysParam) % 不确定变量 delta_pv sdpvar(24, 1, full); delta_wt sdpvar(24, 1, full); % 不确定范围约束 constraints [... -sysParam.delta_pv_max delta_pv sysParam.delta_pv_max, ... -sysParam.delta_wt_max delta_wt sysParam.delta_wt_max]; % 预算约束 constraints [constraints, ... sum(abs(delta_pv) ./ sysParam.delta_pv_max) ... sum(abs(delta_wt) ./ sysParam.delta_wt_max) sysParam.Gamma]; % 对偶变量和线性化后子问题目标 % 具体形式取决于原始内层约束的系数矩阵需要根据你的模型推导 % 这里用简化示意 obj_sp ...; % 对偶目标函数包含x_opt和delta的耦合项 optimize(constraints, -obj_sp, sdpsettings(solver, gurobi, verbose, 0)); u_worst.p_pv sysParam.p_pv_base value(delta_pv); u_worst.p_wt sysParam.p_wt_base value(delta_wt); obj_sp value(obj_sp); end绝对值项在YALMIP里处理需要额外变量实际项目中我用的是引入正负偏差变量的方式把delta_pv拆成delta_pv_pos和delta_pv_neg然后用它们的和来构建预算约束。这样整体保持线性约束的表达形式。4.5 求解器参数和数值预处理YALMIP求解时的sdpsettings参数里有几个值得调。一个是数值容差Gurobi默认的绝对间隙和相对间隙设置比较保守我一般会把相对间隙设成1e-4或者更高避免主问题和子问题在边界上反复震荡。另一个是线程数如果模型规模不大设成单线程反而更稳定多线程在小模型上的收益不明显还容易引入结果不确定性。数值预处理方面建议把所有成本系数和功率量纲统一。我习惯把所有功率用kW电量用kWh成本用元这样数量级差异不会太大。如果有些系数相差超过1e4倍求解器的数值稳定性会明显下降出现莫名其妙的不可行或者无界警告。检查数值尺度是否合理是调试CCG循环的第一步。5. 实测中的坑与调参经验这部分是最干货的部分都是我在代码跑通和项目复现过程中真正踩过的坑。如果你是自己动手写代码大概率会撞上其中几个。5.1 双线性项的处理别偷懒CCG子问题里不确定变量u和对偶变量λ会以乘积形式出现在目标函数中这是max-min结构带来的天然产物。YALMIP可以直接建立两个sdpvar变量的乘积但这样得到的模型是双线性的Gurobi根本解不了。我在处理这个双线性项时的经验是如果是连续变量和连续变量的乘积根据不确定变量的边界值可以通过大M法或者McCormick包络线性近似或者直接利用盒子约束的对偶性质消去变量如果是二进制变量和连续变量的乘积则要先引入中间变量再用大M法线性化。这里有个容易踩的坑McCormick包络是松弛近似得到的最坏场景不一定是最优最坏场景保守性会受影响。所以能通过严格对偶变换直接消去双线性项的情况下尽量用严格方法。5.2 迭代不收敛的排查路线如果CCG迭代超过40轮还在上下界之间震荡不要盲目加大迭代上限这是一条错误的排查路径。更有效的做法是检查以下几点第一检查子问题是否真的求到了最坏场景。把子问题算出来的风光场景和目标函数值记录下来人工验证一下在这个场景下主问题给出的容量方案是否真的会产生这么高的运行成本。如果手工验算成本比子问题给出的低说明子问题的对偶转换或者双线性项处理错了。第二检查主问题是否保留了所有历史场景。我在初版代码里犯过这个错主问题只保留最新场景结果每次迭代都在不同的场景之间跳来跳去上下界永远收敛不到一起。把场景保留改成累积式之后大约第八轮就开始稳步收敛了。第三检查数值容差。上下界差距缩小到一定程度就不动了可能是求解器精度限制。试着把求解器的相对间隙收紧一点或者把目标函数系数统一缩放1000倍排除数值尺度问题。5.3 保守性-经济性平衡的调参Γ的值直接决定了容量配置结果的保守程度。我在项目里做了个实验把Γ从0逐步加到最大预算值观察光伏、风机、储能容量的变化。在没有预算约束时即Γ为0系统完全不考虑不确定性光伏容量大约是600kW储能配置也会比较小总成本最低。随着Γ逐步增大系统为了应对更恶劣的风光场景光伏和风机容量反而可能减小因为不可靠的电源容量在极端场景下形同虚设柴油机和储能的容量则稳步上升。在Γ达到最大值时储能容量可能比Γ为0时高出50%以上柴油机容量也显著增加。这个趋势说明如果某个地区的风光资源波动特别大单纯加大光伏和风机容量并不能提升供电可靠性反而应该配置灵活调节资源如储能和柴油机。在做项目汇报的时候我会把这张Γ-容量曲线直接放进去决策者看到之后更容易理解保守程度对投资的影响。5.4 24时段扩展到8760小时的加速思路容量配置严格来说应该用全年8760小时的时序数据做校核但两阶段鲁棒优化的子问题本身是24时段直接扩展到8760会让变量规模爆炸。我试过直接扩结果主问题单轮求解就要半小时以上迭代根本走不动。工程上我最后采用的方法是典型日法加季节性分解。把一年分成四个季节每个季节选一个代表日然后对四个典型日分别做容量配置取各设备容量的包络值作为最终结果。这样虽然会牺牲一点精度但计算时间能控制在十几分钟内而且结果的可解释性更强。更进阶的做法是用Benders分解和CCG的嵌套结构把8760小时的运行调度问题拆到多个时段块分别求解。这个方向我在后续迭代里验证过速度和精度都比典型日法好但代码复杂度明显上升。如果只是做一个初步的可行性研究典型日法已经够用了。5.5 关于Matlab代码复用和扩展的几个建议代码结构上强烈建议把CCG主循环、主问题构造、子问题构造、参数设置、结果输出拆成独立函数每个函数用清晰注释说明输入输出。我一开始图省事写了个巨型脚本调试时改一个参数要拖半天后来重构成函数化结构后加新设备、改不确定集合都很方便。不确定参数的波动范围设置要基于实际数据。如果项目手头暂时没有历史数据至少也要参考行业报告或仿真数据不要随意设一个对称区间那样算出来的配置方案缺乏实际参考价值。如果你想在这个基础上加新的不确定性来源比如负荷不确定性或者电价不确定性只需要在不确定变量列表里增加对应变量在预算约束里加入对应的偏离项即可整体框架不用大改。CCG算法的兼容性很好这也是它在工程界受欢迎的原因之一。最后提一个我踩过最多次的坑每次迭代开始前记得清空YALMIP的变量缓存和求解器的工作区否则旧模型的变量约束会残留造成模型规模膨胀和求解偏差。YALMIP里sdpvar对象重复定义不会直接报错你以为是在更新模型实际上是在叠加约束这个问题排查起来相当隐蔽。

读完文章,也想定制专属网站?

尧图设计师 24 小时内与您沟通定制方案

免费获取报价 →
↑