电力系统潮流计算详解.docx
- 文档编号:25615622
- 上传时间:2023-06-10
- 格式:DOCX
- 页数:23
- 大小:18.22KB
电力系统潮流计算详解.docx
《电力系统潮流计算详解.docx》由会员分享,可在线阅读,更多相关《电力系统潮流计算详解.docx(23页珍藏版)》请在冰豆网上搜索。
电力系统潮流计算详解
首先声明一下,这些是从网站上转载的,不是本人上编写的
基于MATLAB的电力系统潮流计算
%简单潮流计算的小程序,相关的原始数据数据数据输入格式如下:
%B1是支路参数矩阵,第一列和第二列是节点编号。
节点编号由小到大编写
%对于含有变压器的支路,第一列为低压侧节点编号,第二列为高压侧节点
%编号,将变压器的串联阻抗置于低压侧处理。
%第三列为支路的串列阻抗参数。
%第四列为支路的对地导纳参数。
%第五烈为含变压器支路的变压器的变比
%第六列为变压器是否是否含有变压器的参数,其中“1”为含有变压器,
%“0”为不含有变压器。
%B2为节点参数矩阵,其中第一列为节点注入发电功率参数;第二列为节点
%负荷功率参数;第三列为节点电压参数;第六列为节点类型参数,其中
%“1”为平衡节点,“2”为PQ节点,“3”为PV节点参数。
%X为节点号和对地参数矩阵。
其中第一列为节点编号,第二列为节点对地
%参数。
n=input('请输入节点数:
n=');
n1=input('请输入支路数:
n1=');
isb=input('请输入平衡节点号:
isb=');
pr=input('请输入误差精度:
pr=');
B1=input('请输入支路参数:
B1=');
B2=input('请输入节点参数:
B2=');
X=input('节点号和对地参数:
X=');
Y=zeros(n);
Times=1;%置迭代次数为初始值
%创建节点导纳矩阵
fori=1:
n1
ifB1(i,6)==0%不含变压器的支路
p=B1(i,1);
q=B1(i,2);
Y(p,q)=Y(p,q)-1/B1(i,3);
Y(q,p)=Y(p,q);
Y(p,p)=Y(p,p)+1/B1(i,3)+0.5*B1(i,4);
Y(q,q)=Y(q,q)+1/B1(i,3)+0.5*B1(i,4);
else%含有变压器的支路
p=B1(i,1);
q=B1(i,2);
Y(p,q)=Y(p,q)-1/(B1(i,3)*B1(i,5));
Y(q,p)=Y(p,q);
Y(p,p)=Y(p,p)+1/B1(i,3);
Y(q,q)=Y(q,q)+1/(B1(i,5)^2*B1(i,3));
end
end
Y
OrgS=zeros(2*n-2,1);
DetaS=zeros(2*n-2,1);%将OrgS、DetaS初始化
%创建OrgS,用于存储初始功率参数
h=0;
j=0;
fori=1:
n%对PQ节点的处理
ifi~=isb&B2(i,6)==2
h=h+1;
forj=1:
n
OrgS(2*h-1,1)=OrgS(2*h-1,1)+real(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))+imag(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));
OrgS(2*h,1)=OrgS(2*h,1)+imag(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))-real(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));
end
end
end
fori=1:
n%对PV节点的处理,注意这时不可再将h初始化为0
ifi~=isb&B2(i,6)==3
h=h+1;
forj=1:
n
OrgS(2*h-1,1)=OrgS(2*h-1,1)+real(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))+imag(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));
OrgS(2*h,1)=OrgS(2*h,1)+imag(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))-real(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));
end
end
end
OrgS
%创建PVU用于存储PV节点的初始电压
PVU=zeros(n-h-1,1);
t=0;
fori=1:
n
ifB2(i,6)==3
t=t+1;
PVU(t,1)=B2(i,3);
end
end
PVU
%创建DetaS,用于存储有功功率、无功功率和电压幅值的不平衡量
h=0;
fori=1:
n%对PQ节点的处理
ifi~=isb&B2(i,6)==2
h=h+1;
DetaS(2*h-1,1)=real(B2(i,2))-OrgS(2*h-1,1);
DetaS(2*h,1)=imag(B2(i,2))-OrgS(2*h,1);
end
end
t=0;
fori=1:
n%对PV节点的处理,注意这时不可再将h初始化为0
ifi~=isb&B2(i,6)==3
h=h+1;
t=t+1;
DetaS(2*h-1,1)=real(B2(i,2))-OrgS(2*h-1,1);
DetaS(2*h,1)=real(PVU(t,1))^2+imag(PVU(t,1))^2-real(B2(i,3))^2-imag(B2(i,3))^2;
end
end
DetaS
%创建I,用于存储节点电流参数
i=zeros(n-1,1);
h=0;
fori=1:
n
ifi~=isb
h=h+1;
I(h,1)=(OrgS(2*h-1,1)-OrgS(2*h,1)*sqrt(-1))/conj(B2(i,3));
end
end
I
%创建Jacbi(雅可比矩阵)
Jacbi=zeros(2*n-2);
h=0;
k=0;
fori=1:
n%对PQ节点的处理
ifB2(i,6)==2
h=h+1;
forj=1:
n
ifj~=isb
k=k+1;
ifi==j%对角元素的处理
Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3))+imag(I(h,1));
Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3))+real(I(h,1));
Jacbi(2*h,2*k-1)=-Jacbi(2*h-1,2*k)+2*real(I(h,1));
Jacbi(2*h,2*k)=Jacbi(2*h-1,2*k-1)-2*imag(I(h,1));
else%非对角元素的处理
Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3));
Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3));
Jacbi(2*h,2*k-1)=-Jacbi(2*h-1,2*k);
Jacbi(2*h,2*k)=Jacbi(2*h-1,2*k-1);
end
ifk==(n-1)%将用于内循环的指针置于初始值,以确保雅可比矩阵换行
k=0;
end
end
end
end
end
k=0;
fori=1:
n%对PV节点的处理
ifB2(i,6)==3
h=h+1;
forj=1:
n
ifj~=isb
k=k+1;
ifi==j%对角元素的处理
Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3))+imag(I(h,1));
Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3))+real(I(h,1));
Jacbi(2*h,2*k-1)=2*imag(B2(i,3));
Jacbi(2*h,2*k)=2*real(B2(i,3));
else%非对角元素的处理
Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3));
Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3));
Jacbi(2*h,2*k-1)=0;
Jacbi(2*h,2*k)=0;
end
ifk==(n-1)%将用于内循环的指针置于初始值,以确保雅可比矩阵换行
k=0;
end
end
end
end
end
Jacbi
%求解修正方程,获取节点电压的不平衡量
DetaU=zeros(2*n-2,1);
DetaU=inv(Jacbi)*DetaS;
DetaU
%修正节点电压
j=0;
fori=1:
n%对PQ节点处理
ifB2(i,6)==2
j=j+1;
B2(i,3)=B2(i,3)+DetaU(2*j,1)+DetaU(2*j-1,1)*sqrt(-1);
end
end
fori=1:
n%对PV节点的处理
ifB2(i,6)==3
j=j+1;
B2(i,3)=B2(i,3)+DetaU(2*j,1)+DetaU(2*j-1,1)*sqrt(-1);
end
end
B2
%开始循环**********************************************************************
whileabs(max(DetaU))>pr
OrgS=zeros(2*n-2,1);%!
!
!
初始功率参数在迭代过程中是不累加的,所以在这里必须将其初始化为零矩阵
h=0;
j=0;
fori=1:
n
ifi~=isb&B2(i,6)==2
h=h+1;
forj=1:
n
OrgS(2*h-1,1)=OrgS(2*h-1,1)+real(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))+imag(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));
OrgS(2*h,1)=OrgS(2*h,1)+imag(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))-real(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));
end
end
end
fori=1:
n
ifi~=isb&B2(i,6)==3
h=h+1;
forj=1:
n
OrgS(2*h-1,1)=OrgS(2*h-1,1)+real(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))+imag(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));
OrgS(2*h,1)=OrgS(2*h,1)+imag(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))-real(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));
end
end
end
OrgS
%创建DetaS
h=0;
fori=1:
n
ifi~=isb&B2(i,6)==2
h=h+1;
DetaS(2*h-1,1)=real(B2(i,2))-OrgS(2*h-1,1);
DetaS(2*h,1)=imag(B2(i,2))-OrgS(2*h,1);
end
end
t=0;
fori=1:
n
ifi~=isb&B2(i,6)==3
h=h+1;
t=t+1;
DetaS(2*h-1,1)=real(B2(i,2))-OrgS(2*h-1,1);
DetaS(2*h,1)=real(PVU(t,1))^2+imag(PVU(t,1))^2-real(B2(i,3))^2-imag(B2(i,3))^2;
end
end
DetaS
%创建I
i=zeros(n-1,1);
h=0;
fori=1:
n
ifi~=isb
h=h+1;
I(h,1)=(OrgS(2*h-1,1)-OrgS(2*h,1)*sqrt(-1))/conj(B2(i,3));
end
end
I
%创建Jacbi
Jacbi=zeros(2*n-2);
h=0;
k=0;
fori=1:
n
ifB2(i,6)==2
h=h+1;
forj=1:
n
ifj~=isb
k=k+1;
ifi==j
Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3))+imag(I(h,1));
Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3))+real(I(h,1));
Jacbi(2*h,2*k-1)=-Jacbi(2*h-1,2*k)+2*real(I(h,1));
Jacbi(2*h,2*k)=Jacbi(2*h-1,2*k-1)-2*imag(I(h,1));
else
Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3));
Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3));
Jacbi(2*h,2*k-1)=-Jacbi(2*h-1,2*k);
Jacbi(2*h,2*k)=Jacbi(2*h-1,2*k-1);
end
ifk==(n-1)
k=0;
end
end
end
end
end
k=0;
fori=1:
n
ifB2(i,6)==3
h=h+1;
forj=1:
n
ifj~=isb
k=k+1;
ifi==j
Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3))+imag(I(h,1));
Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3))+real(I(h,1));
Jacbi(2*h,2*k-1)=2*imag(B2(i,3));
Jacbi(2*h,2*k)=2*real(B2(i,3));
else
Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3));
Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3));
Jacbi(2*h,2*k-1)=0;
Jacbi(2*h,2*k)=0;
end
ifk==(n-1)
k=0;
end
end
end
end
end
Jacbi
DetaU=zeros(2*n-2,1);
DetaU=inv(Jacbi)*DetaS;
DetaU
%修正节点电压
j=0;
fori=1:
n
ifB2(i,6)==2
j=j+1;
B2(i,3)=B2(i,3)+DetaU(2*j,1)+DetaU(2*j-1,1)*sqrt(-1);
end
end
fori=1:
n
ifB2(i,6)==3
j=j+1;
B2(i,3)=B2(i,3)+DetaU(2*j,1)+DetaU(2*j-1,1)*sqrt(-1);
end
end
B2
Times=Times+1;%迭代次数加1
end
Times
一个原始数据的例子
节点数5
支路数5
平衡节点编号5
精度pr0.000001
B1(支路参数矩阵)
[120.04+0.25i0.5i10;130.1+0.35i010;230.08+0.30i0.5i10;420.015i01.051;530.03i01.051]
B2(节点参数矩阵)
[0-1.6-0.8i1002;0-2-1i1002;0-3.7-1.3i1002;05+0i1.051.0503;001.051.0501]
X(节点号和对地参数)
[10;20;30;40;50]
程序2:
%本程序采用牛顿-拉夫讯法对一个五节点电力系统进行潮流计算
%n=input('请输入系统节点数目:
n=');
%nl=input('请输入系统支路数目:
nl=');
%ph=input('请输入平衡节点号:
ph=');
%jd=input('请输入误差精度:
jd=');
%B=input('请输入由支路参数形成的矩阵:
B=');
%A=input('请输入各节点参数形成的矩阵:
A=');
clearall
n=5;
nl=5;
ph=5;
jd=1E-6;
%支路参数矩阵
Br=[120.04+0.25*j0.25j1;
130.1+0.35*j01;
230.08+0.30*j0.25*j1;
240.015*j01.05;
350.03*j01.05];
%节点参数矩阵:
2-PQ节点,3-PV节点,4-平衡节点
A=[-1.6-0.8*j102;
-2-1*j102;
-3.7-1.3*j102;
51.051.053;
01.0504];
%雅克比矩阵形成
Y=zeros(n);
e=zeros(1,n);
f=zeros(1,n);
V=zeros(1,n);
fori=1:
nl
p=Br(i,1);
q=Br(i,2);
Y(p,q)=Y(p,q)-1./(Br(i,3)*Br(i,5));
Y(q,p)=Y(p,q);
Y(q,q)=Y(q,q)+1./(Br(i,3))+Br(i,4);
Y(p,p)=Y(p,p)+1./(Br(i,3)*Br(i,5)^2)+Br(i,4);
end
disp('节点导纳矩阵为:
Y=');
disp(Y);
G=real(Y);
B=imag(Y);
fori=1:
n
e(i)=real(A(i,2));
f(i)=imag(A(i,2));
V(i)=A(i,3);
end
fori=1:
n
S(i)=A(i);
end
P=real(S);
Q=imag(S);
%雅克比矩阵求取
Ci=0;
a=1;
NO=2*n;
N=NO-1;
whilea~=0
a=0;
fori=1:
n
ifi~=ph
C(i)=0;
D(i)=0;
forp=1:
n
C(i)=C(i)+G(i,p)*e(p)-B(i,p)*f(p);
D(i)=D(i)+G(i,p)*f(p)+B(i,p)*e(p);
end
P1=e(i)*C(i)+f(i)*D(i);
Q1=f(i)*C(i)-e(i)*D(i);
V2=e(i)^2+f(i)^2;
ifA(i,4)~=3
DP=P(i)-P1;
DQ=Q(i)-Q1;
fork=1:
n
%非平衡节点时,非对角线元素
ifk~=ph&k~=i
X1=-G(i,k)*f(i)+B(i,k)*e(i);
X2=-G(i,k)*e(i)-B(i,k)*f(i);
X3=-X2;
X4=X1;
p=2*i-1;
q=2*k-1;
J(p,q)=X1;
J(p,N)=DP;
m=p+1;
J(m,q)=X3;
J(m,N)=DQ;
q=q+1;
J(p,q)=X2;
J(m,q)=X4;
%非平衡节点时对角线元素
elseifk~=ph&k==i
X1=-D(i)+B(i,i)*e(i)-G(i,i)*f(i);
X2=-C(i)-G(i,i)*e(i)-B(i,i)*f(i);
X3=-C(i)+G(i,i)*e(i)+B(i,i)*f(i);
X4=D(i)+B(i,i)*e(i)-G(i,i)*f(i);
p=2*i-1;
q=2*k-1;
J(p,q)=X1;
J(p,N)=DP;
m=p+1;
J(m,q)=X3;
J(m,N)=DQ;
q=q+1;
J(p,q)=X2;
J(m,q)=X4;
end
end
else
DP=P(i)-P1;
DV=V(i)^2-V2;
fors=1:
n
ifs~=ph&s~=i
X1=-G(i,s)*f(i)+B(i,s)*e(i);
X2=-G(i,s)*e(i)-B(i,s)*f(i);
X5=0;
X6=0;
p=2*i-1;
q=2*s-1;
J(p,q)=X1;
J(p,N)=DP;
m=p+1;
J(m,q)=X5;
J(m,N)=DV;
q=q+1;
J(p,q)=X2;
J(m,q)=X6;
elseifs~=ph&s==i
X1=-D(i)+B(i,i)*e(i)-G(i,i)*f(i);
X2=-C(i)-G(i,s)*e(i)-B(i,s)*f(i);
X5=-2*f(i);
X
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 电力系统 潮流 计算 详解