《电力系统潮流上机》课程设计报告
院 系:电气与电子工程学院
班 级: 电气1108
学 号: 1111550112
学生姓名: 龙日尚
指导教师: 刘宝柱
设计周数: 两周
成 绩:
日期:20##年1月10日
一、课程设计的目的与要求
培养学生的电力系统潮流计算机编程能力,掌握计算机潮流计算的相关知识
二、设计正文(详细内容见附录)
1. 手算
2. 计算机计算
3. 思考题
三、课程设计总结或结论(详细内容见附录)
四、参考文献
1. 《电力系统计算:电子数字计算机的应用》,西安交通大学等合编。北京:水利电力出版社;
2. 《现代电力系统分析》,王锡凡主编,科学出版社;
3. 《电力系统稳态分析》,陈珩,中国电力出版社,1995年,第三版;
附录(设计流程图、程序、表格、数据等)
4. 机算潮流程序及结果
// dierti.cpp : 定义控制台应用程序的入口点。
//
#include "stdafx.h"
struct Line //线路结构体
{
int Num,NumI,NumJ; //线路号 左节点名 右节点名
float R,X,B,K; //电阻 电抗 电纳 变比(K等于1为普通支路, 不等于1为变压器支路的变比)
};
struct Bus //节点结构体
{
int Num ;
float Volt,Phase,GenP,GenQ,LoadP,LoadQ;
int Type;
};
#include"stdio.h"
#include"string.h"
#include"math.h"
#include"stdlib.h"
#define NBUS 4
#define NLINE 4
/* Global variables */
int nL,nB,nVA,nSH;
float X[NBUS];
int L;
double def[2*NBUS];
double mn[50];
void Gauss(double a[50][50],double b[50], int n) /*定义高斯法 */
{
int JS[50];
int i,j,k;
float d,t,x[50];
FILE *fp;
int L=1;
for(i=0;i<50;i++) JS[i]=0;
for(k=0;k<n;k++){
d=0.0;
for(j=k;j<n;j++)
if(fabs(a[k][j])>d){ /*在一行中找到一个最大值赋值d,并用JS[K]记住这个最大值所在的列号*/
d=fabs(a[k][j]);
JS[k]=j;
}
if(fabs(d)<0.000001) /*如果d的数值太小,做为被除数将带来很大的误差 */
L=0;
else {
if(JS[k]!=k)
for(i=0;i<n;i++)
{
t=a[i][k];
a[i][k]=a[i][JS[k]]; /*进行列交换,让最大值始终在对角元上*/
a[i][JS[k]]=t;
}
}
if(L==0)break;
for(j=k+1;j<n;j++)a[k][j]=a[k][j]/a[k][k]; /*对角元上的元素消为1*/
b[k]=b[k]/a[k][k];
for(i=k+1;i<n;i++){
for(j=k+1;j<n;j++)a[i][j]=a[i][j]-a[i][k]*a[k][j]; /*使下三角阵的元素为0*/
b[i]=b[i]-a[i][k]*b[k];
}
}
if(fabs(a[n-1][n-1])>0.00001){ /*用追赶法,解方程组,求未知数x*/
x[n-1]=b[n-1];
for(i=n-2;i>=0;i--){
t=0.0;
for(j=i+1;j<n;j++)t=t+a[i][j]*x[j];
x[i]=(b[i]-t);
}
}
if((fp=fopen("gauss.txt","w"))==NULL) /*将结果写到TXT文件中*/
{printf("err");exit(0);}
for(i=0;i<n;i++){
fprintf(fp,"%lf",x[i]);
mn[i]=x[i];
fprintf(fp,"\n");}
fclose(fp);
if(fp!=NULL) fclose(fp);
}
int _tmain(int argc, _TCHAR* argv[])
{
FILE *fp;
FILE *fpout;
int i,j,k,l,h,n,v;
int i1,i2,i3,kp,kq;
float d1,d2,d3,d4,d5,d6,r,x,g,b,tt,LL,e,ps,qs,shsh,m;
struct Line sL[NLINE];
struct Bus sB[NBUS];
float YG[NBUS+1][NBUS+1],YB[NBUS+1][NBUS+1];
double u[50][2];
i1=i2=i3=0;
d1=d2=d3=d4=d5=d6=ps=qs=0.0;
for(i=0;i<NBUS;i++)
if((fp=fopen("in.txt","r"))==NULL)
{ printf("Can not open the file named 'in.txt' \n");
exit(0);
}
fscanf(fp,"%d,%d,%d",&nB,&nL,&nSH);
for(i=0;i<nB;i++){
sB[i].Num=sB[i].Type=0;sB[i].Volt=1.0;
sB[i].Phase=sB[i].GenP=sB[i].GenQ=sB[i].LoadP=sB[i].LoadQ=0.0;
fscanf(fp,"%d,%f,%f,%f,%f,%f,%f,%d",&i1,&d1,&d2,&d3,&d4,&d5,&d6,&i2);
sB[i].Num=i1;sB[i].Volt=d1;sB[i].Phase=d2;sB[i].GenP=d3;sB[i].GenQ=d4;sB[i].LoadP=d5,sB[i].LoadQ=d6;sB[i].Type=i2;
};
for(i=0;i<nL;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 %f %f %f %f",&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;
}
if(fp!=NULL) fclose(fp);
/*Make Y Matrix*/
for(i=1;i<nB+1;i++)for(j=1;j<nB+1;j++){
YG[i][j]=0.0;
YB[i][j]=0.0;
};
for(l=0; l<nL; 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;
m=sL[l].K;
if(fabs(sL[l].K-1.0)<0.000001) //普通支路
{
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 //变压器支路
{YG[i][i]=YG[i][i]+g/m+g*(m-1)/m;
YG[j][j]=YG[j][j]+g/m+g*(1-m)/m/m;
YB[i][i]=YB[i][i]+b/m+b*(m-1)/m;
YB[j][j]=YB[j][j]+b/m+b*(1-m)/m/m;
YG[i][j]=YG[i][j]-g/m;
YG[j][i]=YG[j][i]-g/m;
YB[i][j]=YB[i][j]-b/m;
YB[j][i]=YB[j][i]-b/m; }
}
/* Check the Y matrix */
if((fp=fopen("GGBB.txt","w"))==NULL){
printf("Can not open the file named 'GGBB.txt' \n");exit(0);}
fprintf(fp,"---Y Matrix---\n");
for(i=1;i<nB+1;i++)for(j=1;j<nB+1;j++)if(fabs(YB[i][j]-0.0)>0.000001)
fprintf(fp,"Y(%3d,%-3d)=(%10.5f,%10.5f)\n",i,j,YG[i][j],YB[i][j]);
if(fp!=NULL) fclose(fp);
/* 节点电压附初值 */
for(i=1;i<nB+1;i++)
{if(sB[i-1].Type==0)
{u[i][0]=0.0;
u[i][1]=1.0;
}
else if(sB[i-1].Type==1)
{u[i][1]=sB[i-1].Volt;
u[i][0]=0.0;
}
else if(sB[i-1].Type==2)
{u[i][1]=sB[i-1].Volt;
u[i][0]= sB[i-1].Phase;
}
}
for(v=1;;v++)/* 迭代次数可以无限大 */
{
/* 节点电压附初值 */
printf("迭代第%d次赋予的电压初值为e+jf:\n",v);
for(i=1;i<nB+1;i++)
printf("%lf,%lf\n",u[i][1],u[i][0]);
printf("\n");
printf("\n");
/* 求偏移量 */
double P_P[10];
double P_Q[10];
double P_UU[10];
for(i=1;i<nB+1;i++)
{if(sB[i-1].Type==2)
{P_P[i]=0.0;
P_Q[i]=0.0;
P_UU[i]=1.05;}
if(sB[i-1].Type==0)
{double tempP=0.0;
double tempQ=0.0;
for(j=1;j<nB+1;j++)
{
tempP+=YG[i][j]*u[j][1]-YB[i][j]*u[j][0];
tempQ+=YG[i][j]*u[j][0]+YB[i][j]*u[j][1];
}
P_P[i]=(sB[i-1].GenP-sB[i-1].LoadP)-tempP*u[i][1]-tempQ*u[i][0];
P_Q[i]=(sB[i-1].GenQ-sB[i-1].LoadQ)-tempP*u[i][0]+tempQ*u[i][1];
P_UU[i]=0.0;
}
if(sB[i-1].Type==1)
{double tempP=0.0;
double tempQ=0.0;
for(j=1;j<nB+1;j++)
{
tempP+=YG[i][j]*u[j][1]-YB[i][j]*u[j][0];
tempQ+=YG[i][j]*u[j][0]+YB[i][j]*u[j][1];
P_P[i]=(sB[i-1].GenP-sB[i-1].LoadP)-tempP*u[i][1]-tempQ*u[i][0];
}
P_UU[i]=sB[i-1].Volt*sB[i-1].Volt-u[i][1]*u[i][1]-u[i][0]*u[i][0];
P_Q[i]=0.0;
}
}
/* 偏移量阵 */
double P_PQ[6];
int a=0;
for(i=1;i<3;i++)
{P_PQ[a]=P_P[i];
a=a+2;
}
a=1;
for(i=1;i<3;i++)
{P_PQ[a]=P_Q[i];
a=a+2;
}
P_PQ[4]=P_P[3];
P_PQ[5]=P_UU[3];
printf("迭代第%d次的偏移量为:\n",v);
for(i=0;i<6;i++)
{printf("%f",P_PQ[i]);
printf("\n");}
printf("\n");
printf("\n");
/* 雅可比矩阵 */
double H[6][6],N[6][6],J[6][6],L[6][6],R[6][6],S[6][6],aa[6],bb[6];
for(i=1;i<5;i++)
{ if(fabs(sB[i-1].Type-2.0)<0.000001)
continue;
else
{for(j=1;j<5;j++)
if(i!=j)
{H[i][j]=-YB[i][j]*u[i][1]+YG[i][j]*u[i][0];
N[i][j]=YG[i][j]*u[i][1]+YB[i][j]*u[i][0];
J[i][j]=-N[i][j];
L[i][j]=H[i][j];
R[i][j]=0;
S[i][j]=0;}
else
{
aa[i]=bb[i]=0.0;
for(n=1;n<5;n++)
{
aa[i]+=YG[i][n]*u[n][1]-YB[i][n]*u[n][0];
bb[i]+=YG[i][n]*u[n][0]+YB[i][n]*u[n][1];
}
H[i][i]=-YB[i][i]*u[i][1]+YG[i][i]*u[i][0]+bb[i];
N[i][i]=YG[i][i]*u[i][1]+YB[i][i]*u[i][0]+aa[i];
J[i][i]=-YG[i][i]*u[i][1]-YB[i][i]*u[i][0]+aa[i];
L[i][i]=YG[i][i]*u[i][0]-YB[i][i]*u[i][1]-bb[i];
R[i][i]=2*u[i][0];
S[i][i]=2*u[i][1];
}
}}
double ss[50][50];
for(i=0;i<6;i++)
for(j=0;j<6;j++)
ss[i][j]=0.0;
for(i=1;i<3;i++)
for(j=1;j<4;j++)
{ss[2*i-2][2*j-2]=H[i][j];
ss[2*i-2][2*j-1]=N[i][j];
ss[2*i-1][2*j-2]=J[i][j];
ss[2*i-1][2*j-1]=L[i][j];
}
i=3;
for(j=1;j<4;j++)
for(j=1;j<4;j++)
{ss[2*i-2][2*j-2]=H[i][j];
ss[2*i-2][2*j-1]=N[i][j];
ss[2*i-1][2*j-2]=R[i][j];
ss[2*i-1][2*j-1]=S[i][j];
}
printf("迭代第%d次的雅可比矩阵为:\n",v);
for(i=0;i<6;i++)
{for(j=0;j<6;j++)
printf("%10f",ss[i][j]);
printf("\n");}
printf("\n");
printf("\n");
Gauss(ss,P_PQ,6);
for(i=1;i<nB;i++)
{u[i][0]=u[i][0]+mn[2*(i-1)];
u[i][1]=u[i][1]+mn[2*i-1];}
double max;
max=fabs(P_PQ[0]);
for(i=0;i<=5;i++)
if (max<fabs(P_PQ[i]))
max=fabs(P_PQ[i]);
if(fabs(max)<0.0001)
{printf("满足精度要求,迭代终止,迭代次数为%d\n",v);
printf("\n");
printf("\n");
break;}
}/* 叠代循环的括号 */
printf("最终求得的节点电压值为e+jf:\n");
for(i=1;i<nB+1;i++)
printf("%lf,%lf\n",u[i][1],u[i][0]);
printf("\n");
printf("\n");
double uu[5],Phase[5];
for(i=1;i<nB+1;i++)
{uu[i]=sqrt(u[i][1]*u[i][1]+u[i][0]*u[i][0]);
Phase[i]=atan(u[i][0]/u[i][1]);}
for(i=1;i<nB+1;i++)
printf("%lf,%lf\n",uu[i],Phase[i]);
*计算线路功率和平衡节点 PV节点功率*/
double P[5],Q[5];
double tempP=0.0;
double tempQ=0.0;
for(i=1;i<nB+1;i++)
{for(j=1;j<nB+1;j++)
{
tempP+=YG[i][j]*u[j][1]-YB[i][j]*u[j][0];
tempQ+=YG[i][j]*u[j][0]+YB[i][j]*u[j][1];
}
P[i]=tempP*u[i][1]+tempQ*u[i][0];
Q[i]=tempP*u[i][0]-tempQ*u[i][1];
tempP=0.0;
tempQ=0.0;
}
for(i=1;i<nB+1;i++)
printf("节点%d注入功率为%lf,%lf\n",i,P[i],Q[i]);
/* 支路功率 */
double V[4][2];
for(i=1;i<5;i++)
for(j=0;j<3;j++)
V[i][j]=u[i][j];
double sP[5][5],sQ[5][5];
double dsq,dsp,dp,sumgen;
for(i=1;i<NBUS+1;i++)
{
for(j=1;j<NBUS+1;j++)
{
sP[i][j]=0.0;
sQ[i][j]=0.0;
}
}
for(l=0; l<nL; 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)
{/*Normal lines or transformers*/
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
{/*abnormal transformer ratio*/
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));
}
}
/* 输电效率 */
dsp=P[4];
dsq=Q[4];
sumgen=P[4];
for(i=0;i<NBUS;i++)
{
dsp+=sB[i].GenP-sB[i].LoadP;
dsq+=sB[i].GenQ-sB[i].LoadQ;
sumgen+=sB[i].GenP;
}
dp=dsp/sumgen*100;
/* 输出功率情况 */
if((fp=fopen("功率情况.txt","w"))==NULL)
{
printf("Can not open the file named '功率情况.txt' \n");
exit(0);
}
fprintf(fp,"---功率情况---\n");
fprintf(fp,"平衡节点功率S=%10.5f+ j%10.5f\n",P[4],Q[4]);
for(i=1;i<NBUS+1;i++)
for(j=1;j<NBUS+1;j++)
if(fabs(sP[i][j]-0.0)>0.000001)
fprintf(fp,"S(%3d,%-3d)=(%10.5f,j%10.5f)\n",i,j,sP[i][j],sQ[i][j]);
fprintf(fp,"网损为%10.5f+j%10.3f,输电效率为%10.3f\n",dsp,dsq,100-dp);
if(fp!=NULL) fclose(fp);
return 0;
}
结果:
1.导纳阵
Y( 1,1 )=( 1.01534, -8.19201)
Y( 1,2 )=( -0.56148, 2.30208)
Y( 1,3 )=( 0.00000, 3.66667)
Y( 1,4 )=( -0.45386, 1.89107)
Y( 2,1 )=( -0.56148, 2.30208)
Y( 2,2 )=( 1.04225, -4.67651)
Y( 2,4 )=( -0.48077, 2.40385)
Y( 3,1 )=( 0.00000, 3.66667)
Y( 3,3 )=( 0.00000, -3.33333)
Y( 4,1 )=( -0.45386, 1.89107)
Y( 4,2 )=( -0.48077, 2.40385)
Y( 4,4 )=( 0.93463, -4.26159)
2.设定电压初值
3.计算功率和电压偏移;
同理可算出
,
,
4.根据求的第一次迭代时雅可比矩阵各元素的公式计算雅可比矩阵各个元素的具体值:
5.求高斯计算后的修正量:
6.计算各节点电压的一次近似值:
返回第三步重新迭代,并校验收敛与否,令。经过 四 次迭代后,收敛条件满足,停止迭代,求出的电压终值:
7.计算出平衡节点4的注入功率。
8.各节点间功率
9.网损为:
10.网损效率为:2.05896%