求解决啊！急!

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=[deltp;deltq];    %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=[h nn;jj l];       %初步形成雅可比矩阵
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

+ h. G: _% b) U o! l; u) s9 e
您的问题，有可能是读取的数据格式，通常会错误，多半原因是格式问题，
, D7 k6 G2 o. Y+ Q) x. c
; p3 C; J2 y$ j6 l+ N; `9 T除非你上传你的文本文件才知道错误在哪里。

账号		自动登录	找回密码
密码			立即加入

[讨论] 求解决啊！急!

马上加入，结交更多好友，共享更多资料，让你轻松玩转电力研学社区！

相关帖子

浏览过的版块