电力系统分析课程设计 潮流计算

1、武汉理工大学电力系统分析课程设计说明书目录摘要11.任务及题目要求22.设计原理2 2.1设计思路 .2 2.2方法原理.33计算步骤64.结果分析8小结10参考文献11附录：源程序12本科生课程设计成绩评定表24摘要潮流计算，就是根据给定电力系统的接线方式和参数等运行条件，确定电力系统稳态运行各母线电压、个支路电流与功率及网损。其是电力系统中应用最广泛、最基本和最重要的一种电气运算。本设计即根据给定3条支路的运行条件，运用牛顿-拉夫逊的潮流计算方法，求出各线路的功率分布及功率损耗，未知节点的电压。同时进行编程，在matlab运行得到结果。 关键词：潮流计算 电力系统 牛顿-拉夫逊法1.任务及

2、题目要求对如下所诉系统编程进行潮流计算：节点数：3 支路数：3 计算精度：0.00010支路1： 0.0300+j0.0900 12支路2： 0.0200+j0.0900 23支路3： 0.0300+j0.0900 31节点1：PQ节点，S（1）=-0.5000-j0.2000节点2：PQ节点，S（2）=-0.6000-j0.2500节点3：平衡节点，U（3）=1.00000.0000经分析可知题目所给的系统为三节点组成的环形回路，且均为线路没有变压器。有两个为PQ节点，这类节点的有功功率P和无功功率Q是给定的，节点电压是待求量。节点3为平衡节点，在潮流分布算出来以前，网络中的功率损失是未知的

3、，因此，至少要有一个节点的有功功率不能给定，这个节点承担了系统的有功功率平衡，故称之为平衡节点。对系统进行潮流计算，要求各线路的功率分布及功率损耗，未知节点的电压等，并进行编程运行得到结果。2. 设计原理2.1 设计思路本设计给出三条支路阻抗，和相关节点数据 节点1：PQ节点，S（1）=-0.5000-j0.2000 节点2：PQ节点，S（2）=-0.6000-j0.2500 节点3：平衡节点，U（3）=1.00000.0000进行潮流计算，求各线路的功率分布及功率损耗，未知节点的电压。常用的计算潮流分析的方法有牛顿拉夫逊法，PQ分解法等。本次设计采用牛顿拉夫逊法进行计算。牛顿拉夫逊法是一种求

4、解非线性方程的数值解法 ，其要点是把非线性方程式的求解过程变成反复地对相应的线性方程式进行求解的过程，即通常所称的逐次线性化过程。由于便于编写程序用计算机求解，应用较广。 2.1使用方法简介牛顿拉夫逊法（NewtonRaphson法）是求解非线性方程代数方程组的有效迭代计算方法，同时是一种由泰勒展开式只取线性项所得到的线性近似。在牛顿拉夫逊法的每一次迭代过程中，对非线性方程通过线性化处理逐步近似。下面以单变量加以说明。 设欲求解的非线性代数方程为 (2-1) 设方程的真实解为x*，则必有f(x*)=0。用牛顿-拉夫逊法求方程真实解x*的步骤如下： 先给出解的近似值它与真解的误差为，取x

5、5;=x°+x°为第一次的修正估值 ，即 (2-2)其中x°为初始估值的增量，即x°=x¹-x°。设函数f(x)具有任意阶导数，即可将上式在x°的邻域展开为泰勒级数，即： （2-3）式中,分别为函数在处的一阶导数，.，n阶导数。如果差值很小，2-3式右端的一切高阶项均可略去。于是，2-3便简化为 0 (2-4) 这是对于变量的修正量的方程式，。解此方程可得修正量 (2-5)用所求的去修正近似解，变得 (2-6)可见，只要f'(x°)0，即可根据上式求出第一次的修正估值x¹，若恰巧有f(x¹

6、;)=0，则方程的真实解即为x*=x¹。若f(x¹)0，则用上述方法由x¹再确定第二次的修正估值x²， 这样的迭代计算可以反复进行下去，迭代计算的通式是 (2-7)迭代过程的收敛判据为 (2-8)或 (2-9) 式中，为预先给定的小正数。 应用牛顿法求解多变量非线性方程组2-1时，假定已给出各变量的初值，. ，令，. 分别为各变量的修正量，使其满足方程 (2-10) 将上式中的n个多元函数在初始值附近展成泰勒级数,并略去含有，的高阶项，便得到. (2-11) 方程式2-11是对于修正量， 的线性方程组,称为牛顿法的修正方程式.利用高斯消去法或三角分解法可

7、以解出修正量，。然后对初始近似值进行修正 (i=1,2,.,n) (2-13)如此反复迭代，在进行k1次迭代时，从求解修正方程式 (2-14)得到修正量，并对各变量进行修正得到： (2-16) 和 (2-17) 式中的X和分别是由n个变量和修正量组成的n维列向量；F(X)是由n个多元函数组成的n维列项量；J是n阶方阵，称为雅可比矩阵，它的第i、j个元素是第n个函数对第j个变量的偏导数；上角标(k)表示阵的每一个元素都在点处取值。迭代过程一直到满足收敛判据 (2-18) 或 (2-19)为止。和为预先给定的小正数。需要指出的是，用牛顿-拉夫逊法求解非线性代数方程时，巧妙地设定一个合适的初始值是十

8、分重要的，若初始值设定的合适，不但可以保证迭代收敛（即向真实值逼近），而且可使迭代次数减少。相反，若初始值设定的不合适，则不但要增加迭代次数，而且还可能迭代发散（即远离真实值）或者循环（即在真实值附近往复变化）。为解决此问题，可在叠代前先粗略画出函数f(x)是曲线，在该曲线与x轴的焦点附近设定合适的初始估值会有裨益。 3计算步骤用牛顿拉夫逊法计算电力系统潮流的基本步骤：(1)形成节点导纳矩阵；(2)将各节点电压设初值 ；(3)将节点初值代入相关求式，求出修正方程式的常数项向量；(4)将节点电压初值代入求式，求出雅可比矩阵元素；(5)求解修正方程，求修正向量；(6)求取节点电压的新值；(7)检查

9、是否收敛，如不收敛，则以各节点电压的新值作为初值自第3步重新开始进行狭义次迭代，否则转入下一步；(8)计算支路功率分布，节点无功功率和平衡节点功率。具体流程图见下页图二。启动输入原始数据形成节点导纳阵给定电压初值e（0）、f（0）k=0根据公式计算用公式计算雅克比矩阵各元素高斯法解修正方程，求是计算平衡节点功率及全部线路功率输出k=k+1否 图2计算流程图4.结果分析 经过用C语言编程运算得到第一次迭代结果为： 节点导纳矩阵: 6.66667+j-20.00000-3.33333+j10.00000-3.33333+j10.00000-3.33333+j10.00000 5.68627+j-2

10、0.58824-2.35294+j10.58824-3.33333+j10.00000-2.35294+j10.58824 5.68627+j-20.58824 Gauss_Seidel迭代结果： 第1次迭代PQ节点e及f的值：0.983500 -0.019500 第1次迭代PQ节点e及f的值：0.973610 -0.033997 第2次迭代PQ节点e及f的值：0.969642 -0.036485 第2次迭代PQ节点e及f的值：0.965875 -0.042673 Jacobi 矩阵: 19.31452 6.80237-9.57480-3.59699-7.5856018.98468 3.596

11、99-9.57480-9.51651-3.6463119.92097 5.76424 3.64631-9.51651-6.9773519.36506 不平衡值：-0.114256-0.054375-0.002280-0.007412 输出 df,de 节点为 1 的电压不平衡值： df=-0.00508 de=-0.00267 节点为 2 的电压不平衡值： df=-0.00006 de=-0.00037 输出迭代过程中的电压值 节点为 1 的电压修正新值： f=-0.04156 e= 0.96698 节点为 2 的电压修正新值： f=-0.04273 e= 0.96551 还不能达到题目中的精

12、度要求，继续进行迭代，经6次迭代后得到结果。 各节点电压为: U1=0.963186 + j-0.042323 U2=0.962333 + j-0.045732 U3=1.000000 + j0.000000 平衡节点功率为: S3=1.118788+j0.518288 输电线路功率如下: 线路1-2的功率: 0.035453 + j-0.004289 线路2-1的功率: -0.035412 + j0.004412 线路2-3的功率: -0.564587 + j-0.254058 线路3-2的功率: 0.572847 + j0.291225 线路3-1的功率: 0.545942 + j0.22

13、7063 线路1-3的功率: -0.535453 + j-0.195598 输电线路上损耗的功率为: 线路1-2上的功率损耗: 0.000041 + j 0.000123 线路2-3上的功率损耗: 0.008259 + j 0.037167 线路3-1上的功率损耗: 0.010488 + j 0.031465 网络总损耗: 0.018789 + j 0.068755 满足精度要求，计算出了潮流计算要求的各节点电压，网络中的功率分布及功率损耗。完成和达到了设计的要求和目的。小结在本次电力系统分析课程实习中，我感觉到考试仅仅体现了学习能力的一部分，而课程设计更能体现了综合的学习实践能力。原本模糊不

14、清的基本概念和实用方法，在操作本次课设的过程中都得到了明确和加强。课程设计让我感受到比考试更大的压力，有了更高的效率和收获，深层次加深了知识的深度层次。望进行更多一些这样的实践性自主行的学习活动，受益良多。在本次课程设计得到了很多老师和同学的支持和帮助，特在此表示感谢，不再一一署名。参考文献：1 何仰赞 温增银 电力系统分析 华中科技大学出版社 20062 于永源 杨绮雯 电力系统分析（第二版）中国电力出版社 20073 大电网最优潮流计算 刘明波、谢敏、 赵维兴 科学出版社 20104 电力系统分析的计算机算法 邱晓燕、 刘天琪 中国电力出版社 20095 薛定宇 控制系统计算机辅助设计 清

15、华大学出版社 1996附录：源程序#include <stdio.h>#include <math.h>#define M 40 /*矩阵阶数*/#define N 10 /*最大迭代次数*/ int i,j,k,a,b,l,z; /* 循环变量 */ int n, /* 节点数 */ m, /* 支路数 */ dd, /*对地支路数*/ pq, /* PQ节点数 */ pv, /* PV节点数 */ byq; /* 变压器数 */ float eps, /* 精度 */ bb, /* 变比 */ max; /*指标标志*/ float yyM; /*中间变量*/ st

16、atic float GMM,BMM,B1MM,B2MM,B3MM; struct jiedian /* 定义节点结构体 */ int num,s; /* num为节点号，s为节点类型*/ float p,q,u,f,v; float dp,dq,df,du; jiedianM; struct zhilu /* 定义支路结构体 */ int num; int p1,p2; /*支路的两个节点*/ float r,x; /*支路的电阻与电抗*/ zhiluM; FILE *fp1,*fp2; /*文件指针*/ void data1() /* 读取数据函数 */ int h,numb,a,b; a

17、=1; fp1=fopen("d:input.txt","r"); if(fp1=NULL) printf(" can not open file !n"); exit(0); fscanf(fp1,"%d,%d,%d,%d,%d,%d,%fn",&n,&m,&dd,&pq,&pv,&byq,&eps); /*输入节点数,支路数,对地支路数,PQ节点数，PV节点数，变压器数和精度*/ j=1;k=pq+1; for(i=1;i<=n;i+) /*输入节点

18、类型的输入功率和节电电压初值*/ fscanf(fp1,"%d,%d",&numb,&h);if(h=1) /*类型h=1是PQ节点*/ fscanf(fp1,",%f,%f,%f,%fn",&jiedianj.p,&jiedianj.q,&jiedianj.u,&jiedianj.f); /*输入PQ节点的数据*/ jiedianj.num=numb; jiedianj.s=h; j+; if(h=2) /*类型h=2是pv节点*/ fscanf(fp1,",%f,%f,%f,%fn"

19、,&jiediank.p,&jiediank.v,&jiediank.u,&jiediank.f); /*输入PV节点的数据*/ jiediank.num=numb; jiediank.s=h; k+; if(h=3) /*类型h=3是平衡节点*/ fscanf(fp1,",%f,%fn",&jiediann.u,&jiediann.f); /*输入平衡节点的数据*/ jiediann.num=numb; jiediann.s=h; for(i=1;i<=m;i+) /*输入支路阻抗*/ fscanf(fp1,"

20、;%d,%d,",&zhilui.num,&h);if(h=0) /*输入非变压器支路阻抗*/fscanf(fp1,"%d,%d,%f,%fn", &zhilui.p1,&zhilui.p2,&zhilui.r,&zhilui.x); if(h=1) /*输入变压器支路阻抗*/ fscanf(fp1,"%f,%d,%d,%f,%fn",&bb,&zhilui.p1,&zhilui.p2,&zhilui.r, &zhilui.x); zhilui.r=zhil

21、ui.r*bb; zhilui.x=zhilui.x*bb; b=m+a; zhilub.num=b; zhilub.p1=zhilui.p1; zhilub.p2=0; zhilub.r=zhilui.r*bb/(1-bb); zhilub.x=zhilui.x*bb/(1-bb); zhilub+1.num=b zhilub+1.p1=zhilui.p2; zhilub+1.p2=0; zhilub+1.r=zhilui.r/(bb-1); zhilub+1.x=zhilui.x/(bb-1); a=a+2; fclose(fp1); if(fp2=fopen("d:output

22、.txt","w")=NULL) printf(" can not open file!n"); exit(0); fprintf(fp2,"nn * 原始数据 *n"); fprintf(fp2," =n") fprintf(fp2," 节点数:%2d 支路数:%2d 对地支路数:%2d 变压器数:%2d n PQ节点数:%2d PV节点数:%2d 精度:%fn",n,m+2*byq,dd,byq,pq,pv,eps); fprintf(fp2," -n"); f

23、or(i=1;i<=pq;i+) fprintf(fp2," PQ节点 节点%2d P%d=%+f Q%d=%+fn", jiediani.num,jiediani.num,jiediani.p,jiediani.num,jiediani.q); for(i=pq+1;i<=pq+pv;i+) fprintf(fp2," PV节点 节点%2d P%d=%+f V%d=%+fn", jiediani.num,jiediani.num,jiediani.p,jiediani.num,jiediani.v); fprintf(fp2,"

24、平衡节点 节点%2d u%d=%+f f%d=%+fn", jiediann.num,jiediann.num,jiediann.u,jiediann.num,jiediann.f); fprintf(fp2," -n"); for(i=1;i<=m+2*byq;i+) fprintf(fp2," 支路%2d 相关节点:%2d,%2d R=%+f X=%+fn", i,zhilui.p1,zhilui.p2,zhilui.r,zhilui.x); fprintf(fp2," =n"); void form_y() /*

25、 形成节点导纳矩阵函数 */ float S; for(i=1;i<=m+2*byq;i+) /*节点导纳矩阵的主对角线上的导纳*/ for(j=0;j<=n;j+) if(zhilui.p1=j)|(zhilui.p2=j) S=zhilui.r*zhilui.r+zhilui.x*zhilui.x; if(S=0) continue; Gjj+=zhilui.r/S; Bjj+=-zhilui.x/S; for(i=1;i<=m+2*byq;i+) /*对地导纳*/ j=zhilui.p1; k=zhilui.p2; if(j=0)|(k=0) S=zhilui.r*zh

26、ilui.r+zhilui.x*zhilui.x; if(S=0) G0i=B0i=Gi0=Bi0=0;continue; if(j=0) G0k+=zhilui.r/S; B0k+=-zhilui.x/S; Gk0=G0k; Bk0=B0k; if(k=0) G0j+=zhilui.r/S; B0j+=-zhilui.x/S; Gj0=G0j; Bj0=B0j; for(k=1;k<=m-dd;k+) /*节点导纳矩阵非主对角线上的导纳*/ i=zhiluk.p1; j=zhiluk.p2;S=zhiluk.r*zhiluk.r+zhiluk.x*zhiluk.x;if(S=0) co

27、ntinue;Gij+=-zhiluk.r/S;Bij+=zhiluk.x/S;Gji=Gij;Bji=Bij; for(i=1;i<=(pq+pv);i+) /*形成B'矩阵*/ for(j=1;j<=(pq+pv);j+) B1ij=Bjiediani.numjiedianj.num; for(i=1;i<=pq;i+) /*形成B''矩阵*/ for(j=1;j<=pq;j+) B2ij=Bjiediani.numjiedianj.num; /*输出节点导纳矩阵、系数矩阵B'、系数矩阵B''*/ fprintf(fp

28、2,"nn * 计算结果 *n"); fprintf(fp2,"=n") fprintf(fp2,"n 节点导纳矩阵为:"); for(i=1;i<=n;i+) fprintf(fp2,"n "); for(j=1;j<=n;j+) fprintf(fp2,"%+8.6f+j%+8.6f ",Gij,Bij); fprintf(fp2,"n -n"); fprintf(fp2,"n 节点导纳矩阵实部G为:"); for(i=1;i<=n;

29、i+) fprintf(fp2,"n "); for(j=1;j<=n;j+) fprintf(fp2,"%+8.6f ",Gij); fprintf(fp2,"n -n"); fprintf(fp2,"n 节点导纳矩阵虚部B为:"); for(i=1;i<=n;i+) fprintf(fp2,"n "); for(j=1;j<=n;j+) fprintf(fp2,"%+8.6f ",Bij); fprintf(fp2,"n -n");

30、fprintf(fp2,"n 系数矩阵B'为:"); for(i=1;i<=pv+pq;i+) fprintf(fp2,"n "); for(j=1;j<=pv+pq;j+) fprintf(fp2,"%+8.6f ",B1ij); fprintf(fp2,"n -n"); fprintf(fp2,"n 系数矩阵B''为:"); for(i=1;i<=pq;i+) fprintf(fp2,"n "); for(j=1;j<=pq

31、;j+) fprintf(fp2,"%+8.6f ",B2ij); fprintf(fp2,"n =n"); void gaoshi(int A,float aaMM,float bbM) /*高斯计算函数*/int i,j,z; int k=0; float a; for(i=k+1;i<A;i+) for(z=0;(z+i)<A;z+) aai+zk=aai+zk/aakk; /*校验过程*/ bbi+z=bbi+z-aai+zk*bbk; for(j=k+1;j<A;j+) aai+zj=aai+zj-aai+zk*aakj; k

32、=k+1; bbA-1=bbA-1/aaA-1A-1; for(i=A-2,k=A-2;i>=0;i-) a=0; for(j=i+1;j<A;j+) a=aaij*bbj+a; bbi=(bbi-a)/aaii; float solve(float max) /*迭代函数*/ float aaM,bbM;/*求节点有功功率不平衡量dp*/ for(i=1;i<=(pq+pv);i+) /*求dp*/ jiediani.dp=jiediani.p; for(j=1;j<=n;j+) jiediani.dp=jiediani.dp-jiediani.u*jiedianj.

33、u*(Gjiediani.numjiedianj.num*cos(jiediani.f-jiedianj.f)+Bjiediani.numjiedianj.num*sin(jiediani.f-jiedianj.f); for(i=1;i<=(pq+pv);i+) /*定义各节点的dP/U*/ if(jiediani.s!=3) yyi=fabs(jiediani.dp/jiediani.u); max=yy1; /*求dP/U的最大值*/ for(i=1;i<=(pq+pv);i+) if(yyi>max)max=yyi; for(i=1;i<=(pq+pv);i+)

34、 fprintf(fp2," 有功功率不平衡量dP%d=%+fn",jiediani.num,jiediani.dp); fprintf(fp2,"n -n"); /*求各节点的相位角f，df*/ a=1; /*装换为高斯函数的标准形式*/ for(i=1;i<=n;i+) if(jiediani.s=3)a=a+1;continue; aai-a=jiediani.dp/jiediani.u; b=1; for(j=1;j<=n;j+) if(jiedianj.s=3)b=b+1;continue; B3i-aj-b=0.0-B1i+1-a

35、j+1-b*jiediani.u; gaoshi(pv+pq),B3,aa); /*代入高斯函数求解方程*/ a=1; for(i=1;i<=n;i+) if(jiediani.s=3)a+;continue; jiediani.df=aai-a; jiediani.f=jiediani.f+jiediani.df; for(i=1;i<=(pq+pv);i+) fprintf(fp2," df%d=%f 电压的相位角f%d=%+fn",jiediani.num,jiediani.df,jiediani.num,jiediani.f); fprintf(fp2,

36、"n -n"); /*求节点的无功功率不平衡量dq*/ for(i=1;i<=n;i+) /*求dq*/ jiediani.dq=jiediani.q; for(j=1;j<=n;j+) jiediani.dq=jiediani.dq-jiediani.u*jiedianj.u*(Gjiediani.numjiedianj.num*sin(jiediani.f-jiedianj.f)-Bjiediani.numjiedianj.num*cos(jiediani.f-jiedianj.f); for(i=1;i<=pq;i+) /*定义各节点的dQ/U*/

37、yyi=fabs(jiediani.dq/jiediani.u); for(i=1;i<=pq;i+) /*求所有dP/U和dQ/U中的最大值*/ if(yyi>max)max=yyi; for(i=1;i<=pq;i+) fprintf(fp2," 无功功率不平衡量dQ%d=%+fn",jiediani.num,jiediani.dq); fprintf(fp2,"n -n"); /*计算各节点电压*/ for(i=1;i<=n;i+) /*装换为高斯函数的标准形式*/ for(j=1;j<=n;j+) B3i-1j-1=

38、0-B2ij; a=1; for(i=1;i<=n;i+) if(jiediani.s!=1)a=a+1;continue; bbi-a=jiediani.dq/jiediani.u; gaoshi(pq,B3,bb); /*代入高斯函数求解方程*/ a=1; for(i=1;i<=n;i+) if(jiediani.s=3)a+;continue; jiediani.du=bbi-a; jiediani.u=jiediani.u+jiediani.du; for(i=1;i<=n;i+) if(jiediani.s=1) fprintf(fp2," dU%d=%f

39、 电压的大小U%d=%+fn",jiediani.num,jiediani.du,jiediani.num,jiediani.u); fprintf(fp2,"n =n"); return(max); /*返回循环标志值*/float mul_Re(x1,y1,x2,y2) /*求实部*/float x1,x2,y1,y2;float x; x=x1*x2-y1*y2; return(x); float mul_Im(x1,y1,x2,y2) /*求虚部*/float x1,x2,y1,y2; float y; y=x1*y2+x2*y1; return(y);

40、void data2() /* 潮流计算结果 */ float sp=0,sq=0,mo1,mo2,mo3,mo4,mo5,mo6,ff,ttM; int i1,j1,i11=0,j11=0; static float PMM,QMM; fprintf(fp2,"nn 各节点电压为:"); for(i=1;i<=n;i+) ttM=jiediani.u;jiediani.u=(jiediani.u)*(cos(jiediani.f); /*由极坐标转化为直角坐标*/ jiediani.f=(ttM)*(sin(jiediani.f); fprintf(fp2,"n U%d=%+8.6f + j %+8.6f", jiediani.num,jiediani.u,jiediani.f); /*输出直角坐标节点电压*/ fprintf(fp2,"n =n"); fprintf(fp2,"nn 平衡节点功率为:"); for(i=1;i<=n;i+) i1=jiediann.num; j1=jiediani.num; sp+=mul_Re(Gi1j1,-Bi1j1,jiediani.u,-jie