《电力系统潮流计算matlab辅助程序设计-潮流计算程序》由会员分享，可在线阅读，更多相关《电力系统潮流计算matlab辅助程序设计-潮流计算程序（21页珍藏版）》请在金锄头文库上搜索。
1、电力系统潮流计算的MATLAB辅助程序设计潮流计算，通常指负荷潮流，是电力系统分析和设计的主要组成部分，对系统规划、安全运行、经济调度和电力公司的功率交换非常重要。此外，潮流计算还是其它电力系统分析的基础，比如暂态稳定，突发事件处理等。现代电力系统潮流计算的方法主要:高斯法、牛顿法、快速解耦法和MATLAB的M语言编写的MATPOWER4.1，这里主要介绍高斯法、牛顿法和快速解耦法。高斯法的程序是lfgauss，其与lfybus、busout和lineflow程序联合使用求解潮流功率。lfybus、busout和lineflow程序也可与牛顿法的lfnewton程序和快速解耦法的decoupl
2、e程序联合使用。(读者可以到MATPOWER主页下载MATPOWER4.1,然后将其解压到MATLAB目录下，即可使用该软件进行潮流计算)一、高斯-赛德尔法潮流计算使用的程序： 高斯-赛德法的具体使用方法读者可参考后面的实例，这里仅介绍各程序的编写格式:lfgauss：该程序是用高斯法对实际电力系统进行潮流计算，需要用到busdata和linedata两个文件。程序设计为输入负荷和发电机的有功MW和无功Mvar，以及节点电压标幺值和相角的角度值。根据所选复功率为基准值将负荷和发电机的功率转换为标幺值。对于PV节点，如发电机节点，要提供一个无功功率限定值。当给定电压过高或过低时，无功功率可能超出
3、功率限定值。在几次迭代之后(高斯-塞德尔迭代为10次)，需要检查一次发电机节点的无功出力，如果接近限定值，电压幅值进行上下5%的调整，使得无功保持在限定值内。lfybus：这个程序需要输入线路参数、变压器参数以及变压器分接头参数。并将这些参数放在名为linedata的文件中。这个程序将阻抗转换为导纳，并得到节点导纳矩阵。busout：该程序以表格形式输出结果，节点输出包括电压幅值和相角，发电机和负荷的有功和无功功率，以及并联电容器或电抗器的有功和无功功率。lineflow：该程序输出线路的相关数据，程序设计输出流入线路终端的有功和无功的功率、线损以及节点功率，还包含整个系统的有功和无功损耗。l
4、fnewton是牛顿-拉夫逊法对实际电力系统潮流计算开发的程序，数据准备和程序格式和高斯-赛德尔法一样，包括程序lfybus，busout和 lineflow。 decouple是快速解耦法对实际电力系统潮流计算开发的程序，同高斯法和牛顿法一样需要用到三个程序：lfybus、busout、lineflow。 二、数据准备为了在MATLAB环境下用高斯法进行潮流计算，必须定义下列变量：基准功率，功率允许误差，加速因子和最大迭代次数。上述变量命名(小写字母)为：basemva、accuracy、accel和maxiter，一般规定为：basemva=100； accuracy=0.001；acce
5、l=1.6；maxiter=80；输入文件准备的第一步是给节点编号，节点号码必须是连续的，但节点数据输入不一定按顺序来编写。此外，还需要下列数据文件：1.节点数据文件busdata：节点信息输入格式为单行输入，输入的数据形成一个矩阵，叫做busdata矩阵。第一列为节点号；第二列为节点类型；第三列和第四列分别为节点电压幅值(标幺值)和相角(单位为度)；第五列和第六列分别为负荷的有功功率和无功功率；第七列到十列分别为发电机的有功功率、无功功率、最小无功出力和最大无功出力；最后一列为并联电容器注入无功功率。第二列的编码用0、1、2来区分PQ节点、平衡节点和PV节点：0表示PQ节点，输入正的有功功率
6、（MW）和无功功率（Mvar），并且要设定节点电压初始估计值，一般幅值和相角分别设为1和0，若已经给定初始值，则用其给定值来代替1和0。1表示平衡节点，且已知该节点的电压幅值和相角。2表示PV节点，要设定该节点的节点电压幅值和发电机的有功功率（MW），并设定发电机的无功最小出力和最大出力（Mvar）。2.线路数据文件linedata 线路数据用节点对的方法来确定，数据包含在称为linedata的矩阵中。第一列和第二列为节点号码，第三列到第五列为线路电阻、电抗及该线路电纳值的一半，以标幺值表示。最后一列为变压器分接头设定值，对线路来说，需要输入1。线路输入为无输入顺序，对变压器来说，左侧的节点号
7、设为分接头端。3.zdata是线路数据输入变量，包括四项，前两项是节点编号，后两项是线路电阻和电抗，均以标幺值表示，函数返回节点导纳矩阵。三、潮流计算的MATLAB程序清单1. lfgauss.m程序清单 % Power flow solution by Gauss-Seidel methodVm=0; delta=0; yload=0; deltad =0;nbus = length(busdata(:,1);kb=;Vm=; delta=; Pd=; Qd=; Pg=; Qg=; Qmin=; Qmax=; Pk=; P=; Qk=; Q=; S=; V=; for k=1:nbusn=b
8、usdata(k,1);kb(n)=busdata(k,2); Vm(n)=busdata(k,3); delta(n)=busdata(k, 4);Pd(n)=busdata(k,5); Qd(n)=busdata(k,6); Pg(n)=busdata(k,7); Qg(n) = busdata(k,8);Qmin(n)=busdata(k, 9); Qmax(n)=busdata(k, 10);Qsh(n)=busdata(k, 11); if Vm(n) = accuracy & iter = maxiteriter=iter+1; for n = 1:nbus; YV = 0+j*0
9、; for L = 1:nbr; if (nl(L) = n & mline(L) = 1), k=nr(L); YV = YV + Ybus(n,k)*V(k); elseif (nr(L) = n & mline(L)=1), k=nl(L); YV = YV + Ybus(n,k)*V(k); end end Sc = conj(V(n)*(Ybus(n,n)*V(n) + YV) ; Sc = conj(Sc); DP(n) = P(n) - real(Sc); DQ(n) = Q(n) - imag(Sc); if kb(n) = 1 S(n) =Sc; P(n) = real(Sc
10、); Q(n) = imag(Sc); DP(n) =0; DQ(n)=0; Vc(n) = V(n); elseif kb(n) = 2 Q(n) = imag(Sc); S(n) = P(n) + j*Q(n); if Qmax(n) = 0 Qgc = Q(n)*basemva + Qd(n) - Qsh(n); if abs(DQ(n) = 10 if DV(n) = 0.045 if Qgc Qmax(n), Vm(n) = Vm(n) - 0.005; DV(n)=DV(n)+.005; end else, end else,end else,end end if kb(n) =
11、1 Vc(n) = (conj(S(n)/conj(V(n) - YV )/ Ybus(n,n); else, end if kb(n) = 0 V(n) = V(n) + accel*(Vc(n)-V(n); elseif kb(n) = 2 VcI = imag(Vc(n); VcR = sqrt(Vm(n)2 - VcI2); Vc(n) = VcR + j*VcI; V(n) = V(n) + accel*(Vc(n) -V(n); end end maxerror=max( max(abs(real(DP), max(abs(imag(DQ) ); if iter = maxiter & maxerror accuracy
