📄 floor.c
字号:
#include <stdio.h>
#include <math.h>
#define M 50 /*矩阵阶数*/
int i,j,k,l,z; /* 循环变量 */
int n, /* 节点数 */
m, /* 支路数 */
dd, /*对地支路数*/
pq, /* PQ节点数 */
pv, /* PV节点数 */
byqn, /* 变压器数*/
number,
ii,
ll=1,
li=1;
float eps, /* 精度 */
max1,max2;
static float G[M][M],B[M][M],
D1[M],D2[M],yue[M][M],xing[M][M];
struct jiedian /* 节点结构体 */
{ int num,s; /* num为节点号，s为节点类型*/
float p,q,e,f,v;
} jiedian[M];
struct zhilu /* 支路结构体 */
{ int num;
int p1,p2; /*支路的两个节点*/
float r,x; /*支路的电阻与电抗*/
} zhilu[M];
struct byq
{ int num,
p1,
p2;
float r,x,s,k;
}byq[M];
FILE *fp1,*fp2;
void data1() /* 读取数据 */
{ int h,numb;
fp1=fopen("d:\\input.txt","r");
if(fp1==NULL)
{ printf(" can not open file !\n");
printf(" 在D盘创建一个文件，并命名为 input.txt. 文件输入内容的格式如下：.\n");
printf("\n ==============================================================\n");
printf(" 第一行：节点数， 支路数， 对地支路数，变压器数，PQ节点数，PV节点数，精度，\n");
printf("换行后输入：第几个节点，节点类型，有功功率，无功功率，电压值，电压角度\n");
printf("(若是PQ节点则类型是1，电压是自己预置的值；\n");
printf("PV节点类型是2，不用输入无功功率，电压是题目给定的值；\n");
printf("平衡节点类型是3，只需输入电压及其角度。)\n");
printf("注意每次输入一个节点应该换行。\n");
printf("\n继续换行输入：第几条支路，支路类型，支路连接的第一个节点号，支路连接的第二个节点号，支路电阻，支路电抗\n ");
printf("注意：支路没有变压器类型是1，有是2。若有变压器则还需要输入其变比和容量\n");
printf("\n ==============================================================\n");
exit(0);
}
fscanf(fp1,"%d,%d,%d,%d,%d,%d,%f\n",&n,&m,&dd,&byqn,&pq,&pv,&eps); /*输入节点数,支路数,对地支路数,PQ节点数PV节点数和精度*/
j=1;k=pq+1;
for(i=1;i<=n;i++) /*输入节点类型的输入功率和节电电压初值*/
{
fscanf(fp1,"%d,%d",&numb,&h);
if(h==1) /*类型h=1是PQ节点*/
{ fscanf(fp1,",%f,%f,%f,%f\n",&jiedian[j].p,&jiedian[j].q,
&jiedian[j].e,&jiedian[j].f);
jiedian[j].num=numb;
jiedian[j].s=h;
j++;
}
if(h==2) /*类型h=2是pv节点*/
{ fscanf(fp1,",%f,%f,%f\n",&jiedian[k].p,&jiedian[k].v,
&jiedian[k].f);
jiedian[k].num=numb;
jiedian[k].e=jiedian[k].v;
jiedian[k].s=h;
k++;
}
if(h==3) /*类型h=3是平衡节点*/
{ fscanf(fp1,",%f,%f\n",&jiedian[n].e,&jiedian[n].f);
jiedian[n].num=numb;
jiedian[n].s=h;
}
}
for(i=1;i<=m+byqn;i++) /*输入支路阻抗*/
{fscanf(fp1,"%d,%d",&number,&ii);
if(ii==1)
{fscanf(fp1,",%d,%d,%f,%f\n", &zhilu[ll].p1,
&zhilu[ll].p2, &zhilu[ll].r, &zhilu[ll].x);
zhilu[ll].num=ll;ll++;}
if(ii==2) /*变压器编号num，所连节点p1，p2，阻抗r，x，变比k，容量s*/
{ fscanf(fp1,",%d,%d,%f,%f,%f,%f\n",&byq[li].p1,
&byq[li].p2,&byq[li].r,&byq[li].x,&byq[li].k,&byq[li].s);
byq[li].num=li;li++;}
}
fclose(fp1);
if((fp2=fopen("d:\\out.txt","w"))==NULL)
{ printf(" can not open file!\n");
exit(0);
}
fprintf(fp2,"\n ************** 原始数据 **************\n");
fprintf(fp2," ===================================================================\n");
fprintf(fp2," 节点数:%2d 支路数:%2d 对地支路数:%2d 变压器数：%2d PQ节点数:%2d PV节点数:%2d 精度:%.5f\n",
n,m,dd,byqn,pq,pv,eps);
fprintf(fp2," -------------------------------------------------------------------\n");
for(i=1;i<=pq;i++)
fprintf(fp2," PQ节点 节点%2d P[%d]=%f Q[%d]=%f\n",
jiedian[i].num,jiedian[i].num,jiedian[i].p,jiedian[i].num,jiedian[i].q);
for(i=pq+1;i<=pq+pv;i++)
fprintf(fp2," PV节点 节点%2d P[%d]=%f V[%d]=%f\n",
jiedian[i].num,jiedian[i].num,jiedian[i].p,jiedian[i].num,jiedian[i].v);
fprintf(fp2," 平衡节点 节点%2d e[%d]=%f f[%d]=%f\n",
jiedian[n].num,jiedian[n].num,jiedian[n].e,jiedian[n].num,jiedian[n].f);
fprintf(fp2," -------------------------------------------------------------------\n");
for(i=1;i<=m;i++)
fprintf(fp2," 支路%2d 相关节点:%2d,%2d R=%f X=%f\n",
i,zhilu[i].p1,zhilu[i].p2,zhilu[i].r,zhilu[i].x);
for(i=1;i<=byqn;i++)
fprintf(fp2," 变压器编号%2d 相关节点:%2d,%2d R=%f X=%f 变比K=%f\n",byq[i].num,
byq[i].p1,byq[i].p2,byq[i].r,byq[i].x,byq[i].k);
fprintf(fp2," ===================================================================\n");
}
void form_y() /* 形成节点导纳矩阵 */
{int bn;
float S,g,b,g1,b1,g2,b2;
for(i=1;i<=m;i++)
for(i=1;i<=m;i++)
/*对地导纳*/
if((zhilu[i].p1==0)||(zhilu[i].p2==0))
{S=zhilu[i].r*zhilu[i].r+zhilu[i].x*zhilu[i].x;
if(S==0) {G[zhilu[i].p1][zhilu[i].p2]=B[zhilu[i].p1][zhilu[i].p2]=G[zhilu[i].p2][zhilu[i].p1]=B[zhilu[i].p2][zhilu[i].p1]=0;continue;}
if(zhilu[i].p1==0)
{G[zhilu[i].p1][zhilu[i].p2]+=zhilu[i].r/S;
B[zhilu[i].p1][zhilu[i].p2]+=-zhilu[i].x/S;
G[zhilu[i].p2][zhilu[i].p1]=G[zhilu[i].p1][zhilu[i].p2];
B[zhilu[i].p2][zhilu[i].p1]=B[zhilu[i].p1][zhilu[i].p2];
}
if(zhilu[i].p2==0)
{G[zhilu[i].p2][zhilu[i].p1]+=zhilu[i].r/S;
B[zhilu[i].p2][zhilu[i].p1]+=-zhilu[i].x/S;
G[zhilu[i].p1][zhilu[i].p2]=G[zhilu[i].p2][zhilu[i].p1];
B[zhilu[i].p1][zhilu[i].p2]=B[zhilu[i].p2][zhilu[i].p1];
}
}
for(i=1;i<=m;i++) /*节点导纳矩阵的主对角线上的导纳*/
for(j=1;j<=n;j++)
if((zhilu[i].p1==j)||(zhilu[i].p2==j))
{ S=zhilu[i].r*zhilu[i].r+zhilu[i].x*zhilu[i].x;
if(S==0) continue;
G[j][j]+=zhilu[i].r/S;
B[j][j]+=-zhilu[i].x/S;
}
for(j=1;j<=n;j++)
for(bn=1;bn<=byqn;bn++)
{if(byq[bn].p1==j)
{ S=byq[bn].r*byq[bn].r+byq[bn].x*byq[bn].x;
g=byq[bn].r/S;
b=-byq[bn].x/S;
g2=g/byq[bn].k;
b2=b/byq[bn].k;
g1=g*((1-byq[bn].k)/(byq[bn].k*byq[bn].k));
b1=b*((1-byq[bn].k)/(byq[bn].k*byq[bn].k));
G[j][j]=G[j][j]+g2+g1;
B[j][j]=B[j][j]+b2+b1;}
if(byq[bn].p2==j)
{S=byq[bn].r*byq[bn].r+byq[bn].x*byq[bn].x;
g=byq[bn].r/S;
b=-byq[bn].x/S;
g2=g/byq[bn].k;
b2=b/byq[bn].k;
g1=g*(byq[bn].k-1)/byq[bn].k;
b1=b*(byq[bn].k-1)/byq[bn].k;
G[j][j]=G[j][j]+g2+g1;
B[j][j]=B[j][j]+b2+b1;}
}
for(k=1;k<=m;k++) /*节点导纳矩阵非主对角线上的导纳*/
{ i=zhilu[k].p1;
j=zhilu[k].p2;
S=zhilu[k].r*zhilu[k].r+zhilu[k].x*zhilu[k].x;
if(S==0) continue;
G[i][j]+=-zhilu[k].r/S;
B[i][j]+=zhilu[k].x/S;
G[j][i]=G[i][j];
B[j][i]=B[i][j];
}
for(bn=1;bn<=byqn;bn++)
{i=byq[bn].p1;j=byq[bn].p2;
S=byq[bn].r*byq[bn].r+byq[bn].x*byq[bn].x;
G[i][j]+=-(byq[bn].r/S)/byq[bn].k;
B[i][j]+=(byq[bn].x/S)/byq[bn].k;
G[j][i]=G[i][j];
B[j][i]=B[i][j];}
/*输出节点导纳矩阵*/
fprintf(fp2,"\n\n ********* 计算结果 *********\n");
fprintf(fp2,"\n 节点导纳矩阵为:");
for(i=1;i<=n;i++)
{ fprintf(fp2,"\n ");
for(j=1;j<=n;j++)
fprintf(fp2,"%8.5f+j%8.5f ",G[i][j],B[i][j]);
}
fprintf(fp2,"\n ==============================================================\n");
fprintf(fp2,"\n B'矩阵为:");
for(i=1;i<=n-1;i++)
{ fprintf(fp2,"\n ");
for(j=1;j<=n-1;j++)
fprintf(fp2,"%8.5f ",B[i][j]);
}
fprintf(fp2,"\n ==============================================================\n");
fprintf(fp2,"\n B''矩阵为:");
for(i=1;i<=pq;i++)
{ fprintf(fp2,"\n ");
for(j=1;j<=pq;j++)
fprintf(fp2,"%8.5f ",B[i][j]);
}
fprintf(fp2,"\n ==============================================================\n");
}
void yueyue()
{
int a ,b;
for(a=1;a<=M;a++)
{ for (b=1;b<=M;b++)
{yue[a][b]=0;xing[a][b]=0;}
}
for(i=1;i<=pq;i++)
{for (j=1;j<=pq;j++)
xing[i][j]=B[i][j];
}
for(i=1;i<=pq+pv;i++)
{for(j=1;j<=pq+pv;j++)
yue[i][j]=B[i][j];
}
}
void xingxing() /* 形成矩阵方程 */
{ float ei,ej,fi,fj;
int i1,j1;
float dP[M];
for(i=1;i<=pq+pv;i++) /*初始化矩阵*/
{ dP[i]=0;
D1[i]=0;
}
for(i=1;i<=pq;i++) /*PQ节点函数残量*/
{ ei=jiedian[i].e;
fi=jiedian[i].f;
for(j=1;j<=n;j++)
{ i1=jiedian[i].num;
j1=jiedian[j].num;
ej=jiedian[j].e;
fj=jiedian[j].f;
dP[i]+=ej*(G[i1][j1]*cos(fi-fj)+B[i1][j1]*sin(fi-fj));
}
dP[i]=jiedian[i].p-dP[i]*ei;
}
for(i=pq+1;i<=pq+pv;i++) /*PV节点函数残量*/
{ ei=jiedian[i].e;
fi=jiedian[i].f;
for(j=1;j<=n;j++)
{ i1=jiedian[i].num;
j1=jiedian[j].num;
ej=jiedian[j].e;
fj=jiedian[j].f;
dP[i]+=ej*(G[i1][j1]*cos(fi-fj)+B[i1][j1]*sin(fi-fj));
}
dP[i]=jiedian[i].p-dP[i]*ei;
}
for(i=1;i<=pq+pv;i++) /*形成函数残量矩阵D[M]*/
D1[i]=dP[i]/jiedian[i].e; /*形成^P */
max1=D1[1];
for(i=2;i<=pq+pv;i++)
if(fabs(D1[i])>fabs(max1))
max1=D1[i];
max1=fabs(max1);
fprintf(fp2,"\n%8.5f\n",max1);
}
void xingxing1() /* 形成矩阵方程 */
{ float ei,ej,fi,fj;
int i1,j1;
float dQ[M];
for(i=1;i<=pq+pv;i++) /*初始化矩阵*/
{
dQ[i]=0;
D2[i]=0;
}
for(i=1;i<=pq;i++) /*PQ节点函数残量*/
{ ei=jiedian[i].e;
fi=jiedian[i].f;
for(j=1;j<=n;j++)
{ i1=jiedian[i].num;
j1=jiedian[j].num;
ej=jiedian[j].e;
fj=jiedian[j].f;
dQ[i]+=ej*(G[i1][j1]*sin(fi-fj)-B[i1][j1]*cos(fi-fj));
}
dQ[i]=jiedian[i].q-dQ[i]*ei;
}
for(i=1;i<=pq;i++)
D2[i]=dQ[i]/jiedian[i].e; /*形成^Q */
max2=D2[1];
for(i=2;i<=pq;i++)
if(fabs(D2[i])>fabs(max2))
max2=D2[i];
max2=fabs(max2);
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -