电力系统分析潮流计算最终完整版

1、电力系统分析潮流计算实验报告姓名：XXXXXX学号：XXXXXXXXX班级：XXXXXXXX一、实验目的掌握潮流计算计算机算法的方法，熟悉 MATLAB勺程序调试方法。二、实验准备根据课程内容，熟悉MATLAB件的使用方法，自行学习 MATLAB?序的基础 语法，并根据所学知识编写潮流计算牛顿拉夫逊法（或PQ分解法）的计算程 序，用相应的算例在 MATLA上进行计算、调试和验证。三、实验要求每人一组，在实验课时内， 调试和修改运行程序， 用算例计算输出潮流结果。四、程序流程启幼iK=Ofl用戎(446)(4斗咼计尊AP叫胡叭虫鬥|是 冷茸裁路功率及平衛节点功率利用式鬥49)(4 50)求雅可

2、比期半-蓉元素输出五、实验程序%本程序的功能是用牛拉法进行潮流计算%原理介绍详见鞠平著电气工程%默认数据为鞠平著电气工程例8.4 所示数据%B1是支路参数矩阵%第一列和第二列是节点编号。节点编号由小到大编写%对于含有变压器的支路，第一列为低压侧节点编号，第二列为高压侧节点编号%第三列为支路的串列阻抗参数，含变压器支路此值为变压器短路电抗%第四列为支路的对地导纳参数，含变压器支路此值不代入计算%第五烈为含变压器支路的变压器的变比，变压器非标准电压比%第六列为变压器是否是否含有变压器的参数，其中“1 ”为含有变压器， “0”为不含有变压器%B2为节点参数矩阵%第一列为节点注入发电功率参数%第二列为

3、节点负荷功率参数%第三列为节点电压参数%第四列%第五列%第六列为节点类型参数，“ 1 ”为平衡节点，“2”为PQ节点，“3”为PV节点参数%X为节点号和对地参数矩阵%第一列为节点编号%第二列为节点对地参数%默认算例%n=4;%n1=4;%isb=4;%pr=0.00001;%B1=1 2 0.1667i 0 0.8864 1;1 3 0.1302+0.2479i 0.0258i 1 0;1 4 0.1736+0.3306i 0.0344i 1 0;3 40.2603+0.4959i 0.0518i 1 0;%B2=0 0 1 0 0 2;0 -0.5-0.3i 1 0 0 2;0.2 0 1.

4、05 0 0 3;0 -0.15-0.1i 1.05 0 0 1;%X=1 0;2 0.05i;3 0;4 0;clear;clc;num=input( 是否采用默认数据？ (1-默认数据 ;2-手动输入 );if num=1n=4;n1=4;isb=4;pr=0.00001;B1=1 2 0.1667i 0 0.8864 1;1 3 0.1302+0.2479i 0.0258i 1 0;1 4 0.1736+0.3306i 0.0344i 1 0;3 40.2603+0.4959i 0.0518i 1 0;B2=0 0 1 0 0 2;0 -0.5-0.3i 1 0 0 2;0.2 0 1.

5、05 0 0 3;0 -0.15-0.1i 1.05 0 0 1;X=1 0;2 0.05i;3 0;4 0;elsen=input( 请输入节点数 :n=); n1=input( 请输入支路数 :n1=); isb=input( 请输入平衡节点号 :isb=); pr=input( 请输入误差精度 :pr=); B1=input( 请输入支路参数 :B1=); B2=input( 请输入节点参数 :B2=); X=input( 节点号和对地参数 :X=); endTimes=1; %迭代次数%创建节点导纳矩阵 Y=zeros(n);for i=1:n1if B1(i,6)=0 %不含变压器的

6、支路 p=B1(i,1);q=B1(i,2); Y(p,q)=Y(p,q)-1/B1(i,3); Y(q,p)=Y(p,q);Y(p,p)=Y(p,p)+1/B1(i,3)+0.5*B1(i,4); Y(q,q)=Y(q,q)+1/B1(i,3)+0.5*B1(i,4);else %含有变压器的支路 p=B1(i,1); q=B1(i,2); Y(p,q)=Y(p,q)-B1(i,5)/B1(i,3); Y(q,p)=Y(p,q);Y(p,p)=Y(p,p)+B1(i,5)/B1(i,3)+(1-B1(i,5)/B1(i,3); Y(q,q)=Y(q,q)+B1(i,5)/B1(i,3)+(B

7、1(i,5)*(B1(i,5)-1)/B1(i,3); endendfor i=1:n1Y(i,i)=Y(i,i)+X(i,2); %计及补偿电容电纳 enddisp(导纳矩阵为：)；disp(Y);%显示导纳矩阵%初始化 OrgS、 DetaS OrgS=zeros(2*n-2,1)； DetaS=zeros(2*n-2,1)；%创建OrgS,用于存储初始功率参数 h=0；j=0；for i=1:n %对 PQ 节点的处理if i=isb&B2(i,6)=2%不是平衡点&是PQ点h=h+1;for j=1:n%公式 8-74%Pi=ei*(Gij*ej-Bij*fj)+fi*(Gij*fj+

8、Bij*ej)%Qi=fi*(Gij*ej-Bij*fj)-ei*(Gij*fj+Bij*ej)OrgS(2*h-1,1)=OrgS(2*h-1,1)+real(B2(i,3)*(real(Y(i,j)*real(B2(j,3)-imag(Y(i,j)*imag(B2(j,3)+imag(B2(i,3)*(real (Y(i,j)*imag(B2(j,3)+imag(Y(i,j)*real(B2(j,3);OrgS(2*h,1) =OrgS(2*h,1)+imag(B2(i,3)*(real(Y(i,j)*real(B2(j,3)-imag(Y(i,j)*imag(B2(j,3)-real(B

9、2(i,3)*(real(Y(i,j)*imag(B2(j ,3)+imag(Y(i,j)*real(B2(j,3);endendendfor i=1:n %对PV节点的处理，注意这时不可再将h初始化为0if i=isb&B2(i,6)=3 %不是平衡点 &是 PV 点h=h+1;for j=1:n%公式 8-75-a%Pi=ei*(Gij*ej-Bij*fj)+fi*(Gij*fj+Bij*ej)%Qi=fi*(Gij*ej-Bij*fj)-ei*(Gij*fj+Bij*ej)OrgS(2*h-1,1)=OrgS(2*h-1,1)+real(B2(i,3)*(real(Y(i,j)*real

10、(B2(j,3)-imag(Y(i,j)*imag(B2(j,3)+imag(B2(i,3)*(real (Y(i,j)*imag(B2(j,3)+imag(Y(i,j)*real(B2(j,3);OrgS(2*h,1)=OrgS(2*h,1)+imag(B2(i,3)*(real(Y(i,j)*real(B2(j,3)-imag(Y(i,j)*imag(B2(j,3)-real(B2(i,3)*(real(Y(i,j )*imag(B2(j,3)+imag(Y(i,j)*real(B2(j,3);endendend%创建PVU用于存储PV节点的初始电压PVU=zeros(n-h-1,1);t

11、=0;for i=1:nif B2(i,6)=3t=t+1;PVU(t,1)=B2(i,3);endend%创建DetaS,用于存储有功功率、无功功率和电压幅值的不平衡量h=0;for i=1:n %对 PQ 节点的处理if i=isb&B2(i,6)=2h=h+1;DetaS(2*h-1,1)=real(B2(i,2)-OrgS(2*h-1,1); %delPiDetaS(2*h,1)=imag(B2(i,2)-OrgS(2*h,1); %delQiendendt=0;for i=1:n %对PV节点的处理，注意这时不可再将h初始化为0if i=isb&B2(i,6)=3h=h+1;t=t+

12、1;DetaS(2*h-1,1)=real(B2(i,1)-OrgS(2*h-1,1); %delPiDetaS(2*h,1)=real(PVU(t,1)A2+imag(PVU(t,1)A2-real(B2(i,3)A2-imag(B2(i,3)A2; %delUiendend% DetaS%创建I,用于存储节点电流参数i=zeros(n-1,1);h=0;for i=1:nif i=isbh=h+1;I(h,1)=(OrgS(2*h-1,1)-OrgS(2*h,1)*sqrt(-1)/conj(B2(i,3);%conj 求共轭endend%创建Jacbi(雅可比矩阵)Jacbi=zeros

13、(2*n-2);h=0;k=0;for i=1:n %对PQ节点的处理if B2(i,6)=2h=h+1;for j=1:nif j=isbk=k+1;if i=j %对角元素的处理Jacbi(2*h-1,2*k-1)=-imag(Y(i,j)*real(B2(i,3)+real(Y(i,j)*imag(B2(i,3)+imag(I(h,1);Jacbi(2*h-1,2*k)=real(Y(i,j)*real(B2(i,3)+imag(Y(i,j)*imag(B2(i,3)+real(I(h,1);Jacbi(2*h,2*k-1)=-Jacbi(2*h-1,2*k)+2*real(I(h,1)

14、;Jacbi(2*h,2*k)=Jacbi(2*h-1,2*k-1)-2*imag(I(h,1);else %非对角元素的处理Jacbi(2*h-1,2*k-1)=-imag(Y(i,j)*real(B2(i,3)+real(Y(i,j)*imag(B2(i,3);Jacbi(2*h-1,2*k)=real(Y(i,j)*real(B2(i,3)+imag(Y(i,j)*imag(B2(i,3);Jacbi(2*h,2*k-1)=-Jacbi(2*h-1,2*k);Jacbi(2*h,2*k)=Jacbi(2*h-1,2*k-1);endif k=(n-1) %将用于内循环的指针置于初始值，以

15、确保雅可比矩阵换行k=0;endendendendendk=0;for i=1:n %对 PV 节点的处理if B2(i,6)=3h=h+1;for j=1:nif j=isbk=k+1;if i=j %对角元素的处理Jacbi(2*h-1,2*k-1)=-imag(Y(i,j)*real(B2(i,3)+real(Y(i,j)*imag(B2(i,3)+imag(I(h,1);Jacbi(2*h-1,2*k)=real(Y(i,j)*real(B2(i,3)+imag(Y(i,j)*imag(B2(i,3)+real(I(h,1);Jacbi(2*h,2*k-1)=2*imag(B2(i,3

16、);Jacbi(2*h,2*k)=2*real(B2(i,3);else %非对角元素的处理Jacbi(2*h-1,2*k-1)=-imag(Y(i,j)*real(B2(i,3)+real(Y(i,j)*imag(B2(i,3);Jacbi(2*h-1,2*k)=real(Y(i,j)*real(B2(i,3)+imag(Y(i,j)*imag(B2(i,3);Jacbi(2*h,2*k-1)=0;Jacbi(2*h,2*k)=0;endif k=(n-1) %将用于内循环的指针置于初始值，以确保雅可比矩阵换行k=0;endendendendenddisp(初始雅可比矩阵为：)；disp(J

17、acbi);%求解修正方程，获取节点电压的不平衡量DetaU=zeros(2*n-2,1)；DetaU=inv(Jacbi)*DetaS; %inv 矩阵求逆% DetaU%修正节点电压 j=0;for i=1:n %对 PQ 节点处理if B2(i,6)=2j=j+1;B2(i,3)=B2(i,3)+DetaU(2*j,1)+DetaU(2*j-1,1)*sqrt(-1);endendfor i=1:n %对 PV 节点的处理if B2(i,6)=3j=j+1;B2(i,3)=B2(i,3)+DetaU(2*j,1)+DetaU(2*j-1,1)*sqrt(-1);endend% B2%开始

18、循环 *while abs(max(DetaU)prOrgS=zeros(2*n-2,1);h=0;j=0;for i=1:nif i=isb&B2(i,6)=2h=h+1;for j=1:nOrgS(2*h-1,1)=OrgS(2*h-1,1)+real(B2(i,3)*(real(Y(i,j)*real(B2(j,3)-imag(Y(i,j)*imag(B2(j,3)+imag(B2(i,3)*(real (Y(i,j)*imag(B2(j,3)+imag(Y(i,j)*real(B2(j,3);OrgS(2*h,1)=OrgS(2*h,1)+imag(B2(i,3)*(real(Y(i,

19、j)*real(B2(j,3)-imag(Y(i,j)*imag(B2(j,3)-real(B2(i,3)*(real(Y(i,j )*imag(B2(j,3)+imag(Y(i,j)*real(B2(j,3);endendendfor i=1:nif i=isb&B2(i,6)=3h=h+1;for j=1:nOrgS(2*h-1,1)=OrgS(2*h-1,1)+real(B2(i,3)*(real(Y(i,j)*real(B2(j,3)-imag(Y(i,j)*imag(B2(j,3)+imag(B2(i,3)*(real (Y(i,j)*imag(B2(j,3)+imag(Y(i,j)

20、*real(B2(j,3);OrgS(2*h,1)=OrgS(2*h,1)+imag(B2(i,3)*(real(Y(i,j)*real(B2(j,3)-imag(Y(i,j)*imag(B2(j,3)-real(B2(i,3)*(real(Y(i,j )*imag(B2(j,3)+imag(Y(i,j)*real(B2(j,3);endendend % OrgS%创建 DetaSh=0;for i=1:nif i=isb&B2(i,6)=2h=h+1;DetaS(2*h-1,1)=real(B2(i,2)-OrgS(2*h-1,1);DetaS(2*h,1)=imag(B2(i,2)-Org

21、S(2*h,1);endendt=0;for i=1:nif i=isb&B2(i,6)=3h=h+1;t=t+1;% DetaS(2*h-1,1)=real(B2(i,2)-OrgS(2*h-1,1); DetaS(2*h-1,1)=real(B2(i,1)-OrgS(2*h-1,1);DetaS(2*h,1)=real(PVU(t,1)A2+imag(PVU(t,1)A2-real(B2(i,3)A2-imag(B2(i,3)A2; endend% DetaS%创建 Ii=zeros(n-1,1);h=0;for i=1:nif i=isbh=h+1;I(h,1)=(OrgS(2*h-1,

22、1)-OrgS(2*h,1)*sqrt(-1)/conj(B2(i,3);endend% I%创建 JacbiJacbi=zeros(2*n-2);h=0;k=0;for i=1:nif B2(i,6)=2h=h+1;for j=1:nif j=isbk=k+1;if i=jJacbi(2*h-1,2*k-1)=-imag(Y(i,j)*real(B2(i,3)+real(Y(i,j)*imag(B2(i,3)+imag(I(h,1);Jacbi(2*h-1,2*k)=real(Y(i,j)*real(B2(i,3)+imag(Y(i,j)*imag(B2(i,3)+real(I(h,1);J

23、acbi(2*h,2*k-1)=-Jacbi(2*h-1,2*k)+2*real(I(h,1);Jacbi(2*h,2*k)=Jacbi(2*h-1,2*k-1)-2*imag(I(h,1);elseJacbi(2*h-1,2*k-1)=-imag(Y(i,j)*real(B2(i,3)+real(Y(i,j)*imag(B2(i,3);Jacbi(2*h-1,2*k)=real(Y(i,j)*real(B2(i,3)+imag(Y(i,j)*imag(B2(i,3);Jacbi(2*h,2*k-1)=-Jacbi(2*h-1,2*k);Jacbi(2*h,2*k)=Jacbi(2*h-1,2

24、*k-1);endif k=(n-1)k=0;endendendendendk=0;for i=1:nif B2(i,6)=3h=h+1;for j=1:nif j=isbk=k+1;if i=jJacbi(2*h-1,2*k-1)=-imag(Y(i,j)*real(B2(i,3)+real(Y(i,j)*imag(B2(i,3)+imag(I(h,1);Jacbi(2*h-1,2*k)=real(Y(i,j)*real(B2(i,3)+imag(Y(i,j)*imag(B2(i,3)+real(I(h,1);Jacbi(2*h,2*k-1)=2*imag(B2(i,3);Jacbi(2*h

25、,2*k)=2*real(B2(i,3);elseJacbi(2*h-1,2*k-1)=-imag(Y(i,j)*real(B2(i,3)+real(Y(i,j)*imag(B2(i,3);Jacbi(2*h-1,2*k)=real(Y(i,j)*real(B2(i,3)+imag(Y(i,j)*imag(B2(i,3);Jacbi(2*h,2*k-1)=0;Jacbi(2*h,2*k)=0;endif k=(n-1)k=0;endendendendend % JacbiDetaU=zeros(2*n-2,1);DetaU=inv(Jacbi)*DetaS;% DetaU%修正节点电压j=0;

26、for i=1:nif B2(i,6)=2j=j+1;B2(i,3)=B2(i,3)+DetaU(2*j,1)+DetaU(2*j-1,1)*sqrt(-1);endendfor i=1:nif B2(i,6)=3j=j+1;B2(i,3)=B2(i,3)+DetaU(2*j,1)+DetaU(2*j-1,1)*sqrt(-1);endend% B2Times=Times+1; %迭代次数加 1enddisp(迭代次数为：);disp(Times);disp(收敛时电压修正量为：)；disp(DetaU);for k=1:nE(k)=B2(k,3)；e(k)=real(E(k)；f(k)=im

27、ag(E(k)；V(k)=sqrt(e(k)A2+f(kF2);sida(k)=atan(f(k)./e(k)*180./pi；end%=计算各输出量 =disp(各节点的实际电压标幺值E为(节点号从小到大排列)：)；disp(E);%显示各节点的实际电压标幺值E用复数表示disp()disp(各节点的电压大小 V为(节点号从小到大排列)：);disp(V);%显示各节点的电压大小V的模值disp();disp(各节点的电压相角 sida为(节点号从小到大排列)：);disp(sida);%显示各节点的电压相for p=1:nC(p)=0;for q=1:nC(p)=C(p)+conj(Y(p

28、,q)*conj(E(q); %计算各节点的注入电流的共轭值 endS(p)=E(p)*C(p);%计算各节点的功率S = 电压 X 注入电流的共轭值enddisp(各节点的功率S为(节点号从小到大排列)：)；disp(S);%显示各节点的注入功率Sline=zeros(n1,5)；disp()；disp(各条支路的首端功率Si为(顺序同您输入B1时一致)：);for i=1:n1p=B1(i,1);q=B1(i,2);Sline(i,1)=B1(i,1);Sline(i,2)=B1(i,2);if B1(i,6)=0 Si(p,q)=E(p)*(conj(E(p)*conj(B1(i,4).

29、/2)+(conj(E(p)*B1(i,5)-conj(E(q)*conj(1./(B1(i,3)*B1(i,5); Siz(i)=Si(p,q);elseSi(p,q)=E(p)*(conj(E(p)*(1-B1(i,5)/B1(i,3)+(conj(E(p)-conj(E(q)*(B1(i,5)/B1(i,3); Siz(i)=Si(p,q);endSSi(p,q)=Si(p,q);Sline(i,3)=Siz(i);ZF=S(,num2str(p),num2str(q),)=,num2str(SSi(p,q);disp(ZF);enddisp();disp(各条支路的末端功率Sj为(顺序

30、同您输入B1时一致)：);for i=1:n1p=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);elseSj(q,p)=E(q)*(conj(E(q)*(B1(i,5)*(B1(i,5)-1)/B1(i,3)+(conj(E(q)-conj(E(p)*(B1(i,5)/B1(i,3); Sjy(i)=Sj(q,p);endSSj(q,p)=Sj(q,p);S

31、line(i,4)=Sjy(i);ZF=S(,num2str(q),num2str(p),)=,num2str(SSj(q,p);disp(ZF);enddisp();disp(各条支路的功率损耗DS为(顺序同您输入B1时一致)：)；for i=1:n1p=B1(i,1)；q=B1(i,2)；DS(i)=Si(p,q)+Sj(q,p)；DDS(i)=DS(i)；Sline(i,5)=DS(i)；ZF=DS(,num2str(p),num2str(q),)=,num2str(DDS(i)；disp(ZF)；enddisp()；线路损disp(各支路首端编号末端编号首端功率末端功率耗)；disp(

32、Sline)；六、运行结果及其分析是否采用默认数据？导纳矩阵为：2.9056 -11.5015i0.0000 + 5.3173i-1.6606 + 3.1617i-1.2450 + 2.3710i初始雅可比矩阵为：（1- 默认数据 ；2- 手动输入 ）10.0000 + 5.3173i0.0000 - 4.6633i0.0000 + 0.0000i0.0000 + 0.0000i-1.6606 + 3.1617i0.0000 + 0.0000i2.4904 - 4.7039i-1.2450 + 2.3710i0.0000 + 0.0000i-0.8298 + 1.5809i-0.8298 +

33、1.5809i 2.0749 - 3.9089i11.1267 2.7603 -5.31730 -3.1617 -1.6606-3.0509 11.87620 -5.3173 1.6606 -3.1617-5.3173 0 5.3173 0000 -5.3173 0 4.009200-3.3198 -1.7436 0 0 4.8217 2.69800 0 0 0 0 2.1000 迭代次数为：4收敛时电压修正量为：1.0e-05 *0.0349-0.2445-0.0101-0.5713-0.0931-0.0073各节点的实际电压标幺值 E为(节点号从小到大排列)：0.9673 - 0.0655i 1.0252 - 0.1666i 1.0495 - 0.0337i 1.0500 + 0.0000i各节点的电压大小 V为（节点号从小到大排列）:0.9695 1.0387 1.0500 1.0500各节点的电压相角 sida 为（节点号从小到大排列 ）：-3.8734 -9.2315 -1.8419 0各节点的功率S为（节点号从小到大排列）：-0.