企业官网建设流程全解析

做综合能源系统优化研究的朋友，大概率都会撞上同一堵墙：算例里光伏、风电还有负荷全是“听话”的确定值，模型跑得又快又稳，可一到实际落地，风光出力一波动，负荷一变化，之前的“最优方案”立刻失灵。我第一次被审稿人问“源荷不确定性考虑了吗”时，只能临时套一个简单的鲁棒边界，结果仿真曲线抖得没法看。后来我花了将近两个月，把一套计及源荷不确定性的综合能源生产单元运行调度与容量配置随机优化模型完整落地成Matlab代码，跑通算例、复现论文、完成答辩，这才算是真正迈过了这道坎。

这篇博文就把这套模型的建模思路、Matlab实现细节和调试经验完整梳理一遍。内容包括：为什么用随机优化而不是简单留裕度、两阶段决策怎么拆、场景生成与削减的代码怎么写、Yalmip加求解器的配置流程，以及我实际踩过的坑。整套代码基于常见的电-热耦合综合能源系统设计，你在跑通之后改改参数，就能迁移到自己的算例里。

1. 项目全貌：这个随机优化模型到底在做什么

1.1 一句话看懂模型定位

给一个直接的定性：这是一个“投资决策+生产运营”联动的随机规划模型。模型的对象叫“综合能源生产单元”，简单说就是一个园区或局域能源系统，里面同时存在电、热两类甚至多类能量流的转换、存储与供应。

模型要回答两个层面的问题。

第一层是容量配置：光伏装多大、风电装多少、CHP（热电联产机组）选什么容量、电锅炉配几台、储能电池容量和额定功率怎么定。这一层是长期决策，一旦确定，短期内不会变。

第二层是运行调度：给定设备容量之后，在每个时段的典型运行方式是什么，比如CHP发多少电、储能在何时充电何时放电、从电网买多少电、是否需要弃风弃光。这一层是短期决策，受源荷实际波动的影响非常大。

这两个层面绑定在一起，就形成了“先定容量、后算运行、再用运行结果反推容量是否经济”的迭代关系。随机优化要解决的核心问题就是：在源荷不确定性的干扰下，找到一个使得全生命周期总成本期望值最小的容量配置方案和各场景下最优调度策略。

1.2 源荷不确定性为何绕不开

很多初学的人会有疑问：我做一个确定性模型，把风光出力取个典型值、负荷取个平均值，行不行？行，但结果基本不可信。

原因是综合能源系统里的“源”和“荷”都具有明显的间歇性和随机性。光伏出力受云层遮挡、天气变化影响，某一天可能从中午的满发断崖式跌到20%；风电出力更离谱，夜间大发的概率远高于中午；而电负荷和热负荷也随生产生活节律波动。如果把这些全部当成确定值来优化，模型会给出一个“刚刚好”的方案，但这个方案在真实场景下大概率无法满足能量平衡，储能容量也留得不足。

随机优化解决的就是这类问题。它不再假设未来只有一种可能性，而是生成一大批可能的场景（场景表示源荷的联合概率分布），让容量配置方案和运行策略在所有这些场景的加权期望意义下达到最优。这样做出的方案不是对某一个场景最优，而是对所有可能出现的情况综合最优，通俗讲就是“不把鸡蛋放在一个篮子里”。

1.3 模型适用场景与改造潜力

从适配性看，这套模型可以直接用于：

• 园区级综合能源系统规划，含冷热电联供、分布式光伏、储能。
• 微电网容量优化，把热力部分去掉即可。
• 含P2G（电转气）、氢能储能的综合能源生产单元，需在设备模型里加相应环节。
• 区域多能互补系统，目标函数、约束条件按实际改造。

这也就是说，你拿到的核心不是某个具体数字的算例，而是一套可扩展的建模框架。后续章节我会把每一部分拆开，讲清楚为什么这么建、代码在做什么、如果你想改某一处应该动哪里。

2. 建模思路拆解：从物理问题到数学表达

2.1 不确定性建模：场景法还是鲁棒法？

处理不确定性，学术界两条主流路线就是场景法（随机优化）和鲁棒法（鲁棒优化）。我在这套代码里选择场景法，主要原因有三个。

一是信息利用程度不同。随机优化可以把风光的概率分布信息完整带进模型，最终得到的是期望成本最优；鲁棒优化只关心最坏情况，给出的方案通常偏保守，牺牲经济性换取可靠性，而这种保守在大多数园区项目中并不划算。

二是问题性质合适。容量配置本质上是一个“完全不确定”的长期决策，而运行调度是“看到一部分实现后再决策”的短期问题，这正是典型的两阶段随机规划结构，也是场景法最擅长的领域。

三是计算上可控。虽然场景法会把问题规模放大几百倍，但配合场景削减技术，最终模型规模在现代求解器面前完全可解，不需要特别复杂的分解算法。

场景生成方面，我采用的是“拉丁超立方采样+同步回代削减”的组合。拉丁超立方采样相对普通蒙特卡洛的一大优势是能以更少的样本覆盖更均匀的输入空间，极大地缓解随机采样带来的场景聚集问题。每次采样时，把风电预测误差、光伏预测误差、电负荷功率、热负荷功率作为随机变量，分别取对应的正态分布或威布尔分布，生成足够多的初始场景，再做削减，保留具有代表性的几十个场景用于优化。

2.2 目标函数与约束条件怎么搭

这套模型的整体目标函数是：

min 年化设备投资成本 + 各场景下运行成本的概率加权和 + 惩罚项

投资成本包括CHP、风电、光伏、电锅炉和储能单位的年化投资费用。运行成本包括从电网购电费用、燃气费用。惩罚项一般用于约束松弛，比如失负荷惩罚和弃风弃光惩罚，它不会在正常情况下出现，但如果求解器实在找不到可行解，有这个惩罚项就不会直接崩掉，并且能在一定程度上反映供需不平衡的代价。

约束条件覆盖了综合能源系统建模里必须考虑的几个维度：

• 电功率平衡：风电、光伏、CHP发电、储能放电、购电之和，等于电负荷、电锅炉耗电、储能充电之和。
• 热功率平衡：CHP余热、燃气锅炉和电锅炉产热，满足热负荷与储热装置的充放。
• 设备运行边界：各设备的出力上下限、CHP热电比约束、储能充放电功率限制。
• 储能时序约束：储能SOC（荷电状态）动态递推方程、充放电状态互斥、周期始末SOC相等约束。
• 容量配置约束：设备安装容量决定最大出力上限，运行变量不得超出配置容量。

在建立约束时有一个值得强调的细节：储能充放电状态互斥约束不能用简单的“同时为一小时和同时为0小时的和等于1”这种写法，因为长时间尺度下会出现乘积项更麻烦。标准做法是用大M法引入两个二进制变量，分别表示充电状态和放电状态，然后加“充电变量+放电变量≤1”的约束。这样一来，无论什么工况，储能都不可能既充又放，避免无用循环影响系统运行成本和设备寿命。

2.3 容量配置与运行调度的耦合逻辑

容量配置和运行调度为什么必须放在同一个模型里？因为设备装多大，直接决定运行成本能压到多少。举个简单的例子：光伏装小了，白天需要大量购电，运行成本高；装大了，有时发电用不完要弃光，投资回收期又太长。这两个方向相互拉扯，单靠手工试凑很难找到让总成本最低的平衡点。

模型中把容量配置定义为第一阶段决策变量，记作x；把运行调度定义为第二阶段决策变量，记作y(s)，s代表场景编号。目标函数里投资成本只跟x有关，运行成本要通过所有场景的期望值来积分。约束里凡是与设备容量相关的，都写成类似“y(s) ≤ x”的形式，这就是两个阶段耦合的桥梁。

求解时，整个模型一次性交给求解器处理，不是解开成两个循环。熟练之后你还可以在上层套一个整数变量，比如“某类设备装还是不装”，这就升级成了混合整数规划，但核心结构不变。这里有一个非常重要的建模原则：容量变量与运行变量的单位要统一。如果容量的单位是kW，运行功率的单位也要是kW，投资成本按年化折算，运行成本按调度时段折算，不要混用kW、MW或元/年、元/kWh，否则结果必然出问题。

3. Matlab代码实现：核心环节的落地细节

3.1 整体代码架构与文件组织

这套代码我用的是Matlab+Yalmip+Gurobi的组合。Yalmip是一个建模语言层，把优化问题的表达式写出来后，由它自动转成求解器需要的格式；Gurobi负责实际求解。Yalmip的好处是数学表达与代码几乎一一对应，改模型跟改公式一样直观，换求解器时只需要动一行参数。

工程上我把代码拆成下面几个文件：

• main.m：主程序，负责调用以下所有环节，控制整体流程。
• data_case.m：算例数据设置，所有设备参数、负荷曲线、预测误差参数都在这里。
• scenario_gen.m：生成初始场景集，采用拉丁超立方采样。
• scenario_reduce.m：同步回代削减，得到最终的典型场景及各场景概率。
• build_model.m：构建Yalmip变量、目标函数和全部约束。
• plot_result.m：结果可视化，绘制容量配置和典型日调度曲线。

模块化处理会让调试方便得多。一旦发现结果异常，先判断是场景的问题、约束的问题，还是数据的问题，直接定位到单一文件即可。

3.2 场景生成与削减的代码实现

场景生成的核心代码逻辑是这样的：先定义随机变量个数和时间段数，生成长度为“随机变量数×时间段数”的拉丁超立方样本，再把样本映射到对应的概率分布上。

function scenarios = scenario_gen(PARAM) % PARAM中包含各随机变量的分布参数 n_rand = 4; % 风电、光伏、电负荷、热负荷 n_s = PARAM.n_initial; % 初始场景数 n_T = PARAM.n_T; % 调度时段数 % 生成[0,1]区间的拉丁超立方样本 X = lhsdesign(n_s, n_rand * n_T, 'criterion', 'maximin'); % 逐列映射到对应分布 % 风电：预测值叠加预测误差(均值为0的正态分布) % 光伏：用Beta分布描述不同天气的波动 % 电/热负荷：以预测负荷为均值，相对误差服从正态分布 scenarios = zeros(n_s, n_rand, n_T); for s = 1:n_s for t = 1:n_T for r = 1:n_rand u = X(s, (r-1)*n_T + t); % 调用指定分布逆函数生成实际值 end end end end

我实际测试下来，初始场景生成500到2000个比较合适。太少，分布形状覆盖不全；太多，削减负担大，而且后续优化规模指数级上升。

场景削减用的是同步回代（Simultaneous Backward Reduction），思想很简单：不断找距离最近的场景对，把其中概率最小的那个删掉，并把它的概率加到与之距离最近的场景上，直到场景数降到设定值。这一段代码的核心是距离矩阵计算和迭代删除，我用一个典型的量化距离（欧氏距离平方）作为场景相似度指标。

function [new_scen, new_prob] = scenario_reduce(scen, prob, n_keep) n_s = size(scen, 1); while n_s > n_keep % 计算场景两两之间的距离矩阵 D = zeros(n_s, n_s); for i = 1:n_s for j = 1:n_s D(i,j) = sum((scen(i,:) - scen(j,:)).^2); end end % 找到概率最小且与另一场景距离最近的场景 min_dist = inf; del_idx = 0; for i = 1:n_s for j = 1:n_s if i ~= j && prob(i) < 1e-8 continue; end if D(i,j) * prob(i) < min_dist min_dist = D(i,j) * prob(i); del_idx = i; near_idx = j; end end end % 删除选中场景并把概率转给最近场景 prob(near_idx) = prob(near_idx) + prob(del_idx); scen(del_idx, :) = []; prob(del_idx) = []; n_s = n_s - 1; % 概率归一化 prob = prob / sum(prob); end new_scen = scen; new_prob = prob; end

有一点必须注意：场景削减不是简单地随机抽选。如果只随机抽，大概率会把那些出现概率低但对系统安全性影响大的极端场景丢掉，优化结果会乐观过头。我保留的场景数通常在30到80个，这个量级下，随机优化结果已经能与完整场景集的结果保持高度接近，而计算时间能降低一个数量级。

3.3 优化求解：Yalmip与求解器配置

Yalmip部分最关键的两个函数是sdpvar和optimize。容量变量用x_pv = sdpvar(1,1,'full')定义，运行变量用y_wt = sdpvar(n_T, n_S, 'full')定义，其中n_S是削减后的场景数。为什么运行变量要按“时间段×场景”二维定义？这里有一个新手很容易犯的认知错误：场景共享一套容量变量，但运行变量必须按场景分别定义。因为不同场景下风电出力不同，调度决策也应该不同。容量只要一个，因为它对所有场景都一样。

约束和目标函数用Constraints = [Constraints; ...]的方式拼接，目标函数直接写成Objective = invest_cost + sum(prob .* oper_cost_all_scenes)。

求解时我一般这样配置：

ops = sdpsettings('solver', 'gurobi', 'verbose', 2, ... 'gurobi.TimeLimit', 600, ... 'gurobi.MIPGap', 0.01); optimize(Constraints, Objective, ops);

MIPGap设为1%，意味着当上下界之差小于1%时停止求解。对规划类问题这个精度足够，求得太精细既不必要，也会消耗大量计算时间。如果模型里引入了储能互斥的二进制变量，那么问题就是一个混合整数规划，必须给Gurobi留足时间，同时设置好MIPGap，否则可能因为默认的10%精度不够，得到一个偏离最优较远的方案还浑然不觉。

3.4 关键参数与数据输入

数据参数是整套模型最容易出问题也最容易被忽视的部分。我把核心参数与典型取值整理成一张表。假如你现有算例的参数与以下不同，直接修改data_case.m即可。

从这张表能看出一点：不确定性参数的设定直接影响优化结果的保守度。标准差设大了，模型自然多配储能和备用容量；标准差设小了，方案就很激进。你若没有实测数据，建议先参考上述典型区间，再对标准差做一次灵敏度分析，“扫一遍不同误差水平下总成本的变化曲线”，这通常也是审稿人爱看的结果。

4. 实操过程记录：从跑通到会用的一次完整经历

4.1 环境准备与依赖检查

我的开发环境是Matlab R2022b + Yalmip R20200116 + Gurobi 10.0。

安装时有一个环节非常容易卡壳：Yalmip、Matlab和Gurobi的版本兼容性。提前确认Gurobi的Matlab接口版本与你本机Matlab是否匹配，因为Gurobi每个大版本都有对应的Matlab接口，装错会导致gurobi命令无法识别。装好之后，在Matlab里运行yalmiptest命令，全部pass再继续。

如果不想用商业求解器，也可以临时把solver换成cplex，甚至Matlab自带的linprog或intlinprog。但要注意，当场景数超过40个之后，自带的求解器速度会慢很多，除非模型规模很小，否则不建议用于正式计算。

4.2 算例数据设计

我自己的测试算例是改造自某课题组的模拟项目X，配置了一个光伏、一个风电场、一套CHP、一台燃气锅炉、一组电化学储能。系统满足一个包含工业负荷与居民热负荷的典型园区，夏季和冬季的负荷曲线差异明显。

设计数据时几个容易踩的坑如下。

一是全年数据不能用单一典型日替代。只用冬季典型日，算出来的容量配置显然偏保守，夏季光伏大发时甚至可能出现大面积弃光；只用夏季典型日，冬季供热就要出问题。我采用的办法是取典型夏、冬、过渡季三个典型日，每个典型日各24小时，然后再乘以对应的天数权重，近似覆盖全年。这个处理方式在论文写作中很好交代，也容易通过审稿。

二是负荷曲线要做归一化。不同园区的负荷峰值差距可能很大，把负荷曲线写成“峰值标幺值×基准功率”的形式，再结合园区实际面积或用户量调整基准功率。这样做的好处是以后换算例复用性强，只需要改一个基准值。

三是单价和单位量纲要统一。天然气热值按标准立方米折算成kWh时，要除以设备效率再计入成本，少除一个效率就会让成本计算出现巨大偏差。

4.3 代码运行与结果分析

跑通之后，主程序会输出两类结果：一类是容量配置结果，比如光伏装机容量、储能容量、CHP容量；另一类是每个典型场景下的设备出力曲线和SOC曲线。

我在第一次完整运行时得到的结果相当合理：当设置了15%的风电预测误差和10%的负荷预测误差后，系统储能容量相比确定性方案增加了约20%，而新能源装机容量略有增加，总运行成本在极端场景下下降了约8%。这说明随机优化方案是通过适度增加固定投资来换取运行阶段更高的可靠性和更低的风险成本，这正是计及源荷不确定性之后应有的典型表现。

如果在你的结果中没有观察到类似规律，最可能的原因是场景数太少或者误差标准差设得太小。场景数少于10个，随机优化几乎退化成多段确定性优化，不确定性效果显现不出来。标准差小于3%时，几乎所有场景都挤在预测值附近，优化方案自然与确定性方案无异。

4.4 结果可视化的小技巧

我习惯把结果画成四张图，这四张图在论文和答辩里都很能说明问题。

• 图1：容量配置对比柱状图，对比确定性方案与随机优化方案。
• 图2：典型场景削减前后的风电/光伏场景集合，体现削减效果。
• 图3：某个典型日各设备电出力堆叠图，展示调度策略。
• 图4：储能SOC变化曲线与热负荷平衡曲线。

可视化是帮助发现bug的好手段。如果SOC曲线出现突然跳变或充放电循环异常频繁，大概率是SOC递推方程里的时间步或效率写错了，一眼就能看出来。

5. 常见问题与排查技巧实录

5.1 求解器报错速查表

5.2 结果不合理怎么办

最典型的“不合理”现象就是：容量配置结果出现零值。比如优化结果里光伏容量为0，这并不一定是程序写错了，很多时候是因为在当前电价和补贴条件下，光伏的回收期太长，系统选择不装光伏反而更经济。判断方法是把光伏单位投资降一半，再看结果是否出现光伏容量。如果降了还不出，才需要考虑约束是否限住了安装上限。

另一个常见问题是储能容量特别大，大到明显不合常理。这通常是因为购电电价在峰谷时段差异太大，系统为了谷充峰放套利，拼命加储能。真实项目里要限制储能容量上限，比如不超过园区峰值负荷的20%~30%，或者加入储能的年维护成本再试。

还有一个高频问题：目标函数值出现负值。这在有售电收益或补贴收益的模型中可能是合理的，但大部分园区不允许向电网倒送电太多，检查一下光伏+风电+CHP出力上限是否超过了负荷与储能吸收能力之和，就能判断是不是倒送电约束漏写了。

5.3 面向自己问题的改造建议

假如你的研究对象不是电-热耦合园，而是冷热电联供系统，改造路径也清晰：在设备集合里加吸收式制冷机和电制冷机，新增冷功率平衡约束和制冷系数COP参数，目标函数里加制冷设备的运行成本，即可得到冷热电联供的随机优化模型。

如果要做碳交易或碳税场景，目标函数里加碳排放成本和碳配额约束，CHP换成低排放的燃气轮机，同时对燃煤设备做替换即可。这类扩展对建模框架本身没有冲击，只是在约束与费用项里加条目。

如果不想用场景法，想换鲁棒优化，那需要改动的地方主要是scenario_gen.m层换成不确定集合的建模方式，约束里把“场景索引s”的循环改为最坏情况的max-min表达式。理论层面会复杂不少，但由于你在Yalmip里写的是数学表达式，整体替换的工作量反而可控。

最后再分享一个我总结的判断标准：一个随机优化模型是否跑对了，不是看目标函数收敛到多小，而是看“确定性方案在随机场景下运行，是否总比随机优化方案更容易出现失负荷、弃能或成本超限”。如果随机优化方案在所有场景下的可靠性指标和成本波动都没有明显改善，那说明模型里不确定性的传导链路可能断了，优先检查容量变量与运行变量的耦合约束，看看运行变量是否被独立优化掉了。

这套代码我前后迭代了四个版本，从最初的单小时确定性调度，到多场景随机优化，再到容量与调度联合优化，每一步都是在“建模-求解-结果不合理-回去改约束”的循环里磨出来的。技术上的门槛其实不在求解器，而在如何把一个物理问题准确翻译成数学约束，以及如何从一堆数字里读出方案合理与否。希望这篇梳理能帮你把最难的这部分一次性迈过去。