// 潮流计算程序.cpp : Defines the entry point for the console application.
//
#include "stdafx.h"
#include
#include
#include
#define N 9 /*总结点数*/
#define M 6 /*PQ结点数*/
#define K 9 /*线路数*/
#define eps 1e-4
#define double float
void guass(int n,int m,double c[],double b[][N],double x[]) /*高斯函数*/
{
double a[N][N],y[N];
int i,j,k;
for(i=0;i {
y[i]=c[i];
for(j=0;j }
for(k=0;k for(i=k+1;i {
for(j=k+1;j y[i]=y[i]-y[k]*a[i][k]/a[k][k];
}
for(i=m-1;i>=0;i--)
{
for(j=i+1;j x[i]=y[i]/a[i][i];
}
}
struct line
{int Lindex;
int Headnode;
int Endnode;
double R;
double X;
double b;
}line[K];
struct line *t;
void main()
{
double r[N]={0.0};
double u[N]={1.04,1.025,1.025,1.0,1.0,1.0,1.0,1.0,1.0};
double p[N]={1,1.63,0.85,0.0,0.0,0.0,-1.25,-0.9,-1.0};
double q[N]={1,1,1,0.0,0.0,0.0,-0.5,-0.3,-0.35};
double g[N][N]={0.0},b[N][N]={0.0};
double h[N][N]={0.0};
double B[N][N]={0.0};
double temp;
double H[N][6];
double lr,lx,lb1,lg,lb;
int i,j;
int ku=0,kr=0,kp=1,kq=1;
void val(double u[N],double g[N][N],double b[N][N],double r[N],int ku, int kr,double h[N][N]); /*函数申明*/
FILE *fp;
if ((fp=fopen("data.txt","r"))==NULL)
{
printf("cannot open this file.\n");
exit(0);
}
//fp=fopen("data.txt","r");
for(i=0;i {
for(j=0;j { fscanf(fp,"%f,",&temp);
H[i][j]=temp;
}
}
fclose(fp);
for(i=0;i {
line[i].Lindex=(int)H[i][0];line[i].Headnode=(int)H[i][1];line[i].Endnode=(int)H[i][2];
line[i].R=H[i][3];line[i].X=H[i][4];line[i].b=H[i][5];}
for(t=line;t {i=t->Headnode-1;
j=t->Endnode-1;
lr=t->R;
lx=t->X;
lb1=t->b;
lg=lr/(lr*lr+lx*lx);
lb=-lx/(lr*lr+lx*lx);
g[i][i]+=lg;
g[j][j]+=lg;
b[i][i]+=lb+lb1;
b[j][j]+=lb+lb1;
h[i][j]=h[j][i]=-lb1;
g[i][j]-=lg;
g[j][i]-=lg;
b[i][j]-=lb;
b[j][i]-=lb;
}
getchar();
printf("\n=====节点导纳矩阵=====\n");
for(i=0;i {
for(j=0;j {
B[i][j]=b[i][j];
printf("%8f,",B[i][j]);
}
printf("\n");
}
printf("\n");
getchar();
printf("\n=====给定初值=====\n");
for(i=0;i printf("u[%d]=%8f p[%d]=%8f q[%d]=%8f\n",i+1,u[i],i+1,p[i],i+1,q[i]);
printf("\n=====迭代求解=====\n");
while(kp==1)
{
double ip,iq,max;
double dp[N],dq[N],dr[N];
double dpu[N],dqu[N];
double y1[N-1],y2[M];
double x1[N-1],x2[M];
for(i=1;i {
ip=0;
for(j=0;j ip=ip+u[j]*(g[i][j]*cos(r[i]-r[j])+b[i][j]*sin(r[i]-r[j]));
dp[i]=p[i]-u[i]*ip;
dpu[i]=dp[i]/u[i];
printf("dp[%d]=%8f\n",i+1,dp[i]);
printf("\n");
}
getchar();
max=fabs(dpu[1]);
for(i=1;i {
if(fabs(dpu[i])>max)
max=fabs(dpu[i]);
}
if (max>=eps)
{
for(i=0;i y1[i]=-dpu[i+1]; /*起值不同,为了对应,故加一*/
guass(1,N-1,y1,B,x1);
for (i=1;i {
dr[i]=x1[i-1]/u[i];
r[i]=r[i]+dr[i];
}
printf("\n===== 第 %d 次迭代后电压 相角初值 =====\n",kr+1);
for(i=1;i printf("r[%d]=%8f\n",i+1,r[i]);
getchar();
printf("\n\n");
kr=kr+1;kq=1;
top:
for(i=N-M;i {
iq=0;
for(j=0;j iq=iq+u[j]*(g[i][j]*sin(r[i]-r[j])-b[i][j]*cos(r[i]-r[j]));
dq[i]=q[i]-u[i]*iq;
dqu[i]=dq[i]/u[i];
printf("dq[%d]=%8f\n",i+1,dq[i]);
printf("\n");
}
max=fabs(dqu[N-M]);
for (i=N-M;i {
if(fabs(dqu[i])>max)
max=fabs(dqu[i]);
}
if(max>=eps)
{
for(i=0;i y2[i]=-dqu[i+N-M]; /*同上,对应加N-M*/
guass(N-M,M,y2,B,x2);
for(i=N-M;i u[i]=u[i]+x2[i-(N-M)];
printf("\n=====第 %d 次迭代电压初值=====\n",ku+1);
for(i=N-M;i printf("u[%d]=%8f\n",i+1,u[i]);
printf("\n\n");
ku=ku+1;kp=1;
}
else
{
kq=0;
if(kp==0)
val(u,g,b,r,ku,kr,h);
}
}
else
{
kp=0;
if(kq==0)
val(u,g,b,r,ku,kr,h);
else
goto top;
}
}
}
void val(double u[N],double g[N][N],double b[N][N],double r[N],int ku, int kr,double h[N][N])
{
double ps=0,pv1=0,pv2=0;
double qs=0,qv1=0,qv2=0;
double p[N][N]={0};
double q[N][N]={0};
double s[N][N]={0};
double dp[N][N]={0};
double dq[N][N]={0};
double ds[N][N]={0};
double dSp=0,dSq=0;
int i,j;
FILE *fp1;
printf("\n=====平衡节点功率 =====\n");
getchar();
for(i=0;i {
ps=ps+u[0]*u[i]*(g[0][i]*cos(r[0]-r[i])+b[0][i]*sin(r[0]-r[i]));
qs=qs+u[0]*u[i]*(g[0][i]*sin(r[0]-r[i])-b[0][i]*cos(r[0]-r[i]));
}
printf("sp=%8f+j(%8f)\n",ps,qs);
printf("\n=====PV 节点功率=====\n");
getchar();
for(i=0;i {
pv1=pv1+u[1]*u[i]*(g[1][i]*cos(r[1]-r[i])+b[1][i]*sin(r[1]-r[i]));
qv1=qv1+u[1]*u[i]*(g[1][i]*sin(r[1]-r[i])-b[1][i]*cos(r[1]-r[i]));
}
printf("sv1=%8f+j(%8f)\n",pv1,qv1);
for(i=0;i {
pv2=pv2+u[2]*u[i]*(g[2][i]*cos(r[2]-r[i])+b[2][i]*sin(r[2]-r[i]));
qv2=qv2+u[2]*u[i]*(g[2][i]*sin(r[2]-r[i])-b[2][i]*cos(r[2]-r[i]));
}
printf("sv2=%8f+j(%8f)\n",pv2,qv2);
for(i=0;i for(j=0;j {
p[i][j]=u[i]*u[i]*(-g[i][j])+u[i]*u[j]*(g[i][j]*cos(r[i]-r[j])+b[i][j]*sin(r[i]-r[j]));
q[i][j]=u[i]*u[i]*(-h[i][j]+b[i][j])+u[i]*u[j]*(g[i][j]*sin(r[i]-r[j])-b[i][j]*cos(r[i]-r[j]));
dp[i][j]=u[i]*u[i]*(-g[i][j])+u[i]*u[j]*(g[i][j]*cos(r[i]-r[j])+b[i][j]*sin(r[i]-r[j]))
+u[j]*u[j]*(-g[j][i])+u[j]*u[i]*(g[j][i]*cos(r[j]-r[i])+b[j][i]*sin(r[j]-r[i]));
dq[i][j]=u[i]*u[i]*(-h[i][j]+b[i][j])+u[i]*u[j]*(g[i][j]*sin(r[i]-r[j])-b[i][j]*cos(r[i]-r[j]))
+u[j]*u[j]*(-h[j][i]+b[j][i])+u[j]*u[i]*(g[j][i]*sin(r[j]-r[i])-b[j][i]*cos(r[j]-r[i]));
}
printf("\n======迭代后电压和相交值 ======\n");
getchar();
for(i=0;i printf("u[%d]=%8f r[%d]=%8f\n",i+1,u[i],i+1,r[i]);
printf("\n=====线路功率=======\n");
for(i=0;i {
getchar();
{ for(j=0;j { printf("s[%d][%d]=%8f+j(%8f)\n",i+1,j+1,p[i][j],q[i][j]);
printf("\n");
}
printf("\n");
}
}
printf("\n");
printf("\n=====线路损耗=====\n");
for(i=0;i {
getchar();
{ for(j=i+1;j { printf("ds[%d][%d]=%8f+j%8f\n",i+1,j+1,dp[i][j],dq[i][j]);
printf("\n");
}
printf("\n");
}
}
printf("\n");
printf("\n=====网络中损耗=====\n");
getchar();
for(i=0;i {
for(j=i+1;j { dSp+=dp[i][j];
dSq+=dq[i][j];
}
}
printf("dS=%8f+j(%8f)\n",dSp,dSq);
printf("\n======迭代次数========\n");
printf("ku=%d\n",ku);
printf("kr=%d\n",kr);
printf("\n=======数据保存=====\n");
fp1=fopen("jieguo.txt","w+");
{
fprintf(fp1,"线路潮流:\n");
for(i=0;i {
for(j=0;j fprintf(fp1,"s[%d][%d]=%8f+j%8f\n",i+1,j+1,p[i][j],q[i][j]);
}
fprintf(fp1,"\n");
fprintf(fp1,"支路损耗:\n");
for(i=0;i {
for(j=i+1;j fprintf(fp1,"ds[%d][%d]=%8f+j%8f\n",i+1,j+1,dp[i][j],dq[i][j]);
}
fprintf(fp1,"\n");
fprintf(fp1,"网络中损耗:\n");
fprintf(fp1,"dS=%8f+j(%8f)\n",dSp,dSq);
fprintf(fp1,"\n");
fprintf(fp1,"PV 节点功率:\n");
fprintf(fp1,"sv1=%8f+j(%8f)\n",pv1,qv1);
fprintf(fp1,"sv2=%8f+j(%8f)\n",pv2,qv2);
fprintf(fp1,"\n");
fprintf(fp1,"平衡节点功率:\n");
fprintf(fp1,"sp=%8f+j(%8f)\n",ps,qs);
fprintf(fp1,"\n");
fprintf(fp1,"节点电压和相交值:\n");
for(i=0;i fprintf(fp1,"u[%d]=%8f r[%d]=%8f\n",i+1,u[i],i+1,r[i]);
}
printf("\n===========结束==============\n");
getchar();
}
