电力网节点导纳矩阵计算例题与程序

电力网节点导纳矩阵计算例题与程序
佘名寰 编写
用计算机解算电力网潮流电压和短路电流问题首先需确定电力网的节点导纳矩阵或节点阻抗矩阵。

本文通过例题介绍用网络拓扑法计算节点导纳矩阵的方法和程序，程序考虑了线路并联电容和变压器支路标么变比不为1时的影响。

程序用MATLAB 语言编写，线路参数均采用标么值。

本文稿用office word 2007 版编写，可供电气专业人员计算相关问题时参考。

1.用网络拓扑计算节点导纳矩阵 1.1网络拓扑矩阵：
【例1.1】 例图1-1是有5 个节点和5条支路的网络，节点5作为基准参考点，1 ,2, 3, 4为独立节点，支路编号和方向图中已标识。

例图1-1
对于具有n 个节点b 条支路的有向图，它的关联矩阵为一个N ×B 的矩阵A a ：
A a =[a ij ]
若支路j 与节点i 相关，且箭头背离节点i ，则a ij =1，若箭头指向节点则a ij =-1，若支路j 与节点i 无关，则a ij =0，图1-1所示的有向图的关联矩阵为
① ② ③ ④ ⑤ 支路编号
A ij =
行编号从上到下为1 2 3 4 5节点编号（5为参考节点） 去掉第5行即为独立节点的关联矩阵。

以下介绍生成网络关联矩阵的M 函数文件 ffm.m ：
% M FUNCTION ffm.m
% Np is number of node point,Nb is number of braches ②Z23
③Z13
①Z21
YC2
YC3
YC1
④Z42
⑤Z53
1:1.05
1.05:1
4
2
3
1
5
Z21=0.04+J0.25 Z23=0.08+J0.30 Z13=0.1+J0.35 Z42=J0.015 Z53=J0.03 YC1=J0.25 YC2=J0.50 YC3=J0.25
% nstart--the start point of branches ,nend -- the end point,
% A -- network incidence matrix
function[A]=ffm(nstart,nend)
global Np Nb
n=length(nstart);
A=zeros(Np,Nb);
for i=1:n
A(nstart(i),i)=1;
A(nend(i),i)=-1;
end
以例图1-1网络为例调用ffm.m文件求其关联矩阵
运算以上程序可得关联矩阵 mm ij如下：
mm =
-1 0 1 0 0
1 1 0 -1 0
0 -1 -1 0 -1
0 0 0 1 0
0 0 0 0 1
Mm ij明显与A ij是相同的。

1.2生成节点导纳矩阵程序：
⑴由网络原始矩阵计算节点导纳矩阵公式
Y=AY s0A t (1-1)
Y ----节点导纳矩阵
A------网络关联矩阵
A t-----A的转置矩阵
Y S0----网络原始导纳矩阵
若网络各支路阻抗为 Z b=[z b1,z b2,……,z bn]
则Z S0=
Y s0=Z s0-1 (1-2)
Y=A Z s0-1A t (1-3)
⑵节点导纳计算程序
以例1-1网络为例，在不计对地电容和变压器变比假定为1条件下，节点导纳矩阵计算程序如下：
clear
global Np Nb
% Np is number of node point,Nb is number of braches,
Np=5;Nb=5;
% nstart--the start point of branches ,nend -- the end point,
% mm -- network incidence matrix
nstart=[2,2,1,4,5];
nend=[1,3,3,2,3];
mm=ffm(nstart,nend);
% zb1,the series impedances of transmission line
% yb1,the series admittances of transmission line
zb1r=[0.04,0.08,0.1,0.0,0.0];
zb1i=[0.25,0.30,0.35,0.015,0.03];
zb1=zb1r+zb1i*j;
yb1=zb1.^(-1);
yb=diag(conj(yb1'));
y=mm*yb*(mm)';
运算程序后可得节点导纳矩阵
y=
1.3787 - 6.5417i -0.6240 + 3.9002i -0.7547 +
2.6415i 0 0
-0.6240 + 3.9002i 1.4539 -73.6789i -0.8299 + 3.1120i 0 +66.6667i 0
-0.7547 +2.6415i -0.8299+ 3.1120i 1.5846 -39.0869i 0 0
+33.3333i
0 0 +66.6667i 0 0 -66.6667i 0
0 0 0 +33.3333i 0 0
-33.3333i
⑶线路对地并联电容的处理
例图1-1网络中，线路①，②两端对地电容导纳均为j0.25,节点1,2,3自导纳相应增加
j0.25,j0.5,j0.25。

程序处理如下：
% yb0 ,the shunt admittances of transmission line
yb0i=[0.25,0.25,0.0,0.0,0.0];
yb0=0+yb0i*j;
y0=diag(conj(yb0'));
yg=mm*y0*(mm)';
yn=diag(diag(yg));
yy=yn-yg;
y=y+yn;
程序中yb0为各线路对地电容导纳，y0为原始导纳矩阵，yn对角矩阵对应各节点自导纳增
加值，yy相应各支路两端对地并联导纳，用于计算支路电流。

y为计及对地电容后的节点导
纳矩阵。

程序运算结果
yn =
0 + 0.2500i 0 0 0 0 0 0 + 0.5000i 0 0 0 0 0 0 + 0.2500i 0 0
0 0 0 0 0
0 0 0 0 0
yy =
0 0 + 0.2500i 0 0 0
0 + 0.2500i 0 0 + 0.2500i 0 0 0 0 + 0.2500i 0 0 0
0 0 0 0 0
0 0 0 0 0
y = 1.3787 - 6.2917i -0.6240 + 3.9002i -0.7547 + 2.6415i 0 0 -0.6240 + 3.9002i 1.4539 -73.1789i -0.8299 + 3.1120i 0 +66.6667i 0 -0.7547 + 2.6415i -0.8299 + 3.1120i 1.5846 -38.8369i 0 0 +33.3333i 0 0 +66.6667i 0 0 -66.6667i 0 0 0 0 +33.3333i 0 0 -33.3333i ⑷ 非标准变压器支路的处理
当支路有变压器并且标么变比不为1时，可将变压器化作П形等值电路，如图1-2所示。

图中理想变压器变比为k:1,接入节点i,高压侧经漏抗Z T 与节点j 相连。

图1-2变压器等值电路
等值电路图中Y T =1/Z T 。

i 侧自导变为y ii =ky t +k(k-1)y t =k 2y t ，改变量为Δy ii =(k 2
-1)y t ; I j 间互导变为y ij =-ky t,Δy ij =-ky t +y t =-(k-1)y t 。

j 侧自导变为y jj =ky t +(1-k)y t =y t,,与K=1时相同;变压器变比不为1时修正导纳矩阵的M 函数文件 fdt1.m
% M FUNCTION fdt1.m
% the transformer turns ratio,with off-nominal tap-setting,fig 1-2 % bt[] two-dimensional array ,line one is point i,line two is point j % p is the transformer turns ratio,k is number of transformers
function[A,B]=fdt1(a,b)
global bt p m k
for i=1:m
for j=1:m
A(i,j)=a(i,j);
B(i,j)=b(i,j);
end
end
for i1=1:k
k1=bt(1,i1);
k2=bt(2,i1);
t=-a(k1,k2);
p1=p(i1);
for i=1:m
for j=1:m
if i==k1
if j==k1
A(i,j)=A(i,j)+(p1*p1-1)*t;
elseif j==k2
A(i,j)=A(i,j)-(p1-1)*t;
B(i,j)=p1*(p1-1)*t;
end
elseif i==k2
if j==k1
A(i,j)=A(i,j)-(p1-1)*t;
B(i,j)=(1-p1)*t;
end
end
end
end
end
程序中bt是2行K列数组，第一行为变压器接入点i的节点编号，第二行为经漏抗后所接j节点编号,k为变压器支路数。

P为一维k列数组，输入各个变压器的变比；
m 为矩阵 a,b,A,B的阶数
a，A是进行变比修正前后的节点导纳矩阵，
b，B是变比修正前后支路对地并联导纳矩阵。

1.3 例图1-1节点导纳和支路对地并联导纳计算程序
例图1-1共有5个节点，5条支路，所以Np=5，Nb=5
支路④，为变压器支路，其 i 节点编号为2，j节点编号为4，变比为1：1/1.05；
支路⑤也为变压器支路，其 i 节点编号为3，j节点编号为5，变比为1：1/1.05；
所以数组bt=[2,3;4,5]，p=[1/1.05,1/1.05]，k=2;
5条支路按箭头方向首末端节点分别为2,2,1,4,5和1,3,3,2,3，所以nstart=[2,2,1,4,5];
nend=[1,3,3,2,3];输入以上数据可得以下程序：
%******* NU11.m , example 1-1, fig 1-1 **************************
% The following Program for node impedances calculation is based on MATLAB2007
clear
% bus 1,2,3 is PQ bus,bus 4 is PV bus,bus 5 is slack bus
global Np Nb bt p m k
% Np is number of node point,Nb is number of braches,
% p the transformer turns ratio,with off-nominal tap-setting,fig 1-1
% bt is two-dimensional array ,line one is point i,line two is point j,
% m is number of node point, k is number of transformers
Np=5;Nb=5;
bt=[2,3;4,5];p=[1/1.05,1/1.05];m=5;k=2;
% nstart--the start point of branches ,nend -- the end point,
% mm -- network incidence matrix
nstart=[2,2,1,4,5];
nend=[1,3,3,2,3];
mm=ffm(nstart,nend);
% zb1,the series impedances of transmission line
% yb1,the series admittances of transmission line
zb1r=[0.04,0.08,0.1,0.0,0.0];
zb1i=[0.25,0.30,0.35,0.015,0.03];
zb1=zb1r+zb1i*j;
yb1=zb1.^(-1);
yb=diag(conj(yb1'));
y=mm*yb*(mm)';
% yb0 ,the shunt admittances of transmission line
yb0i=[0.25,0.25,0.0,0.0,0.0];
yb0=0+yb0i*j;
y0=diag(conj(yb0'));
yg=mm*y0*(mm)';
yn=diag(diag(yg));
yy=yn-yg;
y=y+yn;
[Y,YY]=fdt1(y,yy);
% YY is the shund admittances of the line and transformers at each end,
% Y is node-admittance matrix for network
运行结果为：
Y =
1.3787 - 6.2917i -0.6240 + 3.9002i -0.7547 +
2.6415i 0 0
-0.6240 + 3.9002i 1.4539 -66.9808i -0.8299 + 3.1120i 0 +63.4921i
-0.7547 +2.6415i -0.8299 + 3.1120i 1.5846 -35.7379i 0 0
+31.7460i
0 0 +63.4921i 0 0 -66.6667i 0
0 0 0 +31.7460i 0 0
-33.3333i
>> YY
YY =
0 0 + 0.2500i 0 0 0 0 + 0.2500i 0 0 + 0.2500i 0 + 3.0234i 0
0 0 + 0.2500i 0 0 0 +
1.5117i
0 0 - 3.1746i 0 0 0 0 0 0 - 1.5873i 0 0
参考文献：
1．MATLAB基础与应用简明教程张平等编著北京航空航天大学出版社2001年
2．燃气工程电算方法，冯良张同全惠君编著，同济大学出版社，1998年
3．电力系统计算，西安交通大学清华大学等合编，水利电力出版社，1978年。