您好,欢迎来到意榕旅游网。
搜索
您的当前位置:首页牛顿拉夫逊潮流计算

牛顿拉夫逊潮流计算

来源:意榕旅游网


牛顿—拉夫逊法潮流计算 一、 潮流计算的基本原理 二、实际电力系统中的节点类型实际电力系统中的节点类型 网络中各节点的性质:负荷节点:给定功率P、Q1. 负荷节点:如图中的3、4节点2. 发电机节点:发电机节点负荷节点s4混合节点4512如图中的节点1,可能有两种情况:给定P、Q运行,给定P、V运行3. 负荷发电机混合节点:PQ节点,如图中的节点2负荷节点3s2s3过渡节点4. 过渡节点:PQ为0的给定PQ节点,如图的节点5 5 潮流计算中节点类型划分 1. PQ节点:已知P、Q三、潮流计算中节点类型的划分PQ节点负荷、过渡节点,PQ给定的发电机节点,大部分节点PQ2.PV节点:已知P、V平衡节点s44PV节点152给定PV的发电机节点,具有可调电源的变电所,少量节点PV33. 平衡节点+基准节点:已知V、s2s3PQ节点PQ节点也称为松弛节点,摇摆节点

6

采用极坐标,节点电压表示为

ViViiVicosijsini

节点功率将写成

PiViVj(GijcosijBijsinij)j1 (1) n

QiViVj(GijsinijBijcosij)j1n式中,ijij,是i、j两节点电压的相角差。

方程式把节点功率表示为节点电压的幅值和相角的函数。在有n个节点的系统中,假定第1~m号节点为PQ节点,第m1~n1号节点为PV节点,第n号节点为平衡节点。Vn和n是给定的,PV节点的电压幅值Vm1~Vn1也是给定的。因此,只剩下n1个节点的电压相角1,2,,n1和m个节点的电压幅值

V1,V2,,Vm是未知量。

实际上,对于每一个PQ节点或每一个PV节点都可以列写一个有功功率不平衡量方程式

nPiPisPiPisViVj(GijcosijBijsinij)0j1i1,2,,n1 (2)

而对于每一个PQ节点还可以再列写一个无功功率不平衡量方程式

nQiQisQiQisViVj(GijsinijBijcosij)0j1i1,2,,m (3)

式(2)和式(3)一共包含了n1m个方程式,正好同未知量的数目相等,而比直角坐标形式的方程少了n1m个。

对于方程式(2)和式(3)可以写出修正方程式如下

PHQKNLδ1VVD2 (4)

式中

P1P2;PPn1V1V2;VVmHQ1Q2;QQmV1PijVD2V212n1 (5) Vm是n1n1阶方阵,其元素为HijPiVj;N是n1m阶矩阵,其元

Qij素为NijVj;K是mn1阶矩阵,其元素为KijQiVj;L是mm阶

方阵,其元素为LijVj。

在这里把节点不平衡功率对节点电压幅值的偏导数都乘以该节点电压,相应地把节点电压的修正量都除以该节点的电压幅值,这样,雅可比矩阵元素的表达式就具有比较整齐的形式。

对式(2)和式(3)求偏导数,可以得到雅可比矩阵元素的表达式如下 当ij时

HijViVj(GijsinijBijcosij)NijViVj(GijcosijBijsinij)KijViVj(GijcosijBijsinij)LijViVj(GijsinijBijcosij) (6)

当ij时

2HiiViBiiQi2NiiViGiiPi 2KiiViGiiPi2LiiViBiiQi (7)

二、 程序流程框图

输入原始数据 形成节点导纳阵 给定电压初值e(0)、f(0) k=0 根据公式计算P(k),Q(k)及V2(k) 用公式计算雅克比矩阵各元素 高斯法解修正方程,求e(k),f(k) e(k+1)e(k)e(k),f(k+1)f(k)f(k)修正节点电压 k=k+1 否 max{e(k),f(k) }?是 计算平衡节点功率及全部线路功率 输出

三、 程序清单

% B1为支路参数矩阵:1、支路首端号;2、末端号;3、支路阻抗;4、支路对 地导纳;5、支路的变比;6、支路首端处于K侧为1,1侧为0

% B2为节点参数矩阵:1、该节点发电机功率;2、该节点负荷功率;3、节点电压初始值4、PV节点电压V的给定值;5、节点所接的无功补偿设备的容量 6、节点分类标号:1为平衡节点;2为PQ节点;3为PV节点; clear;

n=4;%input('请输入节点数:n='); nl=4;%input('请输入支路数:nl=');

isb=1;%input('请输入平衡母线节点号:isb='); pr=0.00001;%input('请输入误差精度:pr='); B1=[

1 2 0.1880i -0.6818i 1 1; 1 3 0.1302+0.2479i 0.0129i 1 0; 1 4 0.1736+0.3306i 0.0172i 1 0; 3 4 0.2603+0.4959i 0.0259i 1 0;] B2=[

0 0 1 0 0 2; 0 0.5+0.3i 1 0 0 2; 0.2 0 1.05 1.05 0 3; 0 0.15+0.1i 1.05 1.05 0 1;]

Y=zeros(n);e=zeros(1,n);f=zeros(1,n);V=zeros(1,n);sida=zeros(1,n);S1=zeros(nl); % % %--------------------------------------------------- for i=1:nl

%支路数

if B1(i,6)==0 %左节点处于1侧

p=B1(i,1);q=B1(i,2);

else %左节点处于K侧 p=B1(i,2);q=B1(i,1); end

Y(p,q)=Y(p,q)-1./(B1(i,3)*B1(i,5));

%非对角元

Y(q,p)=Y(p,q); %非对角元

Y(q,q)=Y(q,q)+1./(B1(i,3)*B1(i,5)^2)+B1(i,4)./2; %对角元K侧 Y(p,p)=Y(p,p)+1./B1(i,3)+B1(i,4)./2; end

%对角元1侧

%求导纳矩阵 disp('导纳矩阵 Y='); disp(Y)

%---------------------------------------------------------- G=real(Y);B=imag(Y); for i=1:n

%分解出导纳阵的实部和虚部 %给定各节点初始电压的实部和虚部

e(i)=real(B2(i,3)); f(i)=imag(B2(i,3)); V(i)=B2(i,4); end for i=1:n

%给定各节点注入功率

%PV节点电压给定模值

S(i)=B2(i,1)-B2(i,2); B(i,i)=B(i,i)+B2(i,5); end

%i节点注入功率SG-SL

%i节点无功补偿量

%=========================================================== P=real(S);Q=imag(S); %分解出各节点注入的有功和无功功率 ICT1=0;IT2=1;N0=2*n;N=N0+1;a=0; %迭代次数ICT1、a;不满足收敛要求的节点数IT2

while IT2~=0 % N0=2*n 雅可比矩阵的阶数;N=N0+1扩展列 IT2=0;a=a+1; for i=1:n if i~=isb

%非平衡节点

C(i)=0;D(i)=0; for j1=1:n

C(i)=C(i)+G(i,j1)*e(j1)-B(i,j1)*f(j1);%Σ(Gij*ej-Bij*fj) D(i)=D(i)+G(i,j1)*f(j1)+B(i,j1)*e(j1);%Σ(Gij*fj+Bij*ej) end

P1=C(i)*e(i)+f(i)*D(i);%节点功率P计算eiΣ(Gij*ej-Bij*fj)+fiΣ(Gij*fj+Bij*ej) Q1=C(i)*f(i)-e(i)*D(i);%节点功率Q计算fiΣ(Gij*ej-Bij*fj)-eiΣ(Gij*fj+Bij*ej)

%求i节点有功和无功功率P',Q'的计算值 V2=e(i)^2+f(i)^2;

%电压模平方

%======= 以下针对非PV节点来求取功率差及Jacobi矩阵元素 ========= if B2(i,6)~=3 DP=P(i)-P1;

%非PV节点

%节点有功功率差

DQ=Q(i)-Q1; %节点无功功率差

%================= 求取Jacobi矩阵 =================== for j1=1:n if j1~=isb&j1~=i

%非平衡节点&非对角元

X1=-G(i,j1)*e(i)-B(i,j1)*f(i); % dP/de=-dQ/df X2=B(i,j1)*e(i)-G(i,j1)*f(i); % dP/df=dQ/de X3=X2; % X2=dp/df X3=dQ/de X4=-X1; % X1=dP/de X4=dQ/df p=2*i-1;q=2*j1-1;

J(p,q)=X3;J(p,N)=DQ;m=p+1; % X3=dQ/de J(p,N)=DQ节点无功功率差 J(m,q)=X1;J(m,N)=DP;q=q+1; % X1=dP/de J(m,N)=DP节点有功功率差 J(p,q)=X4;J(m,q)=X2; % X4=dQ/df X2=dp/df elseif j1==i&j1~=isb

%非平衡节点&对角元

X1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i);% dP/de X2=-D(i)+B(i,i)*e(i)-G(i,i)*f(i);% dP/df X3=D(i)+B(i,i)*e(i)-G(i,i)*f(i); % dQ/de X4=-C(i)+G(i,i)*e(i)+B(i,i)*f(i);% dQ/df

p=2*i-1;q=2*j1-1;J(p,q)=X3;J(p,N)=DQ;%扩展列△Q m=p+1;

J(m,q)=X1;q=q+1;J(p,q)=X4;J(m,N)=DP;%扩展列△P J(m,q)=X2; end end else

%========== 下面是针对PV节点来求取Jacobi矩阵的元素 =========== DP=P(i)-P1;

% PV节点有功误差 % PV节点电压误差

DV=V(i)^2-V2; for j1=1:n if j1~=isb&j1~=i

%非平衡节点&非对角元

X1=-G(i,j1)*e(i)-B(i,j1)*f(i); % dP/de X2=B(i,j1)*e(i)-G(i,j1)*f(i); % dP/df X5=0;X6=0;

p=2*i-1;q=2*j1-1;J(p,q)=X5;J(p,N)=DV; % PV节点电压误差 m=p+1;

J(m,q)=X1;J(m,N)=DP;q=q+1;J(p,q)=X6; % PV节点有功误差 J(m,q)=X2;

elseif j1==i&j1~=isb

%非平衡节点&对角元

X1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i);% dP/de X2=-D(i)+B(i,i)*e(i)-G(i,i)*f(i);% dP/df X5=-2*e(i); X6=-2*f(i);

p=2*i-1;q=2*j1-1;J(p,q)=X5;J(p,N)=DV; % PV节点电压误差 m=p+1;

J(m,q)=X1;J(m,N)=DP;q=q+1;J(p,q)=X6; % PV节点有功误差 J(m,q)=X2; end end end end end

%以上为求雅可比矩阵的各个元素及扩展列的功率差或电压差

for k=3:N0 % N0=2*n (从第三行开始,第一、二行是平衡节点) k1=k+1;N1=N; % N=N0+1 即 N=2*n+1扩展列△P、△Q 或 △U

for k2=k1:N1 % 从k+1列的Jacobi元素到扩展列的△P、△Q 或 △U

J(k,k2)=J(k,k2)./J(k,k);% 用K行K列对角元素去除K行K列后的非对角元素进行规格化 end

J(k,k)=1; % 对角元规格化K行K列对角元素赋1 %==================== 回代运算================= if k~=3 % 不是第三行 k > 3 k4=k-1; for k3=3:k4 for k2=k1:N1

% 用k3行从第三行开始到当前行的前一行k4行消去 % k3行后各行上三角元素

J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去运算(当前行k列元素消为0) end %用当前行K2列元素减去当前行k列元素乘以第k行K2列元素 J(k3,k)=0; %当前行第k列元素已消为0 end

if k==N0 %若已到最后一行 break; end

%================== 前代运算=========================== for k3=k1:N0 % 从k+1行到2*n最后一行

for k2=k1:N1 % 从k+1列到扩展列消去k+1行后各行下三角元素 J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去运算

end %用当前行K2列元素减去当前行k列元素乘以第k行K2列元素 J(k3,k)=0; %当前行第k列元素已消为0 end

else %是第三行k=3

%====================== 第三行k=3的前代运算============= for k3=k1:N0 %从第四行到2n行(最后一行) for k2=k1:N1 %从第四列到2n+1列(即扩展列)

J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去运算(当前行3列元素消为0)

end %用当前行K2列元素减去当前行3列元素乘以第三行K2列元素 J(k3,k)=0; %当前行第3列元素已消为0 end end end

%====上面是用线性变换方式高斯消去法将Jacobi矩阵化成单位矩阵===== for k=3:2:N0-1 L=(k+1)./2;

e(L)=e(L)-J(k,N); %修改节点电压实部 k1=k+1;

f(L)=f(L)-J(k1,N); %修改节点电压虚部 end

%------修改节点电压----------- for k=3:N0 DET=abs(J(k,N));

if DET>=pr %电压偏差量是否满足要求 IT2=IT2+1; %不满足要求的节点数加1 end end

ICT2(a)=IT2; %不满足要求的节点数 ICT1=ICT1+1; %迭代次数 end

%用高斯消去法解\"w=-J*V\" disp('迭代次数:'); disp(ICT1);

disp('没有达到精度要求的个数:'); disp(ICT2); for k=1:n

V(k)=sqrt(e(k)^2+f(k)^2); %计算各节点电压的模值

sida(k)=atan(f(k)./e(k))*180./pi; %计算各节点电压的角度 E(k)=e(k)+f(k)*j; %将各节点电压用复数表示 end

%=============== 计算各输出量 =========================== disp('各节点的实际电压标幺值E为(节点号从小到大排列):'); disp(E); %显示各节点的实际电压标幺值E用复数表示 disp('-----------------------------------------------------'); disp('各节点的电压大小V为(节点号从小到大排列):'); disp(V); %显示各节点的电压大小V的模值 disp('-----------------------------------------------------');

disp('各节点的电压相角sida为(节点号从小到大排列):'); disp(sida); %显示各节点的电压相角 for p=1:n C(p)=0; for q=1:n

C(p)=C(p)+conj(Y(p,q))*conj(E(q)); %计算各节点的注入电流的共轭值 end

S(p)=E(p)*C(p); %计算各节点的功率 S = 电压 X 注入电流的共轭值 end

disp('各节点的功率S为(节点号从小到大排列):'); disp(S); %显示各节点的注入功率 disp('-----------------------------------------------------');

disp('各条支路的首端功率Si为(顺序同您输入B1时一致):'); for i=1:nl

p=B1(i,1);q=B1(i,2); if B1(i,6)==0

Si(p,q)=E(p)*(conj(E(p))*conj(B1(i,4)./2)+(conj(E(p)*B1(i,5))... -conj(E(q)))*conj(1./(B1(i,3)*B1(i,5)))); Siz(i)=Si(p,q); else

Si(p,q)=E(p)*(conj(E(p))*conj(B1(i,4)./2)+(conj(E(p)./B1(i,5))... -conj(E(q)))*conj(1./(B1(i,3)*B1(i,5)))); Siz(i)=Si(p,q); end disp(Si(p,q)); SSi(p,q)=Si(p,q);

ZF=['S(',num2str(p),',',num2str(q),')=',num2str(SSi(p,q))]; disp(ZF);

disp('-----------------------------------------------------'); end

disp('各条支路的末端功率Sj为(顺序同您输入B1时一致):'); for i=1:nl

p=B1(i,1);q=B1(i,2); if B1(i,6)==0

Sj(q,p)=E(q)*(conj(E(q))*conj(B1(i,4)./2)+(conj(E(q)./B1(i,5))... -conj(E(p)))*conj(1./(B1(i,3)*B1(i,5)))); Sjy(i)=Sj(q,p); else

Sj(q,p)=E(q)*(conj(E(q))*conj(B1(i,4)./2)+(conj(E(q)*B1(i,5))... -conj(E(p)))*conj(1./(B1(i,3)*B1(i,5)))); Sjy(i)=Sj(q,p); end disp(Sj(q,p)); SSj(q,p)=Sj(q,p);

ZF=['S(',num2str(q),',',num2str(p),')=',num2str(SSj(q,p))]; disp(ZF);

disp('-----------------------------------------------------'); end

disp('各条支路的功率损耗DS为(顺序同您输入B1时一致):'); for i=1:nl

p=B1(i,1);q=B1(i,2); DS(i)=Si(p,q)+Sj(q,p); disp(DS(i)); DDS(i)=DS(i);

ZF=['DS(',num2str(p),',',num2str(q),')=',num2str(DDS(i))]; disp(ZF);

disp('-----------------------------------------------------'); end figure(1); subplot(1,2,1); plot(V);

xlabel('节点号');ylabel('电压标幺值'); grid on; subplot(1,2,2); plot(sida);

xlabel('节点号');ylabel('电压角度'); grid on; figure(2); subplot(2,2,1); P=real(S);Q=imag(S); bar(P);

xlabel('节点号');ylabel('节点注入有功'); grid on; subplot(2,2,2); bar(Q);

xlabel('节点号');ylabel('节点注入无功'); grid on; subplot(2,2,3);

P1=real(Siz);Q1=imag(Siz); bar(P1);

xlabel('支路号');ylabel('支路首端注入有功'); grid on; subplot(2,2,4); bar(Q1);

xlabel('支路号');ylabel('支路首端注入无功'); grid on;

运行结果: 导纳矩阵 Y=

2.9056 -11.1778i 0 + 5.3191i -1.6606 + 3.1617i -1.2450 + 2.3710i 0 + 5.3191i 0 - 5.6600i 0 0 -1.6606 + 3.1617i 0 2.4904 - 4.7233i -0.8298 + 1.5809i -1.2450 + 2.3710i 0 -0.8298 + 1.5809i 2.0749 - 3.9304i

迭代次数: 4

没有达到精度要求的个数: 5 5 2 0

各节点的实际电压标幺值E为(节点号从小到大排列):

1.0000 0.8686 - 0.0940i 1.0500 + 0.0058i 0.9879 - 0.01i

----------------------------------------------------- 各节点的电压大小V为(节点号从小到大排列): 1.0000 0.8736 1.0500 0.9881

-----------------------------------------------------

各节点的电压相角sida为(节点号从小到大排列): 0 -6.1767 0.3150 -1.0945

各节点的功率S为(节点号从小到大排列):

0.4585 + 0.8816i -0.5000 - 0.3000i 0.2000 + 0.2171i -0.1500 - 0.1000i

-----------------------------------------------------

各条支路的首端功率Si为(顺序同您输入B1时一致): 0.5000 + 1.0400i

S(1,2)=0.5+1.04i

----------------------------------------------------- -0.1013 - 0.19i

S(1,3)=-0.10126-0.19i

----------------------------------------------------- 0.0598 - 0.0035i

S(1,4)=0.059786-0.0034655i

----------------------------------------------------- 0.0945 + 0.0678i

S(3,4)=0.094+0.067788i

-----------------------------------------------------

各条支路的末端功率Sj为(顺序同您输入B1时一致): -0.5000 - 0.3000i

S(2,1)=-0.5-0.3i

----------------------------------------------------- 0.1055 + 0.1493i

S(3,1)=0.106+0.14934i

----------------------------------------------------- -0.0592 - 0.0123i

S(4,1)=-0.059161-0.012341i

----------------------------------------------------- -0.0908 - 0.0877i

S(4,3)=-0.090839-0.087659i

-----------------------------------------------------

各条支路的功率损耗DS为(顺序同您输入B1时一致): -0.0000 + 0.7400i

DS(1,2)=-5.5511e-017+0.73997i

----------------------------------------------------- 0.0042 - 0.0056i

DS(1,3)=0.0042042-0.00555i

----------------------------------------------------- 0.0006 - 0.0158i

DS(1,4)=0.00062508-0.015806i

----------------------------------------------------- 0.0037 - 0.0199i

DS(3,4)=0.0037003-0.019872i

-----------------------------------------------------

1.05101-1-2电压标幺值0.95电压角度12节点号34-3-40.9-5-60.85-712节点号34

0.51节点注入有功0节点注入无功123节点号40.50-0.5-0.5123节点号1.0.60.41支路首端注入有功0.2支路首端注入无功123支路号40.500-0.2-0.5123支路号4

因篇幅问题不能全部显示,请点此查看更多更全内容

Copyright © 2019- yrrf.cn 版权所有 赣ICP备2024042794号-2

违法及侵权请联系:TEL:199 1889 7713 E-MAIL:2724546146@qq.com

本站由北京市万商天勤律师事务所王兴未律师提供法律服务