pq_powerflow.txt

PQ分解法潮流程序,Ｃ语言编写

/**********************************************************************/
/*****************PQ分解法潮流程序-BB方案,允许迭代次数2000次******************************/

#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>

  struct Branch_Type
  {
   int i,j;
   double R,X,YK;
  }; 
  struct Compensation_Type
  {
   int i;
   double Y;
  };
  struct Reactance_Type
  {
   int i;
   double X;
  };
  struct Generator_Type
  {
   int i;
   double P,Q;
   double V;
  };
  struct Load_Type
  {
   int i;
   double P,Q;
   double V;
  };
  struct PVNode_Type
  {
   int i;
   double V;
  };
  struct Yii_Type
  {
   double G,B;
  };
  struct Yij_Type
  {
   double G,B;
   int j;
  };
  struct NodalPow
  {
   double P,Q;
  };
  struct NodalVol
  {
    double V,theta;
  };
  struct GeneratorPower
  {
     double P,Q;
  };
  struct U_Type
  {
     double value;
     int j;
  };
  void Datain(int Nb,int Nl,int Nc,int Nr,int Ng,FILE *fp,struct Branch_Type *Branch,struct Compensation_Type *Compensation,struct Reactance_Type *Reactance,struct Load_Type *Load,struct Generator_Type *Generator,struct PVNode_Type *PVNode,int *Npv);/*从文件中读入数据*/
  void AdmittanceMatrix(int N,int Nb,struct Yii_Type *Yii,struct Yii_Type *Yiil,struct Yij_Type *Yij,struct Yij_Type *Yijl,struct Branch_Type *Branch,int *NYseq1,int *NYsum1);/*导纳矩阵形成函数*/
  void AdmittanceMatrixAdd(int Nb,int Nc,int Nr,struct Yii_Type *Yii,struct Yii_Type *Yiil,struct Yij_Type *Yij,struct Yij_Type *Yijl,struct Branch_Type *Branch,struct Compensation_Type *Compensation,struct Reactance_Type *Reactance);/*导纳矩阵追加接地支路函数*/
  void Factorial(int flag,int N,int Npv,struct PVNode_Type *PVNode,int *NUsum,struct Yii_Type *Yii,struct Yij_Type *Yij,int *NYseq,double *D,struct U_Type *U);/*因子表形成函数*/
  void NodePower(int flag,int N,struct NodalVol *NodeVol,struct NodalPow *NodePow,struct Yii_Type *Yii,struct Yij_Type *Yij,int *NYseq);/*节点功率计算函数*/
  void Iteration(int flag,struct Generator_Type *Generator,struct Load_Type *Load,struct PVNode_Type *PVNode,struct NodalVol *NodeVol,struct NodalPow *NodePow,struct GeneratorPower *GenPower,int N,double *DI,double *MaxError,int *ErrNode);/*迭代函数*/
  void FormulaSolution(int flag,struct U_Type *U,double *D,int *NUsum,double *DI,int N,struct NodalVol *NodeVol,double V0);
  void NodeDataOutput(FILE *fp,struct NodalVol *NodeVol,struct Generator_Type *Generator,int N,struct GeneratorPower *GenPower,struct NodalPow *NodePow,struct Load_Type *Load,int Nl);
  void BranchDataOutput(FILE *fp,int Nb,int Nc,int Nr,struct Branch_Type *Branch,struct Compensation_Type *Compensation,struct Reactance_Type *Reactance,struct NodalVol *NodeVol);









main()
{
  int N;
  int Nb;
  int Nc;
  int Nr;
  int Ng;
  int Nl;
  double V0;
  double epsilon;
  struct Branch_Type *Branch;
  struct Compensation_Type *Compensation;
  struct Reactance_Type *Reactance;
  struct Generator_Type *Generator;
  struct Load_Type *Load;
  struct PVNode_Type *PVNode;
  int Npv=0;/*标记PV节点总个数*/
  struct Yii_Type *Yii,*Yiil;
  struct Yij_Type *Yij,*Yijl;
  int *NYseq,*NYsum;
  struct NodalPow *NodePow;
  struct NodalVol *NodeVol;
  struct GeneratorPower *GenPower;
  double *DI1,*DI2;
  double MaxError=0.0,MaxErrorTemp;
  int ErrNode,ErrNodeTemp;
  int Kp=1,Kq=1;


  int n,k,count;
  time_t now;
 



  int i;/*用以标记当前处理行号*/





 

  
  struct U_Type *U1,*U2;
  double *D1,*D2;
  int *NUsum1,*NUsum2;
  

  char FILENAME[20];
  FILE *fp;



 /* if((fp=fopen("README.txt","r"))==NULL)
  {
      printf("部分文件缺失，请确保文件的完整性!!\n");
      printf("按 ENTER 键返回！\n");
      scanf("%c",&ch);
      exit(1);
  }
*/

  /***********打开系统数据文件**********************/
  printf("Please Enter The Filename of The System:");
  gets(FILENAME);
  if((fp=fopen(FILENAME,"r"))==NULL)
  {
    printf("Cannot Find The File:%s\n",FILENAME);
    printf("Press ENTER to Escape!");
    getchar();
    exit(1);
  }



 /******************读系统数据信息****************/
  fscanf(fp,"%d,%d,%d,%d,%d,%d,%lf,%lf",&N,&Nb,&Nc,&Nr,&Ng,&Nl,&V0,&epsilon);


 /******************给支路、发电机、负荷和PV节点分配内存******************/
  Branch=malloc((Nb+1)*sizeof(struct Branch_Type));
  if(Branch==NULL)
  {
    printf("Cannot Get Memory for That Many Values!");
    printf("Press ENTER to Escape!");
    getchar();
    exit(1);
  }
  Compensation=malloc((Nc+1)*sizeof(struct Compensation_Type));
  if(Compensation==NULL)
  {
    printf("Cannot Get Memory for That Many Values!");
    printf("Press ENTER to Escape!");
    getchar();
    exit(1);
  }
  Reactance=malloc((Nr+1)*sizeof(struct Reactance_Type));
  if(Reactance==NULL)
  {
    printf("Cannot Get Memory for That Many Values!");
    printf("Press ENTER to Escape!");
    getchar();
    exit(1);
  }
  Generator=malloc((Ng+1)*sizeof(struct Generator_Type));
  if(Generator==NULL)
  {
    printf("Cannot Get Memory for That Many Values!");
    printf("Press ENTER to Escape!");
    getchar();
    exit(1);
  }
  Load=malloc((Nl+1)*sizeof(struct Load_Type));
  if(Load==NULL)
  {
    printf("Cannot Get Memory for That Many Values!");
    printf("Press ENTER to Escape!");
    getchar();
    exit(1);
  }
  PVNode=malloc(N*sizeof(struct PVNode_Type));
  if(PVNode==NULL)
  {
    printf("Cannot Get Memory for That Many Values!");
    printf("Press ENTER to Escape!");
    getchar();
    exit(1);
  }
/**********内存分配完成***************/





  Datain(Nb,Nl,Nc,Nr,Ng,fp,Branch,Compensation,Reactance,Load,Generator,PVNode,&Npv);/*从文件中读入数据*/




  for(n=0;1;n++)/*打开输出文件*/
  {
    if(FILENAME[n]=='.')
    {
        FILENAME[n]='\0';
        strcat(FILENAME,"out.dat");
        break;
    }
  }
  if((fp=fopen(FILENAME,"w"))==NULL)/*打开数据输出文件*/
  {
    printf("Cannot Find The File:%s\n",FILENAME);
    printf("Press ENTER to Escape!");
    getchar();
    exit(1);
  }
/************为导纳矩阵分配内存,并形成导纳矩阵*******************/

  Yii=malloc((N+1)*sizeof(struct Yii_Type));
  Yiil=malloc((N+1)*sizeof(struct Yii_Type));
  Yij=malloc((N+1)*sizeof(struct Yij_Type));
  Yijl=malloc((N+1)*sizeof(struct Yij_Type));
  NYseq=malloc((N+1)*sizeof(int));
  NYsum=malloc((N+1)*sizeof(int));


  AdmittanceMatrix(N,Nb,Yii,Yiil,Yij,Yijl,Branch,NYseq,NYsum);/*导纳矩阵形成函数*/



/***********分配及导纳矩阵的形成完成*******************/


/******************下面确认导纳矩阵的正确性（仅用于调试）*********************

  fprintf(fp,"The Admittance Matrix is:\n");
  printf("The Admittance Matrix is:\n");

  for(n=1;n<=N;n++)
  {
    fprintf(fp,"Yii[%d]=%lf+j%lf\n",n,Yii[n].G,Yii[n].B);
    fprintf(fp,"Yiil[%d]=%lf+j%lf\n",n,Yiil[n].G,Yiil[n].B);


    printf("Yii[%d]=%lf+j%lf\n",n,Yii[n].G,Yii[n].B);
    printf("Yiil[%d]=%lf+j%lf\n",n,Yiil[n].G,Yiil[n].B);
  }
  for(n=1;n<=N-1;n++)
  {
    for(i=NYseq[n];i<=NYseq[n+1]-1;i++)
    {
        fprintf(fp,"Yij[%d][%d]=%lf+j%lf\n",n,Yij[i].j,Yij[i].G,Yij[i].B);
        fprintf(fp,"Yijl[%d][%d]=%lf+j%lf\n",n,Yijl[i].j,Yijl[i].G,Yijl[i].B);

        printf("Yij[%d][%d]=%lf+j%lf\n",n,Yij[i].j,Yij[i].G,Yij[i].B);
        printf("Yijl[%d][%d]=%lf+j%lf\n",n,Yijl[i].j,Yijl[i].G,Yijl[i].B);
    }
  }

*****************确认导纳矩阵的正确性完成（仅用于调试）*********************/



/*****************下面形成因子表,采用BX方案(在B'中不忽略电阻，在B"中忽略电阻)*****************/
  U1=malloc((N-1)*(N-2)/2*sizeof(struct U_Type));
  U2=malloc((N-1)*(N-2)/2*sizeof(struct U_Type));
  D1=malloc(N*sizeof(double));
  D2=malloc(N*sizeof(double));
  NUsum1=malloc(N*sizeof(int));
  NUsum2=malloc(N*sizeof(int));

  Factorial(1,N,Npv,PVNode,NUsum1,Yii,Yij,NYseq,D1,U1);/*因子表形成函数*/
/**************************

  count=0;
  for(n=1;n<=N-1;n++)
  {
      fprintf(fp,"D[%d]=%lf;\n",n,D1[n]);

     }
  for(i=1;i<=N-1;i++)
  {
      
      for(n=1;n<=NUsum1[i];n++)
      {
          fprintf(fp,"U[%d][%d]=%lf\n",i,U1[count+n].j,U1[count+n].value);

      }
      count+=NUsum1[i];
  }


    count=0;

**************************************/
  AdmittanceMatrixAdd(Nb,Nc,Nr,Yii,Yiil,Yij,Yijl,Branch,Compensation,Reactance);/*导纳矩阵追加接地支路函数*/
  Factorial(2,N,Npv,PVNode,NUsum2,Yii,Yij,NYseq,D2,U2);/*因子表形成函数*/


/********************************

  for(n=1;n<=N-1;n++)
  {
      fprintf(fp,"D[%d]=%lf;\n",n,D2[n]);
 
  }
  for(i=1;i<=N-1;i++)
  {
      
      for(n=1;n<=NUsum2[i];n++)
      {
         fprintf(fp,"U[%d][%d]=%lf\n",i,U2[count+n].j,U2[count+n].value);
      }
      count+=NUsum2[i];
  }



    
 
*************************************/



/****************下面利用所求得的因子表进行迭代求解******************/

   

  DI1=malloc(N*sizeof(double));
  DI2=malloc(N*sizeof(double));
  NodePow=malloc((N+1)*sizeof(struct NodalPow));
  NodeVol=malloc((N+1)*sizeof(struct NodalVol));


  /****先送电压初值******/
  for(i=1;i<=N;i++)
  {
    NodeVol[i].V=V0;
    NodeVol[i].theta=0.0;
  }
  /*for(i=1;i<=Ng;i++)
  {
      n=Generator[i].i;
      NodeVol[n].V=Generator[i].V;
      NodeVol[n].theta=0.0;
  }
   for(i=1;i<=Nl;i++)
  {
      n=Load[i].i;
      NodeVol[n].V=Load[i].V;
      NodeVol[n].theta=0.0;
  }*/

  for(n=1;n<=Npv;n++)
  {
    i=PVNode[n].i;
    NodeVol[i].V=PVNode[n].V;
  }
    
  GenPower=malloc((N+1)*sizeof(struct GeneratorPower)); /*为发电机所发功率数组分配内存*/  
  fprintf(fp,"\t系统潮流计算结果：\n(1)迭代过程纪录：\n迭代次数\t    \t有功迭代\t\t无功迭代\n\t\t   ΔPmax\tP-Node\t ΔQmax\t\tQ-Node\n");
  for(k=0;1;k++)
  {
    fprintf(fp,"   %2d:\t",k);
    if(Kp==1)
    {
        NodePower(1,N,NodeVol,NodePow,Yii,Yij,NYseq);/*节点功率计算函数*/
        Iteration(1,Generator,Load,PVNode,NodeVol,NodePow,GenPower,N,DI1,&MaxError,&ErrNode);/*迭代*/
        fprintf(fp,"\t%10.7lf\t %d\t",MaxError,ErrNode);
        if(MaxError>=epsilon)
            FormulaSolution(1,U1,D1,NUsum1,DI1,N,NodeVol,V0);/*线性方程求解*/
        else
            Kp=0;
    }
    else
        fprintf(fp,"\t\t\t\t");
    if(Kq==1)
    {
        NodePower(2,N,NodeVol,NodePow,Yii,Yij,NYseq);/*节点功率计算函数*/
        Iteration(2,Generator,Load,PVNode,NodeVol,NodePow,GenPower,N,DI2,&MaxError,&ErrNode);/*迭代*/
        fprintf(fp,"%10.7lf\t %d\n",MaxError,ErrNode);
        if(MaxError>=epsilon)
            FormulaSolution(2,U2,D2,NUsum2,DI2,N,NodeVol,V0);/*线性方程求解*/
        else
            Kq=0;

    }
    else
        fprintf(fp,"\n");
    if(Kp==0&&Kq==0)
        break;
    if(k>2000)
    {
        fprintf(fp,"\n迭代次数超过2000次，系统不收敛!\n");
        fprintf(fp,"本系统为%d个节点、%d条支路系统，含有%d个补偿电容和%d个并联电抗\n",N,Nb,Nc,Nr);
        now=time(NULL);
        fprintf(fp,"\t\t处理时间： %s\n",ctime(&now));
        /*fclose(fp);
        exit(1);*/
        goto DATA;
    }

    
    printf("COMPUTING......\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b");
  }


  fprintf(fp,"\n总计迭代次数: %d 次!\n\n",k+1);
  fprintf(fp,"(2)潮流计算结果(节点电压)：\n Node\t\tV\t\tθ\t\tP\t\tQ\n");
  NodeDataOutput(fp,NodeVol,Generator,N,GenPower,NodePow,Load,Nl);
  fprintf(fp,"\n(3)潮流计算结果(支路功率):\n    Branch\t\tP\t\tQ\t\tQc\t\tQr\n");
  BranchDataOutput(fp,Nb,Nc,Nr,Branch,Compensation,Reactance,NodeVol);
  fprintf(fp,"\n本系统为%d个节点、%d条支路系统，含有%d个补偿电容和%d个并联电抗\n",N,Nb,Nc,Nr);