CN102013680B - 一种电力系统快速分解法潮流计算方法 - Google Patents

Abstract

本发明公开了一种电力系统快速分解法潮流计算方法，主要包括以下步骤：原始数据输入和电压初始化；形成导纳矩阵；形成修正方程系数矩阵B′和B″并进行因子表分解；进行P～θ迭代，修正电压相角；进行Q～V迭代，修正电压幅值；判断迭代是否收敛；计算节点功率和支路功率。本发明要求同一次迭代中P～θ迭代和Q～V迭代都收敛，迭代过程才结束，使算法框图更加简单、流程更加清晰。本发明不采用稀疏矩阵技术，矩阵元素存取和计算方便，编程简单；修正方程的系数矩阵按n阶存储，避免了节点换号，降低了编程难度；通过合理的逻辑判断减少计算量，明显提高了计算速度，完全能够满足科研的需要。本发明还可以处理多个平衡节点的电力系统。

Description

一种电力系统快速分解法潮流计算方法

技术领域

本发明涉及一种电力系统潮流计算方法，特别是一种电力系统快速分解法潮流计算方法。

背景技术

潮流计算是电力系统最基本的一种计算项目，也是电力系统其他分析项目如静态安全分析等项目的基础。快速分解法潮流计算方法是目前最常用且最快速的潮流计算方法，该方法采用稀疏矩阵技术和节点优化编号等高级技术，主要流程如下：

A、原始数据输入和电压初始化；

潮流计算的原始数据主要有节点数、支路数、最大迭代次数、收敛精度；支路参数数据，包括输电线路的电阻、电抗、对地电纳或变压器支路的电阻、电抗、变比；节点电压及功率数据，包括节点的类型，发电机的有功功率和无功功率或负荷的有功功率和无功功率、节点给定电压或电压初值；无功补偿数据等。

电力系统潮流计算的节点分为3种类型：PQ节点，此节点的有功功率P和无功功率Q已知，电压幅值V和电压相角θ待求；PV节点，此节点的有功功率P和电压幅值V已知，无功功率Q和电压相角θ待求；平衡节点，此节点的电压幅值V和电压相角θ已知，有功功率P和无功功率Q待求。

电压初始化一般采用平启动，即PV节点和平衡节点的电压幅值取给定值，PQ节点的电压幅值取1.0；所有电压的相角都取0.0。这里电气量的单位采用标幺值。

B、节点优化编号；

潮流计算的修正方程的系数矩阵是稀疏矩阵，可以采用稀疏矩阵技术，节点的编号顺序直接影响系数矩阵消去过程中矩阵的稀疏度，对计算速度影响很大，应该对节点进行优化编号。目前主要有静态节点优化编号、半动态节点优化编号和动态节点优化编号3种编号方法。

C、形成导纳矩阵；

导纳矩阵为：

$Y=\left[\begin{array}{ccccc}{Y}_{11}& Y_{12}& Y_{13}& \·\·\·& Y_{1n}\\ Y_{21}& Y_{22}& Y_{23}& \·\·\·& Y_{2n}\\ Y_{31}& Y_{32}& Y_{33}& \·\·\·& Y_{3n}\\ \·& \·& \·& & \·\\ \·& \·& \·& & \·\\ \·& \·& \·& & \·\\ Y_{n1}& Y_{n2}& Y_{n3}& \·\·\·& Y_{\mathrm{nn}}\end{array}\right]$

导纳矩阵为n阶方阵，式中的n为电力网络的节点数，可以反映出节点间的连接关系，如果两个节点不直接连接，它们的互导纳为0。电力网络中每个节点只与几个节点相连，因此导纳矩阵是一个高度稀疏的矩阵，它的大部分元素为0。

对稀疏矩阵中零元素的存储和运算都是多余的，可以采用稀疏矩阵技术免除对零元素的存储和运算。

D、形成按稀疏矩阵技术存储的修正方程系数矩阵B′和B″，并利用稀疏矩阵运算技术对B′和B″进行因子表分解；

潮流计算的功率方程为：

$P_i=V_i\underset{j\∈i}{\mathrm{\Σ}}{V}_j(G_{\mathrm{ij}}\cos {\mathrm{\θ}}_{\mathrm{ij}}+B_{\mathrm{ij}}\sin {\mathrm{\θ}}_{\mathrm{ij}}),i=1,\·\·\·,n---\left(1\right)$

$Q_i=V_i\underset{j\∈i}{\mathrm{\Σ}}{V}_j(G_{\mathrm{ij}}\sin {\mathrm{\θ}}_{\mathrm{ij}}-B_{\mathrm{ij}}\cos {\mathrm{\θ}}_{\mathrm{ij}}),i=1,\·\·\·,n---\left(2\right)$

式中，V_i和θ_i分别为节点i的电压幅值和电压相角；G和B分别为导纳矩阵的电导部分和电纳部分；j∈i表示节点j与节点i直接相连，n为节点数。

方程(1)和(2)是非线性方程组，通常采用逐次线性化方法迭代求解，先设定电压初值，此电压值一般不可能恰好是方程的解，因此给定功率和计算功率就存在功率偏差，功率偏差方程为：

$\mathrm{\Δ}{P}_i=P_{\mathrm{iS}}-P_i=P_{\mathrm{iS}}-V_i\underset{j\∈i}{\mathrm{\Σ}}{V}_j(G_{\mathrm{ij}}\cos {\mathrm{\θ}}_{\mathrm{ij}}+B_{\mathrm{ij}}\sin {\mathrm{\θ}}_{\mathrm{ij}}),i=1,\·\·\·,n-1---\left(3\right)$

$\mathrm{\Δ}{Q}_i=Q_{\mathrm{iS}}-Q_i=Q_{\mathrm{iS}}-V_i\underset{j\∈i}{\mathrm{\Σ}}{V}_j(G_{\mathrm{ij}}\sin {\mathrm{\θ}}_{\mathrm{ij}}-B_{\mathrm{ij}}\cos {\mathrm{\θ}}_{\mathrm{ij}}),i=1,\·\·\·,m---\left(4\right)$

式中，P_iS和Q_iS分别为节点i的给定有功功率和给定无功功率，m为PQ节点数。

因为平衡节点的有功功率不是给定的，因此式(3)中不包括平衡节点；同理平衡节点和PV节点的无功功率也不是给定的，因此式(4)中不包括平衡节点和PV节点。

对方程(3)和(4)线性化，并结合电力系统的特点进行改进后得到的快速分解法修正方程为：

B′Δθ＝ΔP/V                                   (5)

B″ΔV＝ΔQ/V                                    (6)

式中，ΔP/V和ΔQ/V分别为有功功率和无功功率不平衡量除以电压幅值后的列向量；ΔV和Δθ分别为电压幅值和电压相角修正量列向量；B′为(n-1)阶方阵，包含与除平衡节点外所有节点有关的行和列；B″为m阶方阵，仅包括与PQ节点有关的行和列。

B″为导纳矩阵的虚部，B′为导纳矩阵的虚部，但去掉了支路电阻、对地导纳和非标准变比的影响。

快速分解法潮流计算的修正方程的系数矩阵是稀疏矩阵，可以采用稀疏矩阵技术免除对零元素的存储和运算，为了进一步提高效率，同时采用节点优化编号。但稀疏矩阵技术和节点优化编号编程比较复杂，矩阵元素存取和运算都不直观，不利于对程序维护和修改。

E、设置迭代计数初值；

F、设置收敛标志；

G、计算有功功率不平衡量ΔP，解方程B′Δθ＝ΔP/V，修正电压相角；

求解方程(5)，得到Δθ，对上次迭代的电压相角按下式进行修正：

θ^(k)＝θ^(k-1)-Δθ^(k-1)                             (7)

式中，上标(k)表示第k次迭代的值。

方程(5)和(6)都是线性方程组，可以用高斯消去法求解。

设线性方程一般形式如下：

AX＝B                                 (8)

高斯消去法求解线性方程组包括3部分，分别为：

a)对系数矩阵的消去，消去过程可以采用行消去也可以采用列消去，列消去的公式如下：

$a_{\mathrm{ij}}^{\left(k\right)}=a_{\mathrm{ij}}^{(k-1)}-a_{\mathrm{ik}}^{(k-1)}\·a_{\mathrm{kj}}^{(k-1)}/a_{\mathrm{kk}}^{(k-1)},j=k+1,\·\·\·,n,i=k+1,\·\·\·,n---\left(9\right)$

b)对右端量的前代，公式如下：

$b_i^{\left(k\right)}=b_i^{(k-1)}-a_{\mathrm{ik}}^{(k-1)}\·b_k^{(k-1)}/a_{\mathrm{kk}}^{(k-1)},i=k+1,\·\·\·,n---\left(10\right)$

c)回代求解，公式如下：

$x_i=(b_i^{(i-1)}-\underset{j=i+1}{\overset{n}{\mathrm{\Σ}}}{a}_{\mathrm{ij}}^{(i-1)}\·x_j)/a_{\mathrm{ii}}^{(i-1)},i=n,\·\·\·,1---\left(11\right)$

快速分解法的系数矩阵在迭代过程中是不变的，只是方程右端量是变化的，因此可以先对系数矩阵进行消去，结果保存起来，以后各次迭代中反复使用。所保存的系数矩阵的消去结果也称为因子表。

H、判断有功功率最大不平衡量ΔP_m是否小于收敛精度；

I、计算无功功率不平衡量ΔQ，解方程B″ΔV＝ΔQ/V，修正电压幅值；

求解方程(6)，得到ΔV，对上次迭代的电压幅值按下式进行修正：

V^(k)＝V^(k-1)-ΔV^(k-1)                                         (12)

J、判断无功功率最大不平衡量ΔQ_m是否小于收敛精度；

K、计算平衡节点功率及PV节点的无功功率，计算支路功率。

按式(1)和(2)计算平衡节点功率及PV节点的无功功率。

上述方法虽然能大幅度提高了潮流计算速度、降低内存占用量，但算法采用稀疏矩阵技术和节点优化编号，编程非常麻烦、不易修改和维护。并且快速分解法每次迭代过程包括P～θ迭代，即解方程B′Δθ＝ΔP/V和Q～V迭代即解方程B″ΔV＝ΔQ/V两个半次迭代，只要相邻的两个半次迭代都收敛，迭代过程就结束了。相邻的两个半次迭代都收敛，可以是第k次ΔP_m和ΔQ_m都小于ε，也可以是第k-1次ΔQ_m和第k次ΔP_m都小于ε，只要相邻即可，从图1看，为了考虑这两种情况，程序主流程比较复杂。这样处理不适合科研人员以快速分解法潮流计算为基础进行进一步研究的需要。

如果不采用稀疏矩阵技术，速度很慢；另外B′不包括平衡节点有关的行和列，为(n-1)阶方阵；B″仅包括PQ节点有关的行和列，为m阶方阵；B′的阶数和B″的阶数都与导纳矩阵的阶数不同，也与节点号不对应，编程时需要换号，加大了程序的复杂程度，增加了修改和维护难度，也无法满足科研人员的要求。

发明内容

为解决现有技术存在的上述问题，本发明要设计一种易于编程、修改和调试的电力系统快速分解法潮流计算方法。

为了实现上述目的，本发明的技术方案如下：一种电力系统快速分解法潮流计算方法，包括以下步骤：

A、原始数据输入和电压初始化；

B、形成导纳矩阵；

C、形成n阶的修正方程系数矩阵B′和B″，并进行因子表分解；

D、设置迭代计数初值；

E、计算有功功率不平衡量ΔP，解方程B′Δθ＝ΔP/V，修正电压相角；计算无功功率不平衡量ΔQ，解方程B″ΔV＝ΔQ/V，修正电压幅值；

F、判断P～θ迭代和Q～V迭代是否都收敛；

G、计算平衡节点的功率及PV节点的无功功率，计算支路功率。

步骤C所述的系数矩阵B′和B″都按n阶方阵存储，B′中与平衡节点对应的行和列的元素都为0，B″中与PV节点或平衡节点对应的行和列的元素都为0；形成修正方程系数矩阵B′和B″的步骤如下：

C1、形成修正方程的系数矩阵B′的步骤如下：

C11：不考虑节点的类型，形成系数矩阵B′；

C12：设置当前节点号i＝1；

C13：判断节点i是否为平衡节点，如果不是平衡节点，则转至C15；

C14：把系数矩阵B′第i行和第i列的所有元素置为0；

C15：当前节点号i＝i+1；

C16：判断i是否大于节点数n，如果i大于n，结束，否则转至C13；

C2、形成修正方程的系数矩阵B″的步骤如下：

C21：不考虑节点的类型，取导纳矩阵虚部形成系数矩阵B″；

C22：设置当前节点号i＝1；

C23：判断节点i是否为PQ节点，如果是PQ节点，则转至C25；

C24：把系数矩阵B″第i行和第i列的所有元素置为0；

C25：当前节点号i＝i+1；

C26：判断i是否大于节点数n，如果i大于n，结束，否则转至C23；

C3、形成因子表：消去过程采用列消去，从式(9)可以看出，如果则第k行对第i行的消去就没有必要了，可以通过逻辑判断避免不必要的计算；具体步骤如下：

C31：设置当前行号k＝1；

C32：判断矩阵元素a_kk是否为0，如果a_kk为0，则转至C311；

C33：设置i＝k+1；

C34：判断i是否大于n，如果i大于n，则转至C311；

C35：判断矩阵元素a_ik是否为0，如果a_ik为0，则转至C310；

C36：设置当前列号j＝k+1；

C37：判断j是否大于n，如果j大于n，则转至C310；

C38：根据式a_ij＝a_ij-a_ik a_kj/a_kk对系数矩阵进行消去运算；

C39：令j＝j+1，返回到C37；

C310：令i＝i+1，返回到C34；

C311：令k＝k+1；

C312：判断k是否小于n，如果k小于n，则返回到C32；否则结束。

步骤E所述的修正方程B′Δθ＝ΔP/V和B″ΔV＝ΔQ/V的求解采用因子表技术，利用步骤C形成的因子表对方程的右端向量进行前代和回代就能得到方程的解；本发明中系数矩阵B′和B″都按n阶方阵存储，因此对B′和B″进行因子表分解得到的矩阵中也包含元素都为0的行和列，在前代和回代过程中必须跳过这些行；本发明通过判断矩阵的对角元素是否为0的方法跳过元素都为0的行；前代和回代过程的步骤如下：

E1、前代过程的具体步骤如下：

E11：设置当前行号k＝1；

E12：判断矩阵元素a_kk是否为0，如果a_kk为0，则转至E17；

E13：设置i＝k+1；

E14：判断i是否大于n，如果i大于n，则转至E17；

E15：根据式b_i＝b_i-a_ik b_k/a_kk进行前代运算；

E16：令i＝i+1，返回到E14；

E17：令k＝k+1；

E18：判断k是否小于n，如果k小于n，则返回到E12；否则结束；

E2、回代过程的具体步骤如下：

E21：设置当前行号i＝n；

E22：判断矩阵元素a_ii是否为0，如果a_ii为0，则转至E28；

E23：设置j＝i+1；

E24：判断j是否大于n，如果j大于n，则转至E27；

E25：根据式b_i＝b_i-a_ijb_j进行回代运算；

E26：令j＝j+1，返回到E24；

E27：b_i＝b_i/a_ii；

E28：令i＝i-1；

E29：判断i是否小于1，如果i不小于1，则返回到E22；否则结束。

与现有技术相比，本发明具有以下有益效果：

1、为了简化过程，本发明对现有方法的流程图进行了修改，要求同一次迭代中ΔP_m和ΔQ_m都小于收敛精度ε，迭代过程才结束。与现有的方法相比，可能要多一次Q～V迭代的工作量，但算法框图简单多了，流程更加清晰。并且快速分解法主要计算量在因子表分解上，多一次Q～V迭代的工作量对整个潮流计算的计算时间影响不大。

2、本发明不采用稀疏矩阵技术和节点优化编号技术，简化了编程，使得矩阵元素存取和运算都非常直观，有利于对程序维护和修改。

3、为了降低程序难度，本发明的系数矩阵B′和B″都按n阶方阵存储，B′中与平衡节点对应的行和列的元素都为0，B″中与PV节点或平衡节点对应的行和列的元素都为0，形成因子表时可以通过判断语句跳过这些行。这样程序就可以避免节点换号了，降低了程序的复杂程度。

4、本发明形成系数矩阵B′及其因子表时，并没有限制平衡节点的个数，也不要求对平衡节点特殊编号，因此可以处理多个平衡节点的电力系统。

5、本发明形成因子表时，通过判断语句跳过了元素都为0的行和不需要进行消去运算的行，明显提高了计算速度。

6、本发明提出了一种适合科研人员使用的快速分解法潮流计算方法。本发明通过避免某些不必要的计算来提高潮流计算的速度，效果非常明显。本发明还对快速分解法的主程序结构进行了调整，使得结构更加清晰，更利于理解和维护。本发明具有概念清晰、编程简单、性能良好的特点。本发明方法相比实用算法内存占用量较多，但完全能够接受；计算速度也有所降低，但速度也非常快，完全能够满足科研的需要。

附图说明

本发明共有附图10张，其中：

图1是现有技术的快速分解法潮流计算流程图。

图2是本发明的快速分解法潮流计算流程图。

图3是本发明的修正方程系数矩阵B′的计算流程图。

图4是本发明的修正方程系数矩阵B″的计算流程图。

图5是本发明的因子表的计算流程图。

图6是本发明的修正方程前代过程的流程图。

图7是本发明的修正方程回代过程的流程图。

图8是本发明方法实施例1的网络图。

图9是输电线路的等值电路示意图。

图10是变压器支路的等值电路示意图。

具体实施方式

下面结合附图对本发明的流程及有益效果进行进一步地描述。采用本发明对图8所示的一个简单5节点电力系统算例进行了计算，收敛精度为0.00001，迭代次数为10次。

图8所示的简单电力系统网络，包括2个发电机，3个负荷，2条变压器支路，3条输电线路。该网络共5条母线或节点，其中母线2为平衡节点，此母线的电压幅值和相角给定，母线4为PV节点，此母线的有功功率和电压幅值给定，其它3条母线为PQ节点，这些母线的有功功率和无功功率给定。算例中各电气量的单位除特别说明外均为标幺值。下面以简单5节点网络作为算例来说明本发明的计算过程：

A、原始数据输入和电压初始化

支路参数见表1，母线数据见表2。表2中发电机功率填正数，负荷填负数。

电压初值采用平启动，母线1、母线3、母线5的电压幅值初值为1.0，母线2和母线4的电压幅值初值为1.05，母线1～母线5的电压相角初值都为0.0。

表1 图8实施例的支路参数

表2 图8实施例的母线数据

B、形成导纳矩阵

导纳矩阵如表3：

表3 导纳矩阵元素

C、形成n阶的修正方程系数矩阵B′和B″，并进行因子表分解

1)形成n阶的修正方程系数矩阵B′

首先忽略支路电阻、对地导纳和变压器非标准变比，仅考虑支路电抗形成导纳矩阵的虚部，见表4，然后根据图3所示流程把平衡节点2对应的第2行和第2列的所有元素都置成0，形成修正方程系数矩阵B′，见表5。

表4 仅考虑支路电抗形成的导纳矩阵虚部

表5 修正方程系数矩阵B′

2)对矩阵B′进行因子表分解

根据图5所示流程对表5所示的系数矩阵B′进行因子表分解的过程如下：

首先用第1行对第2～5行进行消去运算。由于a₂₁和a₄₁为0，不必要对第2行和第4行进行消去运算；a₃₁和a₅₁不为0，需要对第3行和第5行进行消去运算。运算结果见表6。注：a₃₁和a₅₁没有必要变成0。

表6 用B′矩阵的第1行进行高斯列消去的结果

第2行的对角元素a₂₂为0，跳过去。

用第3行对第4～5行进行消去运算。由于a₄₃为0，不必要对第4行进行消去运算；a₅₃不为0，需要对第5行进行消去运算。运算结果见表7。

表7 用B′矩阵的第3行进行高斯列消去的结果

用第4行对第5行进行消去运算。由于a₅₄不为0，需要对第5行进行消去运算，运算结果见表8。表8为最终的因子表。

表8 用B′矩阵的第4行进行高斯列消去的结果

3)形成n阶的修正方程系数矩阵B″

首先取得导纳矩阵的虚部，见表9，然后根据图4所示流程把平衡节点2对应的第2行和第2列以及PV节点4对应的第4行和第4列的所有元素都置成0，形成修正方程系数矩阵B″，见表10。

表9 导纳矩阵的虚部

表10 修正方程系数矩阵B″

4)对矩阵B″进行因子表分解

根据图5所示流程对表10所示的系数矩阵B″进行因子表分解的过程如下：

首先用第1行对第2～5行进行消去运算。由于a₂₁和a₄₁为0，不必要对第2行和第4行进行消去运算；a₃₁和a₅₁不为0，需要对第3行和第5行进行消去运算。运算结果见表11。

表11 用B″矩阵的第1行进行高斯列消去的结果

第2行的对角元素a₂₂为0，跳过去。

用第3行对第4～5行进行消去运算。由于a₄₃为0，不必要对第4行进行消去运算；a₅₃不为0，需要对第5行进行消去运算。运算结果见表12。

表12 用B″矩阵的第3行进行高斯列消去的结果

第4行的对角元素a₄₄为0，跳过去。表12为最终的因子表。

D、设置迭代计数初值

设置迭代计数初值k＝0。

E、计算有功功率不平衡量ΔP，解方程B′Δθ＝ΔP/V，修正电压相角；计算无功功率不平衡量ΔQ，解方程B″ΔV＝ΔQ/V，修正电压幅值；

1)首次P～θ迭代

计算有功功率不平衡量ΔP，方程B′Δθ＝ΔP/V前代过程见表13。

表13 首次P～θ迭代计算结果

回代得到方程的解及电压相角修正结果见表14。

表14 电压相角修正结果

2)首次Q～V迭代

计算无功功率不平衡量ΔQ，方程B″ΔV＝ΔQ/V前代过程见表15。

表15 首次Q～V迭代计算结果

回代得到方程的解及电压幅值修正结果见表16。

表16 电压幅值修正结果

F、判断P～θ迭代和Q～V迭代是否收敛

从表13可见有功功率最大不平衡量ΔP_m＝5.0000，从表15可见无功功率最大不平衡量ΔQ_m＝5.5672，均大于收敛精度，不收敛，应该继续迭代。

G、计算平衡节点功率及PV节点的无功功率，计算支路功率

经过10次迭代，潮流计算收敛。节点功率和电压见表17，支路功率计算结果见表18。

表17 图8实施例的节点电压和功率计算结果

表18 图8实施例的支路传输功率计算结果

本方法可以采用任何一种编程语言和编程环境实现，如C语言、C++、FORTRAN、Delphi等。开发环境可以采用Visual C++、Borland C++ Builder、Visual FORTRAN等。

本方法的另一个实施例为一个445节点大型实际电网，潮流计算的收敛精度为0.0001，计算时间为0.156s，计算环境为主频1.10GHz的Intel Pentium的PC机。作为对比，采用4种方法对本算例进行了计算，收敛精度均为0.0001。

4种方法分别如下：

方法1：快速分解法，不采用稀疏矩阵技术；

方法2：快速分解法，不采用稀疏矩阵技术，但通过合理的逻辑判断避免不必要的计算，是本发明方法；

方法3：快速分解法，采用稀疏矩阵技术和节点优化编号技术；

方法4：牛顿法，采用稀疏矩阵技术和节点优化编号技术。

4种方法计算结果如表19：

表19 大型实际电网实施例的计算结果

从表19可见，本发明方法的快速分解法不采用稀疏矩阵技术但通过合理的逻辑判断避免不必要计算，它的计算时间仅为0.156s，与方法4(即采用稀疏矩阵技术的牛顿法)一样快；与不采用稀疏矩阵技术的方法1相比，速度快得多。

从表19可见，不采用稀疏矩阵技术，内存需求量较大，是方法3的90倍、方法4的31倍。即便如此，内存需求量也不到5MB内存，对于目前几个GB的计算机内存和如此大规模的算例来说，内存需求量完全可以接受。

Claims (1)

1.一种电力系统快速分解法潮流计算方法，其特征在于：包括以下步骤：

A、原始数据输入和电压初始化；

B、形成导纳矩阵；

C、形成n阶的修正方程系数矩阵B′和B″，并进行因子表分解；

D、设置迭代计数初值；

E、计算有功功率不平衡量ΔP，解方程B′Δθ＝ΔP/V，修正电压相角；计算无功功率不平衡量ΔQ，解方程B″ΔV＝ΔQ/V，修正电压幅值；

式中：ΔV和Δθ分别为电压幅值和电压相角修正量列向量，V为电压幅值；

F、判断P～θ迭代和Q～V迭代是否都收敛；

P～θ迭代，即解方程B′Δθ＝ΔP/V和Q～V迭代即解方程B″ΔV＝ΔQ/V两个半次迭代，只要相邻的两个半次迭代都收敛，迭代过程就结束了；相邻的两个半次迭代都收敛，第k次ΔPm和ΔQm都小于ε或第k–1次ΔQm和第k次ΔPm都小于ε，只要相邻即可；ΔPm为有功功率最大不平衡量，ΔQm为无功功率最大不平衡量；

G、计算平衡节点的功率及PV节点的无功功率，计算支路功率；

步骤C所述的系数矩阵B′和B″都按n阶方阵存储，B′中与平衡节点对应的行和列的元素都为0，B″中与PV节点或平衡节点对应的行和列的元素都为0；形成修正方程系数矩阵B′和B″的步骤如下：

C1、形成修正方程的系数矩阵B′的步骤如下：

C11：不考虑节点的类型，形成系数矩阵B′；

C12：设置当前节点号i＝1；

C13：判断节点i是否为平衡节点，如果不是平衡节点，则转至C15；

C14：把系数矩阵B′第i行和第i列的所有元素置为0；

C15：当前节点号i＝i+1；

C16：判断i是否大于节点数n，如果i大于n，结束，否则转至C13；

C2、形成修正方程的系数矩阵B″的步骤如下：

C21：不考虑节点的类型，取导纳矩阵虚部形成系数矩阵B″；

C22：设置当前节点号i＝1；

C23：判断节点i是否为PQ节点，如果是PQ节点，则转至C25； 

C24：把系数矩阵B″第i行和第i列的所有元素置为0；

C25：当前节点号i＝i+1；

C26：判断i是否大于节点数n，如果i大于n，结束，否则转至C23；

C3、形成因子表：设线性方程AX＝B，方程系数矩阵A的第i行第j列元素用a_ij表示，方程右端列向量B的第i个元素用b_i表示，形成因子表就是对线性方程AX＝B的系数矩阵A进行高斯消去的过程，消去过程采用列消去，第k行对第i行进行消去的公式为a_ij＝a_ij-a_ika_kj/a_kk，式中j＝k+1,…,n；如果a_ik＝0，则第k行对第i行的消去就没有必要了，因此通过逻辑判断避免不必要的计算；具体步骤如下：

C31：设置当前行号k＝1；

C32：判断矩阵元素a_kk是否为0，如果a_kk为0，则转至C311；

C33：设置i＝k+1；

C34：判断i是否大于n，如果i大于n，则转至C311；

C35：判断矩阵元素a_ik是否为0，如果a_ik为0，则转至C310；

C36：设置当前列号j＝k+1；

C37：判断j是否大于n，如果j大于n，则转至C310；

C38：根据式a_ij=a_ij-a_ika_kj/a_kk对系数矩阵进行消去运算；

C39：令j＝j+1，返回到C37；

C310：令i＝i+1，返回到C34；

C311：令k＝k+1；

C312：判断k是否小于n，如果k小于n，则返回到C32；否则结束;

步骤E所述的修正方程B′Δθ＝ΔP/V和B″ΔV＝ΔQ/V的求解采用因子表技术，利用步骤C形成的因子表对方程的右端向量进行前代和回代就能得到方程的解；系数矩阵B′和B″都按n阶方阵存储，因此对B′和B″进行因子表分解得到的矩阵中也包含元素都为0的行和列，在前代和回代过程中必须跳过这些行；通过判断矩阵的对角元素是否为0的方法跳过元素都为0的行；前代和回代过程的步骤如下：

E1、前代过程的具体步骤如下：

E11：设置当前行号k＝1；

E12：判断矩阵元素a_kk是否为0，如果a_kk为0，则转至E17；

E13：设置i＝k+1； 

E14：判断i是否大于n，如果i大于n，则转至E17；

E15：根据式b_i=b_i-a_ikb_k/a_kk进行前代运算；

E16：令i＝i+1，返回到E14；

E17：令k＝k+1；

E18：判断k是否小于n，如果k小于n，则返回到E12；否则结束；

E2、回代过程的具体步骤如下：

E21：设置当前行号i＝n；

E22：判断矩阵元素a_ii是否为0，如果a_ii为0，则转至E28；

E23：设置j＝i+1；

E24：判断j是否大于n，如果j大于n，则转至E27；

E25：根据式b_i＝b_i–a_ijb_j进行回代运算；

E26：令j=j+1，返回到E24；

E27：b_i=b_i/a_ii；

E28：令i＝i–1；

E29：判断i是否小于1，如果i不小于1，则返回到E22；否则结束。 

Images