求解决啊!急!
clc;clear;
tic;
filename=uigetfile('*.*');
a=textread(filename,'','delimiter',',','emptyvalue',0);%读取TXT文件
c=a(1,:); %平衡节点行数
n=a(1,1); %节点数
b=find(a(:,1)==0);
line=b(2)-b(1)-1; %线路共有行数
knum=b(1)+line; %线路结束行数
lineblock=a(b(1)+1:knum,:); %线路的参数
lineNo=lineblock(1:line,1);
linei=lineblock(1:line,2);
linej=lineblock(1:line,3);
liner=lineblock(1:line,4);
linex=lineblock(1:line,5);
lineb=lineblock(1:line,6);
trans=b(3)-b(2)-1; %变压器共有行数
k1=b(2)+1; %变压器开始行数
k2=b(2)+trans; %结束行数
transblock=a(k1:k2,:); %变压器的参数
transNo=transblock(1:trans,1);
transi=transblock(1:trans,2);
transj=transblock(1:trans,3);
transr=transblock(1:trans,4);
transx=transblock(1:trans,5);
transk=transblock(1:trans,6);
branch=b(4)-b(3)-1; %接地支路共有行数
k3=b(3)+1; %接地支路开始行数
k4=b(3)+1+branch; %结束行数
branchblock=a(k3:k4,:); %接地支路参数
branchNo=branchblock(1:branch,1);
branchi=branchblock(1:branch,2);
branchg=branchblock(1:branch,3);
branchb=branchblock(1:branch,4);
power=b(5)-b(4)-1; %节点功率共有行数
k5=b(4)+1; %节点功率开始行数
k6=b(4)+power; %结束行数
powerblock=a(k5:k6,:); %节点功率参数
powerNo=powerblock(1:power,1);
PQi=powerblock(1:power,2); %powi
PG=powerblock(1:power,3); %powpgi
QG=powerblock(1:power,4); %powqgi
PD=powerblock(1:power,5); %powpdi
QD=powerblock(1:power,6); %powqdi
PV=b(6)-b(5)-1; %PV节点共有行数
k7=b(5)+1; %PV节点开始行数
k8=b(5)+PV; %结束行数
PVblock=a(k7:k8,:); %PV节点参数
PVNo=PVblock(1:PV,1);
PVi=PVblock(1:PV,2);
PVV=PVblock(1:PV,3);
PVqmin=PVblock(1:PV,4);
PVqmax=PVblock(1:PV,5);
clear a;
G=zeros(n,n);
B=zeros(n,n);
for k=(1:line); %线路导纳矩阵
i=linei(k);
j=linej(k);
r=liner(k);
x=linex(k);
b=lineb(k);
GIJ=r/(r*r+x*x);
BIJ=-x/(r*r+x*x);
G(j,i)=-GIJ;
G(i,j)=-GIJ;
B(j,i)=-BIJ;
B(i,j)=-BIJ;
b=lineb(k);
G(i,i)=G(i,i)+GIJ;
B(i,i)=B(i,i)+BIJ+b;
G(j,j)=G(j,j)+GIJ;
B(j,j)=B(j,j)+BIJ+b;
end
for k=(1:trans); %变压器导纳矩阵
i=transi(k);
j=transj(k);
r=transr(k);
x=transx(k);
j=transj(k);
k0=transk(k);
GIJ=r/(r*r+x*x);
BIJ=-x/(r*r+x*x);
G(i,j)=-GIJ/k0;
G(j,i)=G(i,j);
B(i,j)=-BIJ/k0;
B(j,i)=B(i,j);
G(i,i)=G(i,i)+GIJ/k0/k0;
B(i,i)=B(i,i)+BIJ/k0/k0;
G(j,j)=G(j,j)+GIJ;
B(j,j)=B(j,j)+BIJ;
end
for k=(1:branch); %接地支路导纳矩阵
i=branchi(k);
G(i,i)=G(i,i)+branchg(k);
B(i,i)=B(i,i)+branchb(k);
end
Y=sparse(zeros(n,n));
A=sparse(zeros(n,n));
clear j;
X=G+j*B;
X=sparse(X);
Y=(abs(X)); %导纳矩阵完成
A=angle(X); %相角
dp=sparse(PQi,1,(PG-PD)./100,n,n); %节点注入有功功率
dq=sparse(PQi,1,(QG-QD)./100,n,1); %节点注入无功功率
V0=sparse(ones(1,n)'); %节点电压大小初值
VA=sparse(zeros(n,1)); %相角初值
V0(c(2))=c(3); %平衡节点电压
V0(PVi)=PVV; %pv节点电压
ee=1; %精度
k=1; %迭代次数
while(ee>c(4)&&k<20)
ij=repmat(VA,1,n)-repmat(VA',n,1)-angle(X);%角
V=V0.*(cos(VA)+j*sin(VA));%节点电压复数形式
delt0=conj(X*V);
delt1=V.*delt0; %代入节点电压求出的功率
deltp=dp-real(delt1);%有功修正量
deltq=dq-imag(delt1);%无功修正量
deltp(c(2))=0; %有功修正量(处理平衡节点)
deltq(c(2))=0; %无功修正量(处理平衡节点)
deltq(PVi)=0; %无功修正量(处理pv节点)
delt=; %p,q修正量数组
h=diag(diag(V0)*(Y.*sin(ij))*V0)-diag(V0)*(Y.*sin(ij))*diag(V0);
nn=-diag(V0)*(Y.*cos(ij))-diag((Y.*cos(ij))*V0); %n矩阵
jj=-diag(diag(V0)*(Y.*cos(ij))*V0)+diag(V0)*(Y.*cos(ij))*diag(V0);
l=-diag(V0)*(Y.*sin(ij))-diag((Y.*sin(ij))*V0);
%雅克比矩阵完成
nn(:,PVi)=0; %处理pv节点
jj(PVi,:)=0; %处理pv节点
l(PVi,:)=0; %处理pv节点
l(:,PVi)=0; %处理pv节点
l=l+sparse(PVi,PVi,1,n,n);%形成l
ykb=; %初步形成雅可比矩阵
ykb(c(2),:)=0; %处理平衡节点a
ykb(:,c(2))=0; %处理平衡节点a
ykb(c(2),c(2))=1; %处理平衡节点a
ykb(c(2)+n,:)=0; %处理平衡节点v
ykb(:,c(2)+n)=0; %处理平衡节点v
ykb(c(2)+n,c(2)+n)=1;%处理平衡节点v
jie=ykb\delt;
delta=jie(1:n); %相角修正量
deltv=jie((n+1):2*n); %电压大小修正量
VA=VA-delta;
V0=V0-deltv;
ee=max(abs(deltv))
k=k+1;
end
VA=VA.*180./pi %相角
V0 %电压幅值
toc;
求解释,求解答,改了很多次都没用
??? Error using ==> minus
Matrix dimensions must agree.
Error in ==> liweiming2 at 122
deltp=dp-real(delt1);%有功修正量 编程高手解释一下哦! 回复 1# mvpmqc
您的问题,有可能是读取的数据格式,通常会错误,多半原因是格式问题,
除非你上传你的文本文件才知道错误在哪里。
页:
[1]
