CN101404412B - 一种用于静态电压稳定性分析的方法 - Google Patents

Abstract

本发明提出了一种用于静态电压稳定性分析的方法，其利用线性模拟法计算速度快和连续潮流法计算准确的特点，将两者结合起来，可以取得快速和准确地求出某种负荷增长方式下的电网电压稳定裕度。本发明的分析方法既适用于静态电压稳定离线分析也可适用于自动电压控制系统的在线应用。

Description

一种用于静态电压稳定性分析的方法

技术领域

本发明涉及电力系统仿真计算领域，具体涉及一种用于静态电压稳定性分析的方法。 

背景技术

用常规潮流程序可以确定网络中不同点的电压稳定水平。然而这种方法是费力的，而且不能提供作决策非常有用的灵敏度信息。此外，由于雅可比矩阵奇异，常规潮流在崩溃点外无解和在崩溃点附近不能可靠收敛。 

对于调度人员来说，最实用的电压稳定安全指标应是功率裕度指标，这就要求比较准确地确定临界点，指示各负荷节点维持电压稳定性能力的强弱的两个重要参量：负荷点的临界电压和极限功率。为此要解随负荷变化的潮流方程。对于大电力系统来说，要解决三个问题：一是如何规定负荷的变化方式，即各负荷节点负荷如何改变比较接近电力系统实际。二是在临界点(拐点)，常规潮流计算方法(对应于恒功率负荷的牛顿拉夫逊法)是不收敛的。三是在系统负荷较轻时，存在潮流多解，但真正具有实际应用价值的PV曲线，要求其上半段对应于常规潮流解，而不要在各解组之间跳动。 

线性模拟算法，是一种反映系统实际运行状况、计算速度快、适用面广、精度又能满足工程需要的算法。它以良好的潮流初始解为初始条件，利用潮流的线性化形式，以各种非线性条件的变化作为约束条件，同时又以灵敏度指标的变化特征作为临界点位置的判据。 

连续潮流是电压稳定性分析的有力工具，它可以克服接近稳定极限运行状态时的收敛问题。延拓潮流法通过不断更新潮流方程，使得在所有可能的负荷状态下，潮流方程保持为良态，不管在稳定平衡点还是在不稳定平衡点都可有解。延拓法是从初始稳定工作点开始，随着负荷缓慢变化，沿相应PV曲线对下一工作点进行预报、校正，直到勾画出完整的PV曲线连续潮流计算方法。连续潮流分析方法具有鲁棒性和灵活性，它是解决收敛困难的潮流问题的理想方法。然而，这种方法计算非常耗时。 

基于线性模拟法的连续潮流法把线性模拟法和连续潮流法结合起来，可以取得快速和准确的良好效果。算法从基本工况开始逐步增加负荷，用线性模拟方法求得接近的临界点，然后用连续潮流方法求解潮流。通常线性模拟法可以提供直至临界点的解。连续潮流方法只是在需要求临界点和临界点后的准确解时才是必要的。

发明内容

1、基于线性模拟法的改进连续潮流法原理 

设研究系统有n个节点，n_PQ个PQ节点(无功出力越限的发电机组也作为PQ节点处理)，n_PV个PV节点和一个平衡机节点。 

极坐标系统下的常规潮流方程可以写成： 

$\begin{array}{cc}{P}_{\mathrm{Gi}}-P_{\mathrm{Li}}-V_i\underset{j\∈i}{\mathrm{\Σ}}{V}_j(G_{\mathrm{ij}}\cos {\mathrm{\θ}}_{\mathrm{ij}}+B_{\mathrm{ij}}\sin {\mathrm{\θ}}_{\mathrm{ij}})=0& i=\mathrm{1,2},L,n-1\end{array}---\left(1\right)$

$\begin{array}{cc}{Q}_{\mathrm{Gi}}-Q_{\mathrm{Li}}-V_i\underset{j\∈i}{\mathrm{\Σ}}{V}_j(G_{\mathrm{ij}}\sin {\mathrm{\θ}}_{\mathrm{ij}}-B_{\mathrm{ij}}\cos {\mathrm{\θ}}_{\mathrm{ij}})=0& i=\mathrm{1,2},L,n_{\mathrm{PQ}}\end{array}---\left(2\right)$

其中P_Gi，Q_Gi为节点i的发电机出力；V_i，θ_i为节点i的电压幅值和角度； 

P_Li，Q_Li为节点i的负荷；G_ij，B_ij为节点导纳矩阵第(i，j)个元素的实部和虚部。 

若以参数t来表示发电机和负荷的增长，则有： 

P_Gi＝P_Gi0(1+tk_Gi)                          (3) 

P_Li＝P_Li0(1+tk_PLi)                         (4) 

Q_Li＝Q_Li0(1+tk_QLi)                         (5) 

0≤t≤t_cr                                 (6) 

其中参数t＝0，对应基本的发电机出力和负荷水平；t＝t_cr表示临界点的出力和负荷水平。k是指定节点出力或负荷的增长系数。 

将变化的出力和负荷的表达式代入潮流方程，得到新的潮流方程： 

$P_{\mathrm{Gi}0}(1+{\mathrm{tk}}_{\mathrm{Gi}})-P_{\mathrm{Li}0}(1+{\mathrm{tk}}_{\mathrm{PLi}})-V_i\underset{j\∈i}{\mathrm{\Σ}}{V}_j(G_{\mathrm{ij}}\cos {\mathrm{\θ}}_{\mathrm{ij}}+B_{\mathrm{ij}}\sin {\mathrm{\θ}}_{\mathrm{ij}})=0---\left(7\right)$

$Q_{\mathrm{Gi}0}-Q_{\mathrm{Li}0}(1+{\mathrm{tk}}_{\mathrm{QLi}})-V_i\underset{j\∈i}{\mathrm{\Σ}}{V}_j(G_{\mathrm{ij}}\sin {\mathrm{\θ}}_{\mathrm{ij}}-B_{\mathrm{ij}}\cos {\mathrm{\θ}}_{\mathrm{ij}})=0---\left(8\right)$

用矩阵和向量表示，则为： 

H(X，t)＝0 

0≤t≤t_cr                     (9) 

对方程(1)和(2)在工作点处取偏差量，得到矩阵形式的线性化方程式如下： 

$\left[\begin{array}{c}\mathrm{\ΔP}\\ \mathrm{\ΔQ}\end{array}\right]=\left[\begin{array}{cc}{J}_{\mathrm{P\θ}}& J_{\mathrm{PV}}\\ J_{\mathrm{Q\θ}}& J_{\mathrm{QV}}\end{array}\right]\left[\begin{array}{c}\mathrm{\Δ\θ}\\ \mathrm{\ΔV}\end{array}\right]---\left(10\right)$

式中ΔP为节点注入有功的微增列向量，有n-1个元素；ΔQ为节点注入无功的微增列向量，有n_PQ个元素；Δθ为节点电压相角变化列向量，有n-1个元素；ΔV为节点电压幅值变化列向量，有N_PO个元素； 

是极坐标下的系统雅可比矩阵，简单记为J。式(10)即为分析系统静态稳定性的数学模型，也是潮流方程的修正方程式。 

本发明采用中国电科院系统所引进开发的中国版BPA电力系统分析程序工具PSD-BPA进行潮流计算，获得基本工况的潮流解。随着负荷的增长，电力系统从正常运行状态向临界状态的过渡方式可以有单负荷节点增加负荷的方式、多负荷节点增加负荷的方式和全网增加负荷的方式。增加负荷的不同方式所求得的静态电压稳定裕度各有不同的含义。 

从已知的解(O)开始，对于一个规定负荷增长方式至解(A)，在这个过程中假定各负荷节点按照恒定的功率因数增加负荷，则负荷节点有功、无功注入量的改变量可以表示成以下形式： 

其中矩阵 

${\mathrm{\λ}}_i=\mathrm{\Δ}{Q}_{P{Q}_i}/\mathrm{\Δ}{P}_{P{Q}_i}---\left(12\right)$

在系统负荷增长的过渡过程中，要考虑发电机的有功出力调整才符合系统的实际运行情况。在传统的潮流计算中，因为是针对一个固定的运行状态，因此有功的不平衡量是由平衡机承担的，当然这也是一种特殊的发电机有功出力调整方式。不同的发电机有功出力调整方案对应着不同的过渡方式，因此求得的静态电压稳定裕度也是不同的。若不考虑网络的有功损耗，则增加的有功负荷量应与发电机出力增量是相等的。 

假设 $\mathrm{\Δ}{P}_{P{V}_i}={\mathrm{\γ}}_i\underset{j=1}{\overset{n_{\mathrm{PQ}}}{\mathrm{\Σ}}}\mathrm{\Δ}{P}_{P{Q}_j}$

式中 

为第i台发电机组的有功出力增量，γ_i为第i台发电机组的有功出力调整系数， 

为第j个负荷的有功增量。 

写成矩阵形式为： 

$\mathrm{\Δ}{P}_{\mathrm{PV}}=\left[\begin{array}{cccc}{\mathrm{\γ}}_1& L& L& {\mathrm{\γ}}_1\\ {\mathrm{\γ}}_2& L& L& {\mathrm{\γ}}_2\\ & L& L& \\ {\mathrm{\γ}}_{n_{\mathrm{PV}}}& L& L& {\mathrm{\γ}}_{n_{\mathrm{PV}}}\end{array}\right]\mathrm{\Δ}{P}_{\mathrm{PQ}}=H_{\mathrm{\γ}}\mathrm{\Δ}{P}_{\mathrm{PQ}}---\left(13\right)$

γ_i取  $\mathrm{\Δ}{P}_{{\mathrm{PV}}_i}/\underset{j=1}{\overset{n_{\mathrm{PQ}}}{\mathrm{\Σ}}}\mathrm{\Δ}{P}_{{\mathrm{PQ}}_j},$ 把上式代入式(10)， 

$\left[\begin{array}{c}\mathrm{\Δ}{P}_{\mathrm{PQ}}\\ \mathrm{\Δ}{P}_{\mathrm{PV}}\\ M\\ \mathrm{\Δ}{Q}_{\mathrm{PQ}}\end{array}\right]=\left[\begin{array}{cc}{J}_{\mathrm{P\θ}}& J_{\mathrm{PV}}\\ M& M\\ J_{\mathrm{Q\θ}}& J_{\mathrm{QV}}\end{array}\right]\left[\begin{array}{c}\mathrm{\Δ\θ}\\ \mathrm{\ΔV}\end{array}\right]---\left(14\right)$

可以再写成以下形式： 

从上式中可以得到： 

矩阵S_PQ为n_PQ×n_PQ阶，它反映负荷节点有功、无功注入量的变化对负荷电压变化的灵敏度关系，在正常运行状态下，对角线各元素都为负值。当负荷按照一定的过渡方式增加时，矩阵S_PQ元素绝对值向增大的趋势发生变化，  $\frac{d^{2}V}{d{P}_{\mathrm{PQ}}^2}>0.$ 在临界点附近，矩阵中对角线元素绝对值先是剧增(到初始值的几十倍)，然后突然下跌，  $\frac{d^{2}V}{d{P}_{\mathrm{PQ}}^2}\≤0$ 表明已到临界点位置。 

对应某一工作方式下的雅可比矩阵为J_K，给定步长参数为Δt_K，假定多负荷节点增长负荷时按照负荷各自的初始量为比例系数，即有 

ΔP_PQ＝Δt_Kk_PLP_L，            (17) 

则在某一状态下的负荷扰动量由式(10)可以求得状态量的改变量(Δθ_i)_K、(ΔV_i)_K，再由下式求得新的状态量： 

(θ_i)_K+1＝(θ_i)_K+(Δθ_i)_K                    (18)

(V_i)_K+1＝(V_i)_K+(ΔV_i)_K             (19) 

在新的状态量基础上对雅可比矩阵进行修正，得到新的雅可比矩阵J_k+1，然后再由步长参数Δt_K+1推算下一个潮流状态的状态量。每进行一步，都要判断是否达到电压稳定的临界点，即观察矩阵S_PQ对角线元素的值对于有功量的导数是否小于o。若所有对角线元素对于有功量的导数大于零，即 

${\left[\frac{d^{2}V}{d{P}_{\mathrm{PQ}}^2}\right]}_i>0---\left(20\right)$

则判断系统未到达电压稳定临界点，即运行工作点在P-V曲线的上方；若对角线元素对于有功量的导数小于零，则判断已经到达临界点的位置，返回到原来的状态(A)，转入连续潮流法计算。 

连续潮流法是利用如附图所示的包括有预报和校正步的一种迭代方法。如附图所示，预报步从已知的解(A)开始，以一个切线预报来估计对于一个规定负荷增长方式的解(B)；然后校正步，利用常规潮流解(系统负荷值不变)出准确解(C)；负荷进一步增加时，根据新的切线预报电压值。如果新的估计负荷超出了准确解的最大负荷，则以负荷值固定的校正计算就不会收敛，因此采用固定监控点电压的校正计算来求准确解(E)；当接近电压稳定极限时，为了确定准确的最大负荷，在连续预报期间，负荷增量应逐步减少。 

连续潮流的潮流方程的状态变量为[X，t]^T，即[V，θ，t]^T。在预报步，利用线性近似来估计下一步当状态变量(θ，V或t)中的一个变化后的解。对应于初始解状态变量，对(9)式取全微分： 

$d\left[H(X,t)\right]=\frac{\∂H}{\∂X}\mathrm{dX}+\frac{\∂H}{\∂t}\mathrm{dt}=0---\left(21\right)$

即  $\left[\begin{array}{cc}\frac{\∂H}{\∂X}& \frac{\∂H}{\∂t}\end{array}\right]\left[\begin{array}{c}\mathrm{dX}\\ \mathrm{dt}\end{array}\right]=0---\left(22\right)$

其中

为常规潮流方程的雅可比矩阵，  $\left[\begin{array}{c}\mathrm{dX}\\ \mathrm{dt}\end{array}\right]$ 为预报校正要求的切向量。 

注意，由于潮流方程中引入了参数t，因而增加了一个未知状态变量，待求量个数多一个。为了解这个潮流方程需增加一个方程，使切向量有确定的解。这可通过切向量的第 k个分量值dX_k设定为+1或-1来满足。这个分量称连续参数。这样方程(22)变为： 

$\left[\begin{array}{ccc}\frac{\∂H}{\∂X}& & \frac{\∂H}{\∂t}\\ & e_k& \end{array}\right]\left[\begin{array}{c}\mathrm{dX}\\ \mathrm{dt}\end{array}\right]=\left[\begin{array}{c}0\\ \±1\end{array}\right]---\left(23\right)$

$\left[\begin{array}{ccc}\frac{\∂H}{\∂X}& & \frac{\∂H}{\∂t}\\ & e_k& \end{array}\right]$

可称为扩展潮流方程的雅可比矩阵。 

这里e_k表示行向量，除了第k个元素等于1(对应于连续参数)，其他元素都等于零。 

由上可以看出，连续潮流方程类似于常规潮流方程，只是把负荷增加作为参数加进方程中。随着负荷的增加，相应增加发电机的发电量，重新形成新的潮流方程。 

起初，选择负荷参数t为连续参数，相应的切向量分量设定为+1.0。在随后的预报步中，由于知道了切向量，连续参数就被选择为状态变量。在给定解附近，它有最大的变化速率，其斜率符号决定相应切向量分量的符号。当达到最大负荷时，电压通常是具有最大变化的参数。 

一旦切向量找到了，下一步解的预报就可由式(24)给出： 

$\left[\begin{array}{c}X\′\\ t\′\end{array}\right]=\left[\begin{array}{c}{X}^{\left(0\right)}\\ t^{\left(0\right)}\end{array}\right]+\mathrm{\σ}\left[\begin{array}{c}\mathrm{dX}\\ \mathrm{dt}\end{array}\right]---\left(24\right)$

这里，  ${\left[\begin{array}{cc}{X}^{\′}& t^{\′}\end{array}\right]}^T$ 是预测值；[X⁽⁰⁾t⁽⁰⁾]^T是已知潮流解，即预报步开始时的状态变量值；σ是步长，其取值应使下一点的预测值落在其收敛半径内，即在规定的连续参数下潮流解存在。如果对给定的步长，在下一个校正步中得不到潮流解，则要减小步长。 

在校正步中，以预测值代入扩展潮流方程(23)作迭代求解。但这时也会遇到待求量个数比方程个数多1的问题，解决的办法是指定状态变量[X  t]^T的第k个分量x_k的值为预测值，从而增加一个方程，得到新的扩展的潮流方程： 

$\left[\begin{array}{ccc}\frac{\∂H}{\∂X}& & \frac{\∂H}{\∂t}\\ & x_k-\mathrm{\η}& \end{array}\right]\left[\begin{array}{c}\mathrm{dX}\\ \mathrm{dt}\end{array}\right]=\left[0\right]---\left(25\right)$

式(25)中，x_k是待求量[X   t]^T的第k个分量，因其值已被确定，故被称为连续参数；η为预测值[X′  t′]^T的第k个分量

将预测值[X′t′]^T作为初值代入扩展潮流方程(23)，用牛顿拉夫逊法迭代求解，可求得新潮流解。由于引入了规定为x_k的附加方程，使得在临 界运行点处雅可比矩阵非奇异，式(23)和(25)的求解都不会有数值计算的困难，从而能得到完整的P-V曲线。 

如果连续参数是负荷的增加，则校正步将是V-P平面上的垂直线(例如附图中的线段BC)。换句话说，如果电压幅值是一个连续参数，则校正步将是一条水平线(如附图中的线段DE)。 

2、基于线性模拟法的连续潮流法技术步骤 

算法流程图如附图1所示，介绍其技术步骤如下： 

1)输入电网参数； 

2)采用PSD-BPA潮流计算程序获得基本工况的潮流解，设负荷增长方向t＝0； 

3)形成S_PQ矩阵； 

4)计算  $\frac{d^{2}V}{d{P}_{\mathrm{PQ}}^2};$

5)判断是否大于0，若  $\frac{d^{2}V}{d{P}_{\mathrm{PQ}}^2}>0,$ 计算由于负荷扰动量产生的状态量改变量(Δθ_i)_K、(ΔV_i)_K，求得新的状态量(θ_i)_K+1，(V_i)_K+1，设t′＝t+Δt_K，转步骤3)；若  $\frac{d^{2}V}{d{P}_{\mathrm{PQ}}^2}<0,$ 设t＝t′； 

6)预报，由式(23)和(24)求得预测值[X′t′]^T； 

7)校正，将预测值[X′t′]^T带入式(25)，用牛顿法求解得到新潮流解。 

8)保持本次结算结果作为下次计算的初始工作点； 

9)判断dt是否大于0，若dt>0转步骤6)；否则计算终止并输出临界点的计算结果。 

因此，本发明提出一种用于静态电压稳定性分析的方法，其特征在于包括以下步骤： 

1)输入电网参数，包括相关负荷，发电机、线路、变压器等相关参数，利用潮流计算工具进行潮流计算，获得基本工况的潮流解，同时设负荷增长方向t＝0； 

2)形成S_PQ矩阵，该矩阵为n_PQ×n_PQ阶，它反映负荷节点有功、无功注入量的变化对负荷电压变化的灵敏度关系，在正常运行状态下，对角线各元素都为负值； 

3)计算

并判断

是否大于0，若  $\frac{d^{2}V}{d{P}_{\mathrm{PQ}}^2}>0,$ 计算由于负荷扰动量产生的状态量改变量(Δθ_i)_K、(ΔV_i)_K，求得新的状态量(θ_i)_K+1，(V_i)_K+1，设中间变量t′＝t+Δt_K，转步骤2)，若 设t＝t′，转下一步； 

4)预报，由下面两个式子求得预测值[X′t′]^T： 

$\left[\begin{array}{c}\begin{array}{cc}\frac{\&PartialD;H}{\&PartialD;X}& \frac{\&PartialD;H}{\&PartialD;t}\end{array}\\ e_k\end{array}\right]\left[\begin{array}{c}\mathrm{dX}\\ \mathrm{dt}\end{array}\right]=\left[\begin{array}{c}0\\ \&PlusMinus;1\end{array}\right]$

式中 

可称为扩展潮流方程的雅可比矩阵。 

$\left[\begin{array}{c}{X}^{\&prime;}\\ t^{\&prime;}\end{array}\right]=\left[\begin{array}{c}{X}^{\left(0\right)}\\ t^{\left(0\right)}\end{array}\right]+\mathrm{\&sigma;}\left[\begin{array}{c}\mathrm{dX}\\ \mathrm{dt}\end{array}\right]$

式中[X′t′]^T是预测值；[X⁽⁰⁾t⁽⁰⁾]^T是已知潮流解，即预报步开始时的状态变量值；σ是步长； 

5)校正，将预测值[X′t′]^T带入下式，用牛顿法求解得到新潮流解： 

$\left[\begin{array}{c}\begin{array}{cc}\frac{\&PartialD;H}{\&PartialD;X}& \frac{\&PartialD;H}{\&PartialD;t}\end{array}\\ x_k-\mathrm{\&eta;}\end{array}\right]\left[\begin{array}{c}\mathrm{dX}\\ \mathrm{dt}\end{array}\right]=\left[0\right]$

式中，x_k是待求量[X t]^T的第k个分量，η为预测值[X′t′]^T的第k个分量x′_k

6)保持本次计算结果作为下次计算的初始工作点； 

7)判断dt是否大于0，若dt＞0，则转到步骤4)，否则计算终止并输出临界点的计算结果。 

依据本发明的用于静态电压稳定性分析的方法具有以下优点：利用线性模拟法计算速度快和连续潮流法计算准确的特点，将两者结合起来，可以取得快速和准确地求出某种负荷增长方式下的电网电压稳定裕度。本发明的连续潮流法既适用于静态电压稳定离线分析也可适用于自动电压控制系统的在线应用。 

附图说明

图1是依据本发明的用于静态电压稳定性分析的方法的算法流程图； 

图2是依据本发明的用于静态电压稳定性分析的方法的算法示意图。 

具体实施方式

以重庆电网2008年某日实时数据断面为例，该断面共有305条母线，132条线路，91台变压器，负荷增长方式按照全网负荷等功率因数增长，所有发电机也按等功率因数增长参与调节。采用三种算法计算此断面的静态电压稳定裕度，算法1为标准连续潮流法，算法2为线性模拟法，算法3为基于线性模拟法的连续潮流法。计算结果见由表1。 

以基于线性模拟法的连续潮流(算法3)为例，首先输入电网参数，用PSD-BPA潮流计算方法获得基本工况的潮流解，设t＝0，形成S_PQ矩阵并计算

判断

是否大于0，  $\frac{d^{2}V}{d{P}_{\mathrm{PQ}}^2}>0,$ 计算由于负荷扰动量产生的状态量改变量(Δθ_i)_K、(ΔV_i)_K，求得新的状态量(θ_i)_K+1，(V_i)_K+1，设t′＝t+Δt_K，重新计算S_PQ矩阵并计算

判断

是否大于0，先后迭代5次后达到图2中A点，此时  $\frac{d^{2}V}{d{P}_{\mathrm{PQ}}^2}<0,$ 设t＝t′，转入预报步，由式(23)和(24)求得预测值[X′t′]^T；继续转入校正步，将预测值[X′t′]^T带入式(25)，用牛顿法求解得到新潮流解，保持本次结算结果作为下次计算的初始工作点，并判断dt是否大于0，dt>0重新进行预报步计算，迭代6次后，到达崩溃点附近，此时dt<0计算终止并输出临界点的计算结果。 

由表1可以看出，算法1的静态电压稳定裕度为2.352，算法2和算法3的静态电压稳定裕度2.355，计算结果基本一致。算法1计算时间最短，只要1.02秒，算法2计算时间最长，达15.23秒，算法3计算时间介于算法1和算法3之间，达4.16秒，可以满足在线运行的要求。 

表1  静态电压稳定裕度计算结果 

已经根据优选的实施例描述了本发明。显然，在阅读和理解了上述详细说明书后能做出多种修正和替换。本发明意欲的是本申请构建成包括了落入附属的权利要求书或其等同物的范围之内的所有这些修正和替换。

Claims (1)

1.一种用于静态电压稳定性分析的方法，其特征在于包括以下步骤：

1)输入电网参数，包括相关负荷、发电机、线路和变压器的相关参数，利用潮流计算工具进行潮流计算，获得基本工况的潮流解，同时设负荷增长方向t＝0；

2)形成S_PQ矩阵，该矩阵为n_PQ×n_PQ阶，它反映负荷节点有功、无功注入量的变化对负荷电压变化的灵敏度关系，在正常运行状态下，对角线各元素都为负值；

3)计算

并判断

是否大于0，若计算由于负荷扰动量产生的状态量改变量(Δθ_i)_K、(ΔV_i)_K，求得新的状态量(θ_i)_K+1，(V_i)_K+1，设中间变量t′＝t+Δt_K，转步骤2)，若

设t＝t′，转下一步；

为矩阵S_PQ对角线元素对于有功量的导数；(θ_i)_K为节点i的角度状态量；(V_i)_K为节点i的电压幅值状态量；(Δθ_i)_K为节点i的角度状态量的改变量；(ΔV_i)_K为节点i的电压幅值状态量的改变量；(θ_i)_K+1为节点i的状态量改变后新的角度状态量；(V_i)_K+1为节点i的状态量改变后新的电压幅值状态量；Δt_K为给定步长参数；

4)预报，由下面两个式子求得预测值[X′ t′]^T

$\left[\begin{array}{ccc}\frac{\&PartialD;H}{\&PartialD;X}& & \frac{\&PartialD;H}{\&PartialD;t}\\ & e_k& \end{array}\right]\left[\begin{array}{c}\mathrm{dX}\\ \mathrm{dt}\end{array}\right]=\left[\begin{array}{c}0\\ \&PlusMinus;1\end{array}\right]$

式中 $\left[\begin{array}{ccc}\frac{\&PartialD;H}{\&PartialD;X}& & \frac{\&PartialD;H}{\&PartialD;t}\\ & e_k& \end{array}\right]$ 可称为扩展潮流方程的雅可比矩阵，

$\left[\begin{array}{c}{X}^{\&prime;}\\ t^{\&prime;}\end{array}\right]=\left[\begin{array}{c}{X}^{\left(0\right)}\\ t^{\left(0\right)}\end{array}\right]+\mathrm{\&sigma;}\left[\begin{array}{c}\mathrm{dX}\\ \mathrm{dt}\end{array}\right]$

式中[X′ t′]^T是预测值； ${\left[\begin{array}{cc}{X}^{\left(0\right)}& t^{\left(0\right)}\end{array}\right]}^T$ 是已知潮流解，即预报步开始时的状态变量值；σ是步长；

5)校正，将预测值[X′ t′]^T带入下式，用牛顿法求解得到新潮流解：

$\left[\begin{array}{ccc}\frac{\&PartialD;H}{\&PartialD;X}& & \frac{\&PartialD;H}{\&PartialD;t}\\ & x_k-\mathrm{\&eta;}& \end{array}\right]\left[\begin{array}{c}\mathrm{dX}\\ \mathrm{dt}\end{array}\right]=\left[0\right]$

式中，x_k是待求量[X t]^T的第k个分量，η为预测值[X′ t′]^T的第k个分量x′_k；

6)保持本次计算结果作为下次计算的初始工作点；

7)判断dt是否大于0，若dt＞0，则转到步骤4)，否则计算终止并输出临界点的计算结果。

Images