👨‍🎓个人主页

💥💥💞💞欢迎来到本博客❤️❤️💥💥

🏆博主优势：🌞🌞🌞博客内容尽量做到思维缜密，逻辑清晰，为了方便读者。

⛳️座右铭：行百里者，半于九十。

📋📋📋本文目录如下：🎁🎁🎁

目录

💥1 概述

一、统一能路理论的核心框架

二、稳态与动态潮流计算的核心方法

（一）稳态潮流计算

（二）动态潮流计算

三、关键案例与验证

四、技术难点与解决思路

五、应用前景与挑战

结论

📚2 运行结果

2.1 六节点热网动态潮流

2.2  六节点热网稳态潮流

 2.3 七节点气网动态潮流

2.4 七节点气网稳态潮流

🎉3 参考文献

🌈4 Python代码、数据、文章

---

💥1 概述

摘要：潮流计算作为能源网络分析的基础性应用,在各能源网络已形成成熟但不统一的计算模型与方法。为促进不同能源网络研究的学科融合,该文基于统一能路理论,针对天然气网络与供热网络,提出了相适应的潮流计算方法。在此基础上,补充了基值修正的迭代方法以提高潮流计算的精度,并应用了"边值–初值"等效的方法以在动态潮流计算中隐式地给定初始条件。基于统一能路的潮流计算方法不仅统一了不同能源网络的潮流计算,还统一了同一能源网络内的稳态潮流计算与动态潮流计算,奠定了多异质能流在多时间尺度上统一分析的应用基础。此外,相较以有限差分方法为代表的传统动态潮流计算方法,文中所提出的方法在计算性能上实现了显著的提升。

一、统一能路理论的核心框架

统一能路理论借鉴电力系统“场→路”的推演方法论，将电力、天然气、供热等异质能源网络统一映射为“能路模型”，核心是通过数学化简实现跨能源网络的耦合分析：

1. 方法论基础：

   • 从三维偏微分方程描述的电磁场模型出发，抽取能量传播方向（坡印廷矢量方向），忽略正交方向过程，推导出电阻、电感、电容等等效元件模型，形成分布参数时域电路模型 。
   • 通过频域变换（傅里叶变换/相量法）将时域偏微分方程转化为频域常微分方程，再通过二端口等值模型简化为代数方程，实现集总参数频域电路建模 。
   • 针对非正弦激励（如气网），采用傅里叶分解为多频正弦分量，利用线性叠加原理合成响应 。

2. 能路模型的普适性：

   • 不同能源网络（电路、气路、水路、热路）均抽象为由“阻-容-感”能路元件构成的拓扑网络，通过支路特性（如管道阻力、热容）和拓扑约束（基尔霍夫定律）建立统一的网络矩阵与方程 。

---

二、稳态与动态潮流计算的核心方法

（一）稳态潮流计算

1. 数学模型：

   • 天然气网络：基于零频率网络方程，节点分为定压力节点（如气源）和定注入节点（如负荷），通过重排方程求解压力/流量分布 ：
     YgasP=Q
     其中 Ygas​ 为导纳矩阵，P 为节点压力向量，Q 为注入流量向量。
   • 供热网络：
     • 水力稳态：类似电网潮流，求解流量-压力方程；
     • 热力稳态：建立温度传递方程，考虑热损失与混温效应 。

2. 关键技术突破：

   • 基值修正迭代法：因气/热网方程非线性，线性化会引入误差。通过迭代更新流量/压力基值点，逐步逼近真实解（如图1），超线性收敛且误差<0.3% 。
   • 统一求解框架：与电网潮流共享矩阵方程形式，实现多能流协同计算 。

（二）动态潮流计算

1. 数学模型：

   • 动态过程描述为时空偏微分-代数方程组。通过傅里叶变换将任意激励分解为多频分量，对各分量求解频域代数方程，再逆变换叠加得到时域响应 ：
   • “边值-初值”等效法：隐式处理初始条件，避免传统差分法的初值敏感性问题 。

2. 性能优势：

   • 计算效率：传统有限差分法需小时级计算（如6节点热网需363秒），统一能路法仅需1-2秒（300频域分量），加速百倍以上 。
   • 精度保证：频域截断误差可忽略（图2），平均误差<1% 。

   ![[动态激励分解与响应恢复]]
   图：动态边界条件分解与支路温度响应（来源：）

---

三、关键案例与验证

1. 测试网络：

   • 7节点天然气网：动态潮流计算耗时1.89秒（传统法48.72秒），误差0.33% 。
   • 6节点供热网：动态计算1.03秒（传统法363.27秒），误差0.29%；稳态水力计算0.02秒，热力计算<0.01秒 。
   • 实际系统：中国东部气网、北方热网验证了工程适用性 。

2. 结果展示：

   • 显示稳态压力分布、动态温度传播延迟（如热源温度波动67分钟后影响支路末端） 。

---

四、技术难点与解决思路

难点	统一能路理论的解决方案
多能源方程形式差异	统一为频域代数方程，共享矩阵求解框架
非线性方程求解	基值修正迭代法实现超线性收敛
动态计算稳定性差	频域法避免差分步长限制，无稳定性问题
初始条件设置困难	“边值-初值”等效法隐式给定
多尺度耦合	频域分量独立计算，支持并行加速

---

五、应用前景与挑战

1. 三大统一价值：

   • 统一异质能源网络潮流计算模型；
   • 统一稳态/动态计算框架；
   • 统一多时间尺度分析基础 。

2. 现存挑战：

   • 模型扩展：当前仅支持“质调节”供热模式，“量调节”需进一步研究 ；
   • 不确定性处理：可再生能源波动未纳入，需结合鲁棒优化 ；
   • 多尺度耦合：电力（秒级）与气/热网（小时级）动态协同仍需突破 。

3. 未来方向：

   • 高精度设备级能路模型（如储能、电转气）；
   • 基于ML的参数辨识与实时校正；
   • 云边协同的并行计算架构 。

---

结论

统一能路理论通过频域代数化与基值迭代修正，解决了综合能源系统潮流计算的“模型分裂”与“计算低效”两大瓶颈。其动态潮流的百倍加速能力与<1%的误差，为多能流协同优化提供了工程实用工具。未来需进一步攻克多尺度耦合与不确定性分析，推动理论向“全场景-全要素”扩展。

📚2 运行结果

2.1 六节点热网动态潮流

 

2.2  六节点热网稳态潮流

 2.3 七节点气网动态潮流

复现结果：

复现结果：

2.4 七节点气网稳态潮流

 部分代码：

with context('读取数据与处理'):
    pipe_table = pd.read_excel('./7节点气网稳态data.xls', sheet_name='Branch')
    node_table = pd.read_excel('./7节点气网稳态data.xls', sheet_name='Node')
    numpipe = len(pipe_table)  # 支路数
    numnode = len(node_table)  # 节点数
    l = pipe_table['长度(km)'].values * 1e3  # 长度，m
    d = pipe_table['管径(mm)'].values / 1e3  # 管径，m
    lam = pipe_table['粗糙度'].values  # 摩擦系数
    cp = pipe_table['压气机(MPa)'].values * 1e6  # 支路增压，Pa
    c = 340  # 声速
    Apipe = np.pi*d**2/4  # 管道截面积
    v = np.ones(numpipe)*5  # 流速基值


with context('基于能路的潮流计算'):
    MaxIter = 100  # 最大迭代次数
    err = []  # 误差记录
    vb = [v.copy()]  # 基值记录
    for itera in range(MaxIter):
        # 支路参数
        Rg = [lam[i]*v[i]/Apipe[i]/d[i] for i in range(numpipe)]
        Lg = [1/Apipe[i] for i in range(numpipe)]
        Cg = [Apipe[i]/c**2 for i in range(numpipe)]
        Ug = [-lam[i]*v[i]**2/2/c**2/d[i] for i in range(numpipe)]
        
        # 支路导纳矩阵
        Yb1, Yb2, Zb, Ub = [], [], [], []
        f = 0  # 稳态，只取零频率分量
        for i in range(numpipe):
            Z, Y = Rg[i], 0
            za = np.cosh(np.sqrt(Ug[i]**2+4*Z*Y)/2*l[i]) - Ug[i]/np.sqrt(Ug[i]**2+4*Z*Y)*np.sinh(np.sqrt(Ug[i]**2+4*Z*Y)/2*l[i])
            za = za*np.exp(-Ug[i]*l[i]/2)
            zb = -2*Z/np.sqrt(Ug[i]**2+4*Z*Y)*np.sinh(np.sqrt(Ug[i]**2+4*Z*Y)/2*l[i])
            zb = zb*np.exp(-Ug[i]*l[i]/2)
            zc = -2*Y/np.sqrt(Ug[i]**2+4*Z*Y)*np.sinh(np.sqrt(Ug[i]**2+4*Z*Y)/2*l[i])
            zc = zc*np.exp(-Ug[i]*l[i]/2)
            zd = np.cosh(np.sqrt(Ug[i]**2+4*Z*Y)/2*l[i]) + Ug[i]/np.sqrt(Ug[i]**2+4*Z*Y)*np.sinh(np.sqrt(Ug[i]**2+4*Z*Y)/2*l[i])
            zd = zd*np.exp(-Ug[i]*l[i]/2)
            
            Yb1.append((za*zd-zb*zc-za)/zb)  # 稳态计算中，接地支路不起作用
            Yb2.append((1-zd)/zb)  # 稳态计算中，接地支路不起作用
            Zb.append(-zb)
            Ub.append(1-za*zd+zb*zc)
        yb, ub = np.diag(1/np.array(Zb)), np.diag(Ub)
        
        # 节点-支路关联矩阵
        A = np.zeros([numnode, numpipe])
        Ap = np.zeros([numnode, numpipe])
        for row in pipe_table.iterrows():
            A[int(row[1][1])-1, row[0]] = 1
            A[int(row[1][2])-1, row[0]] = -1
            Ap[int(row[1][1])-1, row[0]] = 1
            Ap[int(row[1][2])-1, row[0]] = 0
        
        # 节点导纳矩阵
        Yg_ = np.matmul(np.matmul(A, yb), A.T) - np.matmul(np.matmul(np.matmul(A, yb), ub), Ap.T)
        
        # 节点分类
        fix_G = node_table[node_table['节点类型']=='定注入'].index.values
        fix_p = node_table[node_table['节点类型']=='定压力'].index.values
        Yg_11 = Yg_[fix_G][:,fix_G]
        Yg_12 = Yg_[fix_G][:,fix_p]
        Yg_21 = Yg_[fix_p][:,fix_G]
        Yg_22 = Yg_[fix_p][:,fix_p]
        assert np.linalg.cond(Yg_11)<1e5  #　确认矩阵不奇异
        
        # 形成广义节点注入向量(给定)　与　节点压力向量(给定)
        Gn_1 = node_table[node_table['节点类型']=='定注入']['注入 (kg/s)'].values.reshape([-1,1])  # kg/s
        Gn_1 -= np.matmul(np.matmul(A[fix_G,:], yb), cp.reshape([-1,1]))
        pn2 = node_table[node_table['节点类型']=='定压力']['气压 (MPa)'].values.reshape([-1,1]) * 1e6  # Pa
            
        # 求解零频率网络方程
        pn1 = np.matmul(np.linalg.inv(Yg_11), (Gn_1 - np.matmul(Yg_12, pn2))).real
        Gn_2 = (np.matmul(Yg_21, pn1) + np.matmul(Yg_22, pn2))
        Gn_2 += np.matmul(np.matmul(A[fix_p,:], yb), cp.reshape([-1,1]))
        
        # 计算失配误差
        p = []
        pn1, pn2 = pn1.reshape(-1).tolist(), pn2.reshape(-1).tolist()
        for node in node_table['节点类型'].values:
            p.append(pn1.pop(0) if node=='定注入' else pn2.pop(0))
        p = np.array(p).reshape([-1,1])
        I = np.matmul(A.T, p).reshape(-1) + cp.reshape(-1) - (np.array(Ub).reshape([-1,1])*np.matmul(Ap.T, p)).reshape(-1)
        for i in range(numpipe):
            I[i] = abs(I[i]*np.diag(yb)[i]/Apipe[i]/(np.matmul(Ap.T, p).reshape(-1)[i])*c**2)
        err.append(np.linalg.norm(I-v))
        print('第%d次迭代，失配误差为%.5f'%(itera+1, err[-1]))
       

🎉3 参考文献

部分理论来源于网络，如有侵权请联系删除。

[1] 陈彬彬,孙宏斌,吴文传,等.综合能源系统分析的统一能路理论(三): 稳态与动态潮流计算[J].中国电机工程学报, 2020, 40(15):11.

[2] 陈金富.柔性交流输电系统的潮流计算及控制策略研究[D].华中科技大学,2002.

[3] 贾树森.交直流系统稳态低频减载优化算法和暂态低压减载控制策略研究[D].合肥工业大学,2012.

🌈4 Python代码、数据、文章