mvpmqc 发表于 2011-5-30 19:32:46

求解决啊!急!

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);%有功修正量

mingyu 发表于 2011-5-31 08:35:17

编程高手解释一下哦!

sdad 发表于 2011-6-20 22:19:57

回复 1# mvpmqc



您的问题,有可能是读取的数据格式,通常会错误,多半原因是格式问题,

除非你上传你的文本文件才知道错误在哪里。
页: [1]
查看完整版本: 求解决啊!急!

招聘斑竹