牛顿拉夫逊法潮流计算matlab程序.docx
《牛顿拉夫逊法潮流计算matlab程序.docx》由会员分享，可在线阅读，更多相关《牛顿拉夫逊法潮流计算matlab程序.docx（6页珍藏版）》请在冰豆网上搜索。
牛顿拉夫逊法潮流计算matlab程序
牛顿拉夫逊法潮流计算matlab程序
%电力系统的潮流计算,以下程序参考文献《电力系统毕业设计》中国水利电力出版社
%（该文献用极坐标下的牛顿——拉夫逊方法实现，在此为了与课本一致做了修改）
%为了计算方便将原来的下标做以下修改:
S2S3S4S5U2U3U4U5改为S1S2S3S4U1U2
%U3U4,即原题的平衡点1就变为现在的平衡点5
%1.形成节点导纳矩阵,
yb55=6.250-18.750j;yb51=-5.000+15.000j;yb52=-1.250+3.750j;yb53=0.000-0.000j;yb54=0.000-0.000j;
yb15=-5.000+15.000j;yb11=10.834-32.500j;yb12=-1.667+5.000j;yb13=-1.667+5.000j;yb14=-2.500+7.500j;
yb25=-1.250+3.750j;yb21=-1.667+5.000j;yb22=12.917-38.750j;yb23=-10.000+30.000j;yb24=0.000-0.000j;
yb35=0.000-0.000j;yb31=-1.667+5.000j;yb32=-10.000+30.000j;yb33=12.917-38.750j;yb34=-1.250+3.750j;
yb45=0.000-0.000j;yb41=-2.500+7.500j;yb42=0.000-0.000j;yb43=-1.250+3.750j;yb44=3.750-11.250j;
YB=[yb11yb12yb13yb14yb15;yb21yb22yb23yb24yb25;yb31yb32yb33yb34yb35;yb41yb42yb43yb44yb45;yb51yb52yb53yb54yb55];
%计算各节点功率的不平衡量设U=E+jF；Y=G+Bj；
E
（1）=1.00;E
（2）=1.00;E（3）=1.00;E（4）=1.00;
F
（1）=0;F
（2）=0;F（3）=0;F（4）=0;
G=real（YB）;B=imag（YB）;
%设S=P+Bj；
S
（1）=0.20+0.20i;S
（2）=-0.45-0.15i;S（3）=-0.40-0.05i;S（4）=-0.60-0.10i;
P=real（S）;Q=imag（S）;
k=0;precision=1;
N1=4;
whileprecision>0.00000001
E（5）=1.06;F（5）=0;
form=1:
N1
forn=1:
N1+1
%计算Pi,Qi,设Pi=Pt;Qi=Qt
Pt（n）=（E（m）*（G（m,n）*E（n）-B（m,n）*F（n））+F（m）*（G（m,n）*F（n）+B（m,n）*E（n）））;
Qt（n）=（F（m）*（G（m,n）*E（n）-B（m,n）*F（n））-E（m）*（G（m,n）*F（n）+B（m,n）*E（n）））;
end
%设P,Q的改变量为dP,dQ
dP（m）=P（m）-sum（Pt）;
dQ（m）=Q（m）-sum（Qt）;
end
form=1:
N1
forn=1:
N1+1
%计算HijNijJijLij
H（m,n）=-B（m,n）*E（m）+G（m,n）*F（m）;
N（m,n）=G（m,n）*E（m）+B（m,n）*F（m）;
J（m,n）=-B（m,n）*F（m）-G（m,n）*E（m）;
L（m,n）=G（m,n）*F（m）-B（m,n）*E（m）;
end
end
form=1:
N1
forn=1:
N1+1
Bi（n）=G（m,n）*F（n）+B（m,n）*E（n）;
Ai（n）=G（m,n）*E（n）-B（m,n）*F（n）;
end
%计算Hii,Nii,Jii,Lii,由公式4-44b左侧公式实现,sum（Ai），sum（Bi）用于实现公式中的sigerma从j到n的求和;
H（m,m）=sum（Bi）-（B（m,m）*E（m）+G（m,m）*F（m））+2*G（m,m）*F（m）;
N（m,m）=sum（Ai）-（G（m,m）*E（m）-B（m,m）*F（m））+2*G（m,m）*E（m）;
J（m,m）=-2*B（m,m）*F（m）+sum（Ai）-（G（m,m）*E（m）-B（m,m）*F（m））;
L（m,m）=-2*B（m,m）*E（m）-（sum（Bi）-（B（m,m）*E（m）+G（m,m）*F（m）））;
end
%设雅可比矩阵为JJ，以下语句用来实现雅可比矩阵中对角线上元素HNJL的排列
form=1:
N1
JJ（2*m-1,2*m-1）=H（m,m）;
JJ（2*m-1,2*m）=N（m,m）;
JJ（2*m,2*m-1）=J（m,m）;
JJ（2*m,2*m）=L（m,m）;
end
%以下语句用于实现雅可比矩阵非对角线上元素的排列
form=1:
N1
forn=1:
N1
ifm==n
else
H（m,n）=-B（m,n）*E（m）+G（m,n）*F（m）;
N（m,n）=G（m,n）*E（m）+B（m,n）*F（m）;
J（m,n）=-B（m,n）*F（m）-G（m,n）*E（m）;
L（m,n）=G（m,n）*F（m）-B（m,n）*E（m）;
JJ（2*m-1,2*n-1）=H（m,n）;
JJ（2*m-1,2*n）=N（m,n）;
JJ（2*m,2*n-1）=J（m,n）;
JJ（2*m,2*n）=L（m,n）;
end
end
end
%设由P,Q的改变量组成的8×1矩阵为PQ，由E,F的改变量组成的8×1矩阵为dU
form=1:
N1
PQ（2*m-1）=dP（m）;PQ（2*m）=dQ（m）;
end
dU=inv（JJ）*PQ';
precision=max（abs（dU））;
forn=1:
N1
F（n）=F（n）+dU（2*n-1）;
E（n）=E（n）+dU（2*n）;
end
forn=1:
N1+1
U（n）=E（n）+（F（n））*j;
end
k=k+1;
k-1,dU=dU',PQ,U
end
%计算S（5）,也就是题目中的S1，即平衡节点功率
form=1:
N1+1
I（m）=YB（5,m）*U（m）;
end
S（5）=U（5）*sum（conj（I））
%设网络总损耗为Ss,计算输电效率efficiency
form=1:
N1+1
S0（m）=S（m）;P（m）=real（S（m））;
end
Ss=sum（S0）
efficiency=（abs（P（3）+P（4）+P
（2）））/（P（5）+（P
（1）））*100
%计算线功率S（m,n）,与课本中各元素的相对位置有所不同
form=1:
N1+1
forn=1:
N1+1
S（m,n）=U（m）*（conj（U（m））-conj（U（n）））*conj（-YB（m,n））;
end
end
S
升级会员 