潮流计算上机报告 华电
潮流上机课程设计-华电
课程设计报告( 2011—2012年度第一学期)名称：电力系统潮流上机院系：电气与电子工程学院班级：学号：学生：指导教师：设计周数：两周成绩：日期： 2011年12月19日一、课程设计的目的与要求培养学生的电力系统潮流计算机编程能力，掌握计算机潮流计算的相关知识二、设计正文（详细容见附录）1．手算2．计算机计算3．思考题三、课程设计总结或结论潮流计算是研究电力系统稳态运行的一种基本计算，最初求解电力系统潮流时大多使用手算，但随着电力系统结构的日趋复杂，计算量也越来越大。
复杂电力系统潮流计算中，由于节点数量巨大，所形成的修正方程已经无法通过手算方式解决，尤其是需要迭代次数较多时，手算所需要的时间太长，计算机潮流计算无疑为解决这一问题提供了极大的便利。
计算机潮流计算可以迅速解决复杂网络的潮流计算问题，这是由于无论系统的复杂程度如何，其节点与支路的类型是固定的，所以只需要输入节点与支路的数据，就可以解决任何一个复杂网络的潮流计算问题。
即只需要一次编程，就可以基本上解决所有复杂网络的计算。
需要注意的是，在使用牛顿—拉弗逊发计算潮流时，对于初值要选择比较接近它们的精确解，否则迭代过程可能不收敛。
潮流计算Ｃ语言程序编程过程中需要注意的是，Ｃ语言无法实现复数运算，需要将得到的值的实部与虚部分开储存并计算。
这个过程复杂并且容易出错，编写程序是需要注意。
另外需要注意的一点是：Ｃ语言数组的编号是从零开始的，在程序编写过程中应注意下标的对应。
通过这一次的电力系统潮流计算编程，我不仅对Ｃ语言的编程有了更深刻的理解，也对《电力系统分析》这门课程进行了查漏补缺和巩固，对电力系统的运行也有了更加深入的了解，受益匪浅。
四、参考文献1.《电力系统计算：电子数字计算机的应用》，交通大学等合编。
：水利电力；2.《现代电力系统分析》，王锡凡主编，科学；3.《电力系统稳态分析》，珩，中国电力，2007年，第三版；附录（包括：1.手算过程及结果；2.计算机计算流程图、表格、数据；3.思考题答案）附录一、手算过程及其结果附录二、计算机计算流程图、表格、数据编写潮流计算程序简单系统如下图所示，支路数据如下： 支路14，27，39为变压器支路，参数为100.1,058.0114==K X ，050.1,063.0227==K X 100.1,059.0339==K X其余支路为线路支路，参数为075.02/,072.0019.07878=+=B j Z , 105.02/,101.0012.08989=+=B j Z 153.02/,161.0032.05757=+=B j Z 179.02/,170.0039.06969=+=B j Z 088.02/,085.0010.04545=+=B j Z 079.02/,092.0017.04646=+=B j Z节点数据如下：o U 004.11∠=•025.1,63.122==U P ，025.1,85.033==U P5.025.15j S --= ,3.09.06j S --= ,35.00.18j S --=节点数据线路数据节点导纳矩阵Y B迭代次数：Count_Num=3各节点电压：平衡节点功率：节点编号1: 0.72887-j0.10579 PV节点无功功率：节点编号2: j0.30332PV节点无功功率：节点编号3: j0.17199线路功率及损耗：线路总损耗：ΔS= 0.05887-j0.78048附录三、思考题1、潮流计算的方法有哪些？各有何特点？答：潮流计算分为简单电力网络的手算和复杂电力网络的机算两大类，其中机算又有高斯-赛德尔迭代法，牛顿-拉夫逊迭代法和P-Q分解法。
潮流计算实验
电力系统分析实验报告实验一：潮流计算的计算机算法>> clear;n=10;nl=10;isb=1;pr=0.00001;B1=[120.03512+0.08306i0.13455i10;230.0068+0.18375i0 1.023811;140.05620+0.13289i0.05382i10;450.00811+0.24549i0 1.023811;160.05620+0.13289i0.05382i10;460.04215+0.09967i0.04037i10;670.0068+0.18375i0 1.023811;680.02810+0.06645i0.10764i10;8100.00811+0.24549i011;890.03512+0.08306i0.13455i10] B2=[00 1.1 1.101;001002;00.343+0.21256i1002;001002;00.204+0.12638i1002;001002;00.306+0.18962i1002;001002;0.50 1.1 1.103;00.343+0.21256i1002]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:nlif B1(i,6)==0p=B1(i,1);q=B1(i,2);elsep=B1(i,2);q=B1(i,1);endY(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;Y(p,p)=Y(p,p)+1./B1(i,3)+B1(i,4)./2;enddisp('导纳矩阵Y=');disp(Y)%----------------------------------------------------------G=real(Y);B=imag(Y);for i=1:ne(i)=real(B2(i,3));f(i)=imag(B2(i,3));V(i)=B2(i,4);endfor i=1:nS(i)=B2(i,1)-B2(i,2);B(i,i)=B(i,i)+B2(i,5);endP=real(S);Q=imag(S);ICT1=0;IT2=1;N0=2*n;N=N0+1;a=0;while IT2~=0IT2=0;a=a+1;for i=1:nif i~=isbC(i)=0;D(i)=0;for j1=1:nC(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) endP1=C(i)*e(i)+f(i)*D(i);Q1=C(i)*f(i)-e(i)*D(i);V2=e(i)^2+f(i)^2;if B2(i,6)~=3DP=P(i)-P1;DQ=Q(i)-Q1;for j1=1:nif j1~=isb&j1~=iX1=-G(i,j1)*e(i)-B(i,j1)*f(i);X2=B(i,j1)*e(i)-G(i,j1)*f(i);X3=X2;X4=-X1;p=2*i-1;q=2*j1-1;J(p,q)=X3;J(p,N)=DQ;m=p+1;J(m,q)=X1;J(m,N)=DP;q=q+1;J(p,q)=X4;J(m,q)=X2;elseif j1==i&j1~=isbX1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i);X2=-D(i)+B(i,i)*e(i)-G(i,i)*f(i);X3=D(i)+B(i,i)*e(i)-G(i,i)*f(i);X4=-C(i)+G(i,i)*e(i)+B(i,i)*f(i);p=2*i-1;q=2*j1-1;J(p,q)=X3;J(p,N)=DQ;m=p+1;J(m,q)=X1;q=q+1;J(p,q)=X4;J(m,N)=DP;J(m,q)=X2;endendelseDP=P(i)-P1;DV=V(i)^2-V2;for j1=1:nif j1~=isb&j1~=iX1=-G(i,j1)*e(i)-B(i,j1)*f(i);X2=B(i,j1)*e(i)-G(i,j1)*f(i);X5=0;X6=0;p=2*i-1;q=2*j1-1;J(p,q)=X5;J(p,N)=DV;m=p+1;J(m,q)=X1;J(m,N)=DP;q=q+1;J(p,q)=X6;J(m,q)=X2;elseif j1==i&j1~=isbX1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i);X2=-D(i)+B(i,i)*e(i)-G(i,i)*f(i);X5=-2*e(i);X6=-2*f(i);p=2*i-1;q=2*j1-1;J(p,q)=X5;J(p,N)=DV;m=p+1;J(m,q)=X1;J(m,N)=DP;q=q+1;J(p,q)=X6;J(m,q)=X2;endendendendendfor k=3:N0k1=k+1;N1=N;for k2=k1:N1J(k,k2)=J(k,k2)./J(k,k);endJ(k,k)=1;if k~=3k4=k-1;for k3=3:k4for k2=k1:N1J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);endJ(k3,k)=0;endif k==N0break;endfor k3=k1:N0for k2=k1:N1J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);endJ(k3,k)=0;endelsefor k3=k1:N0for k2=k1:N1J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);endJ(k3,k)=0;endendendfor k=3:2:N0-1L=(k+1)./2;e(L)=e(L)-J(k,N);k1=k+1;f(L)=f(L)-J(k1,N);endfor k=3:N0DET=abs(J(k,N));if DET>=prIT2=IT2+1;endendICT2(a)=IT2;ICT1=ICT1+1;enddisp('迭代次数：');disp(ICT1);disp('没有达到精度要求的个数：');disp(ICT2);for k=1:nV(k)=sqrt(e(k)^2+f(k)^2);sida(k)=atan(f(k)./e(k))*180./pi;E(k)=e(k)+f(k)*j;enddisp('各节点的实际电压标幺值E为(节点号从小到大排列)：');disp(E);disp('-----------------------------------------------------');disp('各节点的电压大小V为(节点号从小到大排列)：');disp(V);disp('-----------------------------------------------------');disp('各节点的电压相角sida为(节点号从小到大排列)：');disp(sida);for p=1:nC(p)=0;for q=1:nC(p)=C(p)+conj(Y(p,q))*conj(E(q));endS(p)=E(p)*C(p);enddisp('各节点的功率S为(节点号从小到大排列)：');disp(S);disp('-----------------------------------------------------');disp('各条支路的首端功率Si为(顺序同输入B1时一致)：');for i=1:nlp=B1(i,1);q=B1(i,2);if B1(i,6)==0Si(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);elseSi(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);enddisp(Si(p,q));SSi(p,q)=Si(p,q);ZF=['S(',num2str(p),',',num2str(q),')=',num2str(SSi(p,q))];disp(ZF);disp('-----------------------------------------------------');enddisp('各条支路的末端功率Sj为(顺序同输入B1时一致)：');for i=1:nlp=B1(i,1);q=B1(i,2);if B1(i,6)==0Sj(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))*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);enddisp(Sj(q,p));SSj(q,p)=Sj(q,p);ZF=['S(',num2str(q),',',num2str(p),')=',num2str(SSj(q,p))];disp(ZF);disp('-----------------------------------------------------');enddisp('各条支路的功率损耗DS为(顺序同输入B1时一致)：');for i=1:nlp=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('-----------------------------------------------------');endfigure(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;B1 =1.00002.0000 0.0351 + 0.0831i 0 + 0.1346i 1.0000 02.00003.0000 0.0068 + 0.1838i 0 1.0238 1.00001.0000 4.0000 0.0562 + 0.1329i 0 + 0.0538i 1.0000 04.00005.0000 0.0081 + 0.2455i 0 1.0238 1.00001.0000 6.0000 0.0562 + 0.1329i 0 + 0.0538i 1.0000 04.0000 6.0000 0.0422 + 0.0997i 0 + 0.0404i 1.0000 06.00007.0000 0.0068 + 0.1838i 0 1.0238 1.00006.0000 8.0000 0.0281 + 0.0664i 0 + 0.1076i 1.0000 08.0000 10.0000 0.0081 + 0.2455i 0 1.0000 1.00008.0000 9.0000 0.0351 + 0.0831i 0 + 0.1346i 1.0000 0 B2 =0 0 1.1000 1.1000 0 1.00000 0 1.0000 0 0 2.00000 0.3430 + 0.2126i 1.0000 0 0 2.00000 0 1.0000 0 0 2.00000 0.2040 + 0.1264i 1.0000 0 0 2.00000 0 1.0000 0 0 2.00000 0.3060 + 0.1896i 1.0000 0 0 2.00000 0 1.0000 0 0 2.0000 0.5000 0 1.1000 1.1000 0 3.00000 0.3430 + 0.2126i 1.0000 0 0 2.0000导纳矩阵Y=Columns 1 through 69.7177 -22.8591i -4.3185 +10.2135i 0 -2.6996 + 6.3834i 0 -2.6996 + 6.3834i-4.3185 +10.2135i 4.5104 -15.3311i -0.1964 + 5.3083i 0 0 00 -0.1964 + 5.3083i 0.2011 - 5.4347i 0 0 0-2.6996 + 6.3834i 0 0 6.4271 -18.7292i -0.1313 + 3.9744i -3.5993 + 8.5110i0 0 0 -0.1313 + 3.9744i 0.1344 - 4.0690i 0-2.6996 + 6.3834i 0 0 -3.5993 + 8.5110i 0 11.8891 -32.7444i0 0 0 0 0 -0.1964 + 5.3083i0 0 0 0 0 -5.3984 +12.7660i0 0 0 0 0 00 0 0 0 0 0Columns 7 through 100 0 0 00 0 0 00 0 0 00 0 0 00 0 0 0-0.1964 + 5.3083i -5.3984 +12.7660i 0 00.2011 - 5.4347i 0 0 00 9.8514 -26.9275i -4.3185 +10.2135i -0.1344 + 4.0690i0 -4.3185 +10.2135i 4.3185 -10.1462i 00 -0.1344 + 4.0690i 0 0.1344 - 4.0690i迭代次数：4没有达到精度要求的个数：17 18 17 0各节点的实际电压标幺值E为(节点号从小到大排列)：Columns 1 through 61.1000 1.0757 - 0.0207i 1.0050 - 0.0780i 1.0772 - 0.0175i 1.0171 - 0.0631i 1.0762 - 0.0152iColumns 7 through 101.0112 - 0.0666i 1.0778 - 0.0051i 1.0996 + 0.0304i 1.0177 - 0.0814i-----------------------------------------------------各节点的电压大小V为(节点号从小到大排列)：1.1000 1.0759 1.0080 1.0773 1.0191 1.0763 1.0134 1.0778 1.1000 1.0209-----------------------------------------------------各节点的电压相角sida为(节点号从小到大排列)：0 -1.1046 -4.4373 -0.9283 -3.5503 -0.8106 -3.7665 -0.27181.5822 -4.5707各节点的功率S为(节点号从小到大排列)：Columns 1 through 60.7165 + 0.2587i 0.0000 + 0.0000i -0.3430 - 0.2126i -0.0000 + 0.0000i -0.2040 - 0.1264i -0.0000 + 0.0000iColumns 7 through 10-0.3060 - 0.1896i -0.0000 + 0.0000i 0.5000 + 0.0089i -0.3430 - 0.2126i-----------------------------------------------------各条支路的首端功率Si为(顺序同输入B1时一致)：0.3485 + 0.0932iS(1,2)=0.3485+0.093157i-----------------------------------------------------0.3441 + 0.2420iS(2,3)=0.34409+0.24201i-----------------------------------------------------0.1904 + 0.0760iS(1,4)=0.19038+0.07599i-----------------------------------------------------0.2044 + 0.1400iS(4,5)=0.20445+0.13999i-----------------------------------------------------0.1777 + 0.0895iS(1,6)=0.17767+0.089525i------------------------------------------------------0.0163 - 0.0055iS(4,6)=-0.016305-0.0054856i-----------------------------------------------------0.3069 + 0.2128iS(6,7)=0.30686+0.21281i------------------------------------------------------0.1477 - 0.0234iS(6,8)=-0.14767-0.02338i-----------------------------------------------------0.3443 + 0.2509iS(8,10)=0.34427+0.25091i------------------------------------------------------0.4925 - 0.1508iS(8,9)=-0.49251-0.15077i-----------------------------------------------------各条支路的末端功率Sj为(顺序同输入B1时一致)：-0.3441 - 0.2420iS(2,1)=-0.34409-0.24201i------------------------------------------------------0.3430 - 0.2126iS(3,2)=-0.343-0.21256i------------------------------------------------------0.1881 - 0.1345iS(4,1)=-0.18815-0.13451i------------------------------------------------------0.2040 - 0.1264iS(5,4)=-0.204-0.12638i------------------------------------------------------0.1755 - 0.1482iS(6,1)=-0.17551-0.14815i-----------------------------------------------------0.0163 - 0.0413iS(6,4)=0.016326-0.041272i------------------------------------------------------0.3060 - 0.1896iS(7,6)=-0.306-0.18962i-----------------------------------------------------0.1482 - 0.1001iS(8,6)=0.14824-0.10014i------------------------------------------------------0.3430 - 0.2126iS(10,8)=-0.343-0.21256i-----------------------------------------------------0.5000 + 0.0089iS(9,8)=0.5+0.0089402i-----------------------------------------------------各条支路的功率损耗DS为(顺序同输入B1时一致)：0.0044 - 0.1488iDS(1,2)=0.0044095-0.14885i-----------------------------------------------------0.0011 + 0.0294iDS(2,3)=0.0010897+0.029445i-----------------------------------------------------0.0022 - 0.0585iDS(1,4)=0.0022306-0.058518i-----------------------------------------------------0.0004 + 0.0136iDS(4,5)=0.00044972+0.013613i-----------------------------------------------------0.0022 - 0.0586iDS(1,6)=0.0021584-0.058629i-----------------------------------------------------0.0000 - 0.0468iDS(4,6)=2.1344e-005-0.046758i-----------------------------------------------------0.0009 + 0.0232iDS(6,7)=0.00085804+0.023186i-----------------------------------------------------0.0006 - 0.1235iDS(6,8)=0.00056584-0.12352i-----------------------------------------------------0.0013 + 0.0384iDS(8,10)=0.001267+0.038353i-----------------------------------------------------0.0075 - 0.1418iDS(8,9)=0.0074931-0.14183i----------------------------------------------------- >>。
电力系统潮流计算实验报告
电力系统潮流上机计算实验报告１1．手算过程已知：节点1：PQ 节点，节点， s(1)= s(1)= -0.5000-j0.3500 节点2：PV 节点，节点， p(2)=0.4000 v(2)=1.0500 p(2)=0.4000 v(2)=1.0500 节点3：平衡节点，：平衡节点，U(3)=1.0000U(3)=1.0000U(3)=1.0000∠∠0.0000 网络的连接图：0.0500+j0.2000 1 0.0500+j0.2000231）计算节点导纳矩阵由2000.00500.012j Z +=Þ71.418.112j y -=; 2000.00500.013j Z +=Þ71.418.113j y -=;\导纳矩阵中的各元素：42.936.271.418.171.418.1131211j j j y y Y -=-+-=+=;71.418.11212jy Y +-=-=; 71.418.11313j y Y +-=-=; =21Y71.418.11212j y Y +-=-=; 71.418.12122j y Y -==; 002323j y Y +=-=;=31Y 71.418.11313j y Y +-=-=; =32Y 002323j y Y +=-=; 71.418.13133j y Y -==;\形成导纳矩阵B Y ：úúúûùêêêëé-++-+-+-+-+--=71.418.10071.418.10071.418.171.418.171.418.171.418.142.936.2j j j j j j j j j Y B 2）计算各PQ PQ、、PV 节点功率的不平衡量，及PV 节点电压的不平衡量：取:000.0000.1)0(1)0(1)0(1j jf e U +=+=000.0000.1)0(2)0(2)0(2j jf e U +=+=节点3是平衡节点，保持000.0000.1333j jf e U +=+=为定值。
华电潮流上机计算实验报告全解
院系：电气与电子工程学院班级：电气1205学号：1121181223学生姓名：王城指导教师：孙英云设计周数：两周成绩：日期：2015年7月7日一、课程设计的目的与要求培养学生的电力系统潮流计算机编程能力，掌握计算机潮流计算的相关知识二、设计正文（详细内容见附录）1．手算2．计算机计算3．思考题三、课程设计总结或结论（详细内容见附录）四、参考文献1.《电力系统计算：电子数字计算机的应用》，西安交通大学等合编。
北京：水利电力出版社；2.《现代电力系统分析》，王锡凡主编，科学出版社；3.《电力系统稳态分析》，陈珩，中国电力出版社，1995年，第三版；附录（设计流程图、程序、表格、数据等）4.机算潮流程序及结果// dierti.cpp : 定义控制台应用程序的入口点。
//#include "stdafx.h"struct Line //线路结构体{int Num,NumI,NumJ; //线路号左节点名右节点名float R,X,B,K; //电阻电抗电纳变比（K等于1为普通支路，不等于1为变压器支路的变比） };struct Bus //节点结构体{int Num ;float Volt,Phase,GenP,GenQ,LoadP,LoadQ;int Type;};#include"stdio.h"#include"string.h"#include"math.h"#include"stdlib.h"#define NBUS 4#define NLINE 4/* Global variables */int nL,nB,nVA,nSH;float X[NBUS];int L;double def[2*NBUS];double mn[50];void Gauss(double a[50][50],double b[50], int n) /*定义高斯法 */{int JS[50];int i,j,k;float d,t,x[50];FILE *fp;int L=1;for(i=0;i<50;i++) JS[i]=0;for(k=0;k<n;k++){d=0.0;for(j=k;j<n;j++)if(fabs(a[k][j])>d){ /*在一行中找到一个最大值赋值d,并用JS[K]记住这个最大值所在的列号*/ d=fabs(a[k][j]);JS[k]=j;}if(fabs(d)<0.000001) /*如果d的数值太小，做为被除数将带来很大的误差 */L=0;else {if(JS[k]!=k)for(i=0;i<n;i++){t=a[i][k];a[i][k]=a[i][JS[k]]; /*进行列交换，让最大值始终在对角元上*/a[i][JS[k]]=t;}}if(L==0)break;for(j=k+1;j<n;j++)a[k][j]=a[k][j]/a[k][k]; /*对角元上的元素消为1*/b[k]=b[k]/a[k][k];for(i=k+1;i<n;i++){for(j=k+1;j<n;j++)a[i][j]=a[i][j]-a[i][k]*a[k][j]; /*使下三角阵的元素为0*/b[i]=b[i]-a[i][k]*b[k];}}if(fabs(a[n-1][n-1])>0.00001){ /*用追赶法，解方程组，求未知数x*/ x[n-1]=b[n-1];for(i=n-2;i>=0;i--){t=0.0;for(j=i+1;j<n;j++)t=t+a[i][j]*x[j];x[i]=(b[i]-t);}}if((fp=fopen("gauss.txt","w"))==NULL) /*将结果写到TXT文件中*/{printf("err");exit(0);}for(i=0;i<n;i++){fprintf(fp,"%lf",x[i]);mn[i]=x[i];fprintf(fp,"\n");}fclose(fp);if(fp!=NULL) fclose(fp);}int _tmain(int argc, _TCHAR* argv[]){FILE *fp;FILE *fpout;int i,j,k,l,h,n,v;int i1,i2,i3,kp,kq;float d1,d2,d3,d4,d5,d6,r,x,g,b,tt,LL,e,ps,qs,shsh,m;struct Line sL[NLINE];struct Bus sB[NBUS];float YG[NBUS+1][NBUS+1],YB[NBUS+1][NBUS+1];double u[50][2];i1=i2=i3=0;d1=d2=d3=d4=d5=d6=ps=qs=0.0;for(i=0;i<NBUS;i++)if((fp=fopen("in.txt","r"))==NULL){ printf("Can not open the file named 'in.txt' \n");exit(0);}fscanf(fp,"%d,%d,%d",&nB,&nL,&nSH);for(i=0;i<nB;i++){sB[i].Num=sB[i].Type=0;sB[i].Volt=1.0;sB[i].Phase=sB[i].GenP=sB[i].GenQ=sB[i].LoadP=sB[i].LoadQ=0.0;fscanf(fp,"%d,%f,%f,%f,%f,%f,%f,%d",&i1,&d1,&d2,&d3,&d4,&d5,&d6,&i2);sB[i].Num=i1;sB[i].Volt=d1;sB[i].Phase=d2;sB[i].GenP=d3;sB[i].GenQ=d4;sB[i].LoadP=d5,sB[i].LoadQ=d6;sB[i].T ype=i2;};for(i=0;i<nL;i++){sL[i].Num=sL[i].NumI=sL[i].NumJ=0;sL[i].R=sL[i].X=sL[i].B=0.0;sL[i].K=1.0;fscanf(fp,"%2d %3d %3d %f %f %f %f",&i1,&i2,&i3,&d1,&d2,&d3,&d4);sL[i].Num=i1;sL[i].NumI=i2;sL[i].NumJ=i3;sL[i].R=d1;sL[i].X=d2;sL[i].B=d3;sL[i].K=d4;}if(fp!=NULL) fclose(fp);/*Make Y Matrix*/for(i=1;i<nB+1;i++)for(j=1;j<nB+1;j++){YG[i][j]=0.0;YB[i][j]=0.0;};for(l=0; l<nL; l++){i=sL[l].NumI;j=sL[l].NumJ;r=sL[l].R;x=sL[l].X;d1=r*r+x*x;g=r/d1;b=-x/d1;m=sL[l].K;if(fabs(sL[l].K-1.0)<0.000001) //普通支路 {YG[i][i]=YG[i][i]+g;YG[j][j]=YG[j][j]+g;YB[i][i]=YB[i][i]+b+sL[l].B;YB[j][j]=YB[j][j]+b+sL[l].B;YG[i][j]=YG[i][j]-g;YG[j][i]=YG[j][i]-g;YB[i][j]=YB[i][j]-b;YB[j][i]=YB[j][i]-b;}else //变压器支路{YG[i][i]=YG[i][i]+g/m+g*(m-1)/m;YG[j][j]=YG[j][j]+g/m+g*(1-m)/m/m;YB[i][i]=YB[i][i]+b/m+b*(m-1)/m;YB[j][j]=YB[j][j]+b/m+b*(1-m)/m/m;YG[i][j]=YG[i][j]-g/m;YG[j][i]=YG[j][i]-g/m;YB[i][j]=YB[i][j]-b/m;YB[j][i]=YB[j][i]-b/m; }}/* Check the Y matrix */if((fp=fopen("GGBB.txt","w"))==NULL){printf("Can not open the file named 'GGBB.txt' \n");exit(0);}fprintf(fp,"---Y Matrix---\n");for(i=1;i<nB+1;i++)for(j=1;j<nB+1;j++)if(fabs(YB[i][j]-0.0)>0.000001) fprintf(fp,"Y(%3d,%-3d)=(%10.5f,%10.5f)\n",i,j,YG[i][j],YB[i][j]);if(fp!=NULL) fclose(fp);/* 节点电压附初值 */for(i=1;i<nB+1;i++){if(sB[i-1].Type==0){u[i][0]=0.0;u[i][1]=1.0;}else if(sB[i-1].Type==1){u[i][1]=sB[i-1].Volt;u[i][0]=0.0;}else if(sB[i-1].Type==2){u[i][1]=sB[i-1].Volt;u[i][0]= sB[i-1].Phase;}}for(v=1;;v++)/* 迭代次数可以无限大 */{/* 节点电压附初值 */printf("迭代第%d次赋予的电压初值为e+jf:\n",v); for(i=1;i<nB+1;i++)printf("%lf,%lf\n",u[i][1],u[i][0]);printf("\n");printf("\n");/* 求偏移量 */double P_P[10];double P_Q[10];double P_UU[10];for(i=1;i<nB+1;i++){if(sB[i-1].Type==2){P_P[i]=0.0;P_Q[i]=0.0;P_UU[i]=1.05;}if(sB[i-1].Type==0){double tempP=0.0;double tempQ=0.0;for(j=1;j<nB+1;j++){tempP+=YG[i][j]*u[j][1]-YB[i][j]*u[j][0];tempQ+=YG[i][j]*u[j][0]+YB[i][j]*u[j][1];}P_P[i]=(sB[i-1].GenP-sB[i-1].LoadP)-tempP*u[i][1]-tempQ*u[i][0]; P_Q[i]=(sB[i-1].GenQ-sB[i-1].LoadQ)-tempP*u[i][0]+tempQ*u[i][1]; P_UU[i]=0.0;}if(sB[i-1].Type==1){double tempP=0.0;double tempQ=0.0;for(j=1;j<nB+1;j++){tempP+=YG[i][j]*u[j][1]-YB[i][j]*u[j][0];tempQ+=YG[i][j]*u[j][0]+YB[i][j]*u[j][1];P_P[i]=(sB[i-1].GenP-sB[i-1].LoadP)-tempP*u[i][1]-tempQ*u[i][0]; }P_UU[i]=sB[i-1].Volt*sB[i-1].Volt-u[i][1]*u[i][1]-u[i][0]*u[i][0]; P_Q[i]=0.0;}}/* 偏移量阵 */double P_PQ[6];int a=0;for(i=1;i<3;i++){P_PQ[a]=P_P[i];a=a+2;}a=1;for(i=1;i<3;i++){P_PQ[a]=P_Q[i];a=a+2;P_PQ[4]=P_P[3];P_PQ[5]=P_UU[3];printf("迭代第%d次的偏移量为：\n",v);for(i=0;i<6;i++){printf("%f",P_PQ[i]);printf("\n");}printf("\n");printf("\n");/* 雅可比矩阵 */double H[6][6],N[6][6],J[6][6],L[6][6],R[6][6],S[6][6],aa[6],bb[6]; for(i=1;i<5;i++){ if(fabs(sB[i-1].Type-2.0)<0.000001)continue;else{for(j=1;j<5;j++)if(i!=j){H[i][j]=-YB[i][j]*u[i][1]+YG[i][j]*u[i][0];N[i][j]=YG[i][j]*u[i][1]+YB[i][j]*u[i][0];J[i][j]=-N[i][j];L[i][j]=H[i][j];R[i][j]=0;S[i][j]=0;}else{aa[i]=bb[i]=0.0;aa[i]+=YG[i][n]*u[n][1]-YB[i][n]*u[n][0];bb[i]+=YG[i][n]*u[n][0]+YB[i][n]*u[n][1];}H[i][i]=-YB[i][i]*u[i][1]+YG[i][i]*u[i][0]+bb[i]; N[i][i]=YG[i][i]*u[i][1]+YB[i][i]*u[i][0]+aa[i]; J[i][i]=-YG[i][i]*u[i][1]-YB[i][i]*u[i][0]+aa[i]; L[i][i]=YG[i][i]*u[i][0]-YB[i][i]*u[i][1]-bb[i]; R[i][i]=2*u[i][0];S[i][i]=2*u[i][1];}}}double ss[50][50];for(i=0;i<6;i++)for(j=0;j<6;j++)ss[i][j]=0.0;for(i=1;i<3;i++)for(j=1;j<4;j++){ss[2*i-2][2*j-2]=H[i][j];ss[2*i-2][2*j-1]=N[i][j];ss[2*i-1][2*j-2]=J[i][j];ss[2*i-1][2*j-1]=L[i][j];}i=3;for(j=1;j<4;j++){ss[2*i-2][2*j-2]=H[i][j];ss[2*i-2][2*j-1]=N[i][j];ss[2*i-1][2*j-2]=R[i][j];ss[2*i-1][2*j-1]=S[i][j];}printf("迭代第%d次的雅可比矩阵为：\n",v);for(i=0;i<6;i++){for(j=0;j<6;j++)printf("%10f",ss[i][j]);printf("\n");}printf("\n");printf("\n");Gauss(ss,P_PQ,6);for(i=1;i<nB;i++){u[i][0]=u[i][0]+mn[2*(i-1)];u[i][1]=u[i][1]+mn[2*i-1];}double max;max=fabs(P_PQ[0]);for(i=0;i<=5;i++)if (max<fabs(P_PQ[i]))max=fabs(P_PQ[i]);if(fabs(max)<0.0001){printf("满足精度要求，迭代终止，迭代次数为%d\n",v); printf("\n");}/* 叠代循环的括号 */printf("最终求得的节点电压值为e+jf:\n");for(i=1;i<nB+1;i++)printf("%lf,%lf\n",u[i][1],u[i][0]);printf("\n");printf("\n");double uu[5],Phase[5];for(i=1;i<nB+1;i++){uu[i]=sqrt(u[i][1]*u[i][1]+u[i][0]*u[i][0]); Phase[i]=atan(u[i][0]/u[i][1]);}for(i=1;i<nB+1;i++)printf("%lf,%lf\n",uu[i],Phase[i]);*计算线路功率和平衡节点 PV节点功率*/double P[5],Q[5];double tempP=0.0;double tempQ=0.0;for(i=1;i<nB+1;i++){for(j=1;j<nB+1;j++){tempP+=YG[i][j]*u[j][1]-YB[i][j]*u[j][0]; tempQ+=YG[i][j]*u[j][0]+YB[i][j]*u[j][1];}P[i]=tempP*u[i][1]+tempQ*u[i][0];Q[i]=tempP*u[i][0]-tempQ*u[i][1];}for(i=1;i<nB+1;i++)printf("节点%d注入功率为%lf,%lf\n",i,P[i],Q[i]);/* 支路功率 */double V[4][2];for(i=1;i<5;i++)for(j=0;j<3;j++)V[i][j]=u[i][j];double sP[5][5],sQ[5][5];double dsq,dsp,dp,sumgen;for(i=1;i<NBUS+1;i++){for(j=1;j<NBUS+1;j++){sP[i][j]=0.0;sQ[i][j]=0.0;}}for(l=0; l<nL; l++){i=sL[l].NumI;j=sL[l].NumJ;r=sL[l].R;x=sL[l].X;d1=r*r+x*x;if(fabs(sL[l].K-1.0)<0.000001){/*Normal lines or transformers*/sP[i][j]=V[i][1]*V[i][1]*g-V[i][1]*V[j][1]*(g*cos(V[i][0]-V[j][0])+b*sin(V[i][0]-V[j][0]));sQ[i][j]=-(V[i][1]*V[i][1]*sL[l].B+V[i][1]*V[i][1]*b+V[i][1]*V[j][1]*(g*sin(V[i][0]-V[j][0])-b*cos(V[i ][0]-V[j][0])));sP[j][i]=V[j][1]*V[j][1]*g-V[i][1]*V[j][1]*(g*cos(V[j][0]-V[i][0])+b*sin(V[j][0]-V[i][0]));sQ[j][i]=-(V[j][1]*V[j][1]*sL[l].B+V[j][1]*V[j][1]*b+V[i][1]*V[j][1]*(g*sin(V[j][0]-V[i][0])-b*cos(V[j ][0]-V[i][0])));}else{/*abnormal transformer ratio*/sP[i][j]=V[i][1]*V[i][1]*g/sL[l].B/sL[l].B-V[i][1]*V[j][1]*(g*cos(V[i][0]-V[j][0])/sL[l].B+b*sin(V[i][ 0]-V[j][0])/sL[l].B);sQ[i][j]=-(V[i][1]*V[i][1]*b/sL[l].B/sL[l].B+V[i][1]*V[j][1]*(g*sin(V[i][0]-V[j][0])/sL[l].B-b*cos(V[i ][0]-V[j][0])/sL[l].B));sP[j][i]=V[j][1]*V[j][1]*g-V[i][1]*V[j][1]*(g*cos(V[j][0]-V[i][0])/sL[l].B+b*sin(V[j][0]-V[i][0])/sL[l ].B);sQ[j][i]=-(V[i][1]*V[i][1]*b+V[i][1]*V[j][1]*(g*sin(V[j][0]-V[i][0])/sL[l].B-b*cos(V[j][0]-V[i][0])/sL [l].B));}}/* 输电效率 */dsp=P[4];sumgen=P[4];for(i=0;i<NBUS;i++){dsp+=sB[i].GenP-sB[i].LoadP;dsq+=sB[i].GenQ-sB[i].LoadQ;sumgen+=sB[i].GenP;}dp=dsp/sumgen*100;/* 输出功率情况 */if((fp=fopen("功率情况.txt","w"))==NULL){printf("Can not open the file named '功率情况.txt' \n");exit(0);}fprintf(fp,"---功率情况---\n");fprintf(fp,"平衡节点功率S=%10.5f+ j%10.5f\n",P[4],Q[4]);for(i=1;i<NBUS+1;i++)for(j=1;j<NBUS+1;j++)if(fabs(sP[i][j]-0.0)>0.000001)fprintf(fp,"S(%3d,%-3d)=(%10.5f,j%10.5f)\n",i,j,sP[i][j],sQ[i][j]); fprintf(fp,"网损为%10.5f+j%10.3f,输电效率为%10.3f\n",dsp,dsq,100-dp);if(fp!=NULL) fclose(fp);return 0;}结果：1．导纳阵Y( 1,1 )=( 1.01534, -8.19201) Y( 1,2 )=( -0.56148, 2.30208) Y( 1,3 )=( 0.00000, 3.66667) Y( 1,4 )=( -0.45386, 1.89107) Y( 2,1 )=( -0.56148, 2.30208) Y( 2,2 )=( 1.04225, -4.67651) Y( 2,4 )=( -0.48077, 2.40385) Y( 3,1 )=( 0.00000, 3.66667) Y( 3,3 )=( 0.00000, -3.33333) Y( 4,1 )=( -0.45386, 1.89107) Y( 4,2 )=( -0.48077, 2.40385) Y( 4,4 )=( 0.93463, -4.26159)2．设定电压初值01.1;01;01)0(3)0(3)0(2)0(2)0(1)0(1j jf e j jf e j jf e +=++=++=+ 3．计算功率和电压偏移；27731.0])()([41)0(11)0(1)0(141)0(11)0(1)0(11)0(11)0(1-=++--=-=∆∑∑==j j j j j j jj s s e B f G ff B eG e P P P P0.05097])()([41)0(11)0(1)0(141)0(11)0(1)0(11)0(11)0(1-=----=-=∆∑∑==jj j j jj jj s s e B f G ef B e G f Q Q Q Q同理可算出52596.0)0(22)0(2-=-=∆P P P s ,0196.0)0(22)0(2=-=∆Q Q Q s 5.0)0(33)0(3=-=∆P P P s ,0.02)0(3232)0(3=-=∆U U U s 4．根据求的第一次迭代时雅可比矩阵各元素的公式计算雅可比矩阵各个元素的具体值：⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤-----------⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡20000.200066667.300003334.40052691.406629.130208.256148.00001821.182612.456148.030208.266667.3030208.256148.006298.803803.1066667.356148.030208.299265.032104.85．求高斯计算后的修正量：⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡∆∆∆∆∆∆=⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡∆∆∆∆∆∆-0.0000000.1276520.023566-0.108546-0.006511-0.007919-2)0(3)0(3)0(2)0(2)0(1)0(11)0()0(3)0(3)0(2)0(2)0(1)0(1U P Q P Q P e f e f e f J 6．计算各节点电压的一次近似值：12765.010855.000792.010000.197643.099349.0)0(3)0(3)1(3)0(2)0(2)1(2)0(1)0(1)1(1)0(3)0(3)1(3)0(2)0(2)1(2)0(1)0(1)1(1=∆+=-=∆+=-=∆+==∆+==∆+==∆+=f f f f f f fffe e e e e e e e e返回第三步重新迭代，并校验收敛与否，令410-=ε。
华北电力大学潮流上机课程设计报告
课程设计报告( 2014—2015年度第一学期) 名称：电力系统潮流上机院系：电气与电子工程学院班级：电气1215 学号：学生姓名：指导教师：王莉丽设计周数：两周成绩：日期： 2015年1月4日一、课程设计的目的与要求培养学生的电力系统潮流计算机编程能力，掌握计算机潮流计算的相关知识二、设计正文（详细内容见附录）1．手算2．计算机计算3．思考题三、课程设计总结或结论四、参考文献1.《电力系统分析基础》，李庚银，机械工业出版社，2011年，第一版；2.《电力系统稳态分析》，陈珩，中国电力出版社，2007年，第三版；附录（包括：1.手算过程及结果；2.计算机计算流程图、表格、数据；3.思考题答案）附录1.手算过程及结果2.计算机计算流程图3.完整的潮流程序头文件（定义变量）：#define Bus_Num 9#define Line_Num 9#define Precision 1e-5struct Bus{int No ;float Voltage,Phase,GenP,GenQ,LoadP,LoadQ;int Type;}gBus[Bus_Num];struct Line{int No,No_I,No_J;float R,X,B,k;}gLine[Line_Num];float gY_G[Bus_Num][Bus_Num],gY_B[Bus_Num][Bus_Num];float gDelta_P[Bus_Num-1],gDelta_Q[Bus_Num-1],gDelta_PQ[2*(Bus_Num-1)]; float gJaccobi[2*(Bus_Num-1)][2*(Bus_Num-1)];float gDelta_f[Bus_Num-1],gDelta_e[Bus_Num-1],gDelta_fe[2*(Bus_Num-1)]; float gf[Bus_Num],ge[Bus_Num];主程序：// flow.cpp: 主项目文件#include"stdafx.h"#include"NEquation.h"#include"math.h"#include"stdio.h"#include"config.h"using namespace System;void test(){NEquation ob1;ob1.SetSize(2);ob1.Data(0,0)=1;ob1.Data(0,1)=2;ob1.Data(1,0)=2;ob1.Data(1,1)=1;ob1.Value(0)=4;ob1.Value(1)=6;ob1.Run();printf("x1=%f\n",ob1.Value(0));printf("x2=%f\n",ob1.Value(1));}void GetData()//Read the data{FILE *fp;int i;fp=fopen("E:\\1121960726\\flow\\data\\data.txt","r");if(fp==NULL){printf("Can not open the file named 'data.txt' \n");return;}for(i=0;i<=Bus_Num-1;i++){fscanf(fp,"%d,%f,%f,%f,%f,%f,%f,%d",&gBus[i].No,&gBus[i].Voltage,&gBus[i].Phase, &gBus[i].GenP,&gBus[i].GenQ,&gBus[i].LoadP,&gBus[i].LoadQ,&gBus[i].Type);}for(i=0;i<=Line_Num-1;i++){fscanf(fp,"%d,%d,%d,%f,%f,%f,%f",&gLine[i].No,&gLine[i].No_I,&gLine[i].No_J, &gLine[i].R,&gLine[i].X,&gLine[i].B,&gLine[i].k);}fclose(fp);}void GetYMatrix(){int i,j,bus1,bus2;float r,x,d,g,b,g1,b1,g2,b2,g3,b3;FILE *fp;for(i=0;i<=Bus_Num-1;i++){for(j=0;j<=Bus_Num-1;j++){gY_G[i][j]=0;gY_B[i][j]=0;}}for(i=0; i<=Line_Num-1; i++){if(gLine[i].k==0){bus1=gLine[i].No_I-1;bus2=gLine[i].No_J-1;r=gLine[i].R;x=gLine[i].X;d=r*r+x*x;g=r/d;b=-x/d;gY_G[bus1][bus1]=gY_G[bus1][bus1]+g;gY_G[bus2][bus2]=gY_G[bus2][bus2]+g;gY_G[bus1][bus2]=gY_G[bus1][bus2]-g;gY_G[bus2][bus1]=gY_G[bus2][bus1]-g;gY_B[bus1][bus1]=gY_B[bus1][bus1]+b+gLine[i].B;gY_B[bus2][bus2]=gY_B[bus2][bus2]+b+gLine[i].B;gY_B[bus1][bus2]=gY_B[bus1][bus2]-b;gY_B[bus2][bus1]=gY_B[bus2][bus1]-b; }else{bus1=gLine[i].No_I-1;bus2=gLine[i].No_J-1;x=gLine[i].X;d=r*r+x*x;g=r/d;b=-x/d;g1=g/gLine[i].k;b1=b/gLine[i].k;g2=g*(1-gLine[i].k)/(gLine[i].k*gLine[i].k);b2=b*(1-gLine[i].k)/(gLine[i].k*gLine[i].k);g3=g*(gLine[i].k-1)/gLine[i].k;b3=b*(gLine[i].k-1)/gLine[i].k;gY_G[bus1][bus1]=gY_G[bus1][bus1]+g1+g2;gY_G[bus2][bus2]=gY_G[bus2][bus2]+g1+g3;gY_G[bus1][bus2]=gY_G[bus1][bus2]-g1;gY_G[bus2][bus1]=gY_G[bus2][bus1]-g1;gY_B[bus1][bus1]=gY_B[bus1][bus1]+b1+b2;gY_B[bus2][bus2]=gY_B[bus2][bus2]+b1+b3;gY_B[bus1][bus2]=gY_B[bus1][bus2]-b1;gY_B[bus2][bus1]=gY_B[bus2][bus1]-b1;}}// output the Y matrixfp=fopen("E:\\1121960726\\flow\\data\\ymatrix.txt","w");if(fp==NULL){printf("Can not open the file named 'ymatrix.txt' \n");return ;}fprintf(fp,"---Y Matrix---\n");for(i=0;i<=Bus_Num-1;i++){for(j=0;j<=Bus_Num-1;j++){fprintf(fp,"Y(%d,%d)=(%10.5f,%10.5f)\n",i+1,j+1,gY_G[i][j],gY_B[i][j]);}}fclose(fp);}void SetInitial(){int i;for(i=0;i<=Bus_Num-1;i++){if(gBus[i].Type==3){gf[i]=gBus[i].Voltage*sin(gBus[i].Phase);ge[i]=gBus[i].Voltage*cos(gBus[i].Phase);}else{gf[i]=0;ge[i]=1;}}}void GetUnbalance(){int i,j;FILE *fp;for(i=0;i<=Bus_Num-2;i++){gDelta_P[i]=gBus[i+1].GenP-gBus[i+1].LoadP;if(gBus[i+1].Type==2) //PV节¨2点ì?gDelta_Q[i]=gBus[i+1].Voltage*gBus[i+1].Voltage-(ge[i+1]*ge[i+1]+gf[i+1]*gf[i+1]);elsegDelta_Q[i]=gBus[i+1].GenQ-gBus[i+1].LoadQ;for(j=0;j<=Bus_Num-1;j++){gDelta_P[i]=gDelta_P[i]-ge[i+1]*(gY_G[i+1][j]*ge[j]-gY_B[i+1][j]*gf[j])-gf[i+1]*(gY_G[i+1][j]*gf[j]+gY_B[i+1][j]*ge[j]);if(gBus[i+1].Type==1) //PQ节¨2点ì?gDelta_Q[i]=gDelta_Q[i]-gf[i+1]*(gY_G[i+1][j]*ge[j]-gY_B[i+1][j]*gf[j])+ge[i+1]*(gY_G[i+1][j]*gf[j]+gY_B[i+1][j]*ge[j]);}}for(i=0;i<=Bus_Num-2;i++) //合?并?é{gDelta_PQ[2*i]=gDelta_P[i];gDelta_PQ[2*i+1]=gDelta_Q[i];}fp=fopen("E:\\1121960726\\flow\\data\\unbalance.txt","w");if(fp==NULL){printf("无法打开文件：”'unbalance.txt' \n");return ;}fprintf(fp,"---Unbalance---\n");for(i=0;i<=2*Bus_Num-3;i++){fprintf(fp,"Unbalance[%d]=%10.5f\n",i+1,gDelta_PQ[i]);}fclose(fp);}void GetJaccobi(){int i,j;float ga[Bus_Num-1],gb[Bus_Num-1];FILE *fp;for(i=0;i<=Bus_Num-2;i++) //计算注入电流{ga[i]=0;gb[i]=0;for(j=0;j<=Bus_Num-1;j++){ga[i]=ga[i]+gY_G[i+1][j]*ge[j]-gY_B[i+1][j]*gf[j];gb[i]=gb[i]+gY_G[i+1][j]*gf[j]+gY_B[i+1][j]*ge[j];}}for(i=0;i<=Bus_Num-2;i++){for(j=0;j<=Bus_Num-2;j++){if(i!=j){gJaccobi[2*i][2*j]=-gY_B[i+1][j+1]*ge[i+1]+gY_G[i+1][j+1]*gf[i+1];gJaccobi[2*i][2*j+1]=gY_G[i+1][j+1]*ge[i+1]+gY_B[i+1][j+1]*gf[i+1];if(gBus[i+1].Type==2) //PV节¨2点ì?{gJaccobi[2*i+1][2*j]=0;gJaccobi[2*i+1][2*j+1]=0;}else//PQ{gJaccobi[2*i+1][2*j]=-gJaccobi[2*i][2*j+1];gJaccobi[2*i+1][2*j+1]=gJaccobi[2*i][2*j];}}else{gJaccobi[2*i][2*j]=-gY_B[i+1][j+1]*ge[i+1]+gY_G[i+1][j+1]*gf[i+1]+gb[i];gJaccobi[2*i][2*j+1]=gY_G[i+1][j+1]*ge[i+1]+gY_B[i+1][j+1]*gf[i+1]+ga[i];if(gBus[i+1].Type==2) //PV节¨2点ì?{gJaccobi[2*i+1][2*j]=2*gf[i+1];gJaccobi[2*i+1][2*j+1]=2*ge[i+1];}else//PQ节点{gJaccobi[2*i+1][2*j]=-gY_G[i+1][j+1]*ge[i+1]-gY_B[i+1][j+1]*gf[i+1]+ga[i];gJaccobi[2*i+1][2*j+1]=-gY_B[i+1][j+1]*ge[i+1]+gY_G[i+1][j+1]*gf[i+1]-gb[i];}}}}fp=fopen(,"w");if(fp==NULL){printf("无法打开文件：” 'jaccobi.txt' \n");return ;}fprintf(fp,"---Jaccobi Matrix---\n");for(i=0;i<=2*Bus_Num-3;i++){for(j=0;j<=2*Bus_Num-3;j++){fprintf(fp,"jaccobi(%d,%d)=%10.5f\n",i+1,j+1,gJaccobi[i][j]);}}fclose(fp);}void GetRevised(){int i,j;FILE *fp;NEquation ob1; //解矩阵方程ob1.SetSize(2*(Bus_Num-1));for(i=0;i<=2*Bus_Num-3;i++)for(j=0;j<=2*Bus_Num-3;j++)ob1.Data(i,j)=gJaccobi[i][j];for(i=0;i<=2*Bus_Num-3;i++)ob1.Value(i)=gDelta_PQ[i];ob1.Run();for(i=0;i<=Bus_Num-2;i++){gDelta_f[i]=ob1.Value(2*i);gDelta_e[i]=ob1.Value(2*i+1);gDelta_fe[2*i]=gDelta_f[i];gDelta_fe[2*i+1]=gDelta_e[i];}fp=fopen("E:\\1121960726\\flow\\data\\revised.txt","w");if(fp==NULL){printf("无法打开文件：” 'revised.txt' \n");return ;}fprintf(fp,"---Revised---\n");for(i=0;i<=2*Bus_Num-3;i++){fprintf(fp,"revised[%d]=%10.5f\n",i+1,gDelta_fe[i]);}fclose(fp);}void GetNewValue(){int i;FILE *fp;for(i=0;i<=Bus_Num-2;i++){gf[i+1]=gf[i+1]+gDelta_f[i];ge[i+1]=ge[i+1]+gDelta_e[i];}fp=fopen("E:\\1121960726\\flow\\data\\newvalue.txt","w");if(fp==NULL){printf("无法打开文件：”'newvalue.txt' \n");return ;}fprintf(fp,"---New Value---\n");for(i=0;i<=Bus_Num-2;i++){fprintf(fp,"f(%d)=%10.5f,e(%d)=%10.5f\n",i+1,gf[i+1],i+1,ge[i+1]); }fclose(fp);}int main(array<System::String ^> ^args){int i,Count_Num;float maxValue;//test();GetData();GetYMatrix();SetInitial();for(Count_Num=0;Count_Num<=100;Count_Num++) {GetUnbalance();GetJaccobi();GetRevised();GetNewValue();maxValue=fabs(gDelta_fe[0]);for(i=1;i<=2*(Bus_Num-1)-1;i++){if(maxValue<fabs(gDelta_fe[i])){maxValue=fabs(gDelta_fe[i]);}}if(maxValue<Precision){break;}}printf("%d\n",Count_Num);for(i=0;i<=Bus_Num-1;i++){printf("%10.5f\n",sqrt(ge[i]*ge[i]+gf[i]*gf[i]));}while(true){}return 0;}计算结果：（1）原始数据：（2）节点导纳矩阵：(3)Jacobi矩阵：（4）NewVoltageValue:(5)不平衡量：（6）修正量：（6）结果：3.回答思考题1.潮流计算的方法有哪些？各有何特点？2. 如果交给你一个任务，请你用已有的潮流计算软件计算北京城市电网的潮流，你应该做哪些工作？（收集哪些数据，如何整理，计算结果如何分析）3. 设计中遇到的问题和解决的办法。
华中科技大学现代电力系统分析潮流计算作业
现代电力系统分析作业------基于Matpower的电力系统潮流计算专业：班级：姓名：学号：目录基于Matpower的电力系统潮流计算 (1)1.本次潮流计算的目的及意义 (1)2.电力系统潮流计算及其意义 (1)3.电力系统潮流计算常规方法 (1)3.1 牛顿-拉夫逊法 (1)3.2 节点电压用直角坐标表示时的牛顿—拉夫逊潮流计算 (3)3.3 牛顿—拉夫逊法潮流计算程序框图 (4)4.选用的潮流计算的系统 (5)5.利用软件matpower计算潮流 (7)5.1 matpower简介 (7)5.2 基态潮流计算 (7)5.21 基态潮流计算条件说明及数据输入 (7)5.22 潮流计算结果 (9)5.23 基态潮流结果计算分析 (10)5.3 最优潮流计算 (11)5.31 最优潮流计算条件说明及数据输入 (12)5.32 最优潮流的理论结果 (12)5.33 最优潮流仿真计算结果 (12)5.23 最优潮流结果计算分析 (13)6.感想与小结 (14)参考文献 (15)附录一、基态潮流计算文件 (15)附录二、最优潮流计算文件 (17)基于Matpower的电力系统潮流计算1.本次潮流计算的目的及意义本次潮流计算的目的及意义主要是了解电力系统潮流计算及其意义，在此基础上，了解电力系统潮流计算的模型以及常规的潮流计算的方法，掌握并熟练使用电力系统潮流计算软件Matpower。
最后，利用Matpower计算美国西部电网WSCC 三机九节点系统的静态潮流及最优潮流并给出分析。
2.电力系统潮流计算及其意义电力系统潮流计算是研究电力系统稳态运行情况的基本电气计算，电力系统潮流计算的任务是根据给定的网络结构及运行条件，求出电网的运行状态，其中包括各母线的电压、各支路的功率分布以及功率损耗等。
潮流计算分为离线计算和在线计算两大类。
离线计算可以用于电力系统调度，并确定系统的运行方式；离线计算的结果还可以用于电力系统规划方案的分析以及优化系统的运行状态；此外离线潮流计算可以作为初值，用于配合系统的故障分析以及稳定性分析。
华北电力大学 电力系统基础 第五章 简单电力系统潮流计算
dU 3ZIZ 3Z ( 3SU22 )
(R jX )( P2 jQ2 ) U2
P2R Q2 X j P2X Q2R U jU
U2
U2
U1 (U2 U)2 U2 arctg U
U2 U
近似计算：
(21.10
j 21.35)
0.0745
j0.0754MVA
S
'
2
(5.62
j3.34) (0.0745
j0.0754)
5.55
j3.26MVA
S
"
1
S'Fra bibliotek2 S1
(5.55
j3.26) (20.15
j13.96)
14.60
j10.70MVA
2 N
1000SN2
QZT
S22
U
2 2
*
U
k
%U
2 N
100SN
PyT
U12
*
P0 1000U
2 N
QyT
U12
*
I0 %SN
100U
2 N
发电厂
PZT
PkU
2 N
S1
2
1000U12 S N2
QZT
U
k
%U
2 N
S1
2
100U
2 2
S
N
PyT
电
电
所
厂
的
的
变
变
压
压
华电潮流上机课程设计报告程序..doc
院系：电气与电子工程学院班级：学号：学生姓名：指导教师：刘宝柱设计周数：成绩：日期：2012年1月5日课程课程设计报告一、课程设计的目的与要求培养学生的电力系统潮流计算机编程能力，掌握计算机潮流计算的相关知识二、设计正文（详细内容见附录）1．手算2．计算机计算3．思考题3.1潮流计算的方法有哪些？各有何特点？答：潮流计算的方法主要有高斯-赛德尔迭代法、牛顿-拉夫逊迭代法和P-Q分解法。
它们各自的特点如下：（1）高斯-赛德尔迭代法分为以节点导纳矩阵为基础的高斯-赛德尔迭代法和以以节点阻抗矩阵为基础的高斯-赛德尔迭代法的原理比较简单，要求的数字计算机的内存量也比较小，但收敛性差，当系统规模变大时，迭代次数急剧上升，往往出现迭代不收敛的情况；而阻抗法改善了电力系统潮流计算导纳法德收敛性问题，在当时获得了广泛的应用，但是，阻抗法的主要缺点是占用计算机的内存很大，每次迭代的计算量很大。
当系统不断扩大时，这些缺点就更加突出。
（2）牛顿-拉夫逊法是数学中求解非线性方程式的典型方法，有较好的收敛性。
只要在迭代过程中尽可能保持方程式系数矩阵的稀疏性，就可以大大提高牛顿潮流计算程序的计算效率，牛顿法在收敛性、内存要求、计算速度方面都超过了阻抗法，成为知道目前仍被广泛采用的方法。
（3）P-Q分解法潮流计算派生于以极坐标表示时的牛顿-拉夫逊法，它根据电力系统的特点，抓住主要矛盾，对纯数学的牛顿法进行了改造。
与牛顿法相比，P-Q分解法的修正方程的系数矩阵B’和B”分别是（n-1）和（m-1）的方阵，替代了原有的（n+m-2）阶系数矩阵J；B’、B”在迭代过程中可以保持不变且为对称的系数矩阵，提高了计算速度，降低了对存储容量的要求。
P-Q分解法在计算速度方面有显著地提高，迅速得到了推广。
3.2如果交给你一个任务，请你用已有的潮流计算软件计算北京城市电网的潮流，你应该做哪些工作？（收集哪些数据，如何整理，计算结果如何分析）答：(1)在进行北京城市电网的潮流计算之前需要了解北京城市电网中所有的节点支路的相关数据，并对节点和支路分类。
华北电力大学潮流上机课程设计报告
课程设计报告( 2011-- 2012年度第一学期)名称：电力系统潮流上机院系：电气与电子工程学院班级：电管0902学号：1091140202 学生：指导教师：麻秀设计周数：两周成绩：日期：2011年1月6日一、课程设计的目的与要求培养学生的电力系统潮流计算机编程能力，掌握计算机潮流计算的相关知识二、设计正文（详细容见附录）1、手算在此我们同样采取计算机计算这种方便快捷的方法。
具体的方法及结果见附录。
2、计算机计算通过计算机编程，我们可以很方便的求解出相关量，其具体计算方法及计算结果见附录。
3、思考题答案见附录。
三、课程设计总结或结论潮流计算是研究电力系统稳态运行的一种基本运算方法，通过对于潮流的计算，实现对潮流的控制，维护电力系统稳定，防止出现超过约束条件的情况，从而使电力系统免遭崩溃的危险。
最初求解电力系统潮流时大多是用手算，这种计算方法不但耗费时间长，精度差，而且对于大容量的电力系统无法做到很好的控制。
随着电力系统结构日趋复杂，计算量也越来越大，仅靠手算难以求解一些复杂网络的潮流。
计算机作为一种处理数据的工具，其计算速度快，准确率高，存储量大，因此它广泛的应用于大规模复杂化的电力系统网络中。
为了能使计算机能进行潮流计算，就必须编制相应的程序，使其能按照人们在程序中设定的计算流程一步一步的执行，直到计算结果达到要求的精度，得到最终想要的结果。
C语言是一种简单实用的语言，利用它进行程序设计可以大大提高电力系统潮流控制的自动化程度。
本次的潮流上机，就是要求我们利用C语言编程，进行计算机潮流计算。
这既是对我们电力系统潮流计算知识的考验，也是对我们C语言编程能力的考验。
是一次各种知识的交汇综合的运用。
通过这次潮流上机，我对电力系统潮流计算有了更深刻的认识，同时，我的C语言编程能力也得到了提高。
我明白了学习不能仅限于课本，还要做到“学以致用”的道理。
在这次潮流设计中，我遇到了由于多问题，其中既包括由于对潮流知识掌握不牢而导致的问题，也包括由于C语言编程能力有限而导致的问题。
电力系统潮流上机计算报告书写要求
电力系统潮流上机计算报告
系别：
班级：
姓名：
学号：
2012年×月×日
一、程序说明
包括：程序设计思想、程序流程图、程序使用说明。
二、给定题目的手算过程（迭代两次）
包括：原题目、节点导纳矩阵、雅克比矩阵、第一次和第二次迭代
结果。
三、给定题目的程序计算结果
包括：原题目、节点导纳矩阵、雅克比矩阵、程序输入和输出文件（误差0.0001）。
四、编程特色与创新
包括：程序能够完成的基本功能；程序能够完成的高级功能（如：是否包括非标准变比变压器支路，是否采用了稀疏矩阵技术，
是否增加了人机对话界面，程序的通用性和实用性如何）。
五、总结
包括：手算结果与程序计算结果的分析比较；本次上机体会，如：独立编程体会、跟踪调试技能的掌握情况、C语言中结构体、
指针、文件输入输出的掌握情况等。
报告要求：
1．报告中除上面的第“三”项外，其他部分必须手写（最好使用黑色水笔）。
2．报告统一采用A4打印纸书写（留出页边距: 1.5~2厘米）。
不使用实验报告纸。
3．封面按上述格式书写。
4．装订统一在左侧1厘米，二个钉。
5．上述五部分内容必须齐全，各部分内容可以扩充。
6．报告书写要求字迹清楚，不得潦草，流程图中的框要用尺子画。
7．报告必须与本人提交程序吻合，否则取消成绩。
8．报告不得有雷同，否则全部取消成绩。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
《电力系统潮流上机》课程设计报告院系：电气与电子工程学院班级：学号：学生姓名：指导教师：刘宝柱设计周数：两周成绩：日期：2012年1月5日一、课程设计的目的与要求培养学生的电力系统潮流计算机编程能力，掌握计算机潮流计算的相关知识 二、设计正文（详细内容见附录）1． 手算： 要求应用牛顿-拉夫逊法或P-Q 分解法手算求解，要求精度为。
节点1为平衡节点，电压︒∠=00.11U &，节点2为PQ 节点，负荷功率6.08.0~2j S +=，节点3是PV 节点，1.1,4.033==U P ，两条支路分别为04.001.013j Z +=，2.005.012j Z +=，对地支路33.030j y =。
2． 计算机计算 1．导纳阵Y( 1,1 )=( , Y( 1,2 )=( , Y( 1,3 )=( , Y( 1,4 )=( , Y( 2,1 )=( , Y( 2,2 )=( , Y( 2,4 )=( , Y( 3,1 )=( , Y( 3,3 )=( , Y( 4,1 )=( , Y( 4,2 )=( , Y( 4,4 )=( , 2．设定电压初值3．计算功率和电压偏移；同理可算出52596.0)0(22)0(2-=-=∆P P P s ,0196.0)0(22)0(2=-=∆Q Q Q s5.0)0(33)0(3=-=∆P P P s ,0.02)0(3232)0(3=-=∆U U U s4．根据求的第一次迭代时雅可比矩阵各元素的公式计算雅可比矩阵各个元素的具体值：5．求高斯计算后的修正量： 6．计算各节点电压的一次近似值： 返回第三步重新迭代，并校验收敛与否，令410-=ε。
经过 四 次迭代后，收敛条件满足，停止迭代，求出的电压终值：7．计算出平衡节点4的注入功率。
8．各节点间功率9．网损为:04487.001787.0~j S-= 10．网损效率为：% 3． 思考题潮流计算的方法有哪些各有何特点答：潮流计算方法主要包括：高斯-赛德尔迭代法、牛顿-拉夫逊迭代法、P-Q 分解法等算法。
各方法特点如下所示：如果交给你一个任务，请你用已有的潮流计算软件计算北京城市电网的潮流，你应该做哪些工作（收集哪些数据，如何整理，计算结果如何分析）答：1.所需要收集的数据：A.电网中所有节点的数据：a．各节点的类型,包括平衡节点、PV 节点、PQ 节点b. 对于平衡节点要了解节点的电压大小相位、及节点所能提供的最大最小有功无功功率c. PV节点要知道节点电压大小注入有功功率及节点所能提供的最大和最小无功功.率d. PQ节点要知道节点的注入有功和无功功率B．电网中所有支路的数据：a.各支路类型，即是否含有变压器b.各支路的电阻、电感、电纳c.各变压器的变比。
2．数据整理：将上述数据资料进行分类整理，并为每个节点及支路编上编号。
将整理的结果写成本实验中所要求的格式（原始数据的 txt 文档），再用本实验所编制的程序进行求解，得到各节点电压、相位，各线路传输功率、损耗，平衡节点注入功率等数值。
3.计算结果分析：可以大致了解整个电网的运行状况并按照电力系统知识依次进行所需的调整。
如调节各节点电压的大小以防止超出规定限度，调整各线路的有功无功分布状况分析整个系统的网损是否达到标准等。
设计中遇到的问题和解决的办法。
1.在设计最开始遇到的问题就是由于以前C语言所学习的知识遗忘较多，对于程序的使用不太熟练，通过老师讲解、看书、询问同学逐渐扫除了这个障碍。
2.设计的时候对于电力系统整体的知识脉络不是很清晰，编写程序时所定义的变量容易弄混淆，编写的程序产生了较多的错误，经过仔细学习和对每一排程序的排查，找出了产生错误的地方。
同时也对潮流计算的过程有了更清楚的认识。
三、课程设计总结或结论该门课程的主要目的是培养学生的电力系统潮流计算机编程能力，掌握计算机潮流计算的相关知识。
潮流计算是电力系统非常重要的分析计算，是电力系统最基本的计算，主要用以研究系统规划和运行中提出的各种问题。
现代电力系统潮流就算过程比较复杂，手算有很大的困难，而计算机作为一种处理数据的工具，具有计算速度快，准确率高的优点，因此现在电力系统的潮流大多也都采用计算机的算法求解。
对于写程序而言，C 语言是一种简单但很实用的语言，利用 C 语言编程实现潮流计算不仅简单易行，而且还能加深对电力系统潮流这一物理感念的理解。
通过这次潮流对 c++的面向对象编程及其过程有了更深入的了解。
对于潮流计算，我更深刻的掌握了牛拉法解潮流方程的方法以及对数学高次多元方程线性化在迭代求解的方法有了更深入的认识。
通过这次课程设计，自己动手编写潮流计算的程序，仔细分析了潮流计算各个步骤，理论与实践相结合,对潮流计算的思想和方法以及过程有了更深入的理解。
在这两周的学习里，我不仅对C语言的学习有了更深层次的理解，还对自己在《电力系统分析》这门课程里所学到的知识进行了查缺补漏，对电力系统的运行分析等有了更加深入的了解，受益匪浅。
四、参考文献1.《电力系统计算：电子数字计算机的应用》，西安交通大学等合编。
北京：水利电力出版社；2.《现代电力系统分析》，王锡凡主编，科学出版社；3.《电力系统稳态分析》，陈珩，中国电力出版社，1995年，第二版；附录（设计流程图、程序、表格、数据等）1.2.程序设计:#include""#include""#include""#include""#define NBUS 5#define NLINE 7/* Global variables */int nL,nSH,nBUS,nVA;double d,t;int L;struct Line{int Num,NumI,NumJ;double R;double X;double B;double K;};struct Bus{int Num ;double Volt,Phase,GenP,GenQ,LoadP,LoadQ; int Type;};struct Shunt{int Num,NumI;double G,B;};void main(){FILE *fp;FILE *fpout;int i,j,k,l,h,c;int i1,i2,i3,kp,kq,LH[100];double d,d1,d2,d3,d4,d5,d6,t,r,x,g,b,tt,LL,e,ps,qs,shsh,sum,w;struct Line sL[NLINE];struct Bus sBus[NBUS];struct Shunt sSH[NBUS];double YG[NBUS][NBUS],YB[NBUS][NBUS];doubleV[100][2],U[100],dPQ[100],PQ[100][2],JJ[100][100]={0},H[100][100],J[100][1 00],N[100][100],L[100][100],sP[NBUS][NBUS]={0},sQ[NBUS][NBUS]={0},dsp,dsq, sumgen,dp;i1=i2=i3=0;d1=d2=d3=d4=d5=d6=ps=qs=;e=;for(i=0;i<NBUS;i++)/*Read the input data*/if((fp=fopen("","r"))==NULL){ printf("Can not open the file named '' \n");exit(0);}fscanf(fp,"%d,%d,%d",&nBUS,&nL,&nSH);for(i=0;i<nBUS;i++){sBus[i].Num=sBus[i].Type=0;sBus[i].Volt=;sBus[i].Phase=sBus[i].GenP=sBus[i].GenQ=sBus[i].LoadP=sBus[i].LoadQ=;fscanf(fp,"%d,%lf,%lf,%lf,%lf,%lf,%lf,%d",&i1,&d1,&d2,&d3,&d4,&d5,&d6,&i2) ;sBus[i].Num=i1;sBus[i].Volt=d1;sBus[i].Phase=d2;sBus[i].GenP=d3;sBus[i].Ge nQ=d4;sBus[i].LoadP=d5,sBus[i].LoadQ=d6;sBus[i].Type=i2;};for(i=0;i<nL;i++){sL[i].Num=sL[i].NumI=sL[i].NumJ=0;sL[i].R=sL[i].X=sL[i].B=;sL[i].K=;fscanf(fp,"%2d %3d %3d %lf %lf %lf %lf",&i1,&i2,&i3,&d1,&d2,&d3,&d4);sL[i].Num=i1;sL[i].NumI=i2;sL[i].NumJ=i3;sL[i].R=d1;sL[i].X=d2;sL[i].B=d3; sL[i].K=d4;}for(i=0;i<nSH;i++){sSH[i].Num=sSH[i].NumI=0;sSH[i].G=sSH[i].B=; fscanf(fp,"%2d %3d %lf",&i1,&i2,&d1);sSH[i].Num=i1;sSH[i].NumI=i2;sSH[i].B=d1; }if(fp!=NULL) fclose(fp);/*Make Y Matrix*/for(i=1;i<nBUS+1;i++)for(j=1;j<nBUS+1;j++){ YG[i][j]=;YB[i][j]=;};for(l=0; l<nL; l++){i=sL[l].NumI;j=sL[l].NumJ;r=sL[l].R;x=sL[l].X;d1=r*r+x*x;g=r/d1;b=-x/d1;if(fabs(sL[l].<{/*Normal lines or transformers*/YG[i][i]=YG[i][i]+g;YG[j][j]=YG[j][j]+g;YB[i][i]=YB[i][i]+b+sL[l].B;YB[j][j]=YB[j][j]+b+sL[l].B;YG[i][j]=YG[i][j]-g;YG[j][i]=YG[j][i]-g;YB[i][j]=YB[i][j]-b;YB[j][i]=YB[j][i]-b;}else{/*abnormal transformer ratio*/if(fabs(sL[l].B)>{YG[i][i]=YG[i][i]+g/sL[l].B/sL[l].B;YB[i][i]=YB[i][i]+b/sL[l].B/sL[l].B;YG[j][j]=YG[j][j]+g;YB[j][j]=YB[j][j]+b;YG[i][j]=YG[i][j]-g/sL[l].B;YG[j][i]=YG[j][i]-g/sL[l].B;YB[i][j]=YB[i][j]-b/sL[l].B;YB[j][i]=YB[j][i]-b/sL[l].B;}elseprintf("节点%2d,%-2d之间变压器变比为0 \n",i,j); }}ype==0){V[i][0]=;V[i][1]=;}for(i=1;i<nBUS+1;i++)if(sBus[i-1].Type==1){V[i][1]=sBus[i-1].Volt;V[i][0]=;}for(i=1;i<nBUS+1;i++)if(sBus[i-1].Type==2){V[i][1]=sBus[i-1].Volt;V[i][0]=sBus[i-1].Phase;}/* 输出电压初值if((fp=fopen("d:\\lx\\电压初值.txt","w"))==NULL){printf("Can not open the file named '电压初值.txt' \n");exit(0);}fprintf(fp,"---电压初值---\n");for(i=1;i<nBUS+1;i++)for(j=1;j<2;j++)fprintf(fp,"Y(%2d )=(%,%\n",i,V[i][0],V[i][1]);if(fp!=NULL) fclose(fp);*/for(c=1;;c++){/* 计算偏移量 */for(i=1;i<nBUS+1;i++){if(sBus[i-1].Type!=2){PQ[i][0]=0;PQ[i][1]=0;for(j=1;j<nBUS+1;j++){PQ[i][0]-=V[i][1]*V[j][1]*(YG[i][j]*cos(V[i][0]-V[j][0])+YB[i][j]*sin( V[i][0]-V[j][0]));PQ[i][1]-=V[i][1]*V[j][1]*(YG[i][j]*sin(V[i][0]-V[j][0])-YB[i][j]*cos( V[i][0]-V[j][0]));}}}h=1;for(i=1;i<nBUS+1;i++){if(sBus[i-1].Type==0){dPQ[h]=PQ[i][0]+sBus[i-1].GenP-sBus[i-1].LoadP;h++;dPQ[h]=PQ[i][1]+sBus[i-1].GenQ-sBus[i-1].LoadQ;h++;}if(sBus[i-1].Type==1){dPQ[h]=PQ[i][0]+sBus[i-1].GenP-sBus[i-1].LoadP;h++;}}/* 输出偏移量if((fp=fopen("d:\\lx\\偏移量 .txt","w"))==NULL){printf("Can not open the file named '偏移量.txt' \n"); exit(0);}fprintf(fp,"--- 偏移量 ---\n");for(i=1;i<h;i++)fprintf(fp,"dPQ(%2d )=(%\n",i,dPQ[i]);if(fp!=NULL) fclose(fp);*//* 计算雅克比矩阵 */for(i=1;i<nBUS+1;i++){for(j=1;j<nBUS+1;j++){if(i==j){H[i][j]=PQ[i][1]-V[i][1]*V[i][1]*YB[i][j];J[i][j]=-PQ[i][0]-V[i][1]*V[i][1]*YG[i][j];N[i][j]=-PQ[i][0]+V[i][1]*V[i][1]*YG[i][j];L[i][j]=-PQ[i][1]-V[i][1]*V[i][1]*YB[i][j];}else{H[i][j]=V[i][1]*V[j][1]*(YG[i][j]*sin(V[i][0]-V[j][0])-YB[i][j]*cos(V[ i][0]-V[j][0]));J[i][j]=-V[i][1]*V[j][1]*(YG[i][j]*cos(V[i][0]-V[j][0])+YB[i][j]*sin(V [i][0]-V[j][0]));N[i][j]=V[i][1]*V[j][1]*(YG[i][j]*cos(V[i][0]-V[j][0])+YB[i][j]*sin(V[ i][0]-V[j][0]));L[i][j]=V[i][1]*V[j][1]*(YG[i][j]*sin(V[i][0]-V[j][0])-YB[i][j]*cos(V[ i][0]-V[j][0]));}}}h=1;for(i=1;i<nBUS+1;i++){l=1;if(sBus[i-1].Type==0){for(j=1;j<nBUS+1;j++){if(sBus[j-1].Type==0){JJ[h][l]=H[i][j];l++;JJ[h][l]=N[i][j];l++;}if(sBus[j-1].Type==1){JJ[h][l]=H[i][j];l++;}}h++;l=1;for(j=1;j<nBUS+1;j++) {if(sBus[j-1].Type==0){JJ[h][l]=J[i][j];l++;JJ[h][l]=L[i][j];l++;}if(sBus[j-1].Type==1){JJ[h][l]=J[i][j];l++;}}h++;}if(sBus[i-1].Type==1){for(j=1;j<nBUS+1;j++){if(sBus[j-1].Type==0){JJ[h][l]=H[i][j];l++;JJ[h][l]=N[i][j];l++;}if(sBus[j-1].Type==1){JJ[h][l]=H[i][j];l++;}}h++;}}/* 输出雅克比矩阵if((fp=fopen("d:\\lx\\雅克比矩阵 .txt","w"))==NULL){printf("Can not open the file named '雅克比矩阵.txt' \n"); exit(0);}fprintf(fp,"--- 雅克比矩阵 ---\n");fprintf(fp," ");for(i=1;i<h;i++){fprintf(fp," %d\t\t",i);}fprintf(fp,"\n");*/for(i=1;i<h;i++){fprintf(fp,"%2d ",i);for(j=1;j<h;j++){fprintf(fp,"%\t",JJ[i][j]);}fprintf(fp,"\n");}if(fp!=NULL) fclose(fp);/*高斯法求解方程组*/l=1;for(i=1;i<h;i++){LH[i]=0;}for(k=1;k<h;k++){d=;for(j=k;j<h;j++){if(fabs(JJ[k][j])>d){d=fabs(JJ[k][j]); /*在一行中找到一个最大值赋值d,并用JS[K]记住这个最大值所在的列号*/LH[k]=j;}}if(fabs(d)< /*如果d的数值太小，做为被除数将带来很大的误差*/{l=0;}else if(LH[k]!=k){for(i=1;i<h;i++){t=JJ[i][k];JJ[i][k]=JJ[i][LH[k]]; /*进行列交换，让最大值始终在对角元上*/JJ[i][LH[k]]=t;}}if(l==0){break;}for(j=k+1;j<h;j++){JJ[k][j]=JJ[k][j]/JJ[k][k]; /*对角元上的元素消为1*/ }dPQ[k]=dPQ[k]/JJ[k][k];for(i=k+1;i<h;i++){for(j=k+1;j<h;j++){JJ[i][j]=JJ[i][j]-JJ[i][k]*JJ[k][j]; /*使下三角阵的元素为0*/}dPQ[i]=dPQ[i]-JJ[i][k]*dPQ[k];}}if(fabs(JJ[h-1][h-1])>{ /*用追赶法，解方程组，求未知数x*/ U[h-1]=dPQ[h-1];for(i=h-2;i>=0;i--){t=;for(j=i+1;j<h;j++){t=t+JJ[i][j]*U[j];}U[i]=(dPQ[i]-t);}}/*输出高斯结果if((fp=fopen("d:\\lx\\高斯.txt","w"))==NULL) {printf("err");exit(0);}for(i=1;i<h;i++){fprintf(fp,"%f",U[i]);fprintf(fp,"\n");}fclose(fp);*//*得到电压值*/h=1;for(i=1;i<nBUS+1;i++){if(sBus[i-1].Type==0){V[i][0]+=U[h];h++;V[i][1]+=U[h];h++;}if(sBus[i-1].Type==1){V[i][0]+=U[h];h++;}}/* 输出电压值if((fp=fopen("d:\\lx\\电压值.txt","w"))==NULL){printf("Can not open the file named '电压值.txt' \n"); exit(0);}fprintf(fp,"---电压初值---\n");for(i=1;i<nBUS+1;i++)for(j=1;j<2;j++)fprintf(fp,"Y(%2d )=(%,%\n",i,V[i][0],V[i][1]);if(fp!=NULL) fclose(fp);*//*求最大变化值*/w=0;for(i=1;i<h;i++){if(U[i]>w){w=U[i];}}if(w<{break;}}/* 输出电压终值if((fp=fopen("d:\\lx\\电压终值.txt","w"))==NULL){printf("Can not open the file named '电压终值.txt' \n"); exit(0);}fprintf(fp,"---电压终值---\n");fprintf(fp,"循环%d次\n",c);for(i=1;i<nBUS+1;i++)for(j=1;j<2;j++)fprintf(fp,"Y(%2d )=(%,%\n",i,V[i][0]*180/,V[i][1]);if(fp!=NULL) fclose(fp);*/ps=0;qs=0;for(i=1;i<nBUS+1;i++)if(sBus[i-1].Type==2){for(j=1;j<nBUS+1;j++){ps+=V[i][1]*V[j][1]*(YG[i][j]*cos(V[i][0]-V[j][0])+YB[i][j]*sin(V[i][0 ]-V[j][0]));qs+=V[i][1]*V[j][1]*(YG[i][j]*sin(V[i][0]-V[j][0])-YB[i][j]*cos(V[i][0 ]-V[j][0]));}}}for(i=1;i<nBUS+1;i++){for(j=1;j<nBUS+1;j++){sP[i][j]=0;sQ[i][j]=0;}}for(l=0; l<nL; l++)i=sL[l].NumI;j=sL[l].NumJ;r=sL[l].R;x=sL[l].X;d1=r*r+x*x;g=r/d1;b=-x/d1;if(fabs(sL[l].<{/*Normal lines or transformers*/sP[i][j]=V[i][1]*V[i][1]*g-V[i][1]*V[j][1]*(g*cos(V[i][0]-V[j][0])+b*s in(V[i][0]-V[j][0]));sQ[i][j]=-(V[i][1]*V[i][1]*sL[l].B+V[i][1]*V[i][1]*b+V[i][1]*V[j][1]*( g*sin(V[i][0]-V[j][0])-b*cos(V[i][0]-V[j][0])));sP[j][i]=V[j][1]*V[j][1]*g-V[i][1]*V[j][1]*(g*cos(V[j][0]-V[i][0])+b*s in(V[j][0]-V[i][0]));sQ[j][i]=-(V[j][1]*V[j][1]*sL[l].B+V[j][1]*V[j][1]*b+V[i][1]*V[j][1]*( g*sin(V[j][0]-V[i][0])-b*cos(V[j][0]-V[i][0])));}else{/*abnormal transformer ratio*/sP[i][j]=V[i][1]*V[i][1]*g/sL[l].B/sL[l].B-V[i][1]*V[j][1]*(g*cos(V[i] [0]-V[j][0])/sL[l].B+b*sin(V[i][0]-V[j][0])/sL[l].B);sQ[i][j]=-(V[i][1]*V[i][1]*b/sL[l].B/sL[l].B+V[i][1]*V[j][1]*(g*sin(V[ i][0]-V[j][0])/sL[l].B-b*cos(V[i][0]-V[j][0])/sL[l].B));sP[j][i]=V[j][1]*V[j][1]*g-V[i][1]*V[j][1]*(g*cos(V[j][0]-V[i][0])/sL[ l].B+b*sin(V[j][0]-V[i][0])/sL[l].B);sQ[j][i]=-(V[i][1]*V[i][1]*b+V[i][1]*V[j][1]*(g*sin(V[j][0]-V[i][0])/s L[l].B-b*cos(V[j][0]-V[i][0])/sL[l].B));}}dsp=ps;dsq=qs;sumgen=ps;for(i=0;i<nBUS;i++){dsp+=sBus[i].GenP-sBus[i].LoadP;dsq+=sBus[i].GenQ-sBus[i].LoadQ;sumgen+=sBus[i].GenP;}dp=dsp/sumgen*100;/* 输出功率情况 */if((fp=fopen("功率情况.txt","w"))==NULL){printf("Can not open the file named '功率情况.txt' \n");exit(0);}fprintf(fp,"---功率情况---\n");fprintf(fp,"平衡节点功率S=%+ j%\n",ps,qs);for(i=1;i<nBUS+1;i++)for(j=1;j<nBUS+1;j++)if(fabs(sP[i][j]>fprintf(fp,"S(%3d,%-3d)=(%,j%\n",i,j,sP[i][j],sQ[i][j]); fprintf(fp,"网损为%+j%,输电效率为%\n",dsp,dsq,100-dp);if(fp!=NULL) fclose(fp);/* 输出最终结果 */if((fp=fopen("最终结果.txt","w"))==NULL){printf("Can not open the file named '最终结果.txt' \n");exit(0);}fprintf(fp,"---最终结果---\n");fprintf(fp,"\n\n\n");fprintf(fp,"一共循环迭代计算%d次\n",c);fprintf(fp,"\n\n\n");fprintf(fp,"---电压终值---\n");for(i=1;i<nBUS+1;i++)for(j=1;j<2;j++)fprintf(fp,"U(%2d )=(%,%\n",i,V[i][0]*180/,V[i][1]);fprintf(fp,"\n\n\n");fprintf(fp,"平衡节点功率S=%+j%\n",ps,qs);fprintf(fp,"\n\n\n");fprintf(fp,"---各节点间功率---\n");for(i=1;i<nBUS+1;i++)for(j=1;j<nBUS+1;j++)if(fabs(sP[i][j]>fprintf(fp,"S(%2d,%-2d)=(%,j%\n",i,j,sP[i][j],sQ[i][j]); fprintf(fp,"\n\n\n");fprintf(fp,"网损为：%+j%。