📄 yangshuo.c
字号:
else /*求i=j时的H、N、J、L*/
{ a=0;b=0;
for(k=1;k<=n;k++)
if(k!=i)
{ k1=jiedian[k].num;
a=a+G[i1][k1]*jiedian[k].e-B[i1][k1]*jiedian[k].f;
b=b+G[i1][k1]*jiedian[k].f+B[i1][k1]*jiedian[k].e;
}
ykb[2*i-1][2*j-1]=2*G[i1][i1]*jiedian[i].f+b; /*H*/
ykb[2*i][2*j]=-2*B[i1][i1]*jiedian[i].e-b; /*L*/
ykb[2*i-1][2*i]=2*G[i1][i1]*jiedian[i].e+a; /*N*/
ykb[2*i][2*i-1]=-2*B[i1][i1]*jiedian[i].f+a; /*J*/
}
}
for(i=pq+1;i<=pq+pv;i++) /*形成pv节点子阵*/
for(j=1;j<n;j++)
{ i1=jiedian[i].num;
j1=jiedian[j].num;
ei=jiedian[i].e;
ej=jiedian[j].e;
fi=jiedian[i].f;
fj=jiedian[j].f;
if(i!=j) /*求i!=j时的H、N*/
{ ykb[2*i-1][2*j-1]=-B[i1][j1]*ei+G[i1][j1]*fi; /*H*/
ykb[2*i-1][2*j]=G[i1][j1]*ei+B[i1][j1]*fi; /*N*/
}/*省略R[i][j]=0及S[i][j]=0*/
else /*求i=j时的H、N、R、S*/
{ a=0;b=0;
for(k=1;k<=n;k++)
if(k!=i)
{ k1=jiedian[k].num;
a+=G[i1][k1]*jiedian[k].e-B[i1][k1]*jiedian[k].f;
b+=G[i1][k1]*jiedian[k].f+B[i1][k1]*jiedian[k].e;
}
ykb[2*i-1][2*j-1]=2*G[i1][i1]*jiedian[i].f+b; /*H*/
ykb[2*i-1][2*i]=2*G[i1][i1]*jiedian[i].e+a; /*N*/
ykb[2*i][2*j-1]=2*fi; /*R*/
ykb[2*i][2*j]=2*ei; /*S*/
}
}
for(i=1;i<=2*(pq+pv);i++) /*把函数残量矩阵编入修正方程*/
ykb[i][2*(pq+pv)+1]=D[i];
}
void print_ykb(int z)
{int i,j;
fprintf(fp,"\n===============================================================\n");
fprintf(fp,"--------------------输出雅可比矩阵和迭代次数---------------------\n");
fprintf(fp,"\n 迭代的次数为 %2d\n", z);
printf("\n输出雅可比矩阵和迭代次数\n");
printf("迭代的次数为 %2d\n", z);
for(i=1;i<=2*(pq+pv);i++)
{printf("\n");fprintf(fp,"\n");
for(j=1;j<=2*(pq+pv);j++)
{ printf("%7.5f ",ykb[i][j]);
fprintf(fp,"%8.5f ",ykb[i][j]);}
}
fprintf(fp," \n=============================================================\n");
}
void solve() /* 求解修正方程组 (用列主元消取法)*/
{ float d,t;
k=1;
do
{
d=ykb[k][k];
l=k;
i=k+1;
do
{
if(fabs(ykb[i][k])>fabs(d)) /*列选主元*/
{ d=ykb[i][k];
l=i;
}
i++;
} while(i<=2*(pq+pv));
if(l!=k)
{ for(j=k;j<=2*(pq+pv)+1;j++)
{ t=ykb[l][j];
ykb[l][j]=ykb[k][j];
ykb[k][j]=t;
}
}
for(j=k+1;j<=2*(pq+pv)+1;j++) /*正消*/
ykb[k][j]/=ykb[k][k];
for(i=k+1;i<=2*(pq+pv);i++)
for(j=k+1;j<=2*(pq+pv)+1;j++)
ykb[i][j]-=ykb[i][k]*ykb[k][j];
k++;
} while(k<=2*(pq+pv));
if(k!=1)
{
for(i=2*(pq+pv);i>=1;i--) /*回代*/
{ t=0;
for(j=i+1;j<=2*(pq+pv);j++)
t+=ykb[i][j]*D[j];
D[i]=ykb[i][2*(pq+pv)+1]-t;
}
}
for(i=1;i<=(pq+pv);i++)
{ jiedian[i].e+=D[2*i];
jiedian[i].f+=D[2*i-1];
}
max=fabs(D[1]);
for(i=1;i<=2*(pq+pv);i++)
if(fabs(D[i])>max)
max=fabs(D[i]);
}
/*输出迭代过程中的电压值*/
void voltage()
{int i,j;
fprintf(fp,"\n=====================================================================\n");
fprintf(fp,"\n--------------------------输出 df,de--------------------------------\n");
for(i=1;i<=(pq+pv);i++)
{fprintf(fp,"\n");
fprintf(fp," 节点为%2d de+jdf=%7.5f+j%7.5f \n",i,D[2*i],D[2*i-1]);}
fprintf(fp,"\n---------------------------------------------------------------------\n");
fprintf(fp,"\n----------------输出迭代过程中的电压值-------------------------------\n");
for(i=1;i<=(pq+pv);i++)
{fprintf(fp,"\n 节点为%2d e+jf=%7.5f+j%7.5f \n",i,jiedian[i].e,jiedian[i].f);}
fprintf(fp,"\n=====================================================================\n");
}
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);
}
void data2() /* 潮流计算结果 */
{ float sp=0,sq=0,mo1,mo2,mo3,mo4,mo5,mo6,ff;
int i1,j1;
static float P[M][M],Q[M][M];
fprintf(fp,"\n ===================================================================\n");
fprintf(fp,"\n\n 各节点电压为:");
for(i=1;i<=n;i++)
fprintf(fp,"\n U[%d]=%7.5f + j %7.5f",
jiedian[i].num,jiedian[i].e,jiedian[i].f);
fprintf(fp,"\n ===================================================================\n");
fprintf(fp,"\n\n 平衡节点功率为:");
for(i=1;i<=n;i++)
{ i1=jiedian[n].num;
j1=jiedian[i].num;
sp+=mul_Re(G[i1][j1],-B[i1][j1],jiedian[i].e,-jiedian[i].f);
sq+=mul_Im(G[i1][j1],-B[i1][j1],jiedian[i].e,-jiedian[i].f);
}
jiedian[n].p=mul_Re(jiedian[n].e,jiedian[n].f,sp,sq);
jiedian[n].q=mul_Im(jiedian[n].e,jiedian[n].f,sp,sq);
fprintf(fp,"\n S[%d]=%7.5f+j %7.5f",jiedian[n].num,jiedian[n].p,jiedian[n].q);
fprintf(fp,"\n ===================================================================\n");
fprintf(fp,"\n\n 线路功率如下:\n");
for(k=1;k<=m;k++)
{ i1=zhilu[k].p1;
j1=zhilu[k].p2;
for(l=1;l<=n;l++)
{ if(jiedian[l].num==i1)
i=l;
if(jiedian[l].num==j1)
j=l;
}
sp=mul_Re(jiedian[i].e-jiedian[j].e,-jiedian[i].f+jiedian[j].f,
-G[i1][j1],B[i1][j1]);
sq=mul_Im(jiedian[i].e-jiedian[j].e,-jiedian[i].f+jiedian[j].f,
-G[i1][j1],B[i1][j1]);
P[i1][j1]=mul_Re(jiedian[i].e,jiedian[i].f,sp,sq);
Q[i1][j1]=mul_Im(jiedian[i].e,jiedian[i].f,sp,sq);
sp=mul_Re(jiedian[j].e-jiedian[i].e,-jiedian[j].f+jiedian[i].f,
-G[j1][i1],B[j1][i1]);
sq=mul_Im(jiedian[j].e-jiedian[i].e,-jiedian[j].f+jiedian[i].f,
-G[j1][i1],B[j1][i1]);
P[j1][i1]=mul_Re(jiedian[j].e,jiedian[j].f,sp,sq);
Q[j1][i1]=mul_Im(jiedian[j].e,jiedian[j].f,sp,sq);
fprintf(fp,"\n 线路%d-%d的功率:",i1,j1);
fprintf(fp," %7.5f + j %7.5f",P[i1][j1],Q[i1][j1]);
fprintf(fp,"\n 线路%d-%d的功率:",j1,i1);
fprintf(fp," %7.5f + j %7.5f",P[j1][i1],Q[j1][i1]);}
for(i=1;i<=(m-dd);i++) /*求平行支路的功率分布（功率按导纳分布）*/
{z=1;
mo3=zhilu[i].r*zhilu[i].r+zhilu[i].x*zhilu[i].x;
mo3=1/sqrt(mo3);
for(j=i+1;j<=(m-dd);j++)
{if((zhilu[i].p1==zhilu[j].p1)&&(zhilu[i].p2==zhilu[j].p2))
{i1=zhilu[i].p1;
j1=zhilu[i].p2;
mo4=zhilu[j].r*zhilu[j].r+zhilu[j].x*zhilu[j].x;
mo4=1/sqrt(mo4);
mo3=mo3+mo4;
z++;}}
if(z>1)
{ fprintf(fp,"\n====================================================================\n");
fprintf(fp,"节点%d-%d之间有%d条平行支路\n",zhilu[i].p1,zhilu[i].p2,z);
fprintf(fp,"功率分布为：\n");
for(k=i;k<z+i;k++)
{ mo5=zhilu[k].r*zhilu[k].r+zhilu[k].x*zhilu[k].x;
mo6=1/sqrt(mo5);
ff=mo6/mo3;
fprintf(fp,"\n 导纳Y=%8.6f的支路",mo6);
fprintf(fp,"\n 支路%d-%d的功率:",zhilu[i].p1,zhilu[i].p2);
fprintf(fp," %8.6f + j %8.6f",ff*P[i1][j1],ff*Q[i1][j1]);
fprintf(fp,"\n 支路%d-%d的功率:",zhilu[i].p2,zhilu[i].p1);
fprintf(fp," %8.6f + j %8.6f\n",ff*P[j1][i1],ff*Q[j1][i1]);
}
i=i+z-1; }
}
fprintf(fp,"\n\n 线路上损耗的功率为:");
for(k=1;k<=m;k++)
{ i=zhilu[k].p1;
j=zhilu[k].p2;
fprintf(fp,"\n 线路%d-%d上的功率损耗:",i,j);
fprintf(fp," %7.5f + j %7.5f",fabs(P[i][j]+P[j][i]),
fabs(Q[i][j]+Q[j][i]));
}
fprintf(fp,"\n\n 网络总损耗:");
sp=0;sq=0;
for(i=1;i<=n;i++)
for(j=1;j<=n;j++)
{ sp+=P[i][j];
sq+=Q[i][j];
}
fprintf(fp," %7.5f + j %7.5f\n",sp,sq);
fprintf(fp,"\n ************************结束************************\n");
fclose(fp);
}
void main() /* 主函数 */
{ z=1;
data1();
form_y();
do
{ form_ykb();
print_ykb(z);
solve();
voltage();
z++;
}
while((z<N)&&(max>=eps));
data2();
printf(" \n!!计算结束，结果存放于输出文件中!!\n");
}
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -