2. 从“单一微网”到“多微网+储能电站”的架构升级
搞过综合能源系统优化的人都有体会,单微网的冷热电联供(CCHP)优化调度已经不算新鲜事了,燃气轮机、余热锅炉、吸收式制冷机这些设备的建模和调度策略,相关研究一抓一大把。但这个项目标题里最值得玩味的是三个关键词:“储能电站服务”、“多微网”、“双层优化配置”。这三个词放在一起,意味着问题的复杂度不是简单叠加,而是质变。
先说“储能电站服务”。单独看储能,无非是电储能、热储能、冷储能各自的充放策略。但一旦储能电站以“服务商”的身份介入多微网系统,它的角色就变了——它不再只是某一个微网的附属设备,而是多个微网之间能量互济的“缓冲池”和“调度杠杆”。储能电站可以同时向多个微网提供电、热、冷三种能量的存储与释放服务,这就把传统的“各自为政”的微网群,变成一个可以跨微网调配能量的协同系统。我在实际项目中体会最深的一点是:储能电站的容量配置和运行策略,必须和微网群的负荷特性、分布式电源出力特性放在同一个优化框架里决策,割裂开来做,经济性和能效都会打折扣。
再说“多微网”。多微网和单微网的本质区别在于“耦合约束”和“利益博弈”。多个微网之间可以通过联络线进行电能交换,储能电站作为第三方服务方介入后,每个微网既可以从储能电站购电购热购冷,也可以把自己富余的能量“存”进储能电站,形成一种动态的供需匹配关系。这种架构下,任何一个微网的调度决策都会影响其他微网的购能成本和储能电站的收益,因此需要统一协调,而不是各自独立求解。
“双层优化配置”则是解决这类问题的主流数学框架。我先说结论:上层解决“怎么建”的问题,下层解决“怎么用”的问题。上层的决策变量是储能电站的容量配置(电储能容量、蓄热罐容量、蓄冷罐容量)和各微网内分布式供能设备的装机容量(燃气轮机、电制冷机、吸收式制冷机等);下层是在给定容量配置方案下,以24小时为调度周期,优化各微网和储能电站的运行策略,包括设备出力、储能充放功率、微网间的交换功率等。上下层通过“容量方案影响运行成本,运行结果反馈回上层目标函数”的方式闭合迭代,最终得出经济性最优的配置方案。
坦白说,如果只做单层优化——比如只优化运行调度,把容量当作已知参数——问题会简单很多,但实际工程中容量和运行是强耦合的:设备容量决定了一次性投资成本,也决定了运行阶段的可调度范围;运行策略又直接决定了全年运行费用。所以两层嵌套寻优是这类问题绕不开的结构。Matlab在这类问题上的优势在于,下层运行优化可以调用成熟的商业求解器(如Cplex、Gurobi)高效求解,上层容量配置可以用智能优化算法(如粒子群、遗传算法)迭代搜索,两者结合的“元启发式+数学规划”混合框架,收敛性和求解速度在实测中都比较理想。
3. 双层优化模型的搭建细节:目标函数、决策变量与耦合机制
3.1 上层的“建设决策层”在优化什么
我在搭建上层模型时,最先关注的是决策变量的完整性和投资成本的计算口径。上层决策变量主要包括:
- 储能电站侧:电储能额定容量、蓄热罐额定容量、蓄冷罐额定容量,以及对应的最大充放功率
- 各微网侧:燃气轮机额定容量、光伏装机容量、吸收式制冷机额定容量、电制冷机额定容量
- 联络线容量:微网间以及微网与储能电站之间的能量交换通道容量
上层的目标函数是全寿命周期内的综合成本最小化,通常表达为等年值投资成本加年运行成本。等年值处理的原因很简单:设备寿命不同,储能电池可能10年,燃气轮机可能20年,如果不折算成年值,直接相加会严重扭曲比选结果。折算公式是:
[ C_{inv} = \sum_{i} \left( \frac{r(1+r)^{N_i}}{(1+r)^{N_i}-1} \right) \cdot c_i \cdot S_i ]
其中r是折现率,一般取8%;Ni是设备i的寿命年限;ci是单位容量投资成本;Si是配置容量。这个公式看起来简单,但有一个容易踩坑的地方——储能电站的容量配置不能只看总容量,还得同时约束“容量”和“功率”两个维度。电池容量决定能量上限,逆变器/换流器功率决定充放电速率上限,两者的比值(即C-rate)直接影响系统能否在负荷高峰时段释放足够的能量。我见过很多初版模型只约束容量不约束功率,结果上层优化出一个“大容量小功率”的储能配置,下层调度时发现功率瓶颈,根本无法满足尖峰负荷的释能需求,整个迭代直接失真。
上层目标函数中的年运行成本来自下层优化返回的最优运行结果。这就是双层模型的耦合机制:下层算出一个最小运行成本,把它代入上层的等年值总成本中,组成了完整的目标值。上层算法根据这个目标值调整容量方案,再次调用下层求解,如此反复迭代,直到达到最大迭代次数或收敛精度。
3.2 下层的“运行调度层”如何建模冷热电联供
下层是在给定容量方案下,以24小时为周期、1小时为步长的运行优化问题。目标函数是日运行成本最小化,主要包括:
- 向电网购电费用
- 天然气购买费用(燃气轮机燃料成本)
- 储能电站充放能量服务费用(注意这里是“服务”模式,微网向储能电站支付费用,储能电站运营商获得收益)
- 设备启停成本
下层模型中的能量流是问题的核心。燃气轮机燃烧天然气发电,发电的同时产生高温烟气,余热通过余热锅炉回收,一部分用于供热,一部分驱动吸收式制冷机制冷——这是典型的“以电定热”和“以热定电”两种运行模式的由来。实际建模中,大多数研究采用“以热定电”模式约束燃气轮机的最小出力,因为如果纯按电负荷定出力,很容易出现大量余热被浪费的情况,系统能效反而不高。这个约束在Matlab中用线性不等式表达即可,但要注意:燃气轮机的热电比范围需要从厂家参数中提取,不能随便拍脑袋定。
各能量流的平衡方程如下:
电力平衡:
P_grid(t) + P_gt(t) + P_pv(t) + P_storage_dis(t) = P_ELoad(t) + P_ec(t) + P_load_proportion(t) + P_storage_ch(t) + P_ws_ch(t)这个式子中容易忽略的是电制冷机的耗电项P_ec(t)——电制冷一开,就是个大耗电设备。夏季工况下电制冷机的耗电可能占到微网总用电的30%以上,如果不精确建模,电力平衡会严重失稳。
热力平衡:
H_gt_recover(t) + H_storage_dis(t) + H_gas_boiler(t) = H_Load(t) + H_storage_ch(t)这里需要换热的效率衔接,余热锅炉的回收效率通常在0.7到0.8之间,吸收式制冷机的COP大约在0.7到1.2之间,而电制冷机的COP在3到5之间。COP的差异直接影响两种制冷方式的经济性比较,这也是下层优化的核心权衡之一:高峰电价时段,电制冷机的运行成本可能高于吸收式制冷机,调度会倾向于用吸收式;低谷电价时段则反过来。
冷力平衡:
C_abs(t) + C_ec(t) + C_storage_dis(t) = C_Load(t) + C_storage_ch(t)3.3 上下层迭代时的数据传递逻辑
我在这里要特别强调一个实操中很容易犯错的点:上下层之间的数据传递不能只传“目标函数值”,还必须把下层的运行状态分量(比如储能电量的时间序列、燃气轮机的出力序列)回传给上层。原因是上层的某些约束——比如储能电站的容量利用率、联络线的最大负载率——需要基于运行结果来校验。如果只传一个孤零零的运行成本数字,上层算法很可能搜索到一个“成本很低但运行不可行”的容量方案,求解器却以为找到了最优解。
我在代码实现中采用的做法是:上层每产生一组新的容量配置方案,就把这些参数写入一个结构体数组,作为全局变量传递给下层模型;下层求解完成后,从求解结果中提取运行成本、各时段储能电量、联络线功率等关键序列,打包返回给上层。这样上层每次评估一个方案时,得到的是一组完整的状态反馈,而不是一个数字。
4. Matlab代码实现:从数学公式到可运行代码的关键拆解
4.1 代码整体架构与文件组织
拿到这个项目,我不建议把代码写成一个大而全的单文件——双层模型本身的逻辑层级本来就复杂,再加上迭代循环,单文件调试起来会非常痛苦。我习惯的分层方式是:
main.m % 主程序,负责参数初始化、上层算法迭代控制、结果输出 upper_level.m % 上层容量优化函数(调用粒子群或遗传算法) lower_level.m % 下层运行优化函数(构建模型并调用求解器) data_parameters.m % 所有基础参数定义(负荷、光伏、电价、设备参数) build_lower_model.m % 在lower_level.m内部使用,构建Yalmip优化模型 plot_results.m % 结果可视化主程序mian.m的核心流程是:
- 加载负荷、光伏、分时电价等基础数据
- 初始化上层优化算法的种群(容量配置方案的初始集合)
- 对种群中每一个个体,调用lower_level.m求解下层的24小时运行优化
- 汇总每个个体的运行成本,计算等年值总成本,作为上层适应度
- 使用粒子群或遗传算法更新种群,重复第3步
- 达到终止条件后,输出最优容量方案与运行结果
4.2 上层算法选型:基于实测建议
上层算法我在不同项目中分别试过遗传算法(GA)、粒子群算法(PSO)和差分进化算法(DE)。实测下来的对比结果是:
- PSO收敛速度最快,但容易陷入局部最优,尤其当决策变量维度较高时,粒子的位置更新容易在某个维度上停滞
- GA的全局搜索能力更强,但收敛速度偏慢,需要较多迭代次数
- DE在平衡收敛速度和全局搜索能力方面比较均衡,且实现起来比GA还简单
对于这个项目,决策变量维度大约在10到15个之间(储能电站三个容量变量加各微网四五个设备容量变量),不算特别高维,用PSO把惯性权重和个体/全局学习因子调好,基本上在50到80次迭代内能稳定收敛。但如果你追求更稳妥的结果,我建议用DE算法,它的变异操作对参数不敏感,对新手更友好。
一个重要的实操建议是:不要在每次迭代中对全部种群个体都重新求解下层模型,因为下层模型是线性或混合整数线性规划,求解时间虽然不长(单次求解在1到3秒量级),但乘以种群规模50再乘以迭代次数100,总时长可能达到一到两个小时。我在工程实践中采用的方法是:先跑一次较小种群(比如30个个体)和较少迭代次数(50次)的预热优化,观察目标函数的收敛趋势,再根据结果决定是否扩大搜索范围。这样可以避免一开始就把计算资源浪费在不必要的精细搜索上。
4.3 下层Yalmip建模的关键代码模式
下层模型我强烈建议用Yalmip工具箱建模,再调用Cplex或Gurobi求解。Yalmip的语法简洁,能大幅减少建模时间,而且可以方便地切换求解器。核心代码如下:
% 定义决策变量(以燃气轮机为例) P_gt = sdpvar(24, 1, 'full'); % 燃气轮机出力 H_gt = sdpvar(24, 1, 'full'); % 燃气轮机余热回收量 u_gt = binvar(24, 1); % 启停状态,0/1变量 % 电储能 P_storage_ch = sdpvar(24, 1, 'full'); % 充电功率 P_storage_dis = sdpvar(24, 1, 'full'); % 放电功率 E_storage = sdpvar(24, 1, 'full'); % 储能电量状态 u_storage_ch = binvar(24, 1); % 充电状态 u_storage_dis = binvar(24, 1); % 放电状态 % 目标函数:购电成本 + 燃料成本 + 储能服务费 objective = sum(P_grid .* price_buy) + ... sum((P_gt ./ eta_gt) .* price_gas) + ... sum(P_storage_dis .* price_storage_service); % 约束条件 Constraints = []; % 电力平衡约束 Constraints = [Constraints, P_grid + P_gt + P_pv + P_storage_dis == ... P_ELoad + P_ec + P_storage_ch + P_ws_ch]; % 燃气轮机出力上下限约束 Constraints = [Constraints, P_gt_min * u_gt <= P_gt <= P_gt_max * u_gt]; % 储能充放电状态互斥约束(最关键的一步) Constraints = [Constraints, u_storage_ch + u_storage_dis <= 1]; Constraints = [Constraints, 0 <= P_storage_ch <= P_storage_ch_max * u_storage_ch]; Constraints = [Constraints, 0 <= P_storage_dis <= P_storage_dis_max * u_storage_dis]; % 储能电量状态转移约束 Constraints = [Constraints, E_storage(1) == E_storage_init + ... P_storage_ch(1) * eta_ch - P_storage_dis(1) / eta_dis]; for t = 2:24 Constraints = [Constraints, E_storage(t) == E_storage(t-1) + ... P_storage_ch(t) * eta_ch - P_storage_dis(t) / eta_dis]; end % 储能电量上下限 Constraints = [Constraints, E_storage_min <= E_storage <= E_storage_max];这里最需要留意的是储能充放电状态互斥约束。很多初学者会忽略u_storage_ch + u_storage_dis <= 1这个约束,导致求解器出现“同时充电和放电”这种物理上不可能的伪最优解。这个约束配合功率上限中的状态变量,才能避免能量凭空消失或产生的bug。另外,储能电量状态转移方程中,充电效率和放电效率不能混用——我见过有人把eta_ch和eta_dis都设成同一个值,结果储能效率虚高,优化结果和经济性分析直接失真。
4.4 双层交互代码的实现细节
双层交互是实现中最容易出bug的部分。我的建议是使用函数嵌套或匿名函数来传递数据,避免全局变量的滥用。具体做法是:
function total_cost = evaluate_individual(capacity_scheme) % capacity_scheme是包含所有容量参数的向量 % 调用下层求解 lower_result = solve_lower_level(capacity_scheme); % 计算等年值投资成本 inv_cost = calculate_investment_cost(capacity_scheme); % 总成本 = 等年值投资成本 + 年运行成本 total_cost = inv_cost + 365 * lower_result.operation_cost; end这里有一个关键的工程决策:年运行成本是否直接乘以365?严格来说,冬季、夏季和过渡季的负荷曲线不同,直接用典型日运行成本乘以365会带来较大误差。更好的做法是分别取冬、夏、过渡季各自的典型日,分别求解下层优化,再按各季节的天数加权平均。但这样会让计算量直接翻三倍,很多研究为了省时间退而求其次。我的折中方案是:如果负荷数据只有一组典型日,那就直接用365倍;如果有三组季节数据,务必用加权平均。后一种方式得到的配置结果往往更接近实际需求——储能容量不会因为只针对夏季典型日优化而配得过大,也不会因为只针对冬季典型日优化而配得过小。
5. 常见问题与排查技巧实录
这个项目我前前后后调试过多轮,踩过的坑比模型本身还多。这里挑几个典型问题,给后来者做个速查表。
5.1 求解器报“infeasible problem”怎么办
下层求解器返回不可行,是双层优化项目里出现频率最高的错误。排查路径按以下顺序去查:
- 先检查储能初末电量约束。很多模型要求一天的调度周期内储能“始末电量相等”(这是为了满足周期性和长期可持续性),但如果初始电量设得不对,比如初始电量低于最小电量约束,直接不可行。
- 再检查燃气轮机的“最小出力”约束。如果某时刻热负荷很小,燃气轮机以最小出力运行时余热已经超出热需求,而此时蓄热罐又是满的——这就是典型的“刚性约束冲突”。遇到这种情况,要么调整燃气轮机的最小出力比例(从0.3降到0.2或者更低),要么允许蓄热罐短时突破容量上限(在物理可行的前提加一个松弛变量)。
- 检查微网间的交换功率约束。多微网系统中,联络线功率上下限如果设得太小,某些时段的功率交换需求无法满足,也会导致不可行。
排查时我强烈建议开启Yalmip的调试模式,在约束构建完成后打印约束数量,并逐个检查是否有NaN或Inf混入:
optimize(Constraints, objective, sdpsettings('debug', 1, 'verbose', 2));debug模式会精确显示是哪一条约束导致了不可行,比人工排查效率高一个量级。
5.2 储能状态的“二进制变量”爆炸问题
如果每个时段都用一个二进制变量表示充放电状态,24小时就是24个0-1变量,再加上启停变量、联络线方向变量等,整个下层的二进制变量可能超过50个。Cplex和Gurobi求解这种规模的MILP问题其实并不吃力,但如果你把时间步长从1小时改成15分钟,时段数就变成96个,二进制变量直接翻四倍,求解时间可能从2秒涨到3分钟以上。
所以我要强调一个细节:在项目早期先用1小时步长把整个框架跑通、把双层交互逻辑理顺,确认可行后再考虑是否细化到15分钟步长。不要一上来就搞高精度时间粒度,否则你会在无休止的求解等待中磨光耐心。
另外,如果遇到求解时间过长的问题,可以尝试用Big-M法将互补约束线性化,或者用启发式规则替代部分二进制变量。比如默认储能电站只允许在电价低谷充电、高峰放电,这样可以在不增加二进制变量的情况下获得接近最优的解。这种方法作为可行性验证阶段的快速近似,实测误差通常控制在5%以内。
5.3 上层算法收敛到“伪最优解”的识别
双层优化的一个特点是:上层算法很容易跑出一个投资和运行的“表面最优解”——总成本很低,但容量配置方案在物理上明显不合理。我最常见到的情形是:储能容量配得特别大,结果全年都用不上几次,但投资成本分摊到每年后仍然不算高,因为折现率低、寿命长。这种解“算出来”经济性不错,实际工程中是完全不可行的。
处理办法是给上层目标函数加一些“工程合理性惩罚项”。具体来说,我在代码中额外统计了储能电站全年的充放电循环次数,如果低于某个阈值(比如150次),就在总成本上增加一个惩罚系数。这个操作不能算数学上很严谨,但工程上非常有效——它能引导算法避开那些“配置过剩但利用不足”的虚假最优解。类似的,我还会对燃气轮机的年利用小时数做约束,避免出现“装了大容量燃气轮机却一年到头不怎么用”的荒唐方案。
5.4 双层嵌套的计算时间优化
前面提到,双层计算的总时间等于上层迭代次数乘以下层单次求解时间。除了调整种群规模和迭代次数这些“软参数”外,还有一个硬性技巧能显著提升效率:为下层模型添加“热启动”。
具体做法是:在双层迭代过程中,把上一次迭代得到的运行变量值(比如燃气轮机出力、储能功率)作为下一次下层求解的初始点,通过sdpsettings('usex0', 1)选项让求解器从初始点出发搜索。实测下来,热启动能使下层单次求解时间下降20%到40%。虽然这个技巧无法保证每次都能加速,但在长时间迭代中,积少成多的效果非常可观。
另一招是复用下层求解结果。如果上层粒子群在迭代过程中生成了与前几代非常接近的容量方案(由于粒子群本身有向历史最优靠拢的特性,这种情况很常见),可以建立一个简单的缓存表:以容量方案的关键参数为键,存储对应的运行成本。当新的个体与缓存中的方案差异小于某个阈值时,直接调用缓存结果,跳过下层求解。这个技巧听起来朴素,但在实际项目中能节省累计时间的三分之一以上。
6. 结果可视化与方案比选:如何让“两层优化”的结论真正可落地
跑完双层优化,拿到一串数字不代表工作结束。从工程交付的角度看,让决策者看懂配置方案、理解运行策略、比较不同方案之间的差异,往往比模型本身更重要。我自己的习惯是至少画三类图:
第一类是“储能电站容量配置与充放电功率的24小时曲线图”。这张图直观展示上下层衔接的核心成果——电价峰谷期储能如何充放电、是否出现放电功率卡在容量边界上的时段、储能电量是否频繁触顶/见底。如果充电功率曲线在低谷时段反复冲到上限,说明储能容量可能仍然偏小;如果电量曲线大部分时间在中间徘徊,说明配置偏大了。
第二类是“微网间联络线功率与能量交换流向图”。多微网系统的核心价值在于能量互济,这张图能直观反映联络线的利用率。如果某条联络线全天的功率始终为零,说明这条线的容量配置成了摆设——要么两个微网的负荷特性过于同质化,要么储能电站的服务定价没有起到引导互济的作用。这类发现是调整模型参数的重要依据。
第三类是“设备年出力热力图”。我习惯用热力图的色块深浅表示各时段各类设备的出力水平,可以快速暴露“某台设备在大部分时间里出力极低但容量却配得很大”的配置浪费问题。这类设备配置调整后往往能显著优化方案的经济性。
可视化用的代码片段不复杂,以储能出力曲线为例:
figure; plot(1:24, P_storage_dis_opt, 'r-o', 'LineWidth', 1.5); hold on; plot(1:24, -P_storage_ch_opt, 'b-s', 'LineWidth', 1.5); stairs(1:24, price_buy_normalized, 'k--', 'LineWidth', 1); legend('放电功率', '充电功率(负值表示)', '分时电价(归一化)'); xlabel('时间/h'); ylabel('功率/kW'); grid on;把充电功率取负值显示,是为了让充放电曲线处在同一个坐标系下,对比电价的峰谷变化时更直观。
我在实际项目中的体会是:双层优化配置的最终价值,不是输出一张容量配置表就结束了,而是要让每个容量决策都能追溯到对应的运行状态和成本构成。比如储能多配100kW,换来的是峰谷套利收益增加多少;燃气轮机放大一档,全年的燃料成本上升多少、购电成本又下降多少。把这种敏感性分析的结果也放进报告里,方案的说服力会强很多。这个项目后续如果需要继续扩展,最值得尝试的方向是引入多场景随机优化——把负荷和光伏出力的不确定性通过典型场景集表达,而不是只用一组确定性曲线,会让配置结果对未来真实运行工况的适应性更强。