电力系统潮流计算

1、Companynumber:0089WT-8898YT-W8CCB-BUUT-202108电力系统课程设计题目：电力系统潮流计算院系名称：电气工程学院专业班级：电气F1206班学生姓名：学号：指导教师：张孝远成绩：目录原始资料指导老师签名-21概述日期:42潮流计算节点介绍4变量的分类5节点的分类53计算方法简介6牛顿拉夫逊法原理6牛顿拉夫逊法概要6牛顿法的框图及求解过程8MATLAB简介94潮流分布计算10系统的一次接线图10参数计算10丰大及枯大下地潮流分布情况14该地区变压器的有功潮流分布数据15重、过载负荷元件统计表175设计心得17参考文献18附录:程序19原始资料一、系统接线图见附

2、件1。二、系统中包含发电厂、变电站、及其间的联络线路。500kV变电站以外的系统以一个等值发电机代替。各元件的参数见附件2。设计任务1、手动画出该系统的电气一次接线图，建立实际网络和模拟网络之间的联系。2、根据已有资料，先手算出各元件的参数，后再用Matlab表格核算出各元件的参数。3、潮流计算1）对两种不同运行方式进行潮流计算，注意110kV电网开环运行。2）注意将电压调整到合理的范围110kV母线电压控制在106kV117kV之间；220kV母线电压控制在220kV242kV之间。附件一：A"一一一:二GEFL：x150H12.5+31.53x40D.2x31.5课程设计地理接线

3、示意图1x31.590+120火电厂110kV线路B20+8水电站1302x10O-2x80水电站272水电站3口24E3水电站火电厂鬥水电站54臼牵引站曰I-220kV线路llOkV变电站_220kV变电站附件二：1、变压器：两个220kV变电站均采用参数一致的三绕组变压器，具体参数如下。220kV变电站参数表咼压侧绕组中压侧绕组低压侧绕组容量12012060(120)电压220110(121)110kV及以下的变电站的变压器省略，即可将负荷直接挂在110kV母线上。而110kV升压变只计及以下参数。110kV变电站参数表序号变电站名容量X1X01A162B283C1205D636E207G

4、8F+9水电站130MW、2*20MVA10水电站272MW、3*11水电站324MW、12水电站54*、2*20MVA13水电站418MW、30MVA14火电50MW厂、+40MVA2、线路：具体参数如下。220kV线路参数表序号线路名称导线牌号线路长度km1ML2x240102ML2x240103MH2x3005110kV线路参数表序号线路名称导线牌号线路长度km1水电站1水电站2150302AB95803BC9514水电站4C150655BD95636CD240607水电站2C240758ED2409水电站3D15010水电站5T节点7011T节点D70512T节点G240513HF15

5、014HD185715水电站2L240/213016LG240103、发电机各发电机的参数如下：XdttXqtt装机容量功率因数水电站130水电站272水电站324水电站530水电站418华鑫电厂5出力情况：水力发电机丰大出力70%,枯大出力20%火力发电机丰大出力80%,枯大出力80%。4、负荷各110kV变电站丰大负荷按该站变电容量的50%估算，枯大负荷按该站变电容量的60%估算。两个220kV变电站的低压侧上各挂10MW的负荷，中压侧各挂20MW负荷。功率因素均为。5、并联电容器两个220kV变电站的低压侧上均装设并联补偿。补偿总量按该站变电容量的20%装设，分组原则以每组电容器的容量不

6、超过10MVar且经济性较好为准。1 概述潮流计算是电力系统最基本最常用的计算。根据系统给定的运行条件，网络接线及元件参数，通过潮流计算可以确定各母线的电压，包括电压的幅值和相角，各元件流过的功率，整个系统的功率损耗等一系列系统中的潮流数据。近几年，对潮流算法的研究仍然是如何改善传统的潮流算法，即高斯-塞德尔法、牛顿法和快速解耦法。牛顿法，由于其在求解非线性潮流方程时采用的是逐次线性化的方法，为了进一步提高算法的收敛性和计算速度，人们考虑采用将泰勒级数的高阶项或非线性项也考虑进来，于是产生了二阶潮流算法。后来又提出了根据直角坐标形式的潮流方程是一个二次代数方程的特点，提出了采用直角坐标的保留非

7、线性快速潮流算法。潮流计算在数学上是多元非线性方程组的求解问题，求解的方法有很多种，牛顿拉夫逊Newton-Raphson法是数学上解非线性方程组的有效方法，有较好的收敛性。将N-R法用于潮流计算是以导纳矩阵为基础的，由于利用了导纳矩阵的对称性，稀疏性及节点编号顺序优划等技巧，使N-R法在收敛性，占用内存，计算速度等方面的优点都超过了阻抗法总结为在电力系统运行方式和规划方案的研究中都需要进行潮流计算以比较运行方式或规划供电方案的可行性、可靠性和经济性。同时为了实时监控电力系统的运行状态也需要进行大量而快速的潮流计算。因此潮流计算是电力系统中应用最广泛、最基本和最重要的一种电气运算。在系统规划设

8、计和安排系统的运行方式时采用离线潮流计算在电力系统运行状态的实时监控中则采用在线潮流计算。2 潮流计算节点介绍常规的电力系统潮流计算中一般具有三种类型的节点：PQ、PV及平衡节点。一个节点有四个变量，即注入有功功率、注入无功功率，电压大小及相角。常规的潮流计算一般给定其中的二个变量:PQ节点（注入有功功率及无功功率），PV节点（注入有功功率及电压的大小），平衡节点（电压的大小及相角）。变量的分类负荷消耗的有功、无功功率一一p、Q、p、QL1L1L2L2电源发出的有功、无功功率一一P、Q、p、QG1G1G2G2母线或节点的电压大小和相位U、U、5、81212在这十二个变量中，负荷消耗的有功和无功

9、功率无法控制，因它们取决于用户，它们就称为不可控变量或是扰动变量。电源发出的有功无功功率是可以控制的自变量，因此它们就称为控制变量。母线或节点电压的大小和相位角一是受控制变量控制的因变量。其中，U、U主要受Q、Q的控制，5、12G1G215主要受p、p的控制。这四个变量就是简单系统的状态变量。2G1G2为了保证系统的正常运行必须满足以下的约束条件：对控制变量对没有电源的节点则为对状态变量U的约束条件则是i对某些状态变量5还有如下的约束条件i节点的分类第一类称PQ节点。等值负荷功率P、Q和等值电源功率P、Q是LiLiGiGi给定的，从而注入功率P、Q是给定的，待求的则是节点电压的大小U和相iii

10、位角5。属于这类节点的有按给定有功、无功率发电的发电厂母线和没有其他i电源的变电所母线。(第二类称PV节点。等值负荷和等值电源的有功功率P、P是给定LiGi的，从而注入有功功率P是给定的。等值负荷的无功功率Q和节点电压的大小iLiU也是给定的。待求的则是等值电源的无功功率Q，从而注入无功功率QiGii和节点电压的相位角5。有一定无功功率储备的发电厂和有一定无功功率电源i的变电所母线都可以作为PV节点；第三类平衡节点。潮流计算时一般只设一个平衡节点。等值负荷功率P、Q是给定的，节点电压的大小和相位也是给定的。担负调整系统频率任LsLs务的发电厂母线往往被选作为平衡节点。3 计算方法简介牛顿拉夫逊

11、法原理牛顿拉夫逊法概要首先对一般的牛顿一拉夫逊法作一简单的说明。已知一个变量X函数为：到此方程时，由适当的近似值X出发，根据：反复进行计算，当X(n)满足适当的收敛条件就是上面方程的根。这样的方法就是所谓的牛顿拉夫逊法这一方法还可以做下面的解释，设第n次迭代得到的解语真值之差，即X(n)的误差为£时，贝,把f(X(n)+£)在X(n)附近对£用泰勒级数展开上式省略去£2以后部分X(n)的误差可以近似由上式计算出来。比较两式，可以看出牛顿一拉夫逊法的休整量和X(n)的误差的一次项相等。用同样的方法考虑，给出n个变量的n个方程：对其近似解X'得修正量

12、AX'可以通过解下边的方程来确定：11式中等号右边的矩阵%都是对于X:X'，,X'的值。这一矩阵称为雅可dx12nn比(JACOBI)矩阵。按上述得到的修正向量AX:AX'，,AX'后，得到如下关12n系这比X:X'，,X'更接近真实值。这一步在收敛到希望的值以前重复进12n行，一般要反复计算满足£为预先规定的小正数，Xn+1是第n次迭代X的近似值。nn牛顿法的框图及求解过程1、用牛顿法计算潮流时，有以下的步骤：(1) 给这各节点电压初始值e(0),f；(2) 将以上电压初始值代入公式，求修正方程的常数项向量AP(0),AQ(0

13、),(AV2)(0)；（3）将电压初始值在带入上述公式，求出修正方程中系数矩阵的各元素。（4）解修正方程式2（o）,曲（o）；（5）修正各节点电压e=e（o）+Ae（o），/=/+Af；（6）将e，/在带入方程式，求出AP（i）,AQ（i）,（AV2）（1）；（7）检验是否收敛，即maxIp（k）,AQ（K）LsII（8）如果收敛，迭代到此结束，进一步计算各线路潮流和平衡节点功率，并打印输出结果。如果不收敛，转回（2）进行下次迭代计算，直到收敛为止。2、程序框图如下MATLAB简介MATLAB是用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境，主要包括MATLAB和

14、Simulink两大部分。是由美国mathworks公司发布的主要面对科学计算、可视化以及交互式程序设计的高科技计算环境。它将数值分析、矩阵计算、科学数据可视化以及非线性动态系统的建模和仿真等诸多强大功能集成在一个易于使用的视窗环境中，为科学研究、工程设计以及必须进行有效数值计算的众多科学领域提供了一种全面的解决方案，并在很大程度上摆脱了传统非交互式程序设计语言（如C、Fortran）的编辑模式，代表了当今国际科学计算软件的先进水平。MATLAB是一种交互式、面向对象的程序设计语言广泛应用于工业界与学术界主要用于矩阵运算同时在数值分析、自动控制模拟、数字信号处理、动态分析、绘图等方面也具有强大

15、的功能。MATLAB程序设计语言结构完整且具有优良的移植性它的基本数据元素是不需要定义的数组。它可以高效率地解决工业计算问题特别是关于矩阵和矢量的计算。MATLAB与C语言和FORTRAN语言相比更容易被掌握。通过M语言可以用类似数学公式的方式来编写算法大大降低了程序所需的难度并节省了时间，从而可把主要的精力集中在算法的构思而不是编程上。目前电子计算机已广泛应用于电力系统的分析计算潮流计算是其基本应用软件之一。现有很多潮流计算方法。对潮流计算方法有五方面的要求（1）计算速度快（2）内存需要少（3）计算结果有良好的可靠性和可信（4）适应性好亦即能处理变压器变比调整、系统元件的不同描述和与其它程序

16、配合的能力强。潮:系统的一次接线图图系统的一次连接图参数计算设定基准值5100MVA，Ub=，则各参数如下。B（1）发电机的次暂态电抗：X=X*Sb/Sn,Zb=Ub2/Sn发电机参数单位（MW）电厂装机容量枯水出力比例丰水出力比例丰大有功丰大无功枯大有功枯大无功短路X*''水电站130水电72电阻：R=Pk*Un2/1000Sn2;电抗：X=UK(%)*UN2/100SN;电导：G=P0/1000UN2;电纳：B=I0(%)*SN/100UN2;式中Un以KV为单位，Sn以MVA为单位，p0、Pk以KW为单位110KV变压器参数变电站名容量XIX0SbR*X*G*B*A16B

17、28C120D63E20GF+水电站130MW、2*20MVA100水电站272MW、3*100水电站324MW、100水电站54*、2*20MVA100火电厂50MW、+40MVA100100220KV三绕组变压器的参数：电阻：R=Pk*Un2/1000Sn2;电抗：X=UK(%)*UN2/100SN;电导：G=P0/1000UN2;电纳：B=I0(%)*SN/100UN2220kV变电站参数表SBR*X*G*B*高压侧绕组中压侧绕组低压侧绕组高压侧绕组(1)100容量120120120中压侧绕组(2)100电压220121低压侧绕组(3)1004)110KV线路参数：r二p/$二=5；X1

18、*=x1*1*SB/。110KV线路参数标么值序号线路名称导线牌号线路长度kmUB(KV)SB(MW)R*X*B*1水电站1水电站2150301151002AB95801151003BC9511151004水电站4C150651151005BD95631151006CD240601151007水电站2C240751151008ED2401151009水电站3电站5T节点7011510011T节点D70511510012T节点G240511510013HD185711510015水电站2L24013011510016LG240101151001

19、7DF150411510018火电D2401115100110KV变电站负荷参数变电站总容曰丰变各站卜卜比变枯电站泳-A瑙偿电容参丰大无T功枯大有T功枯大无功量站匸比变比功功A16B28C120D63E20F44G站别运行方式ABCDEFGH低压侧L低压侧丰大3000枯大30000丰大及枯大下地潮流分布情况电压是衡量电力系统电能质量的标准之一。电压过高或过低，都将对人身及其用电设备产生重大的影响。保证用户的电压接近额定值是电力系统调度的基本任务之一。当系统的电压偏离允许值时，电力系统必须应用电压调节技术调节系统电压的大小，使其维持在允许值范围内。本文经过手算形成了等值电路图，并编写好了程序得出

20、节点电压标幺值，使其满足所要求的调整范围。我们首先对给定的程序输入部分作了简要的分析，程序开始需要我们确定输入节点数、支路数、平衡母线号、支路参数矩阵、节点参数矩阵。（1）为了保证整个系统潮流计算的完整性，我们把凡具有母线及发电机处均选作节点，这样，我们确定电厂一母线上的发电机作为平衡节点，节点号为，其它机组作为PV节点，节点号为，其余节点均为PQ节点，节点号见等值电路图。（2）确定完节点及编号后，各条支路也相应确定了，我们对各支路参数进行了计算。根据所给实际电路图和题中的已知条件，有以下公式计算各输电线路的阻抗和对地支路电容的标幺值和变压器的阻抗标幺值。该地区变压器的有功潮流分布数据（1）该

21、图为丰大潮流模型图：丰大潮流模型图该运行方式的电压合理，负荷分配也均匀，但是有些线路的负载率偏低。比如水电厂2L站的负载率仅仅为%，ML站的线路的负载率也只是%。(2)该图为枯大潮流模型图：)枯大潮流模型图该运行方式的电压合理，负荷分配也均匀，有些线路的负载率偏低但是有些线路的负载率则偏高。比如水电厂2L站的负载率仅仅为,ML站的线路的负载率也只是%D-H站的线路却重载到84%。重、过载负荷元件统计表类型负载率实际是在功率额定容量建议线路限D站负荷到4MW以下5设计心得通过这次课程设计,我发现自己有很多不足的地方,如基础知识掌握不牢固,很多知识点都忘记了,计算速度慢及准确性低,分析问题能力不够

22、全面等等。同时,在设计的过程中遇到很多问题,如怎样使用WORD的工具,计算公式输入,画图等。明白了有些东西看起来很简单,但一旦做起来却需要很多心思,要注意到很多细节问题。要做到能好好理解课本的内容,一定要认认真真做一次计算。因此,完成课程设计使我对课本的内容加深了理解。总体来说,这次的课程设计不单在专业基础方面反映了我的学习还要加倍努力，还在对一些软件的应用需要加强。计算在各种情况下的潮流分布，对于丰大和枯大情况下的潮流分布有了明确的认识。在本次课程设计过程当中，锻炼了自己实际操作分析能力，理论联系实际，对运行中的电力系统，通过潮流计算可以预知各种负荷变化和网络结构的改变会不会危及系统的安全，

23、系统中所有母线的电压是否在允许的范围以内，系统中各种元件（线路、变压器等）是否会出现过负荷，以及可能出现过负荷时应事先采取哪些预防措施等。同时采用MATLAB进行计算机的计算，在计算时采用特殊算法使得潮流计算的过程更快，效率更高。而采用计算机的运算应该是未来的一种趋势，所以我会学习一定的编辑语言如C,C+等，以提高运算准确性和快速性。总体而言，这次的课程设计对我们运用所学知识，发现、提出、分析和解决实际问题、锻炼实践能力的考察，使我们更清楚地知道不足之出，从而提高我们。参考文献1 于永源主编电力系统分析湖南师范大学出版社M.1992年7月2 陈珩编电力系统稳态分析水利电力出版社M.1995年1

24、月第二版3 邱晓燕刘天琪电力系统分析的计算机算法北京中国电力出版2004 李光琦电力系统暂态分析M.北京：水利电力出版社，5 陆敏政主编电力系统习题集水利电力出版社M.199O年附录:程序%潮流计算MATLAB粗略程序？%creatanew_datat=0;s=0;r=0;w=0;number=input('Howmanynodearethere=');%ConvertPqtoanewarrayforii=1:numberifdata(ii,4)=1t=t+1;forjj=1:14new_data1(t,jj)=data(ii,jj);end;a(1,t)=ii;s=s+1;%

25、recordthenumberofthePQnodeend;end;%Convertpvtoanewarrayforii=1:numberifdata(ii,4)=2t=t+1;forjj=1:14new_data1(t,jj)=data(ii,jj);a(1,t)=ii;r=r+1;%recordthenumberofthePVnodeend;end;%Convertset_vtoanewarrayforii=1:numberifdata(ii,4)=3t=t+1;forjj=1:14new_data1(t,jj)=data(ii,jj);end;a(1,t)=ii;w=w+1;end;en

26、d;%creatanew_data2x,y=size(data2)forii=1:xforjj=1:2formm=1:numberifdata2(ii,jj)=a(1,mm)new_data2(ii,jj)=mm;end;end;end;end;forii=1:xforjj=3:14new_data2(ii,jj)=data2(ii,jj);end;end;%creataYY=zeros(number,number);YY=zeros(number,number);yy=zeros(number,number);forii=1:x%forjj=1:14iii=new_data2(ii,1);j

27、jj=new_data2(ii,2);ifnew_data2(ii,5)=2sub=new_data2(ii,6)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6)-new_data2(ii,7)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6)*i;Y(iii,jjj)=-sub./new_data2(ii,14);YY(iii,jjj)=sub./new_data2(ii,14);Y(jjj,iii)=-sub/new_

28、data2(ii,14);YY(jjj,iii)=sub./new_data2(ii,14);yy(iii,jjj)=(ii,14)./(new_data2(ii,14).*new_data2(ii,14).*sub;yy(jjj,iii)=(new_data2(ii,14)-1)./(new_data2(ii,14).*sub;elseY(iii,jjj)=-new_data2(ii,6)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6)+new_data2(ii,7)./(new_data2(ii,7).

29、*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6)*i;YY(iii,jjj)=new_data2(ii,6)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6)-new_data2(ii,7)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6)*i;Y(jjj,iii)=-new_data2(ii,6)./(new_data2(ii,7).*new_data2(ii,7)+new

30、_data2(ii,6).*new_data2(ii,6)+new_data2(ii,7)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6)*i;YY(jjj,iii)=new_data2(ii,6)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6)-new_data2(ii,7)./(new_data2(ii,7).*new_data2(ii,7)+new_data2(ii,6).*new_data2(ii,6)*i;yy

31、(iii,jjj)=new_data2(ii,8)./2.*i;yy(jjj,iii)=new_data2(ii,8)./2.*i;end;%end;end;foriii=1:numberY(iii,iii)=0;end;%forii=1:x%forjj=1:14foriii=1:numberforjj=1:number%ifiii=jjY(iii,iii)=Y(iii,iii)+YY(iii,jj)+yy(iii,jj);%end;end;end;%creatB,Gforii=1:numberforjj=1:numberG(ii,jj)=real(Y(ii,jj);B(ii,jj)=imag

32、(Y(ii,jj);end;end;%creatInitial_PInitial_QInitial_Vforii=1:(s+r)set_P(ii,1)=(new_data1(ii,9)-new_data1(ii,7)./100;end;forii=1:s;set_Q(ii,1)=(new_data1(ii,10)-new_data1(ii,8)./100;end;forii=1:rset_V(ii,1)=new_data1(ii+s,12).*new_data1(ii+s,12);%trytomodifyforsikeofcorrectingend;Initial_p_q_v=set_P;se

33、t_Q;set_V;disp(Initial_p_q_v);%creatInitial_e,Initial_fforii=1:number-1e(ii,1)=1;f(ii,1)=;%changeftotestusedtobeend;e(number,1)=new_data1(number,12);f(number,1)=0;%e(64,1)=;%test118ieee%f(14,1)=0;%e(10,1)=;%e(11,1)=;%e(12,1)=;%e(13,1)=;%StartNEWTOWNCALULATIONfortry_time=1:25%CreateverynodeconsumePQa

34、ndUn=s;m=r;forii=1:(n+m)sum1=0;forjj=1:(n+m+1)sum1=sum1+e(ii,1).*(G(ii,jj).*e(jj,1)-B(ii,jj).*f(jj,1)+f(ii,1).*(G(ii,jj).*f(jj,1)+B(ii,jj).*e(jj,1);end;p(ii,1)=sum1;end;forii=1:nsum2=0;forjj=1:(n+m+1)sum2=sum2+f(ii,1).*(G(ii,jj).*e(jj,1)-B(ii,jj).*f(jj,1)-e(ii,1).*(G(ii,jj).*f(jj,1)+B(ii,jj).*e(jj,1

35、);end;q(ii,1)=sum2;end;disp('q=');disp(q);u=zeros(n+m),1);forii=(n+1):(n+m)u(ii,1)=e(ii,1).*e(ii,1)+f(ii,1).*f(ii,1);end;forii=n+1:(n+m)extra_u(ii-n),1)=u(ii,1);end;disp('extra_u=');disp(extra_u);sum=p;q;extra_u;disp(sum)disp(s);disp(p);%creatJacobiandisp(n);disp(m);forii=1:(n+m)forj

36、j=1:(n+m)PF(ii,jj)=B(ii,jj).*e(ii,1)-G(ii,jj).*f(ii,1);PE(ii,jj)=-G(ii,jj).*e(ii,1)-B(ii,jj).*f(ii,1);elsess=0;qq=0;fornum=1:(n+m+1)ss=ss+G(ii,num).*f(num,1)+B(ii,num).*e(num,1);qq=qq+G(ii,num).*e(num,1)-B(ii,num).*f(num,1);end;PF(ii,jj)=-ss+B(ii,jj).*e(ii,1)-G(ii,jj).*f(ii,1);%TEST+1PE(ii,jj)=-qq-G

37、(ii,jj).*e(ii,1)-B(ii,jj).*f(ii,1);%TEST+1end;end;end;copy=;disp('=copy=')forii=1:nforjj=1:m+nif(ii=jj)QE(ii,jj)=B(ii,jj).*e(ii,1)-G(ii,jj).*f(ii,1);%TEST+1QF(ii,jj)=G(ii,jj).*e(ii,1)+B(ii,jj).*f(ii,1);%TEST+1elsess=0;qq=0;fornum=1:(n+m+1)ss=ss+G(ii,num).*f(num,1)+B(ii,num).*e(num,1);qq=qq+G

38、(ii,num).*e(num,1)-B(ii,num).*f(num,1);end;QF(ii,jj)=-qq+G(ii,jj).*e(ii,1)+B(ii,jj).*f(ii,1);%TEST+1QE(ii,jj)=ss+B(ii,jj).*e(ii,1)-G(ii,jj).*f(ii,1);%TEST+1end;end;end;%disp('QF');%disp(QF);%disp('QE');%disp(QE);UE=zeros(n+m),(n+m);UF=zeros(n+m),(n+m);forii=n+1:n+mforjj=1:(n+m)UE(ii,

39、jj)=0;UF(ii,jj)=0;elsess=0;qq=0;fornum=1:(n+m+1)ss=ss+G(ii,num).*f(num,1)+B(ii,num).*e(num,1);qq=qq+G(ii,num).*e(num,1)-B(ii,num).*f(num,1);end;UF(ii,jj)=-2.*f(ii,1);UE(ii,jj)=-2.*e(ii,1);end;end;end;forii=(n+1):(n+m)forjj=1:(n+m)extra_UE(ii-n),jj)=UE(ii,jj);extra_UF(ii-n),jj)=UF(ii,jj);end;end;%dis

40、p('extra_UE');%disp(extra_UE);%disp('extra_Uf');%disp(extra_UF);Jacobian=PF,PE;QF,QE;extra_UF,extra_UE;%disp('Jacobian=');%disp(Jacobian);%creatsubstractresultsubstract_result=Initial_p_q_v-sum;%disp('substract_result');%disp(substract_result);%calculatedelta_f_edelta

41、_f_e=-inv(Jacobian)*substract_result;%disp(delta_f_e);forii=1:number-1;f(ii,1)=f(ii,1)+delta_f_e(ii,1);e(ii,1)=e(ii,1)+delta_f_e(ii+number-1,1);end;ifmax(substract_result)<1e-4break;end;end;%disp('substract_result');%disp(substract_result);%disp('e=');%disp(e);%disp('f=');

42、%disp(f);forii=1:numberuuu(ii,1)=e(ii,1).*e(ii,1)+f(ii,1).*f(ii,1);U_RESULT(ii,1)=sqrt(uuu(ii,1);end;forii=1:numberforjj=1:numberifii=a(1,jj)Old_Uresult(ii,1)=U_RESULT(jj,1)end;end;end;forii=1:numberOld_Uresult(ii,2)=ii;end;%disp('U_result');%disp(U_RESULT);disp('=');disp('Thelastresultis:'