powerflow.cpp

电力系统牛顿拉夫逊法的潮流程序

/*下面是数据输入格式说明：
为了减少不必要的运算，本潮流程序中将严格规定数据的录入格式.
具体要求如下：
1.数据共分六大组：即基本数据组、支路数据组、补偿电容数据组、并联电抗数据组、发电机节点数据组、负荷节点数据组。

2.在基本数据组内按照以下格式录入以下数据：
  节点个数，支路条数，对地补偿电容支路数，并联电抗数、发电机节点个数，负荷节点个数，系统平均电压，计算精度

3.在支路数据组内按照以下格式录入以下数据：
  对于线路支路：  支路的第一个节点号，支路的第二个节点号，支路的电阻，支路的电抗，支路的导纳
  对于变压器支路：变压器的第一个节点号，变压器的第二个节点号，变压器的电阻，变压器的电抗，变压器的非标准变比
注意：（1）支路两端的节点号应把小号排在前面，大号排在后面；
  （2）各支路按其小节点号的顺序排列；
  （3）对于线路支路，导纳应为Y,而不是计算中的Y/2;
  （4）对于变压器支路，变压器非标准变比所在侧的节点号应为负号。
另外，串联电容也应该在系统中被认为是一条支路，且串联电容的参数应当换算成为相应的电抗值，注意：此值为负。
特别提醒：本程序不支持同节点之间的并联线路，需要先行归化为一条线路，然后求解。

4.在补偿电容数据组内按照以下格式录入以下数据：
  所在节点号、导纳值

5.在并联电抗数据组内按照以下格式录入以下数据：
  所在节点号、电抗值

6.在发电机节点数据组内按照以下格式录入以下数据：
  发电机所在节点号，发电机所发有功，发电机所发无功，发电机节点的电压
注意：（1）对于PV节点，无功应为发电机的无功上限；
  （2）对于PV节点，电压为节点所需维持的电压，应为负值；
  （3）对于PQ节点，电压为计算过程中的节点的电压初值。

7.在负荷节点数据组内按照以下格式录入以下数据：
  负荷所在节点号，负荷有功，负荷无功，负荷节点的电压
注意：（1）负荷节点一般为PQ节点，其电压为计算过程中的节点的电压初值；
  （2）负荷的有功和无功均应写成注入节点功率的形式，即全部为负值；
  （3）特别的，系统的平衡节点应该被排在最后，即作为第N个节点，同时要求该节点为负荷节点；如果该点无负荷，则该点的负荷功率填零。

至于对所求解系统的描述性文字可以放在该文件的末尾，任意长度。
另外，本程序可以自动形成输出文件，其文件名为输入文件名后加"out"，且输出文件名的后缀确定为".dat"。例如，对于程序已给的5节点系统，输入文件名"ps-5"，后缀为".dat"，输出文件名为"ps-5out"，后缀为".dat"。为了保证名称的前后统一，因此建议：将输入文件的名字写成如下格式：****.dat;并且输入文件的文件名的长度应不大于5（不计后缀）。











*/


/**********************************************************************/
/*****************PQ分解法潮流程序-BB方案,允许迭代次数2000次******************************/
/*********
        作者：李常刚，山东大学电气工程学院02级3班，[e-mail]
        源于毕业设计题目，2006年4月3日
**********/
#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;
  }