电力系统分析课程设计

电力系统分析课程设计学院电气工程学院班级学号姓名指导教师时间前言电力系统潮流计算是研究电力系统稳态运营情况的一种计算，它根据给定的运营条件及系统接线情况拟定整个电力系统各部分的运营状态：各母线的电压，各元件中流过的功率，系统的功率损耗等等。在电力系统规划的设计和现有电力系统运营方式的研究中，都需要运用潮流计算来定量地分析比较供电方案或运营方式的合理性。可靠性和经济性。此外，电力系统潮流计算也是计算系统动态稳定和静态稳定的基础。所以潮流计算是研究电力系统的一种很重要和很基础的计算。随着科学技术的发展，电力系统变得越来越复杂，电气工程师掌握一种好的能对电力系统进行仿真的软件是学习和研究的需要。与众多专门的电力系统仿真软件相比，MATLAB软件具有易学、功能强大和开放性好，是电力系统仿真研究的有力工具。目录1设计题目 31.1系统图的拟定 31.2各节点的初值及阻抗参数 42潮流计算 52.1潮流计算概述与发展 52.2复杂电力系统潮流计算 52.3MATLAB概述 62.4牛顿-拉夫逊法原理 72.5牛顿-拉夫逊法解决潮流计算问题 82.6计算机潮流计算的环节 92.7计算机潮流计算流程图 113手算潮流计算 123.1拟定节点类型 123.2求节点导纳矩阵Yb 123.3计算各节点功率的修正方程的初始值（不平衡量） 134计算机算法潮流计算 154.1计算机法潮流计算过程 154.2计算机法潮流计算结果 234.3系统功率分布图 25个人心得 29附录：源程序 30参考文献 401设计题目1.1系统图的拟定选择六节点、环网、两电源和多引出的电力系统，简化电力系统图如图1-1所示，等值导纳图如图1-2所示。运用以直角坐标表达的牛顿-拉夫逊计算如图1-1所示系统中的潮流分布。计算精度规定各节点电压的误差或修正量不大于。图1-1电力系统图图1-2电力系统等值导纳图1.2各节点的初值及阻抗参数该系统中，节点①为平衡节点，保持=1.05+j0为定值，节点⑥为PV节点，其他四个节点都是PQ节点。给定的注入电压标幺值如表1-1、线路阻抗标幺值如表1-2、输出功率标幺值如表1-3。表1-1各节点电压标幺值参数UUUUUU1.051.001.001.001.001.05表1-2线路、变压器阻抗标幺值线路L2L3L4L5T1T2Y/2阻抗0.06+j0.250.04+j0.250.08+j0.300.1+j0.35j0.03j0.015j0.25表1-3节点输出功率节点②③④⑤⑥功率2+j11.8+j0.41.6+j0.83.7+j1.35注：各PQ节点的电压取1是为了方便计算和最后验证程序的对的性。2潮流计算2.1潮流计算概述与发展电力系统潮流计算也分为离线计算和在线计算两种，前者重要用于系统规划设计和安排系统的运营方式，后者则用于正在运营系统的经常监视及实时控制。运用电子数字计算机进行电力系统潮流计算从50年代中期就已经开始。在这2023内，潮流计算曾采用了各种不同的方法，这些方法的发展重要围绕着对潮流计算的一些基本规定进行的。牛顿-拉夫逊法作为一种实用的，有竞争力的电力系统潮流计算方法，是在应用了稀疏矩阵技巧和高斯消去法求修正方程后。牛顿-拉夫逊法是求解非线性代数方程有效的迭代计算。2.2复杂电力系统潮流计算电力系统潮流计算是对复杂电力系统正常和故障条件下稳态运营状态的计算。潮流计算的目的是求取电力系统在给定运营方式下的节点电压和功率分布，用以检查系统各元件是否过负荷、各点电压是否满足规定、功率的分布和分派是否合理以及功率损耗等。对现有电力系统的运营和扩建，对新的电力系统进行规划设计都是以潮流计算为基础。潮流计算结果的用途，例如用于电力系统稳定研究、安全估计或最优潮流等也对潮流计算的模型和方法有直接影响。节点类型：（1）PV节点：柱入有功功率P为给定值，电压也保持在给定数值。（2）PQ节点：诸如有功功率和无功功率是给定的。（3）平衡节点：用来平衡全电网的功率。选一容量足够大的发电机担任平衡全电网功率的职责。平衡节点的电压大小与相位是给定的，通常以它的相角为参考量，即取其电压相角为0。一个独立的电力网中只设一个平衡点。基本环节：（1）形成节点导纳矩阵；（2）将各节点电压设初值U；（3）将节点初值代入相关求式，求出修正方程式的常数项向量；（4）将节点电压初值代入求式，求出雅可比矩阵元素；（5）求解修正方程，求修正向量；（6）求取节点电压的新值；（7）检查是否收敛，如不收敛，则以各节点电压的新值作为初值自第3步重新开始进行狭义次迭代，否则转入下一步；（8）计算支路功率分布，PV节点无功功率和平衡节点柱入功率。2.3MATLAB概述目前电子计算机已广泛应用于电力系统的分析计算，潮流计算是其基本应用之一。现有很多潮流计算方法。对潮流计算方法有五方面的规定：（1）计算速度快；（2）内存需要少；（3）计算结果有良好的可靠性和可信性；（4）适应性好，亦即能解决变压器变比调整、系统元件的不同描述和与其它程序配合的能力强；（5）简朴。MATLAB是一种交互式、面向对象的程序设计语言，广泛应用于工业界与学术界，重要用于矩阵运算，同时在数值分析、自动控制模拟、数字信号解决、动态分析、绘图等方面也具有强大的功能。MATLAB程序设计语言结构完整，且具有优良的移植性，它的基本数据元素是不需要定义的数组。它可以高效率地解决工业计算问题，特别是关于矩阵和矢量的计算。MATLAB与C语言和FORTRAN语言相比更容易被掌握。通过M语言，可以用类似数学公式的方式来编写算法，大大减少了程序所需的难度并节省了时间，从而可把重要的精力集中在算法的构思而不是编程上。2.4牛顿-拉夫逊法原理假设有n个联立的非线性代数方程：假设以给出各变量的初值，，……，，令其分别为个变量的修正量，使满足以上方程，所以：将上式中的n个多元函数在初始值附近分别展开成泰勒级数，并略去具有，，……，的二次及以上阶次的各项，便得：方程可写成：以上方程是对于修正量，，……，的线性方程组，称为牛顿法的修正方程，可解出，，……，。对初始近似解进行修正：（i=1,2,……，n）反复迭代，在进行k+1次迭代时，从求解修正方程式：得到修正量,,……,,对各量进行修正（i=1,2,……，n）迭代过程一直进行到满足收敛判据2.5牛顿-拉夫逊法解决潮流计算问题节点总数为n；PQ节点有m，；PV节点有n-m-1，平衡节点有1个，节点编号按照先PQ节点，再PV节点，最后平衡节点的顺序进行编号，即:1，2，…，m为PQ节点；m+1，m+2，…，n-1为PV节点；n为平衡节点。可形成结点导纳矩阵。导纳矩阵元素可表达为，本文中节点电压以直角坐标形式表达，即。由此下列公式可求出Pi，Qi假设系统中的第1，2，…，m号节点为PQ节点，第i个节点的给定功率为和，对该节点可列方程：假设系统中的第m+1，m+2，…，n-1号节点为PV节点，则对其中每一个节点可列方程:第n号节点为平衡节点，其电压为是给定的，故不参与迭代。修正方程可写成分块矩阵的形式：通过反复求解修正方程，解出各节点的未知量，再通过收敛判据鉴定是否已为真值。从而求得PQ节点的电压V及相角δ的真值，PV节点的Q、δ真值，平衡节点的P、Q真值，以上即为牛顿-拉夫逊迭代法的潮流计算过程，其优点为计算精确，运营速度快。其中的各个环节都可通过MATLAB程序来实现。2.6计算机潮流计算的环节（1）对电力网络的所有参数设初值,涉及电压、相角、有功、无功等。（2）解决非标准变比支路,使其变成标准变比为1的变压器支路。（3）形成节点导纳矩阵Y。（4）计算有功功率的不平衡量ΔPi,从而求出。（5）根据节点的类型形成J。（6）解修正方程式,求各节点的电压的变化量Δei(i=1,2,3...n,i≠s)（7）求各节点相角的新值ei=ei+Δei(i=1,2,3...n,i≠s)（8）计算无功功率的不平衡量ΔQi,从而求出(i=1,2,3...n,i≠s)（9）解修正方程式,求各节点的电压大小的变化量(i=1,2,3...,n,i≠s)。（10）求各节点的电压大小的新值(i=1,2,3...,n,i≠s)。（11）运用个节点的电压的新值自第四步开始下一次迭代。计算平衡节点的功率和线路功率。其中平衡节点的功率的计算公式为线路上的功率为:从而线路上的损耗的功率为:2.7计算机潮流计算流程图选用牛顿-拉夫逊法进行计算机潮流计算。牛顿-拉夫逊法潮流计算流程图如图2-1。图2-1牛顿-拉夫逊法潮流计算流程图3手算潮流计算3.1拟定节点类型节点①为平衡节点节点⑥为PV节点，其他四个节点都是PQ节点。3.2求节点导纳矩阵Yb由图1可知，该系统以串联支路的导纳标幺值和对地并联导纳标幺值得等值电路如图2所示。可得图2的等值节点导纳矩阵。对角线上的元素为：非对角线上的元素为：所以节点导纳矩阵为：Y=3.3计算各节点功率的修正方程的初始值（不平衡量）设各节点电压初始值为：根据上述N-R法的求解过程，将各节点电压的初始值代入（1）式和（2）式，进行潮流计算，得节点功率和节点电压：雅可比矩阵：J=修正各节点电压：4计算机算法潮流计算4.1计算机法潮流计算过程导纳矩阵Y=0-30.2343i0+31.7460i00000+31.7460i14.8252-42.6506i-14.2023+5.9172i0-0.6240+3.9002i00-14.2023+5.9172i15.0311-8.5292i-0.8299+3.1120i0000-0.8299+3.1120i1.5846-5.5035i-0.7547+2.6415i00-0.6240+3.9002i0-0.7547+2.6415i1.3787-72.9583i0+63.4921i00000+63.4921i0-60.4686i功率方程第(1)次差值：Columns1through12004.2619-2.10000.1000-1.8000-0.5500-1.60008.4738-3.700005.0000形成的第(1)次Jacobi矩阵：Columns1through11000000000000000000000000-37.388714.82525.9172-14.2023003.9002-0.6240000-14.8252-47.912614.20235.9172000.62403.90020005.9172-14.2023-8.029215.03113.1120-0.82990000014.20235.9172-15.0311-9.02920.82993.112000000003.1120-0.8299-5.25351.58462.6415-0.7547000000.82993.1120-1.5846-5.75350.75472.64150003.9002-0.6240002.6415-0.7547-63.18451.378763.4921000.62403.9002000.75472.6415-1.3787-82.732100000000000-2.400000000000076.19050Columns12through13000004.26190-2.100000.10000-1.80000-0.55000-1.600008.473863.4921-3.700000-63.49215.0000Jacobi矩阵第(1)次回代运算Columns1through1200-0.14900.1103-0.00240.16570.10150.3342-0.14350.02970-0.0431各个节点电压模01.15431.01600.95861.14391.2023功率方程第(2)次差值：Columns1through1200-0.9196-0.3298-0.2847-0.0127-0.56520.0384-1.51850.1960-0.0019-0.3928形成的第(2)次Jacobi矩阵：Columns1through11000000000000000000000000-47.156723.25685.2333-16.9700004.4126-1.1470000-20.2170-47.588316.97005.2333001.14704.41260003.5786-15.2150-5.660518.19672.9819-1.34740000015.21503.5786-14.7623-6.45821.34742.981900000002.5188-1.7857-3.58994.77952.1211-1.5609000001.78572.5188-1.7466-5.24071.56092.12110004.4413-0.8294002.9982-0.9415-83.48957.152972.6032000.82944.4413000.94152.9982-0.3333-83.28441.88540000000000-2.400000000000-2.737576.1905-1.8854Columns12through1300000-0.91960-0.32980-0.28470-0.01270-0.565200.0384-1.8854-1.518572.60320.1960-0.0862-0.0019-72.6032-0.3928Jacobi矩阵第(2)次回代运算Columns1through12000.04410.01550.0700-0.02100.1790-0.02750.02880.0422-0.0010.0487各个节点电压模01.11200.94350.78211.11701.2023功率方程第(3)次差值：Columns1through1200-0.04170.0017-0.0185-0.0126-0.1340-0.0217-0.03670.0646-0.0024-0.0916形成的第(3)次Jacobi矩阵：Columns1through11000000000000000000000000-44.191723.52304.7526-16.4350004.2309-1.1798000-19.9621-46.331516.43504.7526001.17984.23090003.4621-14.0963-5.087217.05792.7814-1.22400000014.09633.4621-13.4379-6.46741.22402.781400000001.9845-1.5515-1.89914.35041.6691-1.3531000001.55151.9845-1.3054-5.04841.35311.66910004.3026-0.9762002.8902-1.0313-79.882410.076070.7751000.97624.3026001.03132.8902-3.4952-82.57374.56770000000000-2.4019000000000.353576.2524-4.5677Columns12through1300000-0.041700.00170-0.01850-0.01260-0.13400-0.0217-4.5677-0.036770.77510.06460.0111-0.0024-70.7751-0.0916Jacobi矩阵第(3)次回代运算Columns1through12000.00780.00430.0191-0.00630.0768-0.00550.00560.01560.00110.0181各个节点电压模01.10480.92370.70981.11261.2023功率方程第(4)次差值：Columns1through1200-0.00040.0002-0.0002-0.0012-0.0276-0.00580.00060.0016-0.0003-0.0053形成的第(4)次Jacobi矩阵：Columns1through11000000000000000000000000-43.744923.59274.6455-16.3504004.1980-1.1918000-20.0302-45.989216.35044.6455001.19184.19800003.4389-13.7881-4.989916.76802.7273-178813.4389-13.0467-6.42921.18852.727300000001.7502-1.4706-1.12174.24801.4705-1.2806000001.47061.7502-1.1039-4.99841.28061.47050004.2712-1.0337002.8637-1.0684-79.373711.143870.4222001.03374.2712001.06842.8637-4.6945-82.22825.56080000000000-2.3998000000001.503176.1844-5.5608Columns12through1300000-0.000400.00020-0.00020-0.00120-0.02760-0.0058-5.56080.000670.42220.00160.0473-0.0003-70.4222-0.0053Jacobi矩阵第(4)次回代运算Columns1through12000.00210.00140.0056-0.00200.0259-0.00180.00160.00480.00020.0053各个节点电压模01.10290.91780.68571.11141.2023功率方程第(5)次差值：Columns1through1200-0.00000.00000.0000-0.0001-0.0032-0.0007-0.00000.0001-0.0000-0.0005形成的第(5)次Jacobi矩阵：Columns1through11000000000000000000000000-43.630023.62104.6137-16.3286004.1889-1.1958000-20.0561-45.884116.32864.6137001.19584.18890003.4337-13.6963-4.968816.68032.7114-169633.4337-12.9312-6.41371.17762.711400000001.6711-1.4437-0.85644.21621.4034-1.2564000001.44371.6711-1.0344-4.98431.25641.40340004.2621-1.0517002.8559-1.0801-79.237411.491170.3226001.05174.2621001.08012.8559-5.0505-82.12235.86870000000000-2.3993000000001.841976.1690-5.8687Columns12through1300000-0.000000.000000.00000-0.00010-0.00320-0.0007-5.8687-0.000070.32260.00010.0580-0.0000-70.3226-0.0005Jacobi矩阵第(5)次回代运算Columns1through12000.00030.00020.0008-0.00030.0039-0.00030.00020.00070.00000.0008各个节点电压模01.10260.91690.68201.11131.2023功率方程第(6)次差值：1.0e-004*Columns1through1200-0.00260.00070.0032-0.0182-0.7419-0.1647-0.00160.0184-0.0063-0.1036形成的第(6)次Jacobi矩阵：Columns1through11000000000000000000000000-43.612823.62534.6089-16.3254004.1876-1.1964000-20.0599-45.868316.32544.6089001.19644.18760003.4330-13.6825-4.965816.66702.7090-168253.4330-12.9139-6.41131.17602.709000000001.6591-1.4396-0.81614.21141.3932-1.2528000001.43961.6591-1.0238-4.98221.25281.39320004.2607-1.0543002.8548-1.0818-79.217011.543170.3078001.05434.2607001.08182.8548-5.1036-82.10675.91480000000000-2.3993000000001.892476.1670-5.9148Columns12through1300000-0.000000.000000.00000-0.00000-0.00010-0.0000-5.9148-0.000070.30780.00000.0596-0.0000-70.3078-0.0000Jacobi矩阵第(6)次回代运算1.0e-004*Columns1through12000.07640.04960.2048-0.07280.9503-0.06410.05640.17530.00740.1923各个节点电压模01.10260.91690.68191.11131.2023功率方程第(7)次差值：1.0e-007*Columns1through1200-0.00150.00040.0019-0.0106-0.4343-0.0965-0.00100.0108-0.0037-0.06064.2计算机法潮流计算结果迭代次数：6没有达成精度规定的个数：9101010830各节点的电压复数值E为(节点号从小到大排列)：1.20231.0947-0.1316i0.9068-0.1361i0.6128-0.2991i1.1073-0.0932i1.1996-0.0298i各节点的电压模值大小V为(节点号从小到大排列)：1.20231.10260.91690.68191.11131.2023各节点的电压相角sida为(节点号从小到大排列)：0-6.8549-8.5347-26.0187-4.8097-1.4241各节点的功率S为(节点号从小到大排列)：5.0135+1.8332i-2.1000-1.0000i-1.8000-0.4000i-1.6000-0.8000i-3.7000-1.3000i5.0000+2.5555i各条支路的首端功率Si为(顺序同您输入B1时一致)：S(1,2)=5.0135+1.8332iS(2,3)=3.0895+0.48941iS(2,5)=-0.176-0.31075iS(3,4)=0.78738+0.39432iS(4,5)=-0.90639-0.43087iS(6,5)=5+2.5555i各条支路的末端功率Sj为(顺序同您输入B1时一致)：S(2,1)=-5.0135-1.1787iS(3,2)=-2.5874-0.79432iS(5,2)=0.17703-0.29553iS(4,3)=-0.69361-0.36913iS(5,4)=1.123+1.1889iS(5,6)=-5-2.1934i各条支路的功率损耗DS为(顺序同您输入B1时一致)：DS(1,2)=-8.8818e-016+0.65451iDS(2,3)=0.50212-0.3049iDS(2,5)=0.0010207-0.60628iDS(3,4)=0.093766+0.025185iDS(4,5)=0.21659+0.75805iDS(6,5)=0+0.36211i总网损S=0.8135+0.88867i4.3系统功率分布图个人心得在这次课程设计中，我强烈感觉到自己在很多方面的局限性，对别人的依赖性比较强。我想我会在以后的学习中不断去发现自己在很多方面的局限性，并一一改正，希望在以后的工作中不要犯同样的错误。在这次课程设计中，我们尽量按照老师的规定做，但在具体的操作过程中，还是出现了很多的问题。搞完这个课程设计让我感觉电力系统分析是一门很有用的课程。由于我对它的学到的知识比较少。在很多时候我很多东西都不了解。并且走了很多的弯路。并且我感觉自己的知识不够连贯。好些时候都出现了卡壳的情况。这次课程设计后，我一定要重新对电力系统分析这门课程做进一步的了解。对在此过程中遗留下的问题做好好的研究。争取早点对电力系统分析这门课程有个全方位的了解。为在以后的毕业课程设计中多些方案。也为我子以后走上工作岗位，提高自己的专业技能，打下扎实的基础。尚有就是在几天的课程设计中，使我养成了很好的学习习惯，和对学习知识的严谨的态度，同时也养成了积极查阅相关资料的好习惯，好习惯的养成是来之不易的，我相信在以后的学习和工作中，我将继续保持这些良好的习惯，并积极努力的学习。让自己更上一层楼。同时在此也感谢一直指导我的老师，本次课程设计的完毕与老师的指导师分不开的，终在我们的一起努力下，完毕了这门课程设计。在此对指导老师致以衷心的感谢。还要感谢学院图书馆为我提供丰富的参考资料，也感谢班上同学给了我很多宝贵的意见和参考，使我获益很多。附录：源程序>>%本程序的功能是用牛顿——拉夫逊法进行潮流计算%B1矩阵：1、支路首端号；2、末端号；3、支路阻抗；4、线路对地电纳(或变压器导纳)；% 5、支路的变比；6、支路首端处在K侧为1，1侧为0；%7、线路/变压器标记（0/1）变压器参数当支路首端处在K侧标记为1时归算至末端侧，0归算至首端侧%B2矩阵：1、该节点发电机功率；2、该节点负荷功率；%3、PQ节点电压初始值，或平衡节点及PV节点电压的给定值%4、节点所接无功补偿并联电容（感）的电纳%5、节点分类标号：1为平衡节点（应为1号节点）；2为PQ节点；3为PV节点；clear;isb=1;%input('请输入平衡母线节点号：isb=');pr=0.00001;%input('请输入误差精度：pr=');%n=6;%input('请输入节点数：n=');nl=6;%input('请输入支路数：nl=');B1=[120+0.03i01.0511;230.06+0.025i0+0.5i100;250.04+0.25i0+0.5i100;340.08+0.30i0+0.5i100;450.1+0.35i0100;650+0.015i01.0511]B2=[001.2001;02.1+1i1.0002;01.8+0.40i1.0002;01.6+0.8i1.0002;03.7+1.3i1.0002;0-5+0i1.2003]%input('请输入各节点参数形成的矩阵：B2=');%X=[10;20;30;40;50;60]%%n=4;%input('请输入节点数：n=');nl=4;%input('请输入支路数：nl=');%B1=[124+16i0100;134+16i0100;232+8i0100;241.49+48.02i011/11001]%input('请输入由支路参数形成的矩阵：B1=');%B2=[0011501;0011002;020+4i11002;010+6i1002]%input('请输入各节点参数形成的矩阵：B2=');%Y=zeros(n);e=zeros(1,n);f=zeros(1,n);V=zeros(1,n);sida=zeros(1,n);S1=zeros(nl);%%%求导纳矩阵%fori=1:n%ifX(i,2)~=0;%p=X(i,1);%Y(p,p)=1/X(i,2);%end%endfori=1:nl %从1到n1（总支路数） ifB1(i,7)==1%假如是变压器支路 ifB1(i,6)==0 %左节点(首端)处在1侧 p=B1(i,1);q=B1(i,2); else%左节点(首端)处在K侧 p=B1(i,2);q=B1(i,1); end Y(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); %对角元K侧 Y(p,p)=Y(p,p)+1./B1(i,3)+B1(i,4); %对角元1侧+励磁导纳 else%否则为线路支路 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(q,q)=Y(q,q)+1./B1(i,3)+B1(i,4)./2.0000; %对角元j侧+线路电纳的一半 Y(p,p)=Y(p,p)+1./B1(i,3)+B1(i,4)./2.0000; %对角元i侧+线路电纳的一半 endenddisp('导纳矩阵Y=');disp(Y);%给定各节点初始电压及给定各节点注入功率G=real(Y);B=imag(Y); %分解出导纳阵的实部和虚部fori=1:n %给定各节点初始电压的实部和虚部e(i)=real(B2(i,3));f(i)=imag(B2(i,3));V(i)=abs(B2(i,3)); %PV、平衡节点及PQ节点电压模值endfori=1:n %给定各节点注入功率S(i)=B2(i,1)-B2(i,2); %i节点注入功率SG-SLB(i,i)=B(i,i)+B2(i,4); %i节点无功补偿量(电纳值)end%==================用牛顿-拉夫逊法迭代求解非线性代数方程（功率方程）=======================P=real(S);Q=imag(S);%分解出各节点注入的有功和无功功率ICT1=0;IT2=1;N0=2*n;N1=N0+1;a=0;%迭代次数ICT1、a；不满足收敛规定的节点数IT2whileIT2~=0%N0=2*n雅可比矩阵的阶数；N1=N0+1扩展列IT2=0;a=a+1;JZ=['Jacobi矩阵第(',num2str(a),')次消去运算'];JZ1=['Jacobi矩阵第(',num2str(a),')次回代运算'];JZ0=['功率方程第(',num2str(a),')次差值：'];%求取各个节点的功率及功率偏差及PV节点的电压偏差fori=1:n%n个节点2n行(每节点两个方程P和Q或U)p=2*i-1;m=p+1;C(i)=0;D(i)=0;forj1=1:n%第i行共n列(n个节点间互导纳及节点电压相乘即电流)C(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)end%求i节点有功和无功功率P',Q'的计算值P1=C(i)*e(i)+f(i)*D(i);%节点功率P计算eiΣ(Gij*ej-Bij*fj)+fiΣ(Gij*fj+Bij*ej)Q1=C(i)*f(i)-e(i)*D(i);%节点功率Q计算fiΣ(Gij*ej-Bij*fj)-eiΣ(Gij*fj+Bij*ej)V2=e(i)^2+f(i)^2; %电压模平方%===求取功率差及PV节点电压模平方差=========ifi~=isb %非平衡节点(PQ或PV节点)ifB2(i,5)~=3 %非PV节点(只能是PQ节点)J(m,N1)=P(i)-P1; %PQ节点有功功率差J(m,N1)扩展列△PJ(p,N1)=Q(i)-Q1;%PQ节点无功功率差J(p,N1)扩展列△Qelse%PV节点==================J(m,N1)=P(i)-P1; %PV节点有功功率差J(m,N1)扩展列△PJ(p,N1)=V(i)^2-V2; %PV节点电压模平方差J(p,N1)扩展列△Uendend%(ifi~=isb) 非平衡节点(PQ或PV节点)end%(fori=1:n)n个节点2n行(每节点两个方程P和Q或U)form=1:N0JJN1(m)=J(m,N1);enddisp(JZ0);disp(JJN1);%判断功率偏差量及PV节点的电压偏差量是否满足规定fork=3:N0%除去平衡节点1、2号以外的所有节点DET=abs(J(k,N1));ifDET>=pr;%PQ节点的功率偏差量及PV节点的电压偏差量是否满足规定IT2=IT2+1;%不满足规定的节点数加1endendICT2(a)=IT2;%不满足规定的节点数；a为迭代次数ICT1=ICT1+1;%迭代次数ifICT2(a)==0;%当前不满足规定的节点数为零break%退出迭代运算end%以上为求取各个节点的功率及功率偏差及PV节点的电压偏差%=================求取Jacobi矩阵形成修正方程===================fori=2:n%n个节点2n行(每节点两个方程P和Q或U)ifi~=isb %非平衡节点(PQ或PV节点)ifB2(i,5)~=3%下面是针对PQ节点来求取Jacobi矩阵的元素===========C(i)=0;D(i)=0;forj1=1:n%第i行共n列(n个节点间互导纳及节点电压相乘即电流)C(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)endforj1=2:n%第i行共n列(2n个Jacobi矩阵元素dP/de及dP/df或dQ/de及dQ/df)ifj1~=isb&j1~=i %非平衡节点&非对角元X1=-G(i,j1)*e(i)-B(i,j1)*f(i); %X1=dP/de=-dQ/df=-X4X2=B(i,j1)*e(i)-G(i,j1)*f(i); %X2=dP/df=dQ/de=X3X3=X2;%X2=dp/dfX3=dQ/deX4=-X1;%X1=dP/deX4=dQ/dfp=2*i-1;q=2*j1-1;J(p,q)=X3;m=p+1;%X3=dQ/deJ(p,N)=DQ节点无功功率差J(p,N)=DQ;J(m,q)=X1;q=q+1;%X1=dP/deJ(m,N)=DP节点有功功率差J(m,N)=DP;J(p,q)=X4;J(m,q)=X2;%X4=dQ/dfX2=dp/dfelseifj1==i&j1~=isb %非平衡节点&对角元X1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i);%dP/deX2=-D(i)+B(i,i)*e(i)-G(i,i)*f(i);%dP/dfX3=D(i)+B(i,i)*e(i)-G(i,i)*f(i);%dQ/deX4=-C(i)+G(i,i)*e(i)+B(i,i)*f(i);%dQ/dfp=2*i-1;q=2*j1-1;J(p,q)=X3;%扩展列△QJ(p,N)=DQ;m=p+1;J(m,q)=X1;q=q+1;J(p,q)=X4;%扩展列△PJ(m,N)=DP;J(m,q)=X2;endendelse%ifB2(i,5)~=3%否则（即为PV节点）%===============下面是针对PV节点来求取Jacobi矩阵的元素===========forj1=1:nifj1~=isb&j1~=i %非平衡节点&非对角元X1=-G(i,j1)*e(i)-B(i,j1)*f(i);%dP/deX2=B(i,j1)*e(i)-G(i,j1)*f(i);%dP/dfX5=0;X6=0;p=2*i-1;q=2*j1-1;J(p,q)=X5;%PV节点电压误差J(p,N)=DV;m=p+1;J(m,q)=X1;q=q+1;J(p,q)=X6;%PV节点有功误差J(m,N)=DP;J(m,q)=X2;elseifj1==i&j1~=isb %非平衡节点&对角元X1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i);%dP/deX2=-D(i)+B(i,i)*e(i)-G(i,i)*f(i);%dP/dfX5=-2*e(i);X6=-2*f(i);p=2*i-1;q=2*j1-1;J(p,q)=X5;%PV节点电压误差J(p,N)=DV;m=p+1;J(m,q)=X1;q=q+1;J(p,q)=X6;%PV节点有功误差J(m,N)=DP;J(m,q)=X2;endendend%(ifB2(i,5)~=3else)end%(ifi~=isb)end%(fori=1:n)n个节点2n行(每节点两个方程P和Q或U)JZ0=['形成的第(',num2str(a),')次Jacobi矩阵：'];disp(JZ0);disp(J);%===============================以上为形成完整的Jacobi矩阵================================%====下面用高斯消去法对由Jacobi矩阵形成的修正方程进行求解(按列消去、回代)==========fork=3:N0%N0=2*n（从第三行开始，第一、二行是平衡节点）fork1=k+1:N1 %从k+1列的Jacobi元素到扩展列的△P、△Q或△UJ(k,k1)=J(k,k1)./J(k,k);%用K行K列对角元素去除K行K列后的非对角元素进行规格化endJ(k,k)=1;%对角元规格化K行K列对角元素赋1%==================按列消去运算==================================fork2=k+1:N0%从k+1行到2*n最后一行fork3=k+1:N1%从k2+1列到扩展列消去k+1行后各行下三角元素J(k2,k3)=J(k2,k3)-J(k2,k)*J(k,k3);%消去运算end%用当前行K3列元素减去当前行K列元素乘以第k行K3列元素J(k2,k)=0;%当前行第k列元素已消为0endendJZ=['Jacobi矩阵第(',num2str(a),')次消去运算'];JZ1=['Jacobi矩阵第(',num2str(a),')次回代运算'];disp(JZ);disp(J);%====================按列回代运算=======================================fork=N0:-1:3fork1=k-1:-1:3J(k1,N1)=J(k1,N1)-J(