#include<stdio.h>
#include<math.h>
#define eP 0.00001
#define eQ 0.00001
#define Y_(i,j) (*(*(Y+i)-1+j))
#define Yij (*(Yi+j)
#define Yji (*(*(Y+j)-j+i)
#define Pji Yji.G*cos(temp)+Yji.B*sin(temp)
#define Pij Yij.G*cos(temp)+Yij.B*sin(temp)
#define Qij Yij.G*sin(temp)+Yij.B*cos(temp)
#define Qji Yji.G*sin(temp)+Yji.B*cos(temp) /*宏定义*/
struct Nodetype /*节点功率*/
{
float P,Q;
};
struct Linetype/*线路类型*/
{
float G,B,B0,K;
};
void in_node();/*输入节点功率*/
void in_line();/*输入线路类型*/
void B_Form();/*形成导纳矩阵*/
void factor();/*求因子表*/
void solve(float **B,float *X,int N);/*解方程组*/
void PrtNode();/*打印输出节点参数*/
void ErrorMsg(int Flag);/*错误信息*/
int Node;/*节点数*/
int num;/*原始节点序号*/
int NP,NQ=0;/*PQ节点数*/
struct Nodetype *No;/*节点数据*/
struct Linetype **Y;/*线路参数*/
float **BP,**BQ;/*雅克比矩阵*/
float *V;/*节点电压有效值*/
float *Dlta;/*节点电压相角*/
unsigned int count=0;
int i,j,k,m;
float temp;
char *Type;
FILE *in,*out;
int main(void)
{
int kp=1,kq=1;
float *dP,*dPi,*dQ,*dQi;
float Dltai;
struct Linetype *Yi;
struct Nodetype *Noi;
float tP,tQ;
if((in=fopen("data.txt","r"))==NULL)ErrorMsg(1);
if((out=fopen("out.txt","w"))==NULL)ErrorMsg(2);
in_node();
in_line();
B_Form();
factor();
for(i=0;i<NQ;i++)*(V+i)=1;
dP=(float*)malloc(sizeof(float)*NP);
dQ=dP;
loop:
if(kp==0&&kq==0)goto loopEnd;
count++;
if(count==65535)ErrorMsg(99);
kp=0;
for(i=0;i<NP;i++)
{
dPi=dP+i;
Yi=*(Y+i)-i;
Dltai=*(Dlta+i);
*dPi=0;
for(j=0;j<Node;j++)
{
temp=Dltai-*(Dlta+j);
if(i>j)*dPi+=*(V+j)*(Pji);
else *dPi+=*(V+j)*(Pij);
}
*dPi*=*(V+i);
*dPi*=(*(No+i));
if(fabs(*dPi)>0x8fffffff)ErrorMsg(99);
if(fabs(*dPi)>eP)kp=1;
*dPi/=*(V+i);
}
if(kp==0)goto loopQ;
slove(BP,dP,NP);
for(i=0;i<NP;i++)*(Dlta+i)+=(*(dP+i)/(*(V+i)));
loopQ:
if(kp==0&&kq==0)goto loopEnd;
kq=0;
for(i=0;i<NQ;i++)
{
dQi=dQ+i;
Yi=*(Y+i)-i;
Dltai=*(Dlta+i);
*dQi=0;
for(j=0;j<Node;j++)
{
temp=Dltai-*(Dlta+j);
if(i>j)*dQi+=*(V+j)*(Qji);
else *dQi+=*(V+j)*(Qij);
}
*dQi*=*(V+i);
*dQi*=(*(No+i)).Q-*dQi;
if(fabs(*dQi)>0x8fffffff)ErrorMsg(99);
if(fabs(*dQi)>eQ)kq=1;
*dQi/=*(V+i);
}
if(kp==0)goto loop;
slove(BQ,dQ,NQ);
for(i=0;i<NQ;i++)*(V+i)+=*(dQ+i);
goto loop;
loopEnd:
free(dP);
for(i=0;i<NQ;i++)
{
Noi=No+i;
Yi=*(Y+i)-i;
Dltai=*(Dlta+i);
for(j=0;j<Node;j++)
{
temp=Dltai-*(Dlta+j);
if(i>j)(*Noi).Q+=*(V+j)*(Qji);
else (*Noi).Q+=*(V+j)*(Qij);
}
(*Noi).Q*=*(V+i);
}
i=NP;
Dltai=*(Dlta+i);
for(j=0;j<Node;j++)
{
temp=Dltai-*(Dlta+j);
(*Noi).P+=*(V+j)*(Pji);
}
(*Noi).P*=*(V+i);
fprintf(out,"\n\n【潮流计算结果(节点)】\n",count-1);
PrtNode();
fprintf(out,"\n\n【潮流计算结果(线路)】\n");
fprintf(out,"线路 P Q\n");
for(k=0;k<Node;k++)
{
i=*(num+k);
Noi=No+i;
Yi=*(Y+i)-i;
Dltai=*(Dlta+i);
for(m=0;m<Node;m++0
{
j=*(num+m);
if(j==i)continue;
temp=Dltai-*(Dlta+j);
if(j<i)
{
if(Yji.B==0)continue;
tP=*(V+j)*(Pji);
tP=*(V+i)*Yji.G-tP;
tP*=*(V+i);
tQ=*(V+j)*(Qji);
tQ-=*(V+i)*(Yji.B-Yji.B0/Yji.k);
tQ*=*(V+i);
}
else
{
if(Yij.B==0)continue;
tP=*(V+j)*(Pij);
tP=*(V+i)*Yij.G-tP;
tP*=*(V+i);
tQ=*(V+j)*(Qij);
tQ-=*(V+i)*(Yij.B-Yij.B0);
tQ*=*(V+i);
}
fprintf(out,"S[%d,%d]=(%10.6f,%10.6f)\n",k+1,m+1,-tP,-tQ)
}
}
fclose(out);
system("cmd/c,start out.txt");
return(0);
}
*************************************************************************************
void in_node()
{
struct Nodetype *Noi;
fscanf(in,"%d,%d",&Node,&k);
NP=Node-1;
num=(int*)malloc(sizeof(int)*Node);
V=(float*)malloc(sizeof(float)*Node);
Dlta=(float*)malloc(sizeof(float)*Node);
No=(struct Nodetype*)malloc(sizeof(struct Nodetype)*Node);
j=1;
while(k!=0)
{
switch(k)
{
case 1:k=NQ;NQ++;break;
case 2:k=NP-j;j++;break;
case 3:k=NP;break;
default:ErrorMsg(3);
}
Noi=No+k;
fscanf(in,"%d%f%f%f%f",&i,&(Noi).Q,V+k,Dlta+k);
i--;
*(num+i)=k;
fscanf(in,"%d",&k);
}
if(NQ+j!=Node)ErrorMsg(4);
fprintf(out,"【节点参数表】\n");
PrtNode();
fprintf(out,"总节点:%d\n PQ节点%d\n PV节点%d\n",Node,NQ,NP-NQ);
}
*******************************************************************
void in_line()
{
struct Linetype *Yi;
float R,X,k,B;
m=sizeof(struct Linetype);
Y=(struct Linetype**)malloc(m*Node);
for(i=0;i<Node;i++)
{
*(Y+i)=(struct Linetype*)malloc(m*(Node-i));
Yi=*(Y+i)-i;
for(j=i;j<Node;j++)
{Yij.G=Yij.B=Yij.B0=Yij.k;}
}
while(feof(in)==0)
{
fscanf(in,"%d%d%f%f%f%f",&i,&j,&R,&X,&k,&B);
i--;
j--;
i=*(num+i);
j=*(num+j);
(*(*(Y+i))).B+=B;
(*(*(Y+j))).B+=B;
if(k!=1.0)
{
X*=k;R=0;
temp=(1-k)/X;
(*(*(Y+j))).B+=temp;
(*(*(Y+j))).B+=-(temp/k);
B=temp;
k=-k;
}
if(i>j){temp=i;i=j;j=temp;k=1/k;B*=k;}
Yi=*(Y+i)-i;
Yij.B0=B;
Yij.k=k;
temp=R*R+X*X;
R/=temp;
X/=temp;
Yij.G=-R;
Yij.B=X;
(*(*(Y+i))).G+=R;
(*(*(Y+i))).B+=-X;
(*(*(Y+j))).G+=R;
(*(*(Y+j))).B+=X;
}
fclose(in);
fprintf(out,"\n【节点导纳矩阵】\n");
for(k=0;k<Node;k++)
{
i=k;
i=*(num+k);
for(j=0;j<K;j++)
fprintf(out,"\t\t\t");
for(m=k;m<Node;m++)
{
j=*(num+m);
if(i<j)fprintf(out,"(%10.6f,%10.6f)",Y_(i,j).G,Y_(i,j).B);
else fprintf(out,"(%10.6f,%10.6f)",Y_(j,i).G,Y_(j,i).B);
}
fprintf(out,"\n");
}
}
*************************************************************************************
void B_Form()
{
float *BPi,*BQi;
struct Linetype *Yi;
int size=sizeof(float);
BP=(float**)malloc(size*NP);
for(i=0;i<NP;i++)*(BP+i)=(float*)malloc(size*(NP-i));
for(i=0;i<NP;i++)
{
BPi=*(BP+i)-i;
Yi=*(Y+i)-i;
for(j=i;j<NP;j++)
*(BPi+j)=Yij.B;
}
BQ=BP;
}
***************************************************************************************
void factor()
{
float *BPi,*BPk,*BQi;
for(i=0;i<NP;i++)
{
BPi=*(BP+i)-i;
for(k=0;k<i;k++)
{
BPk=*(BP+k)-k;
temp=(*(BPk+i))/(*(BPk+k));
for(j=i;j<NP;j++)(*(BPi+j))-=temp*(*(BPk+j));
}
*(BPi+i)=1/(*(BPi+i));
for(j=i+1;j<NP;j++)*(BPi+j)=*(BPi+i);
}
}
****************************************************************************************
void slove(float **B,float *X,int N)
{
float *Bi,*Xi;
for(i=0;i<N;i++)*(X+i)=-*(X+i);
for(i=0;i<N;i++)
{
Bi=*(B+i)-i;
Xi=X+i;
for(j=i+1;j<N;j++)*(X+j)-=*(Bi+j)**Xi;
*Xi*=*(Bi+i);
}
for(i=N-1;i>=0;i--)
{
Bi=*(B+i)-i;
Xi=X+i;
for(j=N-1;j>i;j--)*Xi-=*(Bi+j)**(X+j);
}
}
****************************************************************************************
void PrtNode()
{
struct Nodetype *Noi;
fprinf(out,"节点类型 P Q V δ\n");
for(i=0;i<Node;i++)
{
j=*(num+i);
Noi=No+j;
if(j<NQ)Type="PQ";else Type="PV";
if(j==NP)Type="BS";
fprinf(out,"%3d%s%10.6f%10.6f%10.6f%10.6f\n",++i,Type,*(Noi).p,*(Noi).Q,*(V+j),*(Dlta+j)/0.017453);
}
}
*****************************************************************************************
void ErrorMsg(int Flag)
{
switch(Flag)
{
case 1:printf("\n\tError(1):Failed to open file\"Data.txt\"!");break;
case 2:printf("\n\tError(2):Failed to creat file\"out.txt\"!");break;
case 3:printf("\n\tError(3):node data error,please check!");break;
case 4:printf("\n\tError(4):lack data node,pleak check!");break;
case 99:printf("k=%d\n\tError(99):it's emanative",count);break;
}
getch();
fclose(out);
exit(Flag);
}