牛顿-拉夫逊潮流计算实例
牛顿—拉夫逊潮流计算实例
摘要：主要根据一个简单的电网潮流计算实例，来说明潮流计算的具体步骤，以及计算过程中出现的一些问题。文章具体从潮流计算原理、功率方程、雅可比矩阵、Matlab 源程序等一些方面说明了潮流计算的一些主要内容。
关键词：节点导纳矩阵、功率方程、雅可比矩阵
一、潮流计算的基本原理
电力系统潮流计算是电力系统运行和规划中最基本和最经常的计算，其任务是在已知某些运行参数的情况下，计算出系统中全部的运行参数，一般来说，各个母线所供负荷的功率是已知的，各个节点电压是未知的（平衡节点除外），可以根据网络结构形成节点导纳矩阵，然后由节点导纳矩阵和网络拓扑结构列写功率方程，由于功率方程里功率是已知的，电压的幅值和相角是未知的，这样潮流计算的问题就转化为求解非线性方程组的问题了。为了便于用迭代法解方程组，需要将上述功率方程改写成功率平衡方程，并对功率平衡方程求偏导，得出对应的雅可比矩阵，给未知节点赋电压初值，一般为额定电压，将初值带入功率平衡方程，得到功率不平衡量，这样由功率不平衡量、雅可比矩阵、节点电压不平衡量（未知的）构成了误差方程，解误差方程，得到节点电压不平衡量，节点电压加上节点电压不平衡量构成新的节点电压初值，将新的初值带入原来的功率平衡方程，并重新形成雅可比矩阵，然后计算新的电压不平衡量，这样不断迭代，不断修正，一般迭代三到五次就能收敛。
二、直角坐标形式的功率方程 ∆P =P 1
1 ∆Q =Q -Q -f ∑(G e -B f ) +e ∑(G f +B e ) =0
∆P =P -P -e ∑(G e -B f ) -f ∑(G f +B e ) =0
∆Q =Q -Q -f ∑(G e -B f ) +e ∑(G f +B e ) =0
∆P =P -P -e ∑(G e -B f ) -f ∑(G f +B e ) =0
∆Q =Q -Q -f ∑(G e -B f ) +e ∑(G f +B e ) =0
∆P =P -P -e ∑(G e -B f ) -f ∑(G f +B e ) =0
∆Q =Q -Q -f ∑(G e -B f ) +e ∑(G f +B e ) =0
∆P =P -P -e ∑(G e -B f ) -f ∑(G f +B e ) =0
∆Q =Q -Q -f ∑(G e -B f ) +e ∑(G f +B e ) =0 g 1-P l 1-e 1∑(G 1j e j -B 1j f j ) -f 1∑(G 1j f j +B 1j e j ) =0j ∈i 1j ∈i 1l 11j j 1j j 1j j 1j j j ∈i j ∈i g 12g 2l 222j j 2j j 22j j 2j j j ∈i 2j ∈i 2g 2l 22j j 2j j 22j j 2j j j ∈i j ∈i 3g 3l 333j j 3j j 33j j 3j j j ∈i 3j ∈i 3g 3l 33j j 3j j 33j j 3j j j ∈i j ∈i 4g 4l 444j j 4j j 44j j 4j j j ∈i 4j ∈i 4g 4l 44j j 4j j 44j j 4j j j ∈i j ∈i 5g 5l 555j j 5j j 55j j 5j j
j ∈i
5j ∈i 5g 5l 55j j 5j j 55j j 5j j
j ∈i j ∈i
我们可以利用它来算初始的功率不平衡量，也可以通 直角坐标形式的功率平衡方程如上，
过对它求导形成雅可比矩阵。另外，由于节点1是平衡节点，故节点1的功率平衡方程不参与误差方程的列写，只作为最后计算平衡节点的功率的公式，又由于所给节点除节点1外均为PQ 节点，故无电压平衡方程
三、雅可比矩阵各元素的表达式
当j ≠i 时：
H ij =∂∆P i =B ij e i -G ij f i ∂f j
N ij =∂∆P i =-G ij e i -B ij f i ∂e j
∂∆Q i =B ij f i +G ij e i ∂f j M ij =
L ij =∂∆Q i =-G ij f i +B ij e i ∂e j
∂∆P i =B ii e i -G ii f i -b i ∂f i
∂∆P i =-G ii e i -B ii f i -a i ∂e i
∂∆Q i =G ii e i +B ii f i -a i ∂f i 当j=i时： H ii =N ii =M ii =
L ii =∂∆Q i =B ii e i -G ii f i +b i ∂e i
a i =∑(G ij e j -B ij f j )
其中，i =1
n n
b i =∑(G ij f j +B ij e j )
i =1
四、心得体会
把潮流计算的任务完成后，有一些体会
对潮流计算也比以前有了进一步的体会，在学习潮流计算时，虽然依次学习了节点导纳矩阵，功率方程、雅可比矩阵，但不能将它们联系起来，更不知道其中的原委，通过程序的编写，知道了其中的联系，也知道了每个方程、矩阵在计算中的作用。
在接到这个课外课题时，由于以前没有用过Matlab ，刚开始感觉有些茫然，不知道从何做起，但后来通过看书，对Matlab 有了一定了解，又觉得和C 语言差不多，没什么难的，真正开始写程序时，发现既不是我们想的那么简单，也不是难得无法下手，说它不是那么简单，是因为Matlab 和C 语言有一定区别，它的数据结构主要是矩阵，另外功率方程里的求和部分也不像WORD 里那样，可以用公式编辑器写出，得用循环实现。说它不是难得无法下手，是因为通过我们查阅资料和自己调试，最后完成了潮流计算的程序并得到了收敛
的结果。
五、程序代码
e=[1.06 1 1 1 1];
f=[0 0 0 0 0];
p2=0.2;
q2=0.2;
p3=0.45;
q3=0.15;
p4=0.4;
q4=0.05;
p5=-0.6;
q5=0.1;
G=[1,-0.02,-0.08,0,0;-0.02,0.18,-0.06,-0.06,-0.04;-0.08,-0.06,0.15,-0.01,0;0,-0.06,-0.01,0.15,-0.08;0,-0.04,0,-0.08,0.12];
B=[0.3,-0.06,-0.24,0,0;-0.06,0.54,-0.18,-0.18,-0.12;-0.24,-0.18,0.45,-0.03,0;0,-0.18,-0.03,0.45,-0.24;0,-0.12,0,-0.24,0.36];
k=0;
maxp=0;
maxq=0;
I=[0,0;0,0;0,0;0,0;0,0];
for v=1:15
for n=1:5
I(1,1)=I(1,1)+G(1,n)*e(n)-B(1,n)*f(n);
I(1,2)=I(1,2)+G(1,n)*f(n)+B(1,n)*e(n);
end
for n=1:5
I(2,1)=I(2,1)+G(2,n)*e(n)-B(2,n)*f(n);
I(2,2)=I(2,2)+G(2,n)*f(n)+B(2,n)*e(n);
end
for n=1:5
I(3,1)=I(3,1)+G(3,n)*e(n)-B(3,n)*f(n);
I(3,2)=I(3,2)+G(3,n)*f(n)+B(3,n)*e(n);
end
for n=1:5
I(4,1)=I(4,1)+G(4,n)*e(n)-B(4,n)*f(n);
I(4,2)=I(4,2)+G(4,n)*f(n)+B(4,n)*e(n);
end
for n=1:5
I(5,1)=I(5,1)+G(5,n)*e(n)-B(5,n)*f(n);
I(5,2)=I(5,2)+G(5,n)*f(n)+B(5,n)*e(n);
end
H=[];
N=[];
M=[];
L=[];
J=[];
P2=0.2-e(2)*I(2,1)-f(2)*I(2,2);
Q2=0.2-f(2)*I(2,1)+e(2)*I(2,2);
P3=-0.45-e(3)*I(3,1)-f(3)*I(3,2);
Q3=-0.15-f(3)*I(3,1)+e(3)*I(3,2);
P4=-0.4-e(4)*I(4,1)-f(4)*I(4,2);
Q4=-0.05-f(4)*I(4,1)+e(4)*I(4,2);
P5=-0.6-e(5)*I(5,1)-f(5)*I(5,2);
Q5=-0.1-f(5)*I(5,1)+e(5)*I(5,2);
for m=2:5
for n=2:5
if(m==n)
H(m,m)=B(m,m)*e(m)-G(m,m)*f(m)-I(m,2);
N(m,m)=-G(m,m)*e(m)-B(m,m)*f(m)-I(m,1);
M(m,m)=G(m,m)*e(m)+B(m,m)*f(m)-I(m,1);
L(m,m)=B(m,m)*e(m)-G(m,m)*f(m)+I(m,2);
else
H(m,n)=B(m,n)*e(m)-G(m,n)*f(m);
N(m,n)=-G(m,n)*e(m)-B(m,n)*f(m);
M(m,n)=-N(m,n);
L(m,n)=H(m,n);
end
end
end
J=[H(2,2),N(2,2),H(2,3),N(2,3),H(2,4),N(2,4),H(2,5),N(2,5);M(2,2),L(2,2),M(2,3),L(2,3),M(2,4),L(2,4),M(2,5),L(2,5);H(3,2),N(3,2),H(3,3),N(3,3),H(3,4),N(3,4),H(3,5),N(3,5);M(3,2),L(3,2),M(3,3),L(3,3),M(3,4),L(3,4),M(3,5),L(3,5);H(4,2),N(4,2),H(4,3),N(4,3),H(4,4),N(4,4),H(4,5),N(4,5);M(4,2),L(4,2),M(4,3),L(4,3),M(4,4),L(4,4),M(4,5),L(4,5);H(5,2),N(5,2),H(5,3),N(5,3),H(5,4),N(5,4),H(5,5),N(5,5);M(5,2),L(5,2),M(5,3),L(5,3),M(5,4),L(5,4),M(5,5),L(5,5)];
A=[];
C=[P2;Q2;P3;Q3;P4;Q4;P5;Q5]
A=-inv(J)*C
maxp=C(1);
maxq=C(2);
for d=1:2:8
if(C(d)>maxp)
maxp=C(d);
end
end
for t=2:2:8
if(C(t)>maxq)
maxq=C(t);
end
end
if(maxp>0.001&&maxq>0.001)
e(2)=e(2)+A(2)
e(3)=e(3)+A(4)
e(4)=e(4)+A(6)
e(5)=e(5)+A(8)
f(2)=f(2)+A(1)
f(3)=f(3)+A(3)
f(4)=f(4)+A(5)
f(5)=f(5)+A(7)
end
end
S1=(1.06+0*i)*[(G(1,1)+i*B(1,1))*(1.06+0*i)+(G(1,2)+i*B(1,2))*(e(2)+i*f(2))+(G(1,3)+i*B(1,3))*(e(3)+i*f(3))]
S12=(1.06+0*i)*[(1.06-e(2)-i*f(2))*(G(1,2)+i*B(1,2))]
S21=(e(2)+i*f(2))*[e(2)+i*f(2)-1.06]*(G(2,1)+i*B(2,1))
S13=(1.06+0*i)*[(1.06-e(3)-i*f(3))*(G(1,3)+i*B(1,3))]
S31=(e(3)+i*f(3))*[e(3)+i*f(3)-1.06]*(G(3,1)+i*B(3,1))
S23=(e(2)+i*f(2))*[e(2)+i*f(2)-e(3)-f(3)*i]*(G(2,3)+i*B(2,3))
S32=(e(3)+i*f(3))*[e(3)+i*f(3)-e(2)-i*f(2)]*(G(3,2)+i*B(3,2))
S24=(e(2)+i*f(2))*[e(2)+i*f(2)-e(4)-f(4)*i]*(G(2,4)+i*B(2,4))
S42=(e(4)+i*f(4))*[e(4)+i*f(4)-e(2)-f(2)*i]*(G(4,2)+i*B(4,2))
S25=(e(2)+i*f(2))*[e(2)+i*f(2)-e(5)-f(5)*i]*(G(2,5)+i*B(2,5))
S52=(e(5)+i*f(5))*[e(5)+i*f(5)-e(2)-f(2)*i]*(G(5,2)+i*B(5,2))
S34=(e(3)+i*f(3))*[e(3)+i*f(3)-e(4)-f(4)*i]*(G(3,4)+i*B(3,4))
S43=(e(4)+i*f(4))*[e(4)+i*f(4)-e(3)-f(3)*i]*(G(4,3)+i*B(4,3))
S45=(e(4)+i*f(4))*[e(4)+i*f(4)-e(5)-f(5)*i]*(G(4,5)+i*B(4,5))
S54=(e(4)+i*f(5))*[e(5)+i*f(5)-e(4)-f(4)*i]*(G(5,4)+i*B(5,4))
参考资料：《电力系统分析》
《电力系统稳态分析》
《数值分析》
《数学建模与实验》
《Matlab 程序设计》