潮流程序matlab

luojinwen · 发表于 2010-10-28 11:30:47

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

您需要登录才可以下载或查看，没有账号？立即加入

×

clear                      %清除工作空间原有的变量
clc
a=textread('IEEE14.dat');       %读数据文件
n=a(1,1);                %节点数
line=a(1,2);             %支路数
blance=a(1,3);             %平衡节点
SB=a(1,4);                %基准功率
kmax=a(1,5);             %最大迭代次数
ab=a(2,1);                %计算精度
t=find(a(:,1)==0);
xianlu=a([3:t(1)-1],:);    %形成线路参数矩阵
jiedi=a([t(1)+1:t(2)-1],:);  %接地支路参数矩阵
bianya=a([t(2)+1:t(3)-1],:); %形成变压器支路参数矩阵
jiedian=a([t(3)+1:t(4)-1],:);%节点功率数据矩阵
pv=a([t(4)+1:t(5)-1],:);    %pv节点数据矩阵
%读线路参数
linei=xianlu(:,2)';
linej=xianlu(:,3)';
r=xianlu(:,4)';
x=xianlu(:,5)';
b0=xianlu(:,6)';

g=r./(r.*r+x.*x);
b=-x./(r.*r+x.*x);
G=sparse(linei,linej,-g,n,n)+sparse(linej,linei,-g,n,n)+sparse(linei,linei,g,n,n)+sparse(linej,linej,g,n,n);
B=sparse(linei,linej,-b,n,n)+sparse(linej,linei,-b,n,n)+sparse(linei,linei,b+b0,n,n)+sparse(linej,linej,b+b0,n,n);
% 读变压器参数
linei=bianya(:,2)';
linej=bianya(:,3)';
r=bianya(:,4)';
x=bianya(:,5)';
k0=bianya(:,6)';
g=r./(r.*r+x.*x);
b=-x./(r.*r+x.*x);
G=G+sparse(linei,linej,-g./k0,n,n)+sparse(linej,linei,-g./k0,n,n)+sparse(linei,linei,g./k0./k0,n,n)+sparse(linej,linej,g,n,n);
B=B+sparse(linei,linej,-b./k0,n,n)+sparse(linej,linei,-b./k0,n,n)+sparse(linei,linei,b./k0./k0,n,n)+sparse(linej,linej,b,n,n);
% 读接地支路参数
jiedidian=jiedi(:,1)';
branchib=jiedi(:,2)';
B=B+sparse(jiedidian,jiedidian,branchib,n,n);
%节点导纳矩阵
Y=sparse(G+i*B);
% 读节点功率数据
pg=jiedian(:,2)/SB;
pd=jiedian(:,4)/SB;
qg=jiedian(:,3)/SB;
qd=jiedian(:,5)/SB;

p=zeros(n,1);
q=zeros(n,1);
p(jiedian(:,1))=pg-pd;
q(jiedian(:,1))=qg-qd;
% 读pv节点数据
pvjiedian=pv(:,1);
vpv=pv(:,2);
v=ones(n,1);                %电压初始值
v(pvjiedian)=vpv; %pv节点(及平衡节点)电压已知
delt=zeros(n,1);    %相角初始值
deltp=zeros(n,1);    % 失配有功功率初值
deltq=zeros(n,1);    % 失配无功功率初值

%迭代求潮流计算
for diedai=1:10             %置最大循环次数
clear i
% 失配功率方程
   V=v.*cos(delt)+i*v.*sin(delt);
   a2=conj(Y*V).*V;
   deltp=p-real(a2);
   deltq=q-imag(a2);

   % 求雅可比矩阵
   for i=1:n
      for j=1:n
            if i~=j
               H(i,j)=-v(i)*v(j)*(G(i,j)*sin(delt(i)-delt(j))-B(i,j)*cos(delt(i)-delt(j)));
            else
               H(i,i)=v(i)*v(i)*B(i,i)+q(i)-deltq(i);
            end
      end
   end

   for i=1:n
      for j=1:n
            if i~=j
               N(i,j)=-v(i)*(G(i,j)*cos(delt(i)-delt(j))+B(i,j)*sin(delt(i)-delt(j)));
            else
               N(i,i)=-v(i)*G(i,i)-(p(i)-deltp(i))/v(i);
            end
      end
   end

   for i=1:n
      for j=1:n
            if i~=j
               J(i,j)=v(i)*v(j)*(G(i,j)*cos(delt(i)-delt(j))+B(i,j)*sin(delt(i)-delt(j)));
            else
               J(i,i)=v(i)*v(i)*G(i,i)-(p(i)-deltp(i));
            end
      end
   end

   for i=1:n
      for j=1:n
            if i~=j
               L(i,j)=-v(i)*(G(i,j)*sin(delt(i)-delt(j))-B(i,j)*cos(delt(i)-delt(j)));
            else
               L(i,i)=v(i)*B(i,i)-(q(i)-deltq(i))/v(i);
            end
      end
   end
   %修正雅克比矩阵
   Ya=[H N;J L];
   Ya([blance blance+n],:)=0;
   Ya(:,[blance blance+n])=0;
   Ya(blance,blance)=1; %平衡节点
   Ya(blance+n,blance+n)=1;
   Ya(pvjiedian+n,:)=0;
   Ya(:,pvjiedian+n)=0;
   for h=1:length(pvjiedian)
         Ya(pvjiedian(h)+n,pvjiedian(h)+n)=1;
      end

   %
   deltp(blance)=0;       %平衡节点deltp=0
   deltq(blance)=0;       %平衡节点deltq=0

   deltq(pvjiedian)=0;    %pv节点deltq=0
   deltpq=[deltp;deltq];
      if max(abs(deltpq))<ab    % 如果结果满足精度
            break                % 退出循环
      end
      R=Ya\(-deltpq);
      delt=delt+R([1:n]);
      deltt=delt*180/pi
      v=v+R([n+1:2*n])
end

数据格式

14	20	1	100	20	0.1
1.00E-06	1
1	1	2	0.01938	0.05917	0.0264
2	1	5	0.05403	0.22304	0.0246
3	2	3	0.04699	0.19797	0.0219
4	2	4	0.05811	0.17632	0.017
5	2	5	0.05695	0.17388	0.0173
6	3	4	0.06701	0.17103	0.0064
7	4	5	0.01335	0.04211	0
11	6	11	0.09498	0.1989	0
12	6	12	0.12291	0.25581	0
13	6	13	0.06615	0.13027	0
14	7	8	0	0.17615	0
15	7	9	0	0.11001	0
16	9	10	0.03181	0.0845	0
19	12	13	0.22092	0.19988	0
20	13	14	0.17093	0.34802	0
4	14	9	0.12711	0.27038	0
5	10	11	0.08205	0.19207	0
0
9	0.19
0
1	4	7	0	0.20912	0.978	0.9
2	4	9	0	0.55618	0.969 , r' c: K! k( P* h0.9	1.1
3	5	6	0	0.25202	0.932	0.9
0
1	60	0	0	0
2	65	42.4	21.7	12.7
3	0	23.39	94.2	19
4	0	0	47.8	-3.9
5	0	0	7.6	1.6
6	85	12.24	11.2	7.5
7	0	0	0	0
8	0	17.36	0	0
9	0	0	29.5	16.6
10	0	0	9	5.8
11	0	0	3.5	1.8
12	0	0	6.1	1.6
13	0	0	13.5	5.8
14	0	0	14.9	5
0
1	1.06	-40	50
2	1.045	-40	50
3	1.01	0	40
6	1.07	-30	40
8	1.09	-30	45
0
1	105	2.45	0.005	50	200
2	44.4	3.51	0.005	20	100
6	40.6	3.89	0.005	20	100
0
0
0

xq2256 · 发表于 2010-10-28 20:56:33

很有用！jixuzheyangyigechengxu a!

luojinwen · 发表于 2010-10-29 11:35:12

希望能有帮助，呵呵

wangwangmoon · 发表于 2010-12-24 11:03:53

修正雅可比矩阵没有看懂

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

潮流程序matlab

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

相关帖子

浏览过的版块