复杂网络 N-R 法潮流分析与计算的设计.doc

目 录摘要21.设计目的与要求31.1 设计目的31.2 设计要求31.3 设计题目42.题目解析42.1 设计思路52.2 基本公式和变量分类52.3电力系统潮流计算计算结果72.3.1节点导纳矩阵72.3.2简化雅可比矩阵82.3.3修正、迭代93.程序设计93.1 MATLAB编程说明93.2潮流计算流程图103.3潮流计算源程序104.小结14参考文献14摘 要关键词： 电力系统 潮流计算 牛顿拉夫逊算法在如今的社会，电力已经成为人们必不可少的需求，而建立结构合理的大型电力系统不仅便于电能生产与消费的集中管理、统一调度和分配，减少总装机容量，节省动力设施投资，且有利于地区能源资源的合理开发利用，更大限度地满足地区国民经济日益增长的用电需要。电力系统建设往往是国家及地区国民经济发展规划的重要组成部分。电力系统的出现，使高效、无污染、使用方便、易于调控的电能得到广泛应用，推动了社会生产各个领域的变化，开创了电力时代，发生了第二次技术革命。电力系统的规模和技术水准已成为一个国家经济发展水平的标志之一。电力系统稳态分析包括潮流计算（或潮流分析)和静态安全分析。潮流计算针对电力革统各正常运行方式，而静态安全分析则要研究各种运行方式下个别系统元件退出运行后系统的状况。其目的是校验系统是否能安全运行，即是否有过负荷的元件或电压过低的母线等。原则上讲，静态安全分析也可U用潮流计算来代替。但是一般静态安全分析需要校验的状态数非常多，用严格的潮流计算来分析这些状态往往计算量过大，因此不得不寻求一些特殊的算法以满足要求。潮流计算在数学上是多元非线性方程组的求解问题，求解的方法有很多种。牛顿拉夫逊法是数学上解非线性方程式的有效方法，有较好的收敛性。将牛顿法用于潮流计算是以导纳矩阵为基础的，由于利用了导纳矩阵的对称性、稀疏性及节点编号顺序优化等技巧，使牛顿法在收敛性、占用内存、计算速度等方面都达到了一定的要求。1.设计目的与要求1.1 设计目的1. 掌握电力系统潮流计算的基本原理；2. 潮流计算是电力系统分析中的一种最基本的计算，他的任务是对给定运行条件确定系统运行状态，如各母线上的电压（幅值及相角）、网络中的功率分布及功率损耗等。具体表现在以下方面： (1)在电网规划阶段,通过潮流计算,合理规划电源容量及接入点,合理规划网架,选择无功补偿 方案,满足规划水平的大、小方式下潮流交换控制、调峰、调相、调压的要求。 (2)在编制年运行方式时,在预计负荷增长及新设备投运基础上,选择典型方式进行潮流计算,发 现电网中薄弱环节,供调度员日常调度控制参考,并对规划、基建部门提出改进网架结构,加快基建进度的建议。 (3)正常检修及特殊运行方式下的潮流计算,用于日运行方式的编制,指导发电厂开机方式,有 功、无功调整方案及负荷调整方案,满足线路、变压器热稳定要求及电压质量要求。 (4)预想事故、设备退出运行对静态安全的影响分析及作出预想的运行方式调整方案。 3. 掌握并能熟练运用一门计算机语言（MATLAB语言或FORTRAN或C语言或C+语言）；4. 采用计算机语言对潮流计算进行计算机编程计算。 1.2设计要求1）根据给定的运行条件，确定图中电力系统潮流计算时各节点的类型、待求量；2）求节点导纳矩阵；3) 赋予各节点电压变量初值后，求解不平衡量；4）形成雅克比矩阵；5）给出潮流方程或功率方程的表达式；6）当用牛顿拉夫逊法计算潮流时，给出修正方程和迭代收敛条件；7）上机编程调试；8）计算分析给定系统潮流分析并与手工计算结果作比较分析。1.3 设计题目2.题目解析2.1设计思路首先写出节点导纳矩阵，并分析各节点的类型，找出待求量。然后，确定潮流方程。最后进行潮流计算，而最后一步，可利用牛顿拉夫逊法潮流分析。其次此电力系统是一个5节点，4支路的电力网络。其中包含3个PQ节点，一个PV节点，和一个平衡节点。由于牛顿拉夫逊法所需解的方程组最少，最后确定采用牛顿拉夫逊法。2.2基本公式和变量分类（1）.节点电压U和节点导纳矩阵Y。（2）.变量分类。在潮流问题中，任何复杂的电力网和电力系统都可以归结为以下元件（参数）组成。 1）.发电机（注入电流或功率）。 2）.负载（负的注入电流或功率）。 3）.输电线支路（电抗、电阻）。 4）.变压器支路（电阻、电抗、变化）。 5）.变压器对地支路（导纳和感纳，本例中忽略）。 6）.母线上的对地支路（阻抗或导纳，本例中忽略）。 7）.线路上的对地支路（一般为线路电容导纳）。（3）.节点功率方程。有节点电压方程展开节点功率方程，此为电力系统的潮流方程的一般形式： （1-1）令，将此式带入（1-1），将实部、虚部分开，得到直角坐标形式的节点功率方程： （1-2） 潮流方程具有的特点是：他能表征电力系统稳态运行特性；其为一组非线性方程，只能用迭代方法求其数值解；方程中的电压U和导纳Y即可表示为直角坐标，又可表示为极坐标。因而潮流方程有多种表达方式极坐标形式、直角坐标形式和混合坐标形式。潮流计算的约束条件，即电压、相角和功率的约束条件。同时，为电压和支路导纳 。 但（1-2）是节点电压的线性方程组，分别为注入有功功率和无功功率。所以只要找出节点电压，这两个物理量，就能算出功率误差，。对于PV节点，电压已知，则： （1-3）所以对于n个节点的网络，其中（m-1）个PQ节点，一个不平衡点，其不平衡方程数目：，i=1，2，n，则有（n-1）个方程。，i=1，2，n，则有（m-1）个方程。，i=m+1，m+2，n，则有（n-m）个方程。 (4)修正方程：形成雅可比矩阵。N-R法的思想是；本例；对F(x)求偏导的式（1-2）、式（1-3），即式（1-2）、式（1-3）中的、是多维变量的函数，对多维变量求偏导（、），并以矩阵的形式表达称为雅可比矩阵。 当j=i时，对角元素为 （1-4） 当时,矩阵非对角元素为： （1-5） 由上式不难看出，雅可比矩阵有以下特点。 雅可比矩阵中的诸元素都是节点电压的函数，因此在迭代过程中，它们将随着节点电压的变化而不断的变化。 雅可比矩阵具有结构对称性，数据不对称。如非对角，。 由式（1-7）可以看出，当导纳矩阵中非对角元素为零时，。雅可比矩阵中相应的元素也为零，即矩阵是非常稀疏的。因此，修正方程的求解同样可以应用稀疏矩阵的求解技巧。正是由于这一点才使N-R法获得广泛的应用。2.3电力系统潮流计算计算结果2.3.1 节点导纳矩阵求得节点导纳矩阵Y = 各节点的导纳值如下：; ； ; ;； ;； ; ; ; ; ; ;.2.3.2 简化雅可比矩阵形成有功迭代和无功迭代的简化雅可比矩阵B/和B/B/= B/= 将B/ 和B/进行三角分解： -0.234654 -0.564073 -0.443749 -0.423068 2.357080 -1.02902 -0.234654 -0.564073 -0.423068 2.3.3 修正、迭代给定PQ节点初值和各节点电压相角初值V1=1.050。 ，V2(0)=V3(0)=1.0，V4=1.12(0)=3(0)=0， 4(0)=01 作第一次有功迭代，按公式计算节点有功功率不平衡量 P2(0)=-0.55-（-0.024037）=-0.525963 P3(0)=-0.30-（-0.022695）=-0.277305 P4(0)=0.500000 P1(0)/V1(0)=0.454545 P2(0)/ V2(0)=-0.525963 P3(0)/V3(0)=-0.2773092做第一次无功迭代，按公式计算无功功率不平衡量，计算时电压相角最新的修正值。 Q2(0)=-0.13-(-0.001550)=-0.039594Q3(0)=-0.18-(-0.14406)=-0.039588Q2(0)/ V2(0)=-0.131553Q3(0)/V3(0)=-0.039588 解修正方程式，可得各节点电压幅值的修正量为 V3(0)）=-0.014855 于是有： V2(1) = V2(0)+V2(1)=0.964776 V3(1) = V3(0)+V3(1)=0.985145 到这里为止，第一轮有功迭代和无功迭代便做完了。3 按公式计算平衡节点功率，得： P1+jQ1=0.367885+j0.264696经过四轮迭代，节点不平衡功率也下降到10-5以下，迭代到此结束。3.程序设计3.1 MATLAB编程说明及元件描述 MATLAB是一种交互式、面向对象的程序设计语言，广泛的应用于工业界和学术界，主要用于矩阵运算，同时数值分析、自动控制模拟、数字信号处理、动态分析，绘图等方面也具有强大的功能。 在MATLAB设计中，原始数据填写是一个很关键的环节，它与程序使用的方便性和灵活性有着直接关系。 3.2潮流计算流程图图3-2 潮流计算流程图5.2 潮流计算源程序据课题题目，本程序把节点1设为平衡节点，节点2、3、4为PQ节点，节点5为PV节点。程序如下：15G(1,1)=10.834;B(1,1)=-32.500;G(1,2)=-1.667;B(1,2)=5.000;G(1,3)=-1.667;B(1,3)=5.000;G(1,4)=-2.500;B(1,4)=7.500;G(1,5)=-5.000;B(1,5)=15.000;G(2,1)=-1.667;B(2,1)=5.000;G(2,2)=12.917;B(2,2)=-38.750;G(2,3)=-10.000;B(2,3)=30.000;G(2,4)=0;B(2,4)=0;G(2,5)=-1.250;B(2,5)=3.750;G(3,1)=-1.667;B(3,1)=5.000;G(3,2)=-10.000;B(3,2)=30.000;G(3,3)=12.917;B(3,3)=-38.750;G(3,4)=-1.250;B(3,4)=3.750;G(3,5)=0;B(3,5)=0;G(4,1)=-2.500;B(4,1)=7.500;G(4,2)=0;B(4,2)=0;G(4,3)=-1.250;B(4,3)=3.750;G(4,4)=3.750;B(4,4)=-11.250;G(4,5)=0;B(4,5)=0;G(5,1)=-5.000;B(5,1)=15.000;G(5,2)=-1.250;B(5,2)=3.750;G(5,3)=0;B(5,3)=0;G(5,4)=0;B(5,4)=0;G(5,5)=6.250;B(5,5)=-18.750;Y=G+j*B;delt(1)=0;delt(2)=0;delt(3)=0;delt(4)=0;u(1)=1.0;u(2)=1.0;u(3)=1.0;u(4)=1.0;p(1)=0.20;q(1)=0.20;p(2)=-0.45;q(2)=-0.15;p(3)=-0.40;q(3)=-0.05;p(4)=-0.60;q(4)=-0.10;k=0;precision=1;N1=4; %the N1 is the amount of the PQ buswhile precision0.00001 delt(5)=0; u(5)=1.06; for m=1:N1 for n=1:N1+1 pt(n)=u(m)*u(n)*(G(m,n)*cos(delt(m)-delt(n)+B(m,n)*sin(delt(m)-delt(n); qt(n)=u(m)*u(n)*(G(m,n)*sin(delt(m)-delt(n)-B(m,n)*cos(delt(m)-delt(n); end pp(m)=p(m)-sum(pt); qq(m)=q(m)-sum(qt); end for m=1:N1 for n=1:N1+1 h0(n)=u(m)*u(n)*(G(m,n)*sin(delt(m)-delt(n)-B(m,n)*cos(delt(m)-delt(n); n0(n)=-u(m)*u(n)*(G(m,n)*cos(delt(m)-delt(n)+B(m,n)*sin(delt(m)-delt(n); j0(n)=-u(m)*u(n)*(G(m,n)*cos(delt(m)-delt(n)+B(m,n)*sin(delt(m)-delt(n); L0(n)=-u(m)*u(n)*(G(m,n)*sin(delt(m)-delt(n)-B(m,n)*cos(delt(m)-delt(n); end H(m,m)=sum(h0)-u(m)2*(G(m,m)*sin(delt(m)-delt(m)-B(m,m)*cos(delt(m)-delt(m); N(m,m)=sum(n0)-2*u(m)2*G(m,m)+u(m)2*(G(m,m)*cos(delt(m)-delt(m)+B(m,m)*sin(delt(m)-delt(m); J(m,m)=sum(j0)+u(m)2*(G(m,m)*cos(delt(m)-delt(m)+B(m,m)*sin(delt(m)-delt(m); L(m,m)=sum(L0)+2*u(m)2*B(m,m)+u(m)2*(G(m,m)*sin(delt(m)-delt(m)-B(m,m)*cos(delt(m)-delt(m); end for m=1:N1 JJ(2*m-1,2*m-1)=H(m,m); JJ(2*m-1,2*m)=N(m,m); JJ(2*m,2*m-1)=J(m,m); JJ(2*m,2*m)=L(m,m); end for m=1:N1 for n=1:N1 if m=n else H(m,n)=-u(m)*u(n)*(G(m,n)*sin(delt(m)-delt(n)-B(m,n)*cos(delt(m)-delt(n); J(m,n)=u(m)*u(n)*(G(m,n)*cos(delt(m)-delt(n)+B(m,n)*sin(delt(m)-delt(n); N(m,n)=-J(m,n); L(m,n)=H(m,n); JJ(2*m-1,2*n-1)=H(m,n); JJ(2*m-1,2*n)=N(m,n); JJ(2*m,2*n-1)=J(m,n); JJ(2*m,2*n)=L(m,n); end end end for m=1:N1 PP(2*m-1)=pp(m); PP(2*m)=qq(m); end uu=-inv(JJ)*PP; precision=max(abs(uu); for n=1:N1 delt(n)=delt(n)+uu(2*n-1); u(n)=u(n)+uu(2*n); end k=k+1;endK=k-1,delt,u%the following program is used to calculate the S5 and S(m,n)for n=1:N1+1 U(n)=u(n)*(cos(delt(n)+j*sin(delt(n);end for m=1:N1+1 I(m)=Y(5,m)*U(m);end S5=U(5)*sum(conj(I) for m=1:N1+1 for n=1:N1+1 S(m,n)=U(m)*(conj(U(m)-conj(U(n)*conj(-Y(m,n); endend4.实验结论该设计课题中，以迭代法思想和牛顿拉夫逊法为基础，通过建Y矩阵.雅可比矩阵.逆矩阵，运用MATLAB编程计算分析，从而实现对复杂网络潮流的计算，大大提高了运算速度。本设计采用的方法简单易懂，适用于任何实际网络。在进行课题设计的过程中，加深了我对潮流计算的认识，尤其是对牛顿拉夫逊潮流计算的求解思路有了比