企业官网建设流程全解析

最近在做一个电网碳排放核算相关的项目，翻了不少EI论文，发现“电力系统碳排放流”这个词出现频率特别高。无论是做源网荷协调优化、碳追踪，还是算负荷侧碳责任，大家都在用这个方法。我前阵子把一个最经典的版本——基于潮流结果的碳排放流计算方法，在IEEE 14节点系统上用Matlab完整复现了一遍。这篇文章就把整个思路、矩阵推导、代码实现和踩坑记录全部写出来，正好给正在看论文但被公式卡住、或者想在自己算例上做碳流分析的读者一份可以直接照做的参考。

1. 碳排放流是个什么东西，为什么大家都在算

1.1 从“发电排放”到“用电责任”

传统意义上我们核算电力系统碳排放，基本都是按发电厂的口径统计：某台火电机组装机多少、发了多少电、烧了多少煤，乘一个排放因子，就得出这家电厂排了多少碳。这种“源头核算”方法本身没问题，但它回答不了几个很现实的问题：电网里流动的电能到底是谁发的？用户用的每一度电，到底对应多少碳排放？西部风电场发的清洁电送到东部，东部负荷的碳排放到底算谁的？

碳排放流就是用来回答这些问题的。它的核心思想很朴素：把发电侧产生的“碳排放”当成一种伴随有功功率流动的物质，顺着输电线路、变压器一级一级往下游传播，最后分摊到每个负荷节点和每条支路上。这样就能算出每个节点的碳势（相当于每个节点的“碳浓度”）、每条支路的碳流率（相当于“碳流量”），以及每个负荷承担了多少碳责任。

这套方法和电力系统里最成熟的潮流计算是天然配套的：潮流算出功率怎么流，碳流就跟着功率怎么走。所以碳排放流本质上是一种基于潮流的后处理分析工具，不会反过头去影响系统运行方式，计算量也小，非常适合做碳排放核算、碳追踪、碳责任分摊这些场景。

1.2 为什么拿IEEE 14节点当试验田

很多刚接触碳排放流的人会问：为什么大家都拿IEEE 14节点系统做复现，而不是直接上IEEE 118或者某省实际电网？

原因有几个。第一，IEEE 14节点规模小但有代表性：14条母线、5台发电机、20条支路（线路加变压器），拓扑不简单到一眼看穿，也没有复杂到手工算不动，刚好能把碳流计算的全流程走通。第二，这个算例是公开标准算例，Matpower里自带case14数据，任何人下载安装Matpower就能直接跑，复现门槛极低。第三，EI论文里做碳排放流验证的，一半以上都用14节点系统作为基础算例，复现它等于和文献里的结果有了一个“公共坐标系”，后续换更大系统或者改算法，对比起来方便。

我在实际做的时候，还特意把case14里默认的零出力机组（节点3、6、8）改成有出力，用来模拟多电源混合输送的场景。这样算出来的碳流分布更有意思，能明显看到清洁电源对下游节点碳势的“稀释”作用。

2. 数学模型：怎么把碳流算出来

2.1 三个基础假设，先立规矩

碳排放流计算不是凭空来的，它建立在几个假设之上。理解这几个假设，后面看公式才不迷糊。

第一个假设：碳排放流只伴随有功功率流动，无功功率不承担碳流传递。这个其实很好理解，碳排放本质上是和电能生产量挂钩的，无功功率不产生实质电能量，所以碳流只沿有功潮流路径传播。

第二个假设：在同一个节点上，来自不同电源的碳排放完全混合，节点上所有流出功率（包括负荷和转出支路）的碳势相同。就好比两杯不同浓度的糖水倒进一个杯子里搅匀，从这杯水里倒出去的任何一杯糖水，浓度都一样。这是碳排放流的“共享分配”原则，也是它和“逐源追踪”类方法的最大区别。

第三个假设：碳排放流在网络中是守恒的。发电机注入多少碳流，最终必须全部分摊到负荷和网络损耗上去，不多不少。基于这个守恒关系，才能在每个节点建立碳流流入等于流出的平衡方程。

这三个假设决定了方法的基本框架。实际工程中还有一些扩展做法，比如把网损按比例分摊到各负荷，或者引入分段碳排放流来处理时序问题，但核心思路都离不开这三点。

2.2 两个核心变量：节点碳势和碳流率

碳排放流体系里的核心变量，说到底是两个。

第一个是节点碳势，记作 e_i，单位一般是 kgCO2/kWh。直观理解就是节点 i 上“每度电对应的碳排放量”。它有点像电力系统里的电压幅值，是整个碳流计算里最关键的中间量。只要所有节点的碳势求出来，后面所有碳流率都唾手可得。

第二个是碳流率，记作 R，单位一般是 kg/h 或 t/h。它表示单位时间内伴随功率流动通过某一断面或注入某节点的碳排放量。具体又分几类：支路碳流率（支路上传播的碳流量）、负荷碳流率（负荷从电网“取走”的碳流量）、发电碳流率（发电机注入电网的碳流量）。碳流率的大小等于对应有功功率乘以该功率携带的碳势，所以它是一个直接反映“碳怎么流”的量。

把两个变量放在一起看：发电机往节点注入功率 PG 的同时，也注入了一个碳流率 PG × e_G（e_G 是发电机的碳排放强度，也就是单位发电量对应的排放）；节点上所有流出功率，无论去负荷还是去支路，都按节点碳势 e_i 携带碳流。

2.3 从两节点推起，看懂矩阵方程来源

矩阵公式 F × e = C 看起来高大上，但它的来源用两节点系统就能推得明明白白。

假设发电机接在节点1，出力 PG1，碳排放强度 e_G1；节点1通过一条线路向节点2送有功 P12；节点2接了一个负荷 PL2。先看节点1：流入节点1的功率只有发电机出力 PG1，没有其他支路注入，所以节点1的碳势就是发电机的碳排放强度，即 e_1 = e_G1。

再看节点2：流入节点2的功率是支路有功 P12，流动过程中携带的碳流率是 P12 × e_1（等于 P12 × e_G1）；流出节点2的功率是负荷 PL2，碳流率是 PL2 × e_2。根据流入等于流出的守恒关系：

P12 × e_G1 = PL2 × e_2

所以：

e_2 = (P12 / PL2) × e_G1

注意，如果 P12 大于 PL2，说明线路上存在网损，那么 e_2 会比 e_G1 大。这个结果其实非常有物理意义：网损相当于白白消耗了一部分碳流，剩下的碳流由更少的负荷功率承担，所以负荷侧的“碳浓度”被抬高了。这也是碳排放流里一个很重要的结论——线损客观上会增加下游用户的碳责任。

把这种“节点流入碳流等于流出碳流”的平衡关系，对每一个节点都写出来，就组成一个线性方程组。把所有节点功率关系整理成矩阵形式，就是后面代码里要解的：

F × e = C

其中 F 是节点有功通量矩阵，对角线元素是节点总注入功率（发电机出力加支路流入），非对角线元素是负的支路流入功率；C 是各节点的发电碳流率向量，等于节点上发电机出力乘对应碳排放强度。解这个方程组，一次性得到所有节点的碳势 e。

3. Matlab复现全过程

3.1 环境准备：MATPOWER装好就成功了一半

碳排放流计算的第一步是拿到准确的潮流结果。自己写牛顿-拉夫逊潮流代码不是不行，但完全没有必要——学术界做电力系统分析事实标准的工具是MATPOWER，开源、免费、精度高，而且自带IEEE 14节点标准数据。

MATPOWER的安装非常简单：去官网下载压缩包，解压后把文件夹路径添加到Matlab路径，或者在Matlab里直接运行文件夹里的install_matpower脚本。装完之后在命令窗口敲 case14，能输出一个14节点的数据体就说明安装成功了。

我建议再顺带跑一下自带测试：

mpc = case14; res = runpf(mpc);

如果 runpf 返回的 res 结果里 c5conv 字段是1（表示收敛），环境就完全OK了。后面所有碳流计算都以 res 里保存的潮流结果作为输入。

3.2 跑潮流和提取数据

MATPOWER的 runpf 返回结果是一个结构体，里面包含 bus、gen、branch 三个最关键的子表。bus 表的第3列是节点有功负荷（MW），gen 表的前两列是发电机所在节点和出力（MW），branch 表的第1、2列是支路首端和末端节点编号，第14列是首端有功潮流 Pf，第16列是末端有功潮流 Pt。

这里的正负号约定特别重要，一定要先弄清楚：Pf 表示从节点 f 注入线路的有功功率，Pt 表示从节点 t 注入线路的有功功率。如果 Pf 大于0，说明实际功率从 f 流向 t；如果 Pf 小于0，说明实际功率从 t 流向 f。Pt 的正负同理，只是相对 t 节点而言。很多复现翻车都翻在这一步，下面代码里会刻意按这个约定处理。

提取基础数据的代码如下：

mpc = case14; res = runpf(mpc); bus = res.bus; gen = res.gen; branch = res.branch; nb = size(bus, 1); % 节点数 nl = size(branch, 1); % 支路数 % 节点发电出力聚合：可能有多个发电机在同一节点 PG = zeros(nb, 1); for k = 1:size(gen, 1) gbus = gen(k, 1); PG(gbus) = PG(gbus) + gen(k, 2); end % 节点有功负荷 PD = bus(:, 3);

这段代码里，给每个节点的所有发电机出力做了聚合，因为case14里一台发电机对应一个节点，但实际系统里一个节点挂多台机很常见，聚合逻辑可以直接复用。

3.3 构造碳流矩阵的核心代码

构造节点有功通量矩阵是碳流计算的核心环节。原理很简单：遍历每条支路，判断功率实际流向，把流入某节点的功率加到该节点的总注入里，同时在矩阵对应位置写上负的流入功率。

关键代码如下：

% 发电机碳排放强度，单位 kgCO2/kWh，按需自行调整 EG_bus = zeros(nb, 1); EG_bus(1) = 0.95; % 节点1：燃煤机组 EG_bus(2) = 0.55; % 节点2：燃气机组 EG_bus(3) = 0; % 节点3：清洁电源 EG_bus(6) = 0; % 节点6：清洁电源 EG_bus(8) = 0; % 节点8：清洁电源 % 构造节点有功通量矩阵 F，以及发电碳流率向量 C F = zeros(nb, nb); for k = 1:nl f = branch(k, 1); t = branch(k, 2); pf = branch(k, 14); pt = branch(k, 16); if pf > 0 % 功率方向 f -> t，流入 t 的实际功率为 -pt inflow_t = -pt; F(t, t) = F(t, t) + inflow_t; F(t, f) = F(t, f) - inflow_t; else % 功率方向 t -> f，流入 f 的实际功率为 -pf inflow_f = -pf; F(f, f) = F(f, f) + inflow_f; F(f, t) = F(f, t) - inflow_f; end end % 加上发电机出力到对应节点对角元 for i = 1:nb F(i, i) = F(i, i) + PG(i); end % 发电碳流率向量 C = PG .* EG_bus; % 求解节点碳势，单位 kgCO2/kWh e_node = F \ C;

这里面有个值得注意的细节：我判断潮流方向只用 pf 的符号，但累加节点流入功率时用的是对端功率 pt（取绝对值）。原因是有功损耗：首端送出的功率经过线路阻抗后会损耗一部分，实际到达末端节点的功率比首端小。如果用 pf 累加流入节点的功率，会导致节点功率不平衡，求出来的节点碳势也会失真。用末端实际流入功率 -pt 来构建矩阵，等价于把网损从下游节点的流入功率里扣除，保证每个节点“流入=流出+损耗”的物理关系成立，数值上更可靠。

解线性方程用的还是 Matlab 反斜杠运算符 F \ C，14阶矩阵求解是毫秒级的事情，完全不用担心性能。

计算完成后，还可以顺手做支路碳流率和负荷碳流率的后处理：

% 支路碳流率，单位 kg/h branch_carbon = zeros(nl, 1); for k = 1:nl f = branch(k, 1); t = branch(k, 2); pf = branch(k, 14); if pf > 0 branch_carbon(k) = pf * 1000 * e_node(f); else branch_carbon(k) = (-pf) * 1000 * e_node(t); end end % 负荷碳流率，单位 kg/h load_carbon = PD .* 1000 .* e_node;

这里功率用MW，碳势用kgCO2/kWh，乘上1000是因为MW等于1000kW，最终算出的碳流率单位是kg/h。如果习惯用t/h，把结果再除以1000即可。

3.4 结果解析与出图

用以上代码跑完case14，会得到每个节点的碳势值。复现时观察到的典型结果大致是这样：电源节点1（燃煤机组）碳势约0.95 kgCO2/kWh，节点2（燃气）约0.55上下，而接有清洁电源的节点3、6、8碳势会明显偏低，清洁电源下游的节点碳势也会被拉低。负荷较重的节点如果上游通道较长、损耗较大，碳势通常会被抬高，这就是前面两节点推导里“网损抬高下游碳浓度”的体现。

论文配图一般画几种类型的图。节点碳势用柱状图最直观：

figure; bar(1:nb, e_node, 'FaceColor', [0.2 0.5 0.8]); xlabel('节点编号'); ylabel('节点碳势 (kgCO2/kWh)'); grid on;

支路碳流率可以用有向图表现，Matlab的graph对象天然支持：

G = digraph(branch(:,1), branch(:,2), branch_carbon); figure; p = plot(G, 'Layout', 'force', 'LineWidth', 2); p.EdgeCData = branch_carbon; colorbar;

线越粗、颜色越深，说明这条支路承载的碳流量越大。这种图放在报告里非常直观，评审和同事看了都能一眼抓住重点。

4. 复现中踩过的坑

4.1 矩阵报奇异的三种可能

我第一次跑通之前，F \ C 直接报矩阵接近奇异，warning刷了一大片。排查下来基本是三个原因。

第一，支路方向判反导致矩阵结构错误。这在新建节点通量矩阵时最容易出现。比如 pf 为负时，如果不判断方向，仍然往 F(t,t) 里累加功率，就会把“流出”当成“流入”，矩阵某一行可能加起来小于等于零。判断方向必须以 pf 符号为准，不能想当然按支路编号从f到t。

第二，节点功率不平衡。如果直接拿 case14 原始数据构建矩阵，而不先跑 runpf 取潮流结果，那么节点注入和流出功率并不严格满足平衡条件，矩阵可能是病态的。记住：构建碳流矩阵的功率数据必须来自潮流计算结果，不是来自算例原始数据。

第三，孤立节点或者零注入节点处理不当。IEEE 14节点系统本身拓扑是连通的，但如果自己改成其他算例，可能存在某些节点没有任何功率流入也没有发电机，那 F 对应行全为零，矩阵必然奇异。标准做法是先检查每个节点的总注入功率，对孤岛节点单独处理或直接删掉。

4.2 支路潮流方向处理的经典误区

这一条我觉得值得单独拿出来讲，因为几乎所有初学碳排放流的人都会在这里出错。Matpower里 pf 为正表示从 f 流向 t，pf 为负表示从 t 流向 f。但很多人图省事，直接取 abs(pf) 然后默认功率是 f 流向 t，结果在环网或者功率倒送的支路上，碳流方向和实际完全相反，节点碳势全是负的。

我建议把方向判断和功率累加写成独立函数，并且用输入pf符号+绝对值双重校验：

function F = add_branch_flow(F, f, t, pf, pt) if pf > 0 inflow = -pt; % 实际流入节点 t 的功率 F(t, t) = F(t, t) + inflow; F(t, f) = F(t, f) - inflow; else inflow = -pf; % 实际流入节点 f 的功率 F(f, f) = F(f, f) + inflow; F(f, t) = F(f, t) - inflow; end end

这样封装好以后，后面换算例、换数据都不用改主逻辑，只改输入就行。

4.3 单位搞错的结果有多离谱

碳排放流里单位混用是最隐蔽的错误。我见过有人把发电出力读出来是标幺值，直接当成有名值参与计算，算出来的碳势比正常值差了100倍。

IEEE 14节点系统基准容量正好是100MVA，所以Matpower里所有功率都是标幺值形式。但 runpf 返回的 bus、gen、branch 表里，功率已经是有名值MW了，不需要再乘基准容量。真正要注意的是自己在计算发电碳流率时，功率用MW，碳势用kgCO2/kWh，最终碳流率单位是kg/h。如果碳势用g/kWh，算出来的结果又会差1000倍。

我习惯做完一步就打印一次中间结果，用两节点手算值去对照。比如先手动构造一个两节点的简单潮流结果，代入代码算一遍，确认结果和手推公式一致，再上IEEE 14节点。

4.4 多发电机节点碳强度聚合

如果某个节点挂了多台不同燃料类型的发电机，不能简单把每台机组的出力乘碳强度再加起来塞进C向量里面事。正确做法是先算节点总发电碳流率，再反推等效碳强度；或者干脆不聚合，把每个发电机作为独立的碳注入源加入方程。

我在代码里用的是先聚合节点发电出力，再用等效碳强度 EG_bus(i) 去乘。这种情况下 EG_bus(i) 应该是该节点所有发电机碳排放流率之和除以总出力：

EF_node = sum(gen_carbon_rate_per_gen) / PG(i);

如果是做EI论文复现，建议把发电机级别的碳流率保留下来，算完碳势后再按节点汇总，这样既能验证全网碳流守恒，又能方便画不同电源贡献占比的图。

5. 实操体会和扩展玩法

5.1 我的几点个人体会

跑通一遍碳排放流计算之后，有几个感受特别深。

第一，碳排放流计算的难度不在于方程本身，而在于潮流结果的正确解读。只要功率方向、单位这些细节处理好了，核心求解就一行 F \ C。很多论文里搞得神乎其神的公式，落到代码上其实很简洁。

第二，IEEE 14节点系统的碳流结果非常适合做“体检”。把各节点碳势排个序，能明显看出清洁电源对局部碳势的拉低作用，也能看到网损较大、供电距离较长的节点碳势被抬高了多少。这种结果拿来写分析报告，比单纯给一个碳排放总量有说服力得多。

第三，参数碳强度的设置对结果影响很大。我在代码里给的0.95和0.55只是常用参考值，不同文献取值差异很大。复现对比时要先确认对方论文用的排放因子，否则数值没有可比性。

5.2 往大系统和其他方向扩展

这套流程从IEEE 14节点换到IEEE 30、39、118节点，代码几乎不用改，只要把case14换成对应算例名称就行。真正需要改的是发电机碳排放强度配置，因为大系统的发电机类型更多、燃料更杂，得逐台查资料或按区域口径统一设置。

除了换算例，还有几个扩展方向我觉得很有意思：一是把碳排放流扩展到时序场景，用96点或者8760小时的潮流序列逐时段算碳流，再做日、月、年累加，得到的是“碳流曲线”；二是把碳排放流和最短路算法结合，做特定电源的碳流溯源，看某个风电场的低碳电到底送去了哪些负荷；三是把节点碳势作为目标函数的一部分，嵌入到最优潮流模型里，做低碳经济调度。

我现在手上在跑的一个项目，就是把14节点这套流程扩展到一个区域实际网架，难点已经不在碳流方程本身，而在数据清洗和发电机碳强度口径统一上。如果你也准备往更大规模系统做，建议先把F矩阵构造这部分好好封装成函数，后面换算例就只是换数据的问题。各位如果在复现过程中遇到其他奇怪的坑，欢迎一起交流。