《潮流计算方法》由会员分享，可在线阅读，更多相关《潮流计算方法（15页珍藏版）》请在金锄头文库上搜索。
1、由于本人参加我们电气学院的电气小课堂，主讲的是计算机算法计算潮流这章，所以潜 心玩了一个星期，下面整理给大家分享下。本人一个星期以来的汗水，弄清楚了计算机算法计算潮流的基础，如果有什么不懂的可 以发信息到邮箱：* 接下来开始弄潮流的优化问题，吼吼！电力系统的潮流计算的计算机算法：以MATLAB 为环境 这里理论不做过多介绍，推荐一本专门讲解电力系统分析的计算机算法的书籍电力系统分析的计算机算法邱晓燕、刘天琪编著。这里以这本书上的例题【2-1】说明计算机算法计算的过程，分别是牛顿拉弗逊算法的 直角坐标和极坐标算法、P-Q分解算法。主要是简单的网络的潮流计算，其实简单网络计算 和大型网络计算并无本
2、质区别，代码里面只需要修改循环迭代的N即可，这里旨在弄清计 算机算法计算潮流的本质。代码均有详细的注释.其中简单的高斯赛德尔迭代法是以我们的电稳教材为例子讲，其实都差不多，只要把导 纳矩阵Y给你，节点的编号和分类给你，就可以进行计算了，不必要找到原始的电气接线 图。理论不多说，直接上代码： 简单的高斯赛德尔迭代法： 这里我们只是迭代算出各个节点的电压值，支路功率并没有计算。S_ij=P_ij+Q_ij=V_i(V_i*- V_j*) * y_ij*可以计算出各个线路的功率在显示最终电压幅角的时候注意在M ATLAB 里面默认的是弧度的形式，需要转化成角 度显示。clear;clc;%电稳书Pa
3、ge 102例题3-5%计算网络的潮流分布 - 高斯-赛德尔算法%其中节点1是平衡节点%节点2、3是PV节点，其余是PQ节点% 如果节点有对地导纳支路%需将对地导纳支路算到自导纳里面%输入原始数据，每条支路的导纳数值，包括自导和互导纳； y=zeros(5,5);y(1,2)=1/(0.0194+0.0592*1i);y(1,5)=1/(0.054+0.223*1i);y(2,3)=1/(0.04699+0.198*1i);y(2,4)=1/(0.0581+0.1763*1i);%由于电路网络的互易性，导纳矩阵为对称的矩阵 for i=1:1:5for j=1:1:5y(j,i)=y(i,j)
4、;endend%节点导纳矩阵的形成Y=zeros(5,5);%求互导纳for i=1:1:5for j=1:1:5if i=jY(i,j)=-y(i,j);endendend%求自导纳for i=1:1:5%这句话是说将y矩阵的第i行的所有元素相加，得到自导纳的值 Y(i,i)=sum(y(i,:);end %上面求得的自导纳不包含该节点的对地导纳数值，需要加上 Y(2,2)=Y(2,2)+0.067*1i;Y(3,3)=Y(3,3)+0.022*1i;Y(4,4)=Y(4,4)+0.0187*1i;Y(5,5)=Y(5,5)+0.0246*1i;% 导纳矩阵的实部和虚部G = real(Y)
5、;B = imag(Y);Qc2=0;Qc3=0;% 原始节点功率 %这里电源功率为正，负荷功率为负S(1)=0;S(2)=-0.217-0.121*1i+Qc2*1i;S(3)=-0.749-0.19*1i+Qc3*1i;S(4)=-0.658+0.039*1i;S(5)=-0.076-0.016*1i;% 节点功率的 P QP = real(S);Q = imag(S);%下面是两个PV节点的无功初始值Q(2) = 0;Q(3) = 0;U=ones(5,1);%1列5行的1矩阵%节点电压初始值U(1)=1.06;U(2)=1.045;U(3)=1.01;U_reg=U;Sum_YUO=O
6、; %中间变量Sum_YUl=O; %中间变量for cont=1:1:6%这里的cont是迭代次数for i=2:1:5forj=1:1:iif i=jSum_YU0 = Sum_YU0 + Y(i,j)*U_reg(j);endendforj=i+1:1:5Sum_YU1 = Sum_YU1 + Y(i,j)*U(j);endU(i)=(P(i)-Q(i)*1i)/conj(U(i)-Sum_YU0-Sum_YU1)/Y(i,i U_reg(i)=U(i);%PV节点计算%下面是把求出的U2、U3只保留其相位，幅值不变 if i=2angle_U2 = angle(U(2); U(2)=1
7、.045*cos(angle_U2)+1.045*sin(angle_U2)*1i;Q(2)=imag( U(2)*( conj(Sum_YU0) + conj(Sum_YU1) + conj(Y(2,2)*U(2) ) );endif i=3angle_U3 = angle(U(3); U(3)=1.01*cos(angle_U3)+1.01*sin(angle_U3)*1i;Q(3)=imag( U(3)*( conj(Sum_YU0) + conj(Sum_YU1) + conj(Y(3,3)*U(3) ) );end% 下面做越界检查%if Q(4)Q_Max% Q(4) = Q_Ma
8、x;%end%if Q(4) 1e-6),cont=cont+1;%for cont=1:1:3% 下面开始计算 delta_P/delta_Q/delta_Ufor i=2:1:5for j=1:1:5Sum_GB1=Sum_GB1 + ( G(i,j)*e(j) - B(i,j)*f(j) ); Sum_GB2=Sum_GB2 + ( G(i,j)*f(j) + B(i,j)*e(j) ); end delta_P(i)=P(i)-e(i)*Sum_GB1-f(i)*Sum_GB2;if i=2 & i=3%不为节点2,3则计算无功delta_Q(i)=Q(i)-f(i)*Sum_GB1+
9、e(i)*Sum_GB2;end if i=2 | i=3%这里计算delta_U的值，始终为零delta_U(i)=U(i)入2-( e(i)入2 + f(i)入2);end Sum_GB1=0;Sum_GB2=0;end%下面计算雅克比矩阵J=zeros(8,8);for ii=2:1:5i=ii-1;forj=1:1:5Sum_GB1=Sum_GB1 +( G(ii,j)*e(j)- B(ii,j)*f(j) );Sum_GB2=Sum_GB2 +( G(ii,j)*f(j)+ B(ii,j)*e(j) );endforjj=2:1:5j=jj-1;if ii=2 &ii=3%PQ节点if ii=jjJ(2*i-1,2*i-1)=-Sum_GB1-G(ii,ii)*e(ii)-B(ii,ii)*f(ii);J(2*i-1,2*i)=-Sum_GB2+B(ii,ii)*e(ii)-G(i
