电力系统稳态分析牛顿拉夫逊法

0引言潮流是配电网络分析的基础，用于电网调度、运行分析、操作模拟和设计规划，同时也是电压优化和网络接线变化所要参考的内容。

潮流计算通过数值仿真的方法把电力系统的详细运行情况呈现给工作人员，从而便于研究系统在给定条件下的稳态运行特点。

随着市场经济的发展，经济利益是企业十分看重的，而线损却是现阶段阻碍企业提高效益的一大因素。

及时、准确的潮流计算结果，可以给出配电网的潮流分布、理论线损及其在网络中的分布，从而为配电网的安全经济运行提供参考。

从数学的角度来看，牛顿-拉夫逊法能有效进行非线性代数方程组的计算且具有二次收敛的特点，具有收敛快、精度高的特点，在输电网中得到广泛应用。

随着现代计算机技术的发展，利用编程和相关软件，可以更好、更快地实现配电网功能，本文就是结合牛顿-拉夫逊法的基本原理，利用C++程序进行潮流计算，计算结果表明该方法具有良好的收敛性、可靠性及正确性。

1牛顿-拉夫逊法基本介绍1.1潮流方程对于N 个节点的电力网络（地作为参考节点不包括在内），如果网络结构和元件参数已知，则网络方程可表示为：=&&YV I （1-1） 式中，Y 为N*N 阶节点导纳矩阵；&V 为N*1维节点电压列向量；&I 为N*1维节点注入电流列向量。

如果不计网络元件的非线性，也不考虑移相变压器，则Y 为对称矩阵。

电力系统计算中，给定的运行变量是节点注入功率，而不是节点注入电流，这两者之间有如下关系：ˆˆ=&&&EI S （1-2）式中，&S 为节点的注入复功率，是N*1维列矢量；ˆ&S 为&S 的共轭；ˆˆi diag ⎡⎤=⎢⎥⎣⎦&&E V 是由节点电压的共轭组成的N*N 阶对角线矩阵。

由（1-1）和（1-2），可得：ˆˆ=&&&S EYV上式就是潮流方程的复数形式，是N 维的非线性复数代数方程组。

将其展开，有：ˆi i i ij jj iP jQ V Y V ∈-=∑&& j=1,2，….，N （1-3）式中，j i ∈表示所有和i 相连的节点j ，包括j i =。

将节点电压用极坐标表示，即令ii i V V θ=∠&，代入式（1-3）中则有： ()i i i i ij ij j j j iP jQ V G jB V θθ∈-=∠-+∠∑()()cos sin i j ij ij ij ij j iV V G jB j θθ∈=+-∑故有：()()cos sin sin cos i i j ij ij ij ij j ii i j ij ij ij ij j i P V V G B Q V V G B θθθθ∈∈⎧=+⎪⎨=-⎪⎩∑∑ i=1,2,…,N （1-4） 式（1-4）是用极坐标表示的潮流方程。

而节点功率误差：(cos sin )θθ∈∆=-+∑SPi ii j ij ij ij ij j iP P V V G B (1-5) (cos sin )θθ∈∆=--∑SP i i i j ij ij ij ij j iQ Q V V G B (1-6)式中：SP i P ，SP i Q 为节点i 给定的有功功率及无功功率。

1.2牛顿-拉夫逊法基本原理1.2.1牛拉法的一般描述牛拉法是把非线性方程式的求解过程变成反复对相应的线性方程式的求解过程，即非线性问题通过线性化逐步近似，这就是牛拉法的核心。

下面以非线性方程式的求解过程来进行说明。

设电力网络的节点功率方程一般形式如下：()SP =y y x （1-7）式中，SP y 为节点注入功率给定值；y 为SP y 对应的物理量和节点电压之间的函数表达式；x 为节点电压。

写成功率偏差的形式：()()SP =-=0f x y y x （1-8）应用牛拉法求解如下。

在给定的初值()0x 处将式（1-8）作一阶泰勒展开：()()()00T x ∂+∆=∂f f x x x定义T∂=∂f J x为潮流方程的雅克比矩阵，0J 为J 在()0x 处的值，则有： ()()10-∆=-x J f x用∆x 修正()0x 就得到x 的新值。

如果用k 表示迭代次数，写成一般的表达式，有：()()()()()()()()11k k k k k k -+⎧∆=-⎪⎨⎪=+∆⎩x J x f x x x x （1-9） 对于潮流收敛的情况，()1k +x 应比()kx 更接近于解点。

收敛条件为：()()max ki f x ε<由简单迭代法收敛性分析的结论知，越接近解点，牛顿-拉夫逊法收敛越快，它具有二阶收敛速度。

由图1.1可以直观地了解牛拉法的步骤：图1.1 牛顿-拉夫逊法的几何解释1.2.2极坐标的牛顿-拉夫逊法在极坐标中，()f x 有如下的形式：()()()()(),,,,nSP SP n r-⎡⎤⎡⎤∆-⎢⎥⎢⎥==⎢⎥⎢⎥∆-⎣⎦⎣⎦P V P P V f x Q V Q Q V θθθθ（1-10） 共2n-r 个方程，状态变量为：1212T T Tn n r V V V θθθ-⎡⎤⎡⎤==⎣⎦⎣⎦L L x V θ 共2n-r 个待求量。

r 个PV 节点的电压幅值给定，不需求解。

潮流雅克比矩阵的维数是（2n-r ）*(2n-r)，结构如下：nT T TTTn r-⎡⎤∂∆∂∆⎢⎥∂∂∂⎢⎥==∂∂∆∂∆⎢⎥⎢⎥∂∂⎣⎦PPfV J x Q Q V θθ 上式右侧的对电压幅值的偏导数项中的电压幅值的阶数减少了1，为使雅克比矩阵的各部分子矩阵具有一致的形式，在实际计算中，常将该项乘以电压幅值，并选取1122////Tn r n r V V V V V V --⎡⎤⎡⎤∆=∆∆∆⎣⎦⎣⎦L V V 作为待求的修正量，则雅克比矩阵可写成：nT T TT Tn r-⎡⎤∂∆∂∆⎢⎥∂∂∂⎢⎥==∂∂∆∂∆⎢⎥⎢⎥∂∂⎣⎦P P V f V J x Q Q V V θθ（1-11）将式（1-10）和（1-11）代入式（1-9）的修正方程即可求得x 的修正量x ∆，用它修正x 直到()()max ki f x ε<为止。

将式（1-11）用下式表示：⎡⎤=⎢⎥⎣⎦H N J ML 其中每个字块的计算公式如下：对角元素： ''2''2''2''2,,,,i i ii i ii i ii ii i i i i ii i i ii i ii ii i i i i ii i ii i ii ii i i i iii i i ii i ii ii i i P Q H V H V H B V P P N V V N V N G V V Q P M V M V M G V Q Q L V V L V L B V V θθ⎧∂∆===+⎪∂⎪⎪∂∆⎪===--∂⎪⎨∂∆⎪===-⎪∂⎪⎪∂∆===-⎪∂⎩（1-12）非对角元素： '''''''''',cos sin ,cos sin ,,i ij i ij j ij ij ij ij ijji ij j i ij j ij ij ij ij ijji ij i ij j ij ijji ij j i ij j ij ijjP H V H V H B G P N V V N V N G B V Q M V M V M N Q L V V L V L H V θθθθθθ⎧∂∆===-⎪∂⎪⎪∂∆⎪===--∂⎪⎨∂∆⎪===-⎪∂⎪⎪∂∆===⎪∂⎩（1-13）2牛顿法潮流计算步骤2.1程序流程图在了解了牛拉法的原理之后，明确程序编写思路，如图 2.1、2.2所示。

其中图2.1中的“计算电压幅值和角度”步骤较多，单独用图2.2表示出来。

图2.1 牛顿法计算潮流的程序框图 图2.2 电压幅值和角度求解步骤框图当不符合收敛的条件“amontk>1”时,即认为计算不收敛。

具体程序见附录。

2.2计算步骤下面讨论的是极坐标形式的牛顿法求解过程，大致分为以下几个步骤： ① 形成节点导纳矩阵；② 给各节点电压设初值（(0)(0),i i V θ）；③根据式（1-12）、（1-13）生成雅克比矩阵（H 、N 、M 、L ）；④ 将节点电压初值代入式(1-5)、式(1-6)，求出修正方程式的常数项向量,i i P Q ∆∆；⑤ 求解修正方程，得到电压幅值和角度；⑥判断是否收敛，若收敛，计算平衡节点和线路功率； ⑦输出结果，并结束。

3算例3.1系统模型本文以图3.1所示电力网络为例，调用基于牛顿-拉夫逊法的C++程序。

1.1:10.3j0.12+0.5j0.1+0.4j0.08+0.4j①②③④图3.1 系统模型其中节点4设为平衡节点，电压标幺值为1.05，计算误差为0.000001。

3.2输入与输出将图3.1所示模型的相关数据放在data.dat文件中图3.2 输入节点和支路数据对各个数字含义的解释如下：网络模型有四个节点，四条支路，编号见图3.1。

第一个零下面三行数为支路参数，分别表示三条支路的起始和终止节点编号，后面的为电阻、电抗和电纳，电导均为0，例如：1 2 0.1 0.4 0.01528。

第二个零下面的为变压器支路，各数字意义同支路参数。

接下去三行均为节点参数，分别表示注入有功功率和无功功率。

调用text.cpp文件，得到运行结果，见图3.3和图3.4。

图3.3 运行结果1图3.4运行结果23.3结果分析将上述仿真结果整理为表格3.1、3.2，其中“+”表示节点i输出功率给节点j，“-”表示节点j输出功率给i(纵向为i，横向为j)。

表3.1节点有功功率输入与输出节点号 1 2 3 41 0 +0.245981 -0.5 -0.0465632 -0.24316 0 0 -0.3129493 0.5 0 0 04 0.0482143 0.319671 0 0表3.2节点无功功率输入与输出节点号 1 2 3 41 0 -0.014708 -0.029001 -0.1361872 0.0110505 0 0 -0.140363 0.097016 0 0 04 0.10464 0.160255 0 0根据表格计算：节点1有功功率：0+0.245981-0.5-0.046563=-0.300582无功功率：0-0.014708-0.029001-0.136187=-0.179896节点2有功功率：-0.24316+0+0-0.312949=-0.556109无功功率：0.0110505+0+0-0.14036=-0.1293095 节点3有功功率：0.5+0+0+0=0.5无功功率：0.097016+0+0+0=0.097016节点4有功功率：0.0482143+0.319671+0+0=0.3678853无功功率：0.10464+0.160255+0+0=0.264895根据已知条件，两个PQ 节点的注入有功、无功分别为：P1=0.3，Q1=0.18； P2=0.55，Q2=0.13 潮流计算误差：10.3005820.3100%0.194%0.3iP P ∆-=⨯=110.180.179896100%0.0578%0.18Q Q ∆-=⨯=220.5561090.55100% 1.11%0.55P P ∆-=⨯=220.130.1293095100%0.531%0.13Q Q ∆-=⨯=可见，误差均在允许范围内。