现代电力系统潮流计算作业

《现代电力系统潮流计算作业》由会员分享，可在线阅读，更多相关《现代电力系统潮流计算作业（22页珍藏版）》请在装配图网上搜索。

1、现代电力系统潮流计算作业0 序章作业要求（A组）：0.1 调用matpower中的runpf函数，分析输入文件中各矩阵定义；0.2 调用某一个算例，输出潮流结果，并分析。0.3 完成0.1和0.2的基础上，分析matpower中牛顿法和快速解耦法，给出流程图，写出newtonpf和fdpf函数每行程序定义。0.4 完成0.3的基础上，制造一个病态潮流算例，并跟踪调试，分析病态原因。1 分析输入文件中各矩阵的定义1.1 MATPOWER的安装MATPOWER工具箱的安装步骤如下：1） 下载matpower压缩包。官方下载网址：http:/www.pserc.cornell.edu/matpowe

2、r/，目前最新版本为6.0b1，稳定版本为5.1，建议下载稳定版本。2） 解压压缩包，得到文件夹matpower5.1，并将文件夹移动到MATLAB所在路径的toolbox文件夹下。我的路径为： C:Program FilesMATLABR2016atoolbox。3） 添加地址到MATLAB路径。打开MATLAB，点击“文件”“设置路径”“添加并包含子文件夹”，找到matpower5.1所在的位置，点击“确定”，再点“保存”“关闭”。4） 测试matpower工具是否安装成功。在MATLAB命令行窗口输入“test_matpower”，出现一系列的测试，均显示“ok”，最后显示“All te

3、sts successful (3256 passed, 682 skipped of 3938)”，则表示安装成功。1.2 矩阵的定义 打开文档“caseformat.m”，或者在MATLAB命令行窗口中输入“help caseformat”，可以得到关于输入矩阵的数据定义。当然，也可以参考docs文件夹下的manual文档，其中对matpower工具箱进行了详细说明。在matpower中，输入矩阵至少包含三种：母线参数矩阵（Bus Data），发电机参数矩阵（Generator Data），支路参数矩阵（Branch Data）。为了进行最优潮流的相关计算，输入矩阵还包含发电机费用参数矩阵

4、（generator cost data）。以下对三种基本的输入参数矩阵数据格式进行详细说明。表1.1 母线参数矩阵主要数据格式说明列名称含义单位1bus number母线编号，从1开始，正整数2bus type节点类型，1为PQ节点，2为PV节点，3为平衡(参考)节点，4为孤立节点3Pd注入负荷的有功功率MW4Qd注入负荷的无功功率MVar5Gs与母线并联的电导MW6Bs与母线并联的电纳MVar7area number电网断面号，一般设置为18Vm母线电压的幅值p.u.9Va母线电压的相角度10baseKV基准电压kV11zone损耗分区号，一般设置为112Vmax工作时，母线的最高电压幅值

5、p.u.13Vmin工作时，母线的最低电压幅值p.u.表1.2 发电机参数矩阵主要数据格式说明列名称含义单位1bus发电机所在母线的编号2Pg接入发电机的有功功率MW3Qg接入发电机的无功功率MVar4Qmax发电机的最大输出无功功率MVar5Qmin发电机的最小输出无功功率MVar6Vg发电机的工作电压幅值p.u.7mBase发电机的功率基准值，默认为baseMVAMVA8status发电机是否投入运行，大于0表示投入运行，小于或等于0表示退出运行9Pmax发电机的最大输出有功功率MW10Pmin发电机的最小输出有功功率MW表1.3 支路参数矩阵主要数据格式说明列名称含义单位1fbus起始母

6、线编号2tbus终止母线编号3r支路电阻p.u.4x支路电抗p.u.5b支路充电电纳p.u.6rateA支路长期运行允许的功率MVA7rateB支路短期运行允许的功率MVA8rateC支路紧急运行允许的功率MVA9ratio支路上变压器的变比，如果支路元件不是变压器，则该值为010angle支路上变压器的转角，如果支路元件不是变压器，则该值为0度11status支路的初始工作状态，1表示投入运行，0表示退出运行12angmin支路最小相角差度13angmax支路最大相角差度1.3 case9数据分析根据以上分析，打开一个算例，比如默认的case9，进行分析。算例case9.m文件包含两个变量和

7、四个矩阵。其中baseMVA=100，表示功率的基准值为100MVA。三个基本的矩阵定义如下。表1.4 算例case9的母线参数矩阵12345678910111213bus_itypePdQdGsBsareaVmVabaseKVzoneVmaxVmin13000011034511.10.922000011034511.10.932000011034511.10.941000011034511.10.95190300011034511.10.961000011034511.10.971100350011034511.10.981000011034511.10.991125500011034511

8、.10.9表1.5 算例case9的发电机参数矩阵123456789101112busPgQgQmaxQminVgmBasestatusPmaxPminPc1Pc2100300-30011001250100021630300-3001100130010003850300-300110012701000131415161718192021Qc1minQc1maxQc2minQc2maxramp_agcramp_10ramp_30ramp_qapf000000000000000000000000000表1.6 算例case9的支路参数矩阵12345678910111213fbustbusrxbra

9、teArateBrateCratioanglestatusangminangmax1400.05760250250250001-360360450.0170.0920.158250250250001-360360560.0390.170.358150150150001-3603603600.05860300300300001-360360670.01190.10080.209150150150001-360360780.00850.0720.149250250250001-3603608200.06250250250250001-360360890.0320.1610.306250250250

10、001-360360940.010.0850.176250250250001-360360根据参数矩阵，可以推测出case9的电力系统单线图，如图1.1所示。该系统是一个环形网络，包含三个带有发电机的母线，其中母线1是平衡节点，母线2和3均为PV节点，其他的母线都是PQ节点。所有的母线电压初始幅值均设置为1（p.u.），相角为0度，电压基准值为345kV。系统包含3个负荷，分别是母线5上的负荷为90+j30（MVA），母线7上的负荷为100+j35（MVA），母线9上的负荷为125+j50（MVA）。支路1-4，3-6，8-2只有电抗值，电阻和电纳均为0，可以推测该支路为变压器支路的等效。图1

11、.1 系统单线图2 计算潮流并分析2.1 调用runpf计算case9系统的潮流在MATLAB命令行窗口输入“runpf”或“runpf(case9)”，或者直接运行“runpf.m”，得到case9系统的潮流计算结果：MATPOWER Version 5.1, 20-Mar-2015 - AC Power Flow (Newton)Newtons method power flow converged in 4 iterations.Converged in 0.03 seconds=| System Summary |=How many? How much? P (MW) Q (MVAr)

12、- - - -Buses 9 Total Gen Capacity 820.0 -900.0 to 900.0Generators 3 On-line Capacity 820.0 -900.0 to 900.0Committed Gens 3 Generation (actual) 320.0 34.9Loads 3 Load 315.0 115.0 Fixed 3 Fixed 315.0 115.0 Dispatchable 0 Dispatchable -0.0 of -0.0 -0.0Shunts 0 Shunt (inj) -0.0 0.0Branches 9 Losses (I2

13、* Z) 4.95 51.31Transformers 0 Branch Charging (inj) - 131.4Inter-ties 0 Total Inter-tie Flow 0.0 0.0Areas 1 Minimum Maximum - -Voltage Magnitude 0.958 p.u. bus 9 1.003 p.u. bus 6 Voltage Angle -4.35 deg bus 9 9.67 deg bus 2 P Losses (I2*R) - 2.46 MW line 8-9Q Losses (I2*X) - 16.74 MVAr line 8-2=| Bu

14、s Data |= Bus Voltage Generation Load # Mag(pu) Ang(deg) P (MW) Q (MVAr) P (MW) Q (MVAr)- - - - - - - 1 1.000 0.000* 71.95 24.07 - - 2 1.000 9.669 163.00 14.46 - - 3 1.000 4.771 85.00 -3.65 - - 4 0.987 -2.407 - - - - 5 0.975 -4.017 - - 90.00 30.00 6 1.003 1.926 - - - - 7 0.986 0.622 - - 100.00 35.00

15、 8 0.996 3.799 - - - - 9 0.958 -4.350 - - 125.00 50.00 - - - - Total: 319.95 34.88 315.00 115.00=| Branch Data |=Brnch From To From Bus Injection To Bus Injection Loss (I2 * Z) # Bus Bus P (MW) Q (MVAr) P (MW) Q (MVAr) P (MW) Q (MVAr)- - - - - - - - - 1 1 4 71.95 24.07 -71.95 -20.75 -0.000 3.32 2 4

16、5 30.73 -0.59 -30.55 -13.69 0.174 0.94 3 5 6 -59.45 -16.31 60.89 -12.43 1.449 6.31 4 3 6 85.00 -3.65 -85.00 7.89 0.000 4.24 5 6 7 24.11 4.54 -24.01 -24.40 0.095 0.81 6 7 8 -75.99 -10.60 76.50 0.26 0.506 4.29 7 8 2 -163.00 2.28 163.00 14.46 0.000 16.74 8 8 9 86.50 -2.53 -84.04 -14.28 2.465 12.40 9 9

17、4 -40.96 -35.72 41.23 21.34 0.266 2.26 - - Total: 4.955 51.312.2 潮流计算结果分析Matpower工具箱的潮流计算结果由四部分组成：程序运行信息，系统概述，母线数据，支路数据。其中，程序运行信号包含潮流计算类型，使用的迭代算法，迭代次数，所用时间。本次潮流计算是交流潮流计算，采用了Newton法，迭代了4次得到了符合精度要求的结果，耗时0.03s。系统概述描述了系统的基本信息。包括系统元件的数量，元件的功率大小，电压和功率损耗的极值。如case9系统包含9个节点，3个发电机，3个负荷，9条支路。总装机容量820MW，在线容量820

18、MW，实际发电 320MW，负荷消耗有功315MW，总网损4.95MW。母线9上电压幅值最小：0.958（p.u.），电压相角也最小：-4.35，母线6上电压幅值最大：1.003（p.u.），母线2上相角最大：9.67，支路8-9上消耗了最多的有功功率：2.46MW。母线数据包含母线电压结果，发电机输出功率，负荷消耗功率，累计功率。表2.1 潮流计算结果母线数据母线编号电压幅值pu电压相角deg发电机有功MW发电机无功MVar负荷消耗有功MW负荷消耗无功MVar11.0000.000*71.9524.07-21.0009.669163.0014.46-31.0004.77185.00-3.65

19、-40.987-2.407-50.975-4.017-90.0030.0061.0031.926-70.9860.622-100.0035.0080.9963.799-90.958-4.350-125.0050.00累计：319.9534.88315.00115.00支路数据包含起始母线注入功率、终止母线注入功率和支路上的功率损耗。表2.1 潮流计算结果支路数据支路编号起始母线终止母线起始母线注入有功MW起始母线注入无功MVar终止母线注入有功MW终止母线注入无功MVar支路有功损耗MW支路无功损耗MVar11471.9524.07-71.95-20.75-0.0003.3224530.73-

20、0.59-30.55-13.690.1740.94356-59.45-16.3160.89-12.431.4496.3143685.00-3.65-85.007.890.0004.2456724.114.54-24.01-24.400.0950.81678-75.99-10.6076.500.260.5064.29782-163.002.28163.0014.460.00016.7488986.50-2.53-84.04-14.282.46512.40994-40.96-35.7241.2321.340.2662.26累计损耗：4.955 51.31通过以上数据可以发现系统的有功功率守恒：发电

21、机发出的总有功功率为319.95MW，负荷消耗的总有功功率为315MW，网损为4.955MW。另外可以得到系统的潮流分布图如图2.1所示。图2.1 系统潮流分布图3 迭代算法分析3.1 牛顿法分析打开newtonpf.m文档，可以看到matpower的牛顿法的介绍和代码。函数的输入参数包含系统的节点导纳矩阵，节点的注入复功率，初始电压，平衡节点、PV节点和PQ节点的标号列向量，以及包含终止误差、最大迭代次数和输出选项的结构体。返回节点电压，收敛标志和迭代次数。通过分析可以得到matpower的牛顿法的程序流程图，如图3.1所示，这和一般的牛顿法潮流计算程序并没有什么区别。图3.1 牛顿法潮流计

22、算程序流程图以下是newtonpf函数的每行程序的定义。function V, converged, i = newtonpf(Ybus, Sbus, V0, ref, pv, pq, mpopt)%NEWTONPF 使用完整的牛顿法求解潮流% V, CONVERGED, I = NEWTONPF(YBUS, SBUS, V0, REF, PV, PQ, MPOPT)% 通过分别给定完整系统的导纳矩阵（针对所有节点），节点的注入% 复功率（针对所有节点），节点电压的初始值，和平衡节点、PV节点和PQ节% 点标号的列向量，求解节点电压。节点电压矢量包含发电机节点（包括平衡% 节点）的设定值和平衡

23、节点的参考角度，以及幅度的大小和角度的初始值。% MPOPT是一个MATPOWER选项结构体，可用于设置终止误差限，最大迭代次数和% 输出选项（有关详细信息，请参阅MPOPTION）。如果未指定此参数，则使用% 默认选项。最终返回节点电压相量，收敛标志以及迭代次数。% 参考RUNPF.% 缺省参数设置if nargin 1 % 如果该标志位大于1，则保存进度信息到文档 fprintf(n it max P & Q mismatch (p.u.); fprintf(n- -); fprintf(n%3d %10.3e, i, normF);end% 将进度信息输出到文档中if normF 1 %

24、 如果该标志位大于1，则保存进度信息到文档 fprintf(nConverged!n);% 将该字符串输出到指定文档中 endend% 进行牛顿迭代while (converged & i Vm Va = angle(V); % 所有节点电压相位-Va % 计算修正方程式的常数项 mis = V .* conj(Ybus * V) - Sbus; % 计算误差 F = real(mis(pv); % PV节点的有功功率误差 real(mis(pq); % PQ节点的有功功率误差 imag(mis(pq) ; % PQ节点的无功功率误差 % 校验收敛性 normF = norm(F, inf);

25、 % F的无穷范数（等效于取最大值） if mpopt.verbose 1 % 如果该标志位大于1，则保存进度信息到文档 fprintf(n%3d %10.3e, i, normF); % 将进度信息存到文档中 end if normF tol % 误差是否小于误差限 converged = 1; % 收敛标志置1 if mpopt.verbose % 如果该标志位不等于0，则保存最终信息到文档 fprintf(nNewtons method power flow converged in %d iterations.n, i); end % 将该字符串输出到指定文档中 endendif mp

26、opt.verbose % 如果该标志位不等于0，则保存最终信息到文档 if converged % 如果收敛标志为0，即超出了迭代次数，则输出不收敛的信息到文档 fprintf(nNewtons method power flow did not converge in %d iterations.n, i); endend3.2 快速解耦法分析在matpower中，快速解耦法分为XB型和BX型。调用方式如下：runpf(case9,mpoption(pf.alg,FDXB); % XB型runpf(case9,mpoption(pf.alg,FDBX); % BX型快速解耦法、牛顿法和高斯

27、-塞德尔法的比较如表3.1所示，对于不同节点的算例，牛顿法和快速解耦法都具有很强的收敛性，迭代次数并没有明显的增加，而高斯-塞德尔法收敛性不佳，迭代次数增加明显，甚至超过了预设的最大迭代次数。表3.1 各种潮流计算方法的比较算例节点数牛顿法FDXBFDBX高斯-塞德尔法948-77-72122448-89-85355737-79-851811838-77-6NC30059-99-8NC接下来对快速解耦法进行分析，得到流程图如图3.2所示。图3.2 快速解耦法潮流计算程序流程图以下是fdpf函数的每行程序的定义。function V, converged, i = fdpf(Ybus, Sbus

28、, V0, Bp, Bpp, ref, pv, pq, mpopt)%FDPF 使用快速解耦法求解电力系统潮流.% 通过分别给定完整系统的导纳矩阵（针对所有节点），节点的注入% 复功率（针对所有节点），节点电压的初始值，和平衡节点、PV节点和PQ节% 点标号的列向量，求解节点电压。节点电压矢量包含发电机节点（包括平衡% 节点）的设定值和平衡节点的参考角度，以及幅度的大小和角度的初始值。% MPOPT是一个MATPOWER选项结构体，可用于设置终止误差限，最大迭代次数和% 输出选项（有关详细信息，请参阅MPOPTION）。如果未指定此参数，则使用% 默认选项。最终返回节点电压相量，收敛标志以及迭

29、代次数。% 参考RUNPF.% 缺省参数设置if nargin 7 % 如果输入参数少于7项 mpopt = mpoption; % 则设置mpopt的缺省值为mpoptionend% 求解选项tol = mpopt.pf.tol; % 终止误差限max_it = mpopt.pf.fd.max_it; % 最大迭代次数if have_fcn(matlab) & have_fcn(matlab, vnum) 1 % 如果该标志位大于1，则保存进度信息到文档 fprintf(niteration max mismatch (p.u.) ); fprintf(ntype # P Q ); fpri

30、ntf(n- - - -); fprintf(n - %3d %10.3e %10.3e, i, normP, normQ);end % 将进度信息输出到文档中if normP tol & normQ 1 % 如果该标志位大于1，则保存进度信息到文档 fprintf(nConverged!n); % 将该字符串输出到指定文档中 endend% 降低B矩阵维度Bp = Bp(pv; pq, pv; pq); % 通过索引，降低Bp的维度Bpp = Bpp(pq, pq); % 通过索引，降低Bpp的维度% 矩阵B的LU分解if lu_vec % 选择系统支持的lu语句 Lp, Up, pp, q

31、p = lu(Bp, vector); % 矩阵Bp的LU分解 Lpp, Upp, ppp, qpp = lu(Bpp, vector);% 矩阵Bpp的LU分解 junk, iqp = sort(qp);% 行向量qp的排序 junk, iqpp = sort(qpp); % 行向量qpp的排序 % , iqp = sort(qp); % junk没有用到，用代替 % , iqpp = sort(qpp);% junk没有用到，用代替else Lp, Up, Pp = lu(Bp); % 矩阵Bp的LU分解 Lpp, Upp, Ppp = lu(Bpp);% 矩阵Bpp的LU分解end% 进

32、行P和Q迭代while (converged & i 1 % 如果该标志位大于1，则保存进度信息到文档 fprintf(n P %3d %10.3e %10.3e, i, normP, normQ); end % 将进度信息输出到文档中 if normP tol & normQ 1 % 如果该标志位大于1，则保存进度信息到文档 fprintf(n Q %3d %10.3e %10.3e, i, normP, normQ); end % 将进度信息输出到文档中 if normP tol & normQ tol % 误差是否小于误差限 converged = 1; % 收敛标志置1 if mpop

33、t.verbose % 如果该标志位不等于0，则保存最终信息到文档 fprintf(nFast-decoupled power flow converged in %d P-iterations and %d Q-iterations.n, i, i); end % 将该字符串输出到指定文档中 break; % 跳出循环，结束迭代 endendif mpopt.verbose % 如果该标志位不等于0，则保存最终信息到文档 if converged % 如果收敛标志为0，即超出了迭代次数，则输出不收敛的信息到文档 fprintf(nFast-decoupled power flow did n

34、ot converge in %d iterations.n, i); endend4 病态潮流分析4.1 重负荷导致的病态潮流如图4.1所示，将原9节点的系统9#母线上的负荷由原来的125+j50改为125+j300，便构成了一个病态问题。如表4.1所示，逐渐增大9#母线上的负荷，可以清楚地发现节点电压的幅值迅速减小，相位逐渐滞后，牛顿法迭代次数迅速增大，直到问题不收敛，雅克比矩阵的条件数也迅速增大。经过仿真分析，可以发现，雅克比矩阵的条件数最大达到了，所以修正方程是一个病态方程。该问题在给定的运行条件下，潮流方程无解，如果想让系统有解，则必须调整系统的参数，使得系统回到有解区域。表4.1

35、9节点系统9#母线负荷对系统潮流计算的影响负荷（MVA）节点电压(p.u.)牛顿法迭代次数雅克比矩阵最大条件数125+j500.958-4.3504101125+j1000.907-4.4084102125+j1500.850-4.5514104125+j2000.780-4.8555107125+j2500.684-5.5795111125+j2800.574-7.0107165125+j300-迭代200次不收敛图4.1 重负荷导致的病态问题系统单线图通过设置断点、单步调试和查看变量等方式，可以清楚地看到9#节点重负荷情况下，系统潮流的计算过程。如表4.2所示。第7次迭代时，雅克比矩阵接近

36、奇异，条件数急剧增大，节点电压也变得非常大，问题无解。表4.2 9节点系统9#母线电压随迭代次数的变化迭代次数9#节点电压（p.u.）雅克比矩阵的条件数010-10.754-2.96599.0520.607-5.548108.2830.491-7.970116.6240.787-1.644472.3450.628-4.950107.1760.515-7.425114.6974.75259.7388464.382.55830.6951596.991.416173.7279753.4100.626158.7832367.7110.12573.611276170121.151-120.7457697.4130.52417