完整版电力系统潮流程序设计.docx
《完整版电力系统潮流程序设计.docx》由会员分享，可在线阅读，更多相关《完整版电力系统潮流程序设计.docx（17页珍藏版）》请在冰豆网上搜索。
完整版电力系统潮流程序设计
潮流程序的设计与开发
1数据结构的设计
要求：
将设备铭牌参数和有名值作为原始输入数据，潮流结果以有名值的形式输出。
支路数据与节点数据分别以一个矩阵的形式输入，矩阵的每一行表示每一个节点或每一条支路，矩阵的每一列表示不同的参数数据。
2变量命名设计
变量名称
程序变量表示
变量名称
程序变量表示
节点导纳矩阵
无功功率
电导
视在功率
电钠
雅克比矩阵
电压幅值
平衡节点编号
电压相角
有功不平衡量
有功功率
无功不平衡量
电压相角修正量
电压幅值修正量
3程序流程图
4程序代码
4.1主程序
n=input（'请输入节点数:
n='）;
l=input（'请输入支路数:
l='）;
%支路数不要包括三绕组变压器
sw=input（'请输入平衡节点号:
sw='）;
ac=input（'请输入误差精度:
ac='）;
SB=input（'请输入基准功率:
SB='）;
B1=input（'请输入支路参数:
B1='）
%支路参数节点参数和对地支路参数均以矩阵形式储存
%第一列储存支路编号
%第二列与第三列分别储存支路的两个端点，分别为p,q
%第四列储存支路阻抗
%第五列储存支路对地导纳,注意对地导纳不要除以2
%第六列储存该支路是否含有变压器，有为1，无为0
%第七列储存变压器变比k,p指向q的变压器变比为k:
1，且k大于等于1
%第八列储存变压器短路损耗
%第九列储存变压器短路电压百分数
%第十列储存变压器空载损耗
%第十一列储存变压器空载电流百分数
%第十二列储存变压器低压侧电压
%第十三列储存变压器额定功率
%第十四列储存归算所取基准电压
%注意，将三绕组变压器转换为双绕组变压器输入
A1=input（'请输入节点参数:
A1='）;
%第一列为节点编号
%第二列为注入发电功率
%第三列为负荷功率
%第四列为节点电压幅值,为方便起见，以标幺值形式表示
%第五列为节点电压相角
%第六列储存节点对地导纳
%第七列为节点的类型，1为PQ节点，2为PV节点，3为平衡节点
%%%首先求解双绕组变压器参数
fori=1:
l
ifB1（i,6）==1
Zt（i）=B1（i,8）*B1（i,12）^2/（1000*B1（i,13）^2）;
Xt（i）=B1（i,9）*B1（i,12）^2/（100*B1（i,13））;
Gt（i）=B1（i,10）/（1000*B1（i,12）^2）;
Bt（i）=B1（i,11）*B1（i,13）/（100*B1（i,12）^2）;
end
end
B2（:
1:
7）=B1（:
1:
7）;
fori=1:
l
B2（i,8）=Zt（i）;
B2（i,9）=Xt（i）;
end
fori=1:
l
ifB1（i,6）~=1
ZB（i）=B1（i,14）^2/SB;
YB（i）=1/ZB（i）;
B2（i,4）=B1（i,4）/ZB（i）;
B2（i,5）=B1（i,5）/YB（i）;
end
end
fori=1:
l
ifB1（i,6）==1
ZB（i）=B1（i,14）^2/SB;
YB（i）=1/ZB（i）;
B2（i,8）=B2（i,8）/ZB（i）;
B2（i,9）=B2（i,9）/ZB（i）;
Gt（i）=Gt（i）/YB（i）;
Bt（i）=Bt（i）/YB（i）;
A1（B1（i,3）,6）=Gt（i）-Bt（i）*（1i）;
end
end
%%%下面求解节点导纳矩阵
Y=zeros（n）;
fori=1:
n
ifA1（i,6）~=0
Y（i,i）=A1（i,6）;
end
end
B3=B2;
fori=1:
l
p=B3（i,2）;
q=B3（i,3）;
ifB3（i,6）==1%%%含有变压器支路支路
Y（p,p）=Y（p,p）+1./（（B3（i,8）+B3（i,9）*（1i））*B3（i,7）^2）;
Y（p,q）=Y（p,q）-1./（（B3（i,8）+B3（i,9）*（1i））*B3（i,7））;
Y（q,p）=Y（p,q）;
Y（q,q）=Y（q,q）+1./（（B3（i,8）+B3（i,9）*（1i）））;
else%%%无变压器支路
Y（p,p）=Y（p,p）+1./B3（i,4）+B3（i,5）/2;
Y（p,q）=Y（p,q）-1./B3（i,4）;
Y（q,p）=Y（p,q）;
Y（q,q）=Y（q,q）+1./B3（i,4）+B3（i,5）/2;
end
end
Y
A2=A1;
A2（:
2）=A1（:
2）/SB;
A2（:
3）=A1（:
3）/SB;%功率参数标幺化
G=real（Y）;
B=imag（Y）;
S=A2（:
2）-A2（:
3）;
P=real（S）;
Q=imag（S）;
V=A2（:
4）;
delta=A2（:
5）;
DeltaS=bphl（n,sw,A2,V,G,delta,B,P,Q）;
J=jcb（G,B,V,delta,n,B3,A2,sw）;
[ddelta,dV]=xzl（n,J,DeltaS,A2）;
e=1;
while（max（ddelta）>ac|max（dV）>ac）
a=0;
b=0;
fori=1:
n
switchA2（i,7）
case1
delta（i）=delta（i）+ddelta（i-b）;
V（i）=V（i）+dV（i-a-b）;
case2
delta（i）=delta（i）+ddelta（i-b）;
V（i）=V（i）;
a=a+1;
case3
delta（i）=delta（i）;
V（i）=V（i）;
b=b+1;
end
end
DeltaS=bphl（n,sw,A2,V,G,delta,B,P,Q）;
J=jcb（G,B,V,delta,n,B3,A2,sw）;
[ddelta,dV]=xzl（n,J,DeltaS,A2）;
e=e+1;
end
V
delta
e
%%%%下面求平衡节点功率
v=V.*cos（delta）+V.*sin（delta）*（1i）;
forj=1:
n
yu（j）=conj（Y（sw,j））*conj（v（j））;
end
S（sw）=sum（yu）*v（sw）;
input（'平衡节点的功率为'）;
S（sw）
B2（sw,1）=S（sw）;
%%%%下面求解线路功率
fori=1:
n
forj=1:
n
Sl（i,j）=v（i）*（conj（v（i））*conj（A2（i,6））+conj（v（i）-v（j））*conj（-Y（i,j）））;
end
end
input（'线路功率为'）;
Sl
%%%%线路上损耗的功率
fori=1:
n
forj=1:
n
DertaS1（i,j）=（Sl（i,j）+Sl（j,i））/2;
end
end
input（'线路上损耗的功率为'）;
DertaSz=sum（sum（DertaS1））
4.2计算功率不平衡量程序
functionDeltaS=bphl（n,sw,A2,V,G,delta,B,P,Q）
%计算功率不平衡量
fori=1:
n
ifA2（i,7）~=sw
EP（i）=0;
EQ（i）=0;
forj=1:
n
EP（i）=EP（i）+V（i）*V（j）*（G（i,j）*cos（delta（i）-delta（j））+B（i,j）*sin（delta（i）-delta（j）））;
EQ（i）=EQ（i）+V（i）*V（j）*（G（i,j）*sin（delta（i）-delta（j））-B（i,j）*cos（delta（i）-delta（j）））;
end
P1=EP（i）;
Q1=EQ（i）;
ifA2（i,7）==1%PQ节点
p=2*i-1;
DeltaS（p）=P（i）-P1;
p=p+1;
DeltaS（p）=Q（i）-Q1;
else
p=2*i-1;
DeltaS（p）=P（i）-P1;
p=p+1;
DeltaS（p）=0;
end
end
end
DeltaS（2*sw-1）=[];
DeltaS（2*sw-1）=[];
fori=1:
n
ifA2（i,7）==2
DeltaS（2*i）=[];
end
end
DeltaS=DeltaS';
End
4.3计算雅克比矩阵程序
function[J]=jcb（G,B,V,delta,n,B3,A2,sw）
%计算雅克比矩阵
fori=1:
n
ifA2（i,7）==1%PQ节点
forj=1:
n
ifj~=i&j~=sw
H=V（i）*V（j）*（G（i,j）*sin（delta（i）-delta（j））-B（i,j）*cos（delta（i）-delta（j）））;
J1=-V（i）*V（j）*（G（i,j）*cos（delta（i）-delta（j））+B（i,j）*sin（delta（i）-delta（j）））;
N=V（i）*V（j）*（G（i,j）*cos（delta（i）-delta（j））+B（i,j）*sin（delta（i）-delta（j）））;
L=V（i）*V（j）*（G（i,j）*sin（delta（i）-delta（j））-B（i,j）*cos（delta（i）-delta（j）））;
p=2*i-1;
q=2*j-1;
J（p,q）=H;
m=p+1;
J（m,q）=J1;
q=q+1;
J（p,q）=N;
J（m,q）=L;
elseifj==i&j~=sw
H1=0;
forh=1:
n
H1=H1+（-V（i）*V（h）*（G（i,h）*sin（delta（i）-delta（h））-B（i,h）*cos（delta（i）-delta（h））））;
end
H=H1+V（i）*V（i）*（-B（i,i））;
J2=0;
forh=1:
n
J2=J2+（V（i）*V（h）*（G（i,h）*cos（delta（i）-delta（h））+B（i,h）*sin（delta（i）-delta（h））））;
end
J1=J2-V（i）*V（i）*G（i,i）;
N1=0;
forh=1:
n
N1=N1+（V（i）*V（h）*（G（i,h）*cos（delta（i）-delta（h））+B（i,h）*sin（delta（i）-delta（h））））;
end
N=N1-V（i）*V（i）*G（i,i）+2*V（i）*V（i）*G（i,i）;
L1=0;
forh=1:
n
L1=L1+（V（i）*V（h）*（G（i,h）*sin（delta（i）-delta（h））-B（i,h）*cos（delta（i）-delta（h））））;
end
L=L1-V（i）*V（i）*（-B（i,i））-2*V（i）*V（i）*B（i,i）;
p=2*i-1;
q=2
升级会员 