做了好几年综合能源系统优化调度我越来越觉得“有功-无功协同”这五个字的分量远比你想象中重。尤其在“碳中和”这个目标下电力网络和天然气网络深度耦合你要是还像以前那样只盯有功经济调度、无功靠就地补偿硬扛那电压越限、气压波动、碳排放超标这些问题迟早会一起找上门。这篇文章要讲的电气互联系统有功-无功协同优化模型本质上是把电网的潮流约束、气网的稳态流量约束和耦合设备的运行特性放进同一个优化问题里同时求解。目标只有一个在满足安全边界的前提下让运行成本和碳排放都尽量低。模型本身的规模并不大难点全在数学表达和Matlab代码实现上。1. 项目定位与建模思路为什么必须把有功和无功放进一个模型1.1 传统调度方式到底缺了什么以前做电力系统经济调度最常用的做法是“先有功、后无功”先用最优潮流算出各台机组的有功出力再单独做一次无功优化调整发电机机端电压、无功补偿装置和变压器分接头把电压拉回合格范围。这套思路在火电占比高、负荷曲线平稳的年代没什么大问题因为电源分布相对固定无功储备也充足。但碳中和目标带来两个根本性变化一是新能源大规模接入风机和光伏的出力随机波动导致节点电压频繁越限二是电气互联系统里多了燃气轮机和电转气P2G这类耦合设备电网和气网之间不再是“各管各的”而是存在双向影响。燃气机组多发电气网负荷就增加气压跟着下降P2G设备启动电网有功负荷上升局部电压又会变化。这时候如果你还沿袭“先有功、后无功”的两步走策略经常出现的情况是有功调度结果算出来很经济但到了无功优化阶段发现某个节点电压怎么调都调不回来只能回来修改有功计划一来一回浪费大量时间甚至可能陷入迭代不收敛的死循环。1.2 电气互联系统里到底有哪些设备在“互联”所谓电气互联系统在我的理解里包含三块电网、气网、耦合设备。电网大家很熟悉有母线、线路、变压器、发电机还有无功补偿装置气网则包括气源、输气管道、压缩机和气负荷。耦合设备是电气互联系统的灵魂最容易理解的有两类一类是燃气轮机组它消耗天然气、输出电能是“气转电”的通道也是碳排放的主要来源另一类是P2G也就是电转气它利用电力电解水制氢、再通过甲烷化反应合成天然气是“电转气”的通道能够消纳多余的新能源电力同时为气网补充气量。这两类设备的存在让电网的有功潮流量、无功电压水平、气网的压力分布、气源开度全都搅在了一起。燃气机组出力变化会改变气网节点平衡P2G耗电又会改变电网节点负荷任何单一网络的独立优化都很难得到全局最优解。1.3 这个模型定位是什么简单说这是一个扩展的最优潮流问题叠加了天然气网的稳态调度。它适合三类人参考一是做电力系统优化调度研究的学生需要把电-气耦合建模搞扎实二是刚进综合能源方向、想在Matlab里快速复现一个完整优化算例的工程师三是做园区级多能互补规划的人想评估燃气机组和P2G的联合运行经济效益。模型的目标函数可以做成经济成本最小也可以做成碳排放最小或者用碳交易机制把碳排放折算成成本放进目标函数里。我后文统一采用“运行成本碳交易成本”最小这样既贴近碳中和语境又能让优化结果体现出减排信号。2. 核心数学模型拆解电网潮流、气网流量与耦合约束2.1 电网侧建模交流潮流和无功设备约束电网侧采用极坐标形式的交流潮流方程。对每个节点i有功和无功需要满足P_g,i - P_d,i V_i ∑ V_j (G_ij cosθ_ij B_ij sinθ_ij)Q_g,i Q_c,i - Q_d,i V_i ∑ V_j (G_ij sinθ_ij - B_ij cosθ_ij)其中P_g和Q_g是发电机出力P_d和Q_d是负荷Q_c是无功补偿装置出力V_i和θ_i是节点电压幅值和相角G_ij和B_ij是节点导纳矩阵的实部和虚部。无功优化必须配套以下约束发电机无功出力上下限、无功补偿装置容量上下限、节点电压幅值上下限以及支路潮流或电流限制。注意支路限制不能省否则优化模型可能为了降网损把某些线路压得特别狠计算结果虽然“数值好看”实际运行完全没法用。2.2 气网侧建模Weimouth方程引入非凸性气网稳态模型里管道流量通常用Weimouth方程描述F_mn sign(π_m - π_n) · C_mn · sqrt(|π_m - π_n|)这里π p²采用气压平方而不是直接使用压力是为了让方程形式上更接近线性方便后续松弛处理。C_mn是管道常数由管径、长度、温度等因素决定。每个气网节点需要满足气流平衡流入节点的管道流量加上气源注入量等于节点气负荷加上流出节点的管道流量还要扣掉压缩机消耗的气量。此外还有节点气压上下限、气源产气量上下限。Weimouth方程里带符号函数和开方是非凸非线性约束直接丢给求解器大概率会得到“Infeasible problem”。常用的处理办法是做二阶锥松弛把原来的等式约束松弛为F_mn² ≤ C_mn² · (π_m - π_n)并且要求π_m ≥ π_n这一步看着只是把等号变成小于等于但性质完全不同。松弛后问题变成二阶锥规划可以用Gurobi、Cplex、Mosek这类商业求解器高效求解。至于松弛得到的解是否满足原等式需要在求解后回代检查一般只要目标函数往经济性方向压管道流量会自然贴近上限松弛间隙通常在可接受范围内。2.3 耦合设备与碳排放模型燃气轮机组耗气量和有功出力之间可以做线性化处理H_g a_g · P_g b_g这里的H_g表示耗气量a_g是边际气耗率b_g是空载气耗。额定工况附近线性化误差不大工程上完全够用。P2G设备则简单很多输入电力、输出天然气G_p2g η_p2g · P_p2gη_p2g是电转气综合效率一般取50%到70%。效率太低时P2G在经济上不占优势但碳中和模型里它承担了消纳新能源电力和降低碳排放的功能所以不能只看电-气转换效率。碳排放计算的思路是燃气机组按耗气量乘排放因子外购电力也要计及等效碳排放。为了体现碳中和目标我采用碳交易机制给系统分配一个免费碳排放配额超过配额的部分需要按碳价购买。这样碳排放就从硬约束变成了目标函数里的惩罚项优化器会在“多发电多花钱”和“多买碳配额”之间自动权衡。2.4 “协同”到底协同在哪里传统解耦模型里有功潮流方程和无功潮流方程分开处理气网和电网更是各算各的。协同优化模型则把所有决策变量放在同一个求解框架里有功出力P_g、无功出力Q_g、电压幅值V、相角θ、无功补偿Q_c、气源产气量G_s、管道流量F、气压平方π、P2G耗电和产气量全部一次性求出。这样一来燃气机组出力一旦调整不仅影响电网潮流的P方程还通过耗气量H_g影响气网的流量平衡P2G耗电则同时改变电网节点负荷和气网气源注入。优化器能够看到全局所有耦合关系得出的解才是真正意义上的全局最优。3. Matlab实现全流程环境配置、变量定义与求解命令3.1 环境准备YALMIP加Gurobi别在这步省时间Matlab代码实现我建议用YALMIP建模配Gurobi求解。YALMIP是MATLAB下的建模工具箱语法接近数学表达式能把二阶锥约束、线性约束混在一起写非常省心。Gurobi是商业求解器学术用户可以免费申请许可证安装过程也不复杂。先说工具箱路径设置。下载YALMIP后解压到任意目录在MATLAB里执行“addpath(genpath(‘YALMIP文件夹路径’))”再用“savepath”保存。Gurobi装完后在MATLAB里运行其自带的setup命令就能自动配置路径。这里要特别提醒一个很多人踩过的坑MATLAB 2023版本之后如果默认编码不是UTF-8写中文注释会出现乱码处理不好甚至会在运行脚本时报错。我的建议是在“预设项-编辑器/调试器-语言”里把编码统一改为UTF-8然后重启MATLAB。代码注释尽量用简洁中文关键约束行加英文关键词方便排查。3.2 数据准备用结构体把电网和气网数据整理干净数据管理建议用结构体不要散成几十个全局变量。电网数据至少包括母线编号、负荷有功无功、发电机接入母线、发电机有功/无功上下限、支路首末端母线、支路电阻电抗、对地电纳、变压器变比等。气网数据包括气源节点、气源产气上下限、管道首末端节点、管道常数C、节点气压上下限、气负荷分布。耦合数据单独建一个结构体燃气机组接在哪个电网节点、对应消耗哪个气网节点的气、气耗系数是多少P2G接在哪个电网节点、输出到哪个气网节点、转换效率是多少。我在调试不同算例时发现把数据集中放Excel或CSV里统一读取比直接写在m文件里更灵活。换算例时只要替换数据文件代码主体一行都不用改。3.3 变量定义与约束构建的核心片段决策变量通过sdpvar定义关键代码如下% 电网决策变量 P_g sdpvar(n_gen, 1, full); Q_g sdpvar(n_gen, 1, full); V sdpvar(n_bus, 1, full); theta sdpvar(n_bus, 1, full); Q_c sdpvar(n_shunt, 1, full); % 气网决策变量 Pi sdpvar(n_gasbus, 1, full); G_s sdpvar(n_gassource, 1, full); F sdpvar(n_pipe, 1, full); % 耦合设备决策变量 H_g sdpvar(n_gt, 1, full); P_p2g sdpvar(n_p2g, 1, full); G_p2g sdpvar(n_p2g, 1, full);潮流约束写成循环或矩阵形式。气网管道约束用二阶锥形式for k 1:n_pipe fro pipe_fr(k); too pipe_to(k); Constraints [Constraints, F(k)^2 C(k)^2 * (Pi(fro) - Pi(too))]; Constraints [Constraints, Pi(fro) Pi(too)]; end耦合约束则写为% 燃气机组耗气 for k 1:n_gt Constraints [Constraints, H_g(k) a_g(k) * P_g(gen_idx(k)) b_g(k)]; end % P2G电转气 for k 1:n_p2g Constraints [Constraints, G_p2g(k) eta_p2g(k) * P_p2g(k)]; end目标函数按“购电成本购气成本碳交易成本”构建最后调用ops sdpsettings(solver, gurobi, verbose, 2); sol optimize(Constraints, Objective, ops); if sol.problem 0 P_g_opt value(P_g); V_opt value(V); Pi_opt value(Pi); else disp(求解失败请检查约束); end这里有个细节值得强调YALMIP默认变量范围是[-∞, ∞]所有电气量一定要显式添加上下界。以电压为例不约束上下界的话求解器可能跑来一个极其离谱的数值但目标函数竟然还更小那就是彻底的错误解。3.4 求解器的数值稳定性调参Gurobi整体很稳但遇到数值病态问题也可能给出“Numerical trouble”之类的警告。最常见的诱因是量纲不一致电网用有名值气网用MPa结果相差好几个数量级导致约束矩阵条件数很大。我的做法是电网侧统一用标幺值气网侧气压平方π也做归算让数量级都在0.01到100之间。如果模型还是难收敛可以试试给目标函数加一个极小值罚项比如对机组出力加0.0001倍的二次罚项能让求解器在多个几乎相同的最优解里快速稳定下来。4. 仿真结果对比与物理量解读4.1 测试系统怎么搭为了快速验证模型我建议先搭一个小系统比如6节点电网配合7节点气网。电网由6条母线和3台发电机构成其中1台是燃气机组气网由7个节点、6条管道、2个气源构成。P2G设备接在第5个电网节点同时向第6个气网节点注气。我调试时把负荷水平和气负荷水平都做了归一化处理保证系统在没有任何优化的情况下电压和气压都处于合理范围内。这样后续优化结果的提升量才可信不会一上来就是越限状态。4.2 解耦优化与协同优化结果对比我分别跑了传统“解耦”模型和“有功-无功协同”模型结果整理成表格。这里的数据是示意性的但趋势完全符合物理直觉指标解耦优化有功-无功协同优化总运行成本/万元28.627.4网损率/%3.83.1节点电压最大偏差/p.u.0.0620.031气网压力最大偏差/MPa0.180.11系统碳排放/t12601215电压偏差是衡量无功优化的核心指标之一。协同优化把燃气机组有功出力和P2G耗电放在一起优化让有功分布和无功分布彼此照应电压被压在更靠近1.0 p.u.的水平网损率也随之下降。碳排放的下降也很直观协同优化会主动调整燃气机组和普通火电的出力比例更容易安排P2G在负荷低谷期运行把燃气机组的部分出力替换掉碳排放自然降下来。4.3 碳价变化如何影响运行决策碳交易价格是碳中和目标里最有意思的一个参数。我连续把碳价从50元/吨调到300元/吨观察系统运行方式的变化。碳价低的时候系统倾向于多开燃气机组、少用P2G因为气价相对便宜碳价升高之后燃气机组的运行成本明显上升P2G的利用率开始增加。再往上走系统还会出现一个拐点部分燃气机组直接减出力外购电力增加。但外购电也会带来等效碳排放如果电网侧的清洁比例不高可能出现“本地减排、外调增加”的碳排放转移现象。这个点特别值得做研究的人关注它说明碳中和不是一个纯技术优化问题还需要配套的碳排放核算边界设计。5. 常见报错排查与调试技巧实录5.1 求解器不可用或“No suitable solver”出现这条提示十有八九是YALMIP没找到Gurobi。先检查Gurobi能否在MATLAB里单独运行再检查路径是否正确。如果手头没有商业求解器我建议先换用免费的SDPT3或SeDuMi它们也支持二阶锥约束只是求解速度慢一些。但对小规模测试系统完全够用。5.2 模型返回“Infeasible problem”不可行是所有优化模型里最让人头大的问题。我调试时一般按照“先电网、后气网、再耦合”的次序逐段放开。具体做法是先把气网约束全部注释掉把燃气机组当成普通发电机组只跑电网OPF确认电网自身约束没问题再把电网约束锁死单独调度气网模型确认气网无约束冲突最后把耦合约束逐条加回去。如果加上某一类约束后立刻不可行问题通常出在耦合设备的“数量级”或“流向”上。比如燃气机组气耗系数写得太大气网气源根本供不上优化器就会告诉你无解。这种情况要么调整气源上限要么修正气耗系数而不是盲目去放松电压限值。5.3 电压初值导致收敛异常交流潮流方程里有V×V×cosθ这类乘积项本质是非凸问题Gurobi处理的是松弛后的模型求解稳定性整体不错但初值太离谱仍可能导致中断。我在调试时会把电压初值固定为1.0相角初值固定为0也就是常说的平启动。只需要在sdpvar定义后给变量加个低权重罚项或者初始化时用“assign(V, ones(n_bus,1))”。5.4 中文注释乱码问题MATLAB中文注释乱码最直接的原因是文件编码和编辑器编码不匹配。我建议在保存m文件时把编码改成UTF-8如果是旧版本遗留的GBK编码文件直接另存为UTF-8即可。实际操作中还有个更隐蔽的问题代码里出现了中文引号、中文括号YALMIP解析时会直接报错。排查办法很简单把报错行整行删掉重新手打一遍不要用编辑器从网页复制粘贴。5.5 松弛不紧怎么处理Weimouth方程做二阶锥松弛后如果发现管道流量平方远小于右侧压力差说明松弛并不紧解出来可能偏向某些低压节点。我在这种情况会采取两个办法一个是提高目标函数里碳价的权重引导气网和电网的耦合设备更积极运行另一个是给关键管道增加微小的流量罚项让优化器倾向于把管道流量推向约束边界。如果仍然不紧就要考虑是不是管道常数C或者节点气压上下限设置有问题很多人把节点气压下限设得过低导致压力差足够大管道流量反而游刃有余松弛自然紧不起来。5.6 调试心态与节奏做这类电气互联系统优化最容易犯的错误就是想一步到位。我自己的实践体会是第一次跑通模型时目标函数越简单越好甚至可以先把气网方程里的压缩机消耗忽略只保留管道流量平衡和压力约束等整个框架跑通之后再把压缩机、碳交易、P2G等扩展模块逐步加回来。这个“减法到加法”的过程看起来慢实际是最快的。因为每加一个模块你都能清楚地知道它对目标函数和约束的影响一旦结果异常你也能快速定位到底哪个模块在捣乱。最后再分享一个我每次复现这类模型都会固定使用的小技巧在输出模块里把电网节点电压、气网节点压力、机组出力和碳排放量同时画在一张多子图里。不要只在数据层面抠数值直接看趋势图往往比看一堆数字更能发现模型问题。多试几次你会对“电网和气网到底是如何耦合”有非常直观的感觉这东西比任何论文里的公式都来得深刻。