电力系统静态安全分析

电力系统稳态分析

电力系统稳态分析是研究电力系统运行和规划方案最重要和最基本的手段，其任务是根据给定的发电运行方式及系统接线方式求解电力系统的稳态运行情况，包括各母线定的电压、各元件中通过的功率等等。在电力系统运行方式和规划方案研究中，都需要通过潮流计算来分析比较运行方式或规划方案的可行性、可靠性和经济性。牛顿-拉夫逊法是电力系统潮流计算的常用算法之一，它收敛性好，迭代次数少，在电力系统潮流计算中得到了广泛的应用。

本文介绍了电力系统潮流计算方法中的牛顿-拉夫逊法的相关知识及其基本原理，并用MATLAB编写程序，最后通过算例来验证该程序的正确性。

关键词：电力系统潮流计算；牛顿-拉夫逊法；MATLAB

1 绪论

电力系统稳态分析是研究电力系统运行和规划方案最重要和最基本的手段，其任务是根据给定的发电运行方式及系统接线方式求解电力系统的稳态运行情况，包括各母线定的电压、各元件中通过的功率等等。在电力系统运行方式和规划方案研究中，都需要通过潮流计算来分析比较运行方式或规划方案的可行性、可靠性和经济性。

1.1潮流计算概述

电力系统潮流计算也分为离线计算和在线计算两种，前者主要用于系统规划设计和安排系统的运行方式，后者则用于正在运行系统的经常监视及实时控制。

利用电子数字计算机进行电力系统潮流计算从50年代中期就已经开始。在这20年内，潮流计算曾采用了各种不同的方法，这些方法的发展主要围绕着对潮流计算的一些基本要求进行的。对潮流计算的要求可以归纳为下面几点：

（1）计算方法的可靠性或收敛性；

（2）对计算机内存量的要求；

（3）计算速度；

（4）计算的方便性和灵活性。

电力系统潮流计算问题在数学上是一组多元非线性方程式求解问题，其解法都离不开迭代。因此，对潮流计算方法，首先要求它能可靠地收敛，并给出正确答案。由于电力系统结构及参数的一些特点，并且随着电力系统不断扩大，潮流计算的方程式阶数也越来越高，对这样的方程式并不是任何数学方法都能保证给出正确答案的。这种情况成为促使电力系统计算人员不断寻求新的更可靠方法的重要因素。

1.2 潮流计算的意义及其发展

电力系统潮流计算是电力系统分析中的一种最基本的计算，是对复杂电力系统正常和故障条件下稳态运行状态的计算。潮流计算的目标是求取电力系统在给定运行状态的计算。即节点电压和功率分布，用以检查系统各元件是否过负荷。各点电压是否满足要求，功率的分布和分配是否合理以及功率损耗等。对现有电力系统的运行和扩建，对新的电力系统进行规划设计以及对电力系统进行静态和暂态稳定分析都是以潮流计算为基础。潮流计算结果可用如电力系统稳态研究，安全估计或最优潮流等对潮流计算的模型和方法有直接影响。实际电力系统的潮流技术那主要采用牛顿-拉夫逊法。

在运行方式管理中，潮流是确定电网运行方式的基本出发点；在规划领域，需要进行潮流分析验证规划方案的合理性；在实时运行环境，调度员潮流提供了多个在预想操作情况下电网的潮流分布以校验运行可靠性。在电力系统调度运行的多个领域都涉及到电网潮流计算。潮流是确定电力网络运行状态的基本因素，潮流问题是研究电力系统稳态问题的基础和前提。

在用数字计算机解电力系统潮流问题的开始阶段，普遍采取以节点导纳矩阵为基础的逐次代入法。这个方法的原理比较简单，要求的数字计算机内存量比较下，适应50年代电子计算机制造水平和当时电力系统理论水平。但它的收敛性较差，当系统规模变大时，迭代次数急剧上升，在计算中往往出现迭代不收敛的情况。这就迫使电力系统计算人员转向以阻抗矩阵为基础的逐次代入法。阻抗法改善了系统潮流计算问题的收敛性，解决了导纳法无法求解的一些系统的潮流计算，在60年代获得了广泛的应用。阻抗法的主要缺点是占用计算机内存大，每次迭代的计算量大。当系统不断扩大时，这些缺点就更加突出。为了克服阻抗法在内存和速度方面的缺点，60年代中期发展了以阻抗矩阵为基础的分块阻抗法。这个方法把一个大系统分割为几个小的地区系统，在计算机内只需要存储各个地区系统的阻抗矩阵及它们之间联络线的阻抗，这样不仅大幅度地节省了内存容量，同时也提高了计算速度。

克服阻抗法缺点的另一途径是采用牛顿-拉夫逊法。这是数学中解决非线性方程式的典型方法，有较好的收敛性。在解决电力系统潮流计算问题时，是以导纳矩阵为基础的，因此，只要我们能在迭代过程中尽可能保持方程式系数矩阵的稀疏性，就可以大大提高牛顿法潮流程序的效率。自从60年代中期，在牛顿法中利用了最佳顺序消去法以后，牛顿法在收敛性。内存要求。速度方面都超过了阻抗法，成为60年代末期以后广泛采用的优秀方法。

1.3 本文的主要工作

本文介绍了电力系统潮流计算方法中的牛顿-拉夫逊法的相关知识及其基本原理，并用MATLAB编写程序，最后通过算例来验证该程序的正确性。

2 牛顿-拉夫逊法潮流计算基本理论

牛顿法（又称牛顿-拉夫逊法）是求解非线性方程式的典型方法。该方法有较好的收敛性，迭代次数少，在电力系统潮流计算中也得到应用。目前，牛顿法潮流计算是最为广泛、效果最好的一种潮流计算方法。

2.1 牛顿法基本原理

牛顿法是把非线性方程式的求解过程变成反复对相应的线性方程式的求解过程，即非线性问题通过线性化逐步近似，这就是牛顿法的核心。下面以非线性方程式的求解过程来进行说明。

设有非线性方程式

                           （1.1）

设为该方程式的初值，为初值的修正量。如果求得，则就可以得到真解

                        （1.2）

为此，将式

                       （1.3）

按泰勒级数展开

（1.4）

如果接近真值，则相对来说是足够小，所以可以略去所有的高次项。因此

                  （1.5）

可得

                        （1.6）

由于式（1.5）是式（1.4）的简化结果，所以由式（1.6）得到后，还不能得到方程式（1.1）的真解。实际上，用对修正后得到的：

                   （1.7）

只是向真正解更逼近一些。现在如果再以作为初值，求解

                                  （1.8）

就能得到更趋近真正解的

                      （1.9）

这样反复下去，就构成了不断求解非线性方程式的逐次线性化过程。第t次迭代时的参数方程为

                 （1.10）

或                        （1.11）

上式左端可以看成是近似解引起的误差，当时，就满足了原方程式（1.1），因而就成为该方程的解。式中是函数在点的一次导数，也就是曲线在点的斜率，如图1.1所示，

                     （1.12）

修正量则是由点的切线与横轴的交点来确定，由图1.1可以直观的看出牛顿法的求解过程。

图1.1 牛顿法的几何解释

现在把牛顿法推广到多变量非线性方程组的情况。设有变量的非线性联立方程组：

                        (1.13)

给定各变量初值，假设为其修正量，并使其满足

       （1.14）

对以上个方程式分别按泰勒级数展开，当忽略所组成的二次项和高次项时，可以得到

   （1.15）

式中：为函数对自变量的偏导数在点（）处的值。

把上式写成矩阵形式：

          （1.16）

这是变量的线性方程组，称为牛顿法的修正方程，通过它可以解出，并可以进一步求得

                       (1.17)

式中向真正解逼近了一步，如果再以它们作为初值重复解式（1.16）修正方程式，等到更接近真解的，如此迭代下去，并按式（1.17）进行修正，直到满足收敛要求为止并停止迭代计算，这就构成了牛顿法的迭代过程。

一般第t次迭代式的修正方程为

       (1.18)

上式可以简写为

                     (1.19)

式中

               (1.20)

为第次迭代时函数的误差相量；

                   (1.21)

称为的第次迭代时的雅克比矩阵；

                         (1.22)

为第t次迭代时的修正量相量。

同样，也可以写出类似（1.17）的算式

                     (1.23)

这样反复交替的解式（1.19）及式（1.23）就可以使逐步趋近方程式的真正解。为了判断收敛情况，可采用一下两个不等式中的一个：

                        (1.24)

                         (1.25)

式中，为预先给定的很小正数。

2.2 牛顿法潮流求解

在潮流计算中，节点极坐标功率方程式

                (1.26)

                (1.27)

式中：表示号后的节点都直接与节点相连，并且包括的情况。

节点功率误差

               (1.28)

              (1.29)

式中：，为节点给定的有功功率及无功功率。 

当采用的是极坐标模型时，待求量是各节点电压的幅值和角度。对于PV节点来说，节点电压幅值是给定的，不再作为变量。同时，该点不能预先给定无功功率，这样，方程式中也就失去了约束作用。因此，在迭代过程中应该取消与PV节点有关的无功功率方程式。只有当迭代结束后，即各节点电压向量求得之后，才利用这些方程式来求各PV节点应维持的无功功率。同理，由于平衡节点电压幅值和相角都是给定量，因此与平衡节点有关的方程式也不参与迭代过程。迭代结束后，利用平衡节点的功率方程式来确定其有功功率和无功功率。

设系统节点总数为，PV节点共个，PQ节点共个，平衡节点1个（在潮流计算中，该类节点一般在系统中只设一个）。为了求解方便，把平衡节点排在最后，即设为第个节点，则潮流要解的方程式应包括

，，共个方程式。

，，共个方程式。

以上方程式的待求量为各节点电压的角度及电压幅值，其有个。由于中不包括PV节点的电压幅值，所以共有个。这样未知量共有个，恰好可由以上个方程式求出。

将式(1.28)、式(1.29)按泰勒级数展开，略去高次项，可得

(1.30)

  或                                     (1.31)

式中，电压幅值的修正量采用的形式没有特殊意义，只是为了使雅克比矩阵中各元素具有比较相似的表达式。

利用简单的微分运算对式（1.28）、式（1.29）取偏导，可以得到雅克比矩阵中各元素的表达式

非对角元素（）

                 (1.32)

              (1.33)

               (1.34)

              (1.35)

对角元素

      (1.36)

    (1.37)

    (1.38)

     (1.39)

2.3 牛顿法潮流计算步骤

下面讨论的是用极坐标形式的牛顿法求解过程，大致分为以下几个步骤：

形成节点导纳矩阵；

给各节点电压设初值（）；

将节点电压初值代入式(1.28)、式(1.29)，求出修正方程式的常数项向量；

将节点电压初值代入式(1.32)~式(1.39)，求出雅可比矩阵元素；

求解修正方程式(1.30)，求出变量的修正向量；

根据式(1.17)或式(1.23)，求取节点电压的新值；

图1.2 牛顿法计算潮流的程序框图

检查是否收敛，由式(1.24)和式(1.25)可知，若电压趋近于真解时，功率偏移量将趋于零.如不收敛，则以各节点电压的新值作为初值自第步重新开始下一次迭代，否则转入下一步。

计算支路功率分布，PV节点无功功率和平衡节点注入功率。

输出结果，并结束。

图1.2为牛顿法计算潮流的程序框图。图中为事先给定的最大迭代次数，当实际迭代次数时，即认为计算不收敛。下面就据此流程图编写程序。

3 牛顿-拉夫逊法潮流计算编程

本文讨论用极坐标形式的牛顿法潮流的求解。下面就几个子函数作简要说明。

3.1 潮流计算数据文件建立

在一般的潮流计算中，需输入的原始数据包括节点（母线）数据；支路（包括接地支路，变压器支路等）数据。

3.1.1 节点数据

节点数据按一下的格式进行输入

格式：[节点编号，节点电压，节点相角，节点注入有功，节点注入无功，节点类型]

其中，节点类型取值为1、2、3，分别对应PQ节点、PV节点和平衡节点；节点注入功率如果是负荷，取负值；否则取正值；节点编号的原则是先PQ节点，然后PV节点，最后平衡节点。

3.1.2 支路数据

对于支路(包括对地支路和变压器），原始数据格式定义如下

格式：[节点I，节点J,，线路电阻R，线路电抗X，线路电导G，线路电纳B，变压器变比K]

其中，K=0表示普通线路；K>0，表示变压器非标准变比在j侧；K<0，表示变压器非标准变比在i侧。

于是，可以通过以上的格式进行节点和线路数据输入。在实际中，我们通过把节点和线路数据存在文件中，通过子程序1将数据文件的打开进行数据输入，程序代码详见附录A。

3.2 节点重新排序

在前面讲到的电力系统的潮流计算中，总是假设节点编号依据PQ节点、PV节点、平衡节点的次序进行编号的。但是，实际电力系统的节点编号具有一定的随意性，并不总是按照这个原则形成原始数据的。

因此，当原始节点数据不按照PQ节点、PV节点、平衡节点的次序进行编号时，需要对节点重新编号。该程序见子程序2，程序代码详见附录A。

在潮流计算完成后还需要重新将计算结果重新还原为原来的节点编号。该程序见子程序9，程序代码详见附录A。

3.3 计算节点导纳矩阵

由电力系统分析的知识可知，节点导纳矩阵的计算可以归纳如下：

1) 节点导纳矩阵是方阵，其阶数就等于网络中除去参考节点外的节点数。参考节点一般取大地，编号为零。

2) 节点导纳矩阵是稀疏矩阵，其各行非对角元素中非零元素的个数等于对应节点所连接的不接地支路数。

3) 节点导纳矩阵的对角元素，即各该节点的自导纳等于相应节点所连支路的导纳之和。

                                       (3.1)

式中，为节点与节点间支路阻抗的倒数；当节点有接地支路时，还应包括=0的情况。

4) 节点导纳矩阵的非对角元素等于节点与节点之间的导纳的负值。

                                      (3.2)

5) 节点导纳矩阵一般是对称矩阵，这是网络的互易特性所决定的。从而，一般只要求求取这个矩阵的上三角或下三角部分。

按照以上的算式，对于实际网络可以根据给定的支路参数和连接情况，直观而简单地求出导纳矩阵。对于电力系统中的支路参数，可以用一下的算式来表示

                        (3.3)

设导纳矩阵为Y，开始时设Y=0。

3.3.1 普通线路

对于普通的线路来说，，等值电路如图3.1所示。由图可知

Y(I,I)=Y(I,I)+Yt+Ym；

Y(J,J)=Y(J,J)+Yt+Ym；

Y(I,J)=Y(I,J)-Yt；

Y(J,I)=Y(I,J)；

图3.1 普通线路等值电路

3.3.2对地支路

对于对地支路来说，，等值电路如图3.2所示。由图可知

Y(I,I)=Y(I,I)+Ym；

图2.2 对地支路等值电路

3.3.3变压器支路

对于对地支路来说，因非标准变比究竟在哪一侧的不同，有下面两种情况：

1）变压器支路（>0，非标准变比在j侧），等值电路如图3.3所示。由图可知

Y(I,I)=Y(I,I)+Yt+Ym；

Y(J,J)=Y(J,J)+K*K*Yt；

Y(I,J)=Y(I,J)-K*Yt；

Y(J,I)=Y(I,J)；

图2.3 变压器等值电路（>0，非标准变比在j侧）

2）变压器支路（k<0，非标准变比在i侧），等值电路如图3.4所示。由图可知

Y(I,I)=Y(I,I)+Yt+Ym;

Y(J,J)=Y(J,J)+K*K*Yt;

Y(I,J)=Y(I,J)-K*Yt;

Y(J,I)=Y(I,J); 

图2.4 变压器等值电路（<0，非标准变比在i侧）

由上述的基本原理和计算公式来编写程序，可以得到形成节点导纳矩阵的子程序3，程序代码详见附录A。

3.4 功率偏差的计算

按式(1.28)和式(1.29)进行计算，就可很容易地得到功率偏差。该程序见子程序3，程序代码详见附录A。

3.5 计算雅克比矩阵

按式(1.32)和式(1.39)进行计算，就可得到雅克比矩阵。该程序见子程序4，程序代码详见附录A。

3.6 列主元消去法求解方程组

在求解式(1.30)或式(1.31)时，运用列主元消去法进行求解。对于求解线性方程组等价于求解

                                    (3.4)

由上述的基本原理和计算公式来编写程序，可以得到求解方程组的子程序5，程序代码详见附录A。

3.7 计算节点功率注入

按式(1.26)和式(1.27)进行计算，就可很容易地得到PV节点的无功和平衡节点的功率。该程序见子程序7，程序代码详见附录A。

3.8 计算支路潮流

潮流计算的目的不仅仅是为了得到节点信息，更重要的是希望得到系统中每条线路的功率流动情况。计算每条线路流过的功率的计算公式如下：

3.8.1 普通线路

由普通线路的等值电路图3.1，可得

                (3.5)

                        (3.6)

                    (3.7)

3.8.2 变压器支路

对于变压器支路等值电路如图3.3和图3.4所示，它的计算同普通线路。

3.8.3 对地支路

由对地支路的等值电路图3.2，可得

                                        (3.8)

                                (3.9)

由上述的基本原理和计算公式来编写程序，可以得到计算支路潮流的子程序8，程序代码详见附录A。

3.9 牛顿-拉夫逊法潮流计算程序

根据图1.2 所示的牛顿法计算潮流的程序框图编写程序，可以得到牛顿法计算潮流的主程序，程序代码详见附录A。

此外，潮流结果的输出。该程序见子程序10，程序代码详见附录A。

4 算例

4.1系统模型

电力系统模型如图4.1所示，系统的参数见图。利用该模型进行牛顿-拉夫逊潮流计算。

图4.1 电力系统模型

在程序求解潮流时，数据的输入采用M文件的输入方式。节点和支路数据保存在lw_system.m中，如图4.2所示。

图4.2 原始节点和支路数据

4.2 计算结果

在MATLAB的Command Window中输入

跳出如下对话框

打开lw_system.m，计算完成后打开Result.m。

图4.3 牛顿法潮流计算结果

如图4.3所示，得到了牛顿法潮流计算的结果。每次迭代时，产生的雅克比矩阵、节点功率偏差和节点电压值见附录B。

4.3 结果分析

在前面提到的潮流计算中，总是假设节点编号依据PQ节点、PV节点、平衡节点的次序进行编号的。但是，有些时候需要改变节点类型，如改变原来的平衡节点3为PV节点，将原来的PV节点4变为平衡节点等。因此，当原始节点数据不按照PQ节点、PV节点、平衡节点的次序进行编号时，需要编写一个节点重新编号子程序。然而，该节点编号方法并不是一种高效的编号方法。当考虑雅可比矩阵和节点导纳矩阵的稀疏性时，为了减少消去过程或分解过程中产生的新注入非零元素，可以考虑节点优化编号方法。在减少计算量、提高计算速度等方面，该方法更优。

此外，雅克比矩阵中各个元素都是节点电压向量的函数，因此在迭代过程中，他们将随各节点电压向量的变化而不断变化。这一点是影响了牛顿法潮流计算效率的重要因数，因为不仅每次迭代都要重新计算雅克比矩阵，而且还需要重新进行三角分解。参考文献[2]建议采用部分更新雅克比元素以减少运算量，对牛顿法潮流程序进行改进。本文采用的是用极坐标表示，程序中对PV节点处理比较方便。但是，当采用直角坐标时在迭代过程中避免了三角函数的运算，因而每次迭代速度快些。在实际中，我们可以多方面综合考虑来选取采用哪种方法。

在求解线性方程组时，采用的是列主元消去法，该方法能满足一般的实际计算。但是，在实际中经常遇到要对方程组反复多次求解每次改变其常数项，而系数矩阵是不变的。为了提高计算速度，可以利用因子表求解。

5 结论

通过用MATLAB编写牛顿-拉夫逊潮流程序对5节点系统的计算，与参考文献[2]所给结果相同，即可证明该潮流程序的正确性。但是，此程序很难能运用到实际的中去。上面仅仅是对5节点的系统进行计算，实际的电力系统的节点数成百上千万。那么采用上面的程序计算，计算量相当大，会导致计算的速度较慢，同时占用的内存量较大，往往达不到性能指标的要求。总之，该潮流程序对一些节点较少的系统还是适用的，不过在许多方面还有待完善和优化。

参考文献

[1] 张伯明，陈寿孙，严正. 高等电力网络分析[M]. 北京：清华大学出版社，1996.

[2] 王锡凡. 现代电力系统分析[M]. 北京：科学出版社，2003.

[3] 韩祯祥. 电力系统分析[M]. 杭州：浙江大学出版社，2009.

[4] 郗忠梅，李有安，赵法起，等. 基于Matlab的电力系统潮流计算[J]. 山东农业大学学报，2010，41(2):291-294.

[5] 徐劲松，宁玉琳，杨永峰. 基于Matlab的电力系统PQ分解法潮流计算研究[J]. 电气传动自动化，2011，33(2):10-18.

[6] 王忠礼，段慧达，高玉峰. MATLAB应用技术[M]. 北京：清华大学出版社，2007.

附    录A

主程序 Powerflow_NR.m

子程序1  OpDF_.m

子程序2  Num_.m

子程序3  y_.m

子程序4  dPQ_.m

子程序5  Jac_.m

子程序6  Gauss_.m

子程序7  Lineflow_.m

子程序8  Busflow_.m

子程序9  ReNum_.m

子程序10  Printf_.m

附    录B