潮流计算.docx
- 文档编号:9244902
- 上传时间:2023-05-17
- 格式:DOCX
- 页数:19
- 大小:17.91KB
潮流计算.docx
《潮流计算.docx》由会员分享,可在线阅读,更多相关《潮流计算.docx(19页珍藏版)》请在冰点文库上搜索。
潮流计算
#include"stdio.h"
#include"stdlib.h"
#include"string.h"
#include"math.h"
#defineNBUS5
#defineNLINE7
/*Globalvariables*/
intnL,nSH,nBUS,nVA;
doubled,t;
intL;
structLine
{
intNum,NumI,NumJ;
doubleR;
doubleX;
doubleB;
doubleK;
};
structBus
{
intNum;
doubleVolt,Phase,GenP,GenQ,LoadP,LoadQ;
intType;
};
structShunt
{
intNum,NumI;
doubleG,B;
};
voidmain()
{FILE*fp;
FILE*fpout;
inti,j,k,l,h,c;
inti1,i2,i3,kp,kq,LH[100];
doubled,d1,d2,d3,d4,d5,d6,t,r,x,g,b,tt,LL,e,ps,qs,shsh,sum,w;
structLinesL[NLINE];
structBussBus[NBUS];
structShuntsSH[NBUS];
doubleYG[NBUS][NBUS],YB[NBUS][NBUS];
doubleV[100][2],U[100],dPQ[100],PQ[100][2],JJ[100][100]={0},H[100][100],J[100][100],N[100][100],L[100][100],sP[NBUS][NBUS]={0},sQ[NBUS][NBUS]={0},dsp,dsq,sumgen,dp;
i1=i2=i3=0;d1=d2=d3=d4=d5=d6=ps=qs=0.0;
e=0.00001;
for(i=0;i /*Readtheinputdata*/ if((fp=fopen("data.txt","r"))==NULL) {printf("Cannotopenthefilenamed'in.txt'\n"); exit(0); } fscanf(fp,"%d,%d,%d",&nBUS,&nL,&nSH); for(i=0;i sBus[i].Num=sBus[i].Type=0;sBus[i].Volt=1.0; sBus[i].Phase=sBus[i].GenP=sBus[i].GenQ=sBus[i].LoadP=sBus[i].LoadQ=0.0; fscanf(fp,"%d,%lf,%lf,%lf,%lf,%lf,%lf,%d",&i1,&d1,&d2,&d3,&d4,&d5,&d6,&i2); sBus[i].Num=i1;sBus[i].Volt=d1;sBus[i].Phase=d2;sBus[i].GenP=d3;sBus[i].GenQ=d4;sBus[i].LoadP=d5,sBus[i].LoadQ=d6;sBus[i].Type=i2; }; for(i=0;i sL[i].Num=sL[i].NumI=sL[i].NumJ=0; sL[i].R=sL[i].X=sL[i].B=0.0;sL[i].K=1.0; fscanf(fp,"%2d%3d%3d%lf%lf%lf%lf",&i1,&i2,&i3,&d1,&d2,&d3,&d4);sL[i].Num=i1;sL[i].NumI=i2;sL[i].NumJ=i3;sL[i].R=d1;sL[i].X=d2;sL[i].B=d3;sL[i].K=d4; } for(i=0;i sSH[i].Num=sSH[i].NumI=0;sSH[i].G=sSH[i].B=0.0; fscanf(fp,"%2d%3d%lf",&i1,&i2,&d1); sSH[i].Num=i1;sSH[i].NumI=i2;sSH[i].B=d1;} if(fp! =NULL)fclose(fp); /*MakeYMatrix*/ for(i=1;i YG[i][j]=0.0;YB[i][j]=0.0;}; for(l=0;l { i=sL[l].NumI; j=sL[l].NumJ; r=sL[l].R; x=sL[l].X; d1=r*r+x*x; g=r/d1; b=-x/d1; if(fabs(sL[l].K-1.0)<0.000001) {/*Normallinesortransformers*/ YG[i][i]=YG[i][i]+g; YG[j][j]=YG[j][j]+g; YB[i][i]=YB[i][i]+b+sL[l].B; YB[j][j]=YB[j][j]+b+sL[l].B; YG[i][j]=YG[i][j]-g; YG[j][i]=YG[j][i]-g; YB[i][j]=YB[i][j]-b; YB[j][i]=YB[j][i]-b; } else {/*abnormaltransformerratio*/ if(fabs(sL[l].B)>0.000001) { YG[i][i]=YG[i][i]+g/sL[l].B/sL[l].B; YB[i][i]=YB[i][i]+b/sL[l].B/sL[l].B; YG[j][j]=YG[j][j]+g; YB[j][j]=YB[j][j]+b; YG[i][j]=YG[i][j]-g/sL[l].B; YG[j][i]=YG[j][i]-g/sL[l].B; YB[i][j]=YB[i][j]-b/sL[l].B; YB[j][i]=YB[j][i]-b/sL[l].B; } else printf("节点%2d,%-2d之间变压器变比为0\n",i,j); } } //ChecktheYmatrix if((fp=fopen("GGBB.txt","w"))==NULL){ printf("Cannotopenthefilenamed'GGBB.txt'\n");exit(0);} fprintf(fp,"---YMatrix---\n");for(i=1;i if(fp! =NULL)fclose(fp); /*设定电压初值*/ for(i=1;i if(sBus[i-1].Type==0) { V[i][0]=0.0; V[i][1]=1.0; } for(i=1;i if(sBus[i-1].Type==1) { V[i][1]=sBus[i-1].Volt; V[i][0]=0.0; } for(i=1;i if(sBus[i-1].Type==2) {V[i][1]=sBus[i-1].Volt; V[i][0]=sBus[i-1].Phase; } /*输出电压初值*/ if((fp=fopen("d: \\lx\\电压初值.txt","w"))==NULL) { printf("Cannotopenthefilenamed'电压初值.txt'\n");exit(0);} fprintf(fp,"---电压初值---\n"); for(i=1;i for(j=1;j<2;j++) fprintf(fp,"Y(%2d)=(%10.5f,%10.5f)\n",i,V[i][0],V[i][1]);if(fp! =NULL)fclose(fp);for(c=1;;c++) /*计算偏移量*/ for(i=1;i { if(sBus[i-1].Type! =2) { PQ[i][0]=0; PQ[i][1]=0; for(j=1;j PQ[i][0]-=V[i][1]*V[j][1]*(YG[i][j]*cos(V[i][0]-V[j][0])+YB[i][j]*sin(V[i][0]-V[j][0]));PQ[i][1]-=V[i][1]*V[j][1]*(YG[i][j]*sin(V[i][0]-V[j][0])-YB[i][j]*cos(V[i][0]-V[j][0])); } } } h=1; for(i=1;i { if(sBus[i-1].Type==0) { dPQ[h]=PQ[i][0]+sBus[i-1].GenP-sBus[i-1].LoadP; h++; dPQ[h]=PQ[i][1]+sBus[i-1].GenQ-sBus[i-1].LoadQ; h++; } if(sBus[i-1].Type==1) { dPQ[h]=PQ[i][0]+sBus[i-1].GenP-sBus[i-1].LoadP; h++; } } /*输出偏移量 if((fp=fopen("d: \\lx\\偏移量.txt","w"))==NULL) { printf("Cannotopenthefilenamed'偏移量.txt'\n");exit(0);} fprintf(fp,"---偏移量---\n"); for(i=1;i fprintf(fp,"dPQ(%2d)=(%10.5f)\n",i,dPQ[i]); if(fp! =NULL)fclose(fp);*/ /*计算雅克比矩阵*/ for(i=1;i { for(j=1;j { if(i==j) { H[i][j]=PQ[i][1]-V[i][1]*V[i][1]*YB[i][j];J[i][j]=-PQ[i][0]-V[i][1]*V[i][1]*YG[i][j]; N[i][j]=-PQ[i][0]+V[i][1]*V[i][1]*YG[i][j];L[i][j]=-PQ[i][1]-V[i][1]*V[i][1]*YB[i][j]; } else{H[i][j]=V[i][1]*V[j][1]*(YG[i][j]*sin(V[i][0]-V[j][0])-YB[i][j]*cos(V[i][0]-V[j][0]));J[i][j]=-V[i][1]*V[j][1]*(YG[i][j]*cos(V[i][0]-V[j][0])+YB[i][j]*sin(V[i][0]-V[j][0]));N[i][j]=V[i][1]*V[j][1]*(YG[i][j]*cos(V[i][0]-V[j][0])+YB[i][j]*sin(V[i][0]-V[j][0]));L[i][j]=V[i][1]*V[j][1]*(YG[i][j]*sin(V[i][0]-V[j][0])-YB[i][j]*cos(V[i][0]-V[j][0])); } } } h=1; for(i=1;i { l=1; if(sBus[i-1].Type==0) { for(j=1;j { if(sBus[j-1].Type==0) { JJ[h][l]=H[i][j]; l++; JJ[h][l]=N[i][j]; l++; } if(sBus[j-1].Type==1) { JJ[h][l]=H[i][j]; l++; } } h++; l=1; for(j=1;j { if(sBus[j-1].Type==0) { JJ[h][l]=J[i][j]; l++; JJ[h][l]=L[i][j]; l++; } if(sBus[j-1].Type==1) { JJ[h][l]=J[i][j]; l++; } } h++; } if(sBus[i-1].Type==1) { for(j=1;j { if(sBus[j-1].Type==0) { JJ[h][l]=H[i][j]; l++; JJ[h][l]=N[i][j]; l++; } if(sBus[j-1].Type==1) { JJ[h][l]=H[i][j]; l++; } } h++; } } /*输出雅克比矩阵 if((fp=fopen("d: \\lx\\雅克比矩阵.txt","w"))==NULL) { printf("Cannotopenthefilenamed'雅克比矩阵.txt'\n"); exit(0); } fprintf(fp,"---雅克比矩阵---\n");fprintf(fp,""); for(i=1;i { fprintf(fp,"%d\t\t",i); } fprintf(fp,"\n");*/ for(i=1;i { fprintf(fp,"%2d",i); for(j=1;j { fprintf(fp,"%10.5f\t",JJ[i][j]); } fprintf(fp,"\n"); } if(fp! =NULL)fclose(fp); /*高斯法求解方程组*/ l=1; for(i=1;i {LH[i]=0;} for(k=1;k {d=0.0; for(j=k;j { if(fabs(JJ[k][j])>d) { d=fabs(JJ[k][j]); /*在一行中找到一个最大值赋值d,并用JS[K]记住这个最大值所在的列号*/LH[k]=j; } } if(fabs(d)<0.000001) /*如果d的数值太小,做为被除数将带来很大的误差*/ {l=0;} elseif(LH[k]! =k) { for(i=1;i { t=JJ[i][k]; JJ[i][k]=JJ[i][LH[k]]; /*进行列交换,让最大值始终在对角元上*/ JJ[i][LH[k]]=t; } } if(l==0) {break;} for(j=k+1;j { JJ[k][j]=JJ[k][j]/JJ[k][k]; /*对角元上的元素消为1*/ } dPQ[k]=dPQ[k]/JJ[k][k]; for(i=k+1;i { for(j=k+1;j { JJ[i][j]=JJ[i][j]-JJ[i][k]*JJ[k][j]; /*使下三角阵的元素为0*/ } dPQ[i]=dPQ[i]-JJ[i][k]*dPQ[k]; } } if(fabs(JJ[h-1][h-1])>0.00001) { /*用追赶法,解方程组,求未知数x*/U[h-1]=dPQ[h-1]; for(i=h-2;i>=0;i--) { t=0.0; for(j=i+1;j { t=t+JJ[i][j]*U[j]; } U[i]=(dPQ[i]-t); } } /*输出高斯结果 if((fp=fopen("d: \\lx\\高斯.txt","w"))==NULL){printf("err");exit(0);} for(i=1;i { fprintf(fp,"%f",U[i]);fprintf(fp,"\n");} fclose(fp);*/ /*得到电压值*/ h=1; for(i=1;i { if(sBus[i-1].Type==0) { V[i][0]+=U[h]; h++; V[i][1]+=U[h]; h++; } if(sBus[i-1].Type==1) { V[i][0]+=U[h]; h++; } } /*输出电压值 if((fp=fopen("d: \\lx\\电压值.txt","w"))==NULL) { printf("Cannotopenthefilenamed'电压值.txt'\n");exit(0);} fprintf(fp,"---电压初值---\n"); for(i=1;i for(j=1;j<2;j++) fprintf(fp,"Y(%2d)=(%10.5f,%10.5f)\n",i,V[i][0],V[i][1]); if(fp! =NULL)fclose(fp);*/ /*求最大变化值*/ w=0; for(i=1;i { if(U[i]>w) { w=U[i]; } } if(w<0.00001) { break; } } /*输出电压终值 if((fp=fopen("d: \\lx\\电压终值.txt","w"))==NULL) { printf("Cannotopenthefilenamed'电压终值.txt'\n");exit(0); } fprintf(fp,"---电压终值---\n"); fprintf(fp,"循环%d次\n",c); for(i=1;i for(j=1;j<2;j++) printf(fp,"Y(%2d)=(%10.5f,%10.5f)\n",i,V[i][0]*180/3.1415926,V[i][1]); if(fp! =NULL)fclose(fp);*/ ps=0;qs=0; for(i=1;i { if(sBus[i-1].Type==2) { for(j=1;j { ps+=V[i][1]*V[j][1]*(YG[i][j]*cos(V[i][0]-V[j][0])+YB[i][j]*sin(V[i][0]-V[j][0]));qs+=V[i][1]*V[j][1]*(YG[i][j]*sin(V[i][0]-V[j][0])-YB[i][j]*cos(V[i][0]-V[j][0])); } } } for(i=1;i { for(j=1;j {sP[i][j]=0;sQ[i][j]=0;} } for(l=0;l {i=sL[l].NumI; j=sL[l].NumJ; r=sL[l].R; x=sL[l].X; d1=r*r+x*x; g=r/d1; b=-x/d1; if(fabs(sL[l].K-1.0)<0.000001) {/*Normallinesortransformers*/sP[i][j]=V[i][1]*V[i][1]*g-V[i][1]*V[j][1]*(g*cos(V[i][0]-V[j][0])+b*sin(V[i][0]-V[j][0]));sQ[i][j]=-(V[i][1]*V[i][1]*sL[l].B+V[i][1]*V[i][1]*b+V[i][1]*V[j][1]*(g*sin(V[i][0]-V[j][0])-b*cos(V[i][0]-V[j][0])));sP[j][i]=V[j][1]*V[j][1]*g-V[i][1]*V[j][1]*(g*cos(V[j][0]-V[i][0])+b*sin(V[j][0]-V[i][0]));sQ[j][i]=-(V[j][1]*V[j][1]*sL[l].B+V[j][1]*V[j][1]*b+V[i][1]*V[j][1]*(g*sin(V[j][0]-V[i][0])-b*cos(V[j][0]-V[i][0]))); } else {/*abnormaltransformerratio*/sP[i][j]=V[i][1]*V[i][1]*g/sL[l].B/sL[l].B-V[i][1]*V[j][1]*(g*cos(V[i][0]-V[j][0])/sL[l].B+b*sin(V[i][0]-V[j][0])/sL[l].B);sQ[i][j]=-(V[i][1]*V[i][1]*b/sL[l].B/sL[l].B+V[i][1]*V[j][1]*(g*sin(V[i][0]-V[j][0])/sL[l].B-b*cos(V[i][0]-V[j][0])/sL[l].B));sP[j][i]=V[j][1]*V[j][1]*g-V[i][1]*V[j][1]*(g*cos(V[j][0]-V[i][0])/sL[l].B+b*sin(V[j][0]-V[i][0])/sL[l].B);sQ[j][i]=-(V[i][1]*V[i][1]*b+V[i][1]*V[j][1]*(g*sin(V[j][0]-V[i][0])/sL[l].B-b*cos(V[j][0]-V[i][0])/sL[l].B));
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 潮流 计算