对角元素按节点编号顺序存放在对角元素结构体数组中：
structYii_Type
{
doubleG,B;
}Yii[300];
其中G存放对角元素的实部，B存放对角元素的虚部，每个数组的元素个数与系统节点数相等。
由于上三角矩阵中非对角元素和系统中不接地支路一一对应，非对角元素的个数等于网络中不接地支路数Nb。
为了节约内存和提高计算速度，在计算机内存中只储存非零元素，把非零非对角元素“挤实”在一起。
为了识别非对角元素的行号和列号，我们在每个元素后存放相应的列下标。
非对角元素结构体数组定义为：
structYij_Type
{
doubleG,B;
intj;
}Yij[400];
按照这样的排列，取一个互导纳，就可以同时把该元素的列号取出来，为了判断该元素的行号，需要借助于数组NYseq[300]。
数组NYseq[300]按导纳矩阵行号的顺序存放各行非对角元素的首地址（事实上，存放的是各行第一个非对角元素在导纳矩阵非零非对角元素中的顺序号）。
本程序中还定义导纳矩阵中各行非对角元素的个数NYsum[300]。
一般我们有
NYsum[i]=NYseq[i+1]-NYseq[i]（3-7）
由于非对角元素是逐行向下排列的，所以就很容易判断出各非对角元素的行号。
按照对支路原始数据处理的要求，可以得出支路排列的顺序和导纳矩阵非对角元素的排列顺序完全一样。
因此只要顺序取出支路数据，按照式（3-1）求出倒数取负号之后，连同该支路的大节点号（即列下标）顺序送入Yij数组，就形成了导纳矩阵的上三角部分。
3.3导纳矩阵形成过程及框图
在本程序中，适应P-Q分解法的需要，导纳矩阵分为两步形成。
第一步只用不接地支路构成导纳矩阵，不考虑接地支路（包括变压器非标准变比）的影响。
这里同时形成两个导纳矩阵，即通常意义上的系数矩阵Y和不考虑输电线路电阻的系统导纳矩阵Y
，以适应BX法的要求，这两个导纳矩阵实际上只是半成品。
Y主要用来形成BX法所要求的第一个因子表，当该因子表形成后，就在半成品的基础上把接地支路及变压器非标准变比的影响加进去，形成完整的系统导纳矩阵。
其中Y
用来形成BX法第二个因子表，而Y将在整个迭带求解过程及线路潮流计算过程中发挥作用。
只考虑不接地支路构成导纳矩阵的程序框图如图3-2所示。
整个形成过程需要把不接地支路扫描一遍，对每条不接地支路作两方面的工作。
首先把阻抗求倒数并取负号后连同大节点号送到导纳矩阵非对角元素数组Yij（对应于Y）和Yij1（对应于Y
）中形成非对角元素，然后把阻抗的倒数累加到该支路两端节点的自导纳上去[见式（3-2）]。
为了累加形成对角元素，在计算开始时应对数组Yii和Yii1清零[见图3-2中⑴、⑵框]。
⑵框中NYsum为临时工作数组，定义为NYsum[300]，在其中累计导纳矩阵各行非对角元素的个数，因此也需要预先清零。
由于两个导纳矩阵的结构是相同的，共用一个NYsum数组。
图3-2形成不接地支路的导纳矩阵
对不接地支路的扫描用一个循环语句来完成[图中⑶框]。
⑷框把支路的有关数据送进中间工作单元，因为支路为变压器节点号可能为负，所以在这里对节点号取绝对值。
在⑸框中，把阻抗的倒数即支路导纳放到中间工作单元Gij、Bij中，而把支路电抗的倒数放到中间工作单元b_ij中。
⑺框判断支路是否为变压器支路。
若为变压器支路，则导纳需除以变压器非标准变比后再取负号送到导纳矩阵非对角元素数组Yij和Yij1中；否则直接取负号送到导纳矩阵非对角元素数组中[⑻框]。
⑼框向Yij和Yij1数组送列号。
这样就把支路阻抗数据变成了导纳矩阵的上三角部分。
在⑽～⑿框中根据式（3-2）累计有关节点的自导纳。
在⒀框中统计小节点号的不接地支路数目，这也就是导纳矩阵上三角部分每行非对角元素的个数。
至此，完成了一条不接地支路的处理；当循环由1做到Nb时形成了只考虑不接地支路的导纳矩阵。
⒁～⒃框是由NYsum数组根据式（3-7）形成NYseq数组。
3.4追加接地支路的程序框图
图3-3追加对地支路框图
追加接地支路包括两部分内容，即追加对地电容支路和考虑变压器非标准变比的影响，其程序框图如图3-3所示。
整个计算过程需要对支路数据再进行依次扫描，扫描是由一个循环语句来控制[图中⑴框]。
在⑵框中把支路有关数据送入中间工作单元。
⑵框根据节点号i、j的符号判断所取的支路是输电线路还是变压器支路。
当i、j中任一个为负时，为变压器支路，否则为输电线路。
当所取的支路为输电线路是，转入⑿～⒁框，向相应的节点累计自导纳部分。
当所取支路为变压器支路时，转入⑷～⑾框。
在⑷框中判断非标准变比设在支路的哪一侧。
如1章中所述，当i<0时，非标准变比就设在i侧，否则设在j侧。
由图3-1（b）及式（3-5）可知，在非标准变比侧自导纳应累计
yij，但在形成不接地支路的导纳矩阵时，该点自导纳累计了
，因此需要再追加累计（
-1）
=（1-
）Yij。
图中⑹～⑼框就是完成这些运算的。
在非标准变比侧自导纳应累计yij，在形成不接地支路的导纳矩阵时，该点自导纳同样累计了
，因此需要再追加累计（K-1）
=（1-K）Yij，在⑽框和⑾框中完成这些运算。
这样，顺次把支路数据扫描、处理一遍，就形成了描述网络的完整的导纳矩阵。
4.稀疏系数矩阵线性方程式的求解
4.1修正方程式的解法及计算公式
在P-Q分解法潮流计算的迭代过程中，需要反复求解修正方程式
△P/V=B
△θV（4-1）
△Q/V=B
△V（4-2）
如前所述，这两个方程的系数矩阵（B
和B
）在迭代过程中保持不变，只要求对不断变化的常数项（即误差项△P/V、△Q/V）求解出相应的修正量△θV及△V。
在这种情况下，可以先将系数矩阵进行三角分解，然后只要用分解出的三角矩阵（或因子表）对不同常数项进行前代及回代的运算，即可得到要求的修正量。
系数矩阵的三角分解可以利用递推公式求得，也可以利用高斯消去法求得，这两种方法在运算量及内存量上都是等效的。
本程序利用高斯消去法对系数矩阵进行三角分解并形成因子表的计算方法。
在P-Q分解法潮流程序中，为了在迭代过程中轮流求解式（4-1）及式（4-2）需要形成两个因子表,为此，可以把式（4-1）及式（4-2）统一为如下的形式：
B△X=△I（4-3）
在P-θ迭代时，式中B即B
，△X为△θV，△I为△P/V；在Q-V迭代时，式中B为B
，△X为△V，△I为△Q/V。
以下简单归纳以下形成因子表及常数项进行前代及回代运算的有关公式，对于式（4-3）的系数矩阵进行B三角分解以后，可以得到以下形式的因子表：
（4-4）
因子表中上三角部分元素组成了上三角矩阵U：
U=
其中元素
Uij=Bij（i）（i因子表中对角元素组成了对角矩阵D：
D=
其中元素
Dii=
（4-6）
式（4-5）中Bij（i）为因子表中上三角部分i行j列元素，其上标（i）表示此元素由原来系数矩阵元素Bij经过i次运算得来。
这i次运算中包括i-1次消去运算及一次规格化运算：
Bij（k）=Bij（k-1）-Bik（k-1）Bkj（k）=Bij（k-1）-Bik（k-1）Ukj
（k=1,2,…,i-1;j=k+1,k+2,…,n-1）（4-7）
Bij（i）=
（j=i+1,i+2,…,n-1）（4-8）
由于系数矩阵B为对称矩阵，在因子表中不需要保留其下三角部分。
在形成因子表的过程中需要用到下三角部分的元素Bik（k-1）（kBik（k-1）=Bki（k）Bkk（k-1）=Uki
（4-9）
利用式（4-6）～（4-9）就可以由对称系数矩阵B有规律的计算出因子表的所有元素。
对式（4-3）来说，当系数矩阵为对称矩阵时，其具体计算公式为：
△Ij（i）=△Ij（j-1）-Uij△Ii（i-1）（j=2,3,…,n-1;i=1,2,…,j-1）（4-10）
△Ij（j）=Djj△Ij（j-1）（4-11）
显然，顺次取因子表各元素参加一次运算就完成了前代过程的计算。
回代过程计算公式为：
△xi=△Ii（i）-
（j=n-1,n-2,…,1）（4-12）
4.2因子表形成程序框图
由于式（4-3）中的系数矩阵和导纳矩阵B具有相同的结构和稀疏性，因而在一般情况下因子表的上三角矩阵U也是稀疏的。
但由于在消去过程中可能会增加一些新的非零元素，所以U和B的结构不完全相同。
为了节约内存和提高计算速度，在因子表中将不保留零元素，而把上三角矩阵中非零元素紧凑的排列在一起，对每个上三角矩阵U的非零元素都附一个列下标，并对其中每一行都给出非零元素的数目。
以下框图（图4-1）来说明形成因子表的具体过程。
首先对框图中符号作一说明。
flag：
是一个标识变量，当flag=1时，该框图所表示的程序形成第一个因子表（与B
相对应）。
当flag=2时形成第二个因子表（与B
相对应）。
flag在该子程序调用前根据主程序
图4-1形成因子表框图
要求置1或置2，作为参数传个子程序。
n_pv：
为PV节点数组的计数变量，在形成第二个因子表时，用它可以判断系数矩阵中应该去掉哪些列和哪些行。
i：
为因子表正在形成的行号变量。
j：
为列下标变量。
n_u：
是因子表上三角矩阵元素计数变量。
i_above：
为消去行号计数变量（是被消行号计数变量，因子表按行消去过程中，依次取行上面的到行）。
i_pv：
为PV节点的节点号变量。
Btemp：
为临时变量。
count:
为临时计数变量。
形成因子表所需要的原始数据可由以下数组取得，这些数组的定义及内容详见1章和3章：
Load：
负荷功率数组。
PVNode：
PV节点数组。
Yij：
导纳矩阵非对角元素数组。
NYseq：
导纳矩阵各行非对角元素的首地址数组。
Yii：
导纳矩阵对角元素数组。
形成的因子表将放在以下数组中：
NUsum：
存放因子表上三角矩阵各行非对角元素数。
D：
存放因子表的对角元素。
U：
存放因子表上三角矩阵元素，定义为结构体数组：
structU_Type
{
doublevalue;
intj;
}U1[300],U2[300];
其中存放该元素的数值，存放该元素的列下标。
最后，在框图中还有一个很重要的工作数组B[N]，形成因子表的运算主要在这个数组中进行。
现在分别介绍框图中各个部分的作用。
图4-1中A、B、C三部分的作用是“传递”，通过这三部分的工作把系数矩阵中待消行（i行）的元素按其下标稀疏的排列在工作数组B中。
D框的作用是把工作数组B中的待消行元素按照式（4-7）进行消去运算。
最后通过E框的工作把数组B中的元素按式（4-8）进行格式化，并把数组B中非零元素搜集起来，紧密的排列到U数组中去。
这样，从第一行（i=1）做到第N-1行（i=N-1），就形成了完整的因子表。
整个计算过程是一行号为循环变量的。
以下将详细讨论图中各细框的工作情况。
首先介绍当flag=1时即形成第一因子表时的工作情况。
A部分把工作数组B从i+1到N-1单元全部充零，并把导纳矩阵对角元素的虚部送到B数组的第i个单元。
B部分把导纳矩阵非对角元素的虚部按其列下标送到B数组的相应单元中去。
在flag=1时，程序不执行C部分中的运算，直接转入到D部分。
至此，系数矩阵的第i行元素已稀疏的按其列下标排列在B数组中，以下将在D部分中按行对工作数组（即待消行i）中的元素进行消去运算。
一般的说，在形成因子表第i行各元素时，工作数组应该与i-1行以前已形成的各行因子表元素进行消去运算。
因此，在D部分中安排了消去行号i_above=1到i_above=i-1的循环过程。
该框开始时，对n_u赋值1，为顺序取用因子表上三角矩阵已形成各元素作好准备。
由于因子表及系数矩阵B都具有稀疏的特性，消去过程比较复杂，以下用图所示的例子来说明。
图4-2中因子表的第一行及第二行已经形成。
为了清楚起见，在图4-1中已将这两行因子表元素展开排列，实际上在数组U中它们是密集排列的（见图4-3）。
因子表第一行
因子表第二行
待消的第三行
图4-2形成因子表时的消去过程
在图4-2所示的例子中，第一行有三个非对角元素，第二行有三个非对角元素。
由式（4-5）可知，图中各元素的具体意义为
U13=B13
（1）,U15=B15
（1）,…,U28=B28
（2）
图4-3因子表上三角矩阵的存放形式
工作数组中B33为系数矩阵第三行的对角元素，B36及B38为非对角元素。
首先讨论因子表第一行与工作数组的消去过程。
根据式（4-7），可以写出
B3j
（1）=B3j-B31
B1j
（1）
式中：
下标应该由待消行号3开始，即j=3,4…。
因为系数矩阵是对称矩阵，如上所述，在因子表中不必保留下三角部分，因此上式中B31还必须利用式（4-9）求出：
B31=
B13
（1）=
U13
得到后，就可以进行消去运算。
先从，即对角元素起
B33
（1）=B33-B31
B13
（1）=B33-B31
U13（4-13）
当j=4时：
B34
（1）=B34-B31
B14
（1）=B34-B31
U14（4-14）
但在图4-2所示的例子中是零元素，因此上式变为
B34
（1）=B34
这样，式（4-14）所表示的运算对工作数组（即待消行）中的元素没有任何影响。
因此，在第一行与第三行进行消去运算时，只需要从因子表第一行中列下标为3的元素开始，顺次取以后的元素（见图4-3），按照其列下标与工作数组中相应的元素进行消去运算。
在本例中，除了按式（4-13）进行计算以外，还应进行以下两次运算：
B35
（1）=B35-B31
U15
（4-15）
B36
（1）=B36-B31
U16
式（4-15）中，B35为零（见图4-2），即系数矩阵中B35为零元素，但与第一行进行消去运算后变成了非零元素，因而在工作数组（即待消行）中出现了一个注入元素。
现在讨论第二行对第三行的消去运算。
根据式（4-7），消去过程的计算公式为
B3j
（2）=B3k
（1）-B32
（1）
B2j
（2）（4-16）
式中：
B32
（1）=
B23
（2）=
U23
由图4-2可知，在该例中U23=0，因此B32
（1）也等于零，这样式（4-16）变为
B3j
（2）=B3j
（1）
也就是说，由于U23是零元素，因此第二行不必对第三行进行消去运算。
在一般情况下，当工作数组中待消行号为i，而第i
（i
行对i行的消去运算过程。
在图4-1的D部分中，为了判断i
行是否需要对工作数组中的待消行进行消去运算，必须对因子表i
行中元素逐个检查其列下标是否等于待消行的行号i[D部分中⑴]。
当查到列下标等于i的元素时，首先按照式（4-9）求出Bi,i
（i
-1）[D部分中⑵]，然后用该行剩余元素（包括列下标为i的
升级会员 