readital.c

来自「电力系统分析计算 学习调试程序 UNIX / LINUX / CYGWIN 系统」· C语言 代码 · 共 1,096 行 · 第 1/3 页

      if (ACptr->N==0) {
        Nac++;
        ACptr->Num=ACptr->N=Nac;
        KV=Vlevels[KVl];
        sprintf(ACptr->Name,"%5s %6.0lf",Code,KV);
        Vlmax=Vmax[KVl];
        Vlmin=Vmin[KVl];
        if (KV==0){
          fCustomPrint(stderr,"Input Line-> %d\n%s",LineNum,Line);
          ErrorHalt("Base voltage is zero.");
          KV=1;
        }
        else if (Vlmax>Vlmin && Vlmax>=KV && Vlmin<=KV) {
          ACptr->Vlmax=Vlmax/KV;
          ACptr->Vlmin=Vlmin/KV;
        }
        ACptr->KV=KV;
      }
      ACptr->Pg+=-GetValue(Line,14,6,0)/Sn;
      ACptr->PgMax+=-GetValue(Line,25,5,0)/Sn;
      if (ACptr->PgMax==0) ACptr->PgMax=99999999.;
      ACptr->Smax+=-GetValue(Line,30,5,0)/Sn;
      if (ACptr->Smax==0) ACptr->Smax=99999999.;
      if (ACptr->PgMax>ACptr->Smax) {
        ACptr->PgMax=ACptr->Smax;
        fCustomPrint(stderr,"***Warning: Bus %d %s has PgMax > Smax.\n",ACptr->N,ACptr->Name);
        fCustomPrint(stderr,"            PgMax will be set to Smax.\n");
      }
      if (ACptr->Pg>ACptr->PgMax) {
        ACptr->Pg=ACptr->PgMax;
        fCustomPrint(stderr,"***Warning: Bus %d %s has Pg > PgMax.\n",ACptr->N,ACptr->Name);
        fCustomPrint(stderr,"            Pg will be set to PgMax.\n");
      }
      ACptr->Qg+=-GetValue(Line,36,6,0)/Sn;
      ACptr->Qmax+=-GetValue(Line,47,5,0)/Sn;
      ACptr->Qmin+=-GetValue(Line,43,4,0)/Sn;
      if (ACptr->Qmax<=ACptr->Qmin) {
        fCustomPrint(stderr,"Input Line-> %d\n%s",LineNum,Line);
        ErrorHalt("AC bus Q limits are wrong: Qmin >= Qmax.");
      }
    }

    /* Synchronous Condenser */
    else if (Line[79]=='S') {
      N=GetInt(Line,1,3);
      GetStr(Line,4,1,1,Zone);
      KVl=GetInt(Line,5,1);
      sprintf(Code,"%3d%1s%1d",N,Zone,KVl);
      ACptr=ACbusInList2(0,Code,Nac,1,TRUE);
      if (ACptr->N==0) {
        Nac++;
        ACptr->Num=ACptr->N=Nac;
        KV=Vlevels[KVl];
        sprintf(ACptr->Name,"%5s %6.0lf",Code,KV);
        Vlmax=Vmax[KVl];
        Vlmin=Vmin[KVl];
        if (KV==0){
          fCustomPrint(stderr,"Input Line-> %d\n%s",LineNum,Line);
          ErrorHalt("Base voltage is zero.");
          KV=1;
        }
        else if (Vlmax>Vlmin && Vlmax>=KV && Vlmin<=KV) {
          ACptr->Vlmax=Vlmax/KV;
          ACptr->Vlmin=Vlmin/KV;
        }
        ACptr->KV=KV;
      }
      ACptr->Qg+=-GetValue(Line,36,6,0)/Sn;
      ACptr->Qmax+=-GetValue(Line,47,5,0)/Sn;
      ACptr->Qmin+=-GetValue(Line,43,4,0)/Sn;
      if (ACptr->Qmax<=ACptr->Qmin) {
        fCustomPrint(stderr,"Input Line-> %d\n%s",LineNum,Line);
        ErrorHalt("AC bus Q limits are wrong: Qmin >= Qmax.");
      }
    }

    /* Shunt Capacitor */
    else if (Line[79]=='Q') {
      N=GetInt(Line,1,3);
      GetStr(Line,4,1,1,Zone);
      KVl=GetInt(Line,5,1);
      sprintf(Code,"%3d%1s%1d",N,Zone,KVl);
      ACptr=ACbusInList2(0,Code,Nac,1,TRUE);
      if (ACptr->N==0) {
        Nac++;
        ACptr->Num=ACptr->N=Nac;
        KV=Vlevels[KVl];
        sprintf(ACptr->Name,"%5s %6.0lf",Code,KV);
        Vlmax=Vmax[KVl];
        Vlmin=Vmin[KVl];
        if (KV==0){
          fCustomPrint(stderr,"Input Line-> %d\n%s",LineNum,Line);
          ErrorHalt("Base voltage is zero.");
          KV=1;
        }
        else if (Vlmax>Vlmin && Vlmax>=KV && Vlmin<=KV) {
          ACptr->Vlmax=Vlmax/KV;
          ACptr->Vlmin=Vlmin/KV;
        }
        ACptr->KV=KV;
      } else KV=ACptr->KV;
      KV=GetValue(Line,43,4,0)/KV;
      if (KV!=0)  ACptr->B+=(-GetValue(Line,47,5,0)/Sn)/(KV*KV);
    }

  /* --------------- AC element data -------------------------------- */
    else if (Line[79]=='L' || Line[79]=='T') {
      N=GetInt(Line,1,3);
      GetStr(Line,4,1,1,Zone);
      KVl=GetInt(Line,5,1);
      sprintf(Code,"%3d%1s%1d",N,Zone,KVl);
      ACptr=ACbusInList2(0,Code,Nac,1,TRUE);
      if (ACptr->N==0) {
        Nac++;
        ACptr->Num=ACptr->N=Nac;
        KV=Vlevels[KVl];
        sprintf(ACptr->Name,"%5s %6.0lf",Code,KV);
        Vlmax=Vmax[KVl];
        Vlmin=Vmin[KVl];
        if (KV==0){
          fCustomPrint(stderr,"Input Line-> %d\n%s",LineNum,Line);
          ErrorHalt("Base voltage is zero.");
          KV=1;
        }
        else if (Vlmax>Vlmin && Vlmax>=KV && Vlmin<=KV) {
          ACptr->Vlmax=Vlmax/KV;
          ACptr->Vlmin=Vlmin/KV;
        }
        ACptr->KV=KV;
      } else KV=ACptr->KV;
      N=GetInt(Line,7,3);
      GetStr(Line,10,1,1,Zone);
      KVl=GetInt(Line,11,1);
      sprintf(Code,"%3d%1s%1d",N,Zone,KVl);
      ACptrp=ACbusInList2(0,Code,Nac,1,TRUE);
      if (ACptrp->N==0) {
        Nac++;
        ACptrp->Num=ACptrp->N=Nac;
        KVp=Vlevels[KVl];
        sprintf(ACptrp->Name,"%5s %6.0lf",Code,KVp);
        Vlmax=Vmax[KVl];
        Vlmin=Vmin[KVl];
        if (KVp==0){
          fCustomPrint(stderr,"Input Line-> %d\n%s",LineNum,Line);
          ErrorHalt("Base voltage is zero.");
          KVp=1;
        }
        else if (Vlmax>Vlmin && Vlmax>=KVp && Vlmin<=KVp) {
          ACptrp->Vlmax=Vlmax/KVp;
          ACptrp->Vlmin=Vlmin/KVp;
        }
        ACptrp->KV=KVp;
      } else KVp=ACptrp->KV;
      if (ACptr==ACptrp){
        fCustomPrint(stderr,"Input Line-> %d\n%s",LineNum,Line);
        ErrorHalt("Both AC element buses are the same.");
      }
      if (EstimatePl) {
        if (Line[20]=='T') {
          fCustomPrint(stderr,"***Warning: The P load estimate at bus %d %s may be incorrect\n",ACptrEstimate->N,ACptrEstimate->Name);
          fCustomPrint(stderr,"            due to lack of P flow information in the corresponding element data:\n");
          fCustomPrint(stderr,"Input Line-> %d\n%s",LineNum,Line);
        }
        else {
          if (ACptrEstimate==ACptr) ACptrEstimate->Pl+=-GetValue(Line,14,6,0)/Sn;
          else ACptrEstimate->Pl+=GetValue(Line,14,6,0)/Sn;
        }
      }
      if (EstimateQl) {
        if (Line[42]=='T') {
          fCustomPrint(stderr,"***Warning: The Q load estimate at bus %d %s may be incorrect\n",ACptrEstimate->N,ACptrEstimate->Name);
          fCustomPrint(stderr,"            due to lack of Q flow information in the corresponding element data:\n");
          fCustomPrint(stderr,"Input Line-> %d\n%s",LineNum,Line);
        }
        else {
          if (ACptrEstimate==ACptr) ACptrEstimate->Ql+=-GetValue(Line,36,6,0)/Sn;
          else ACptrEstimate->Ql+=GetValue(Line,36,6,0)/Sn;
        }
      }
      if (KV>KVp) {
        KVs=KV; KV=KVp; KVp=KVs;
        ACptrs=ACptr; ACptr=ACptrp; ACptrp=ACptrs;
      }
      if (Line[79]=='L') strcpy(Type,"L");
      else {
        Tmax=GetValue(Line,67,5,2);
        Tmin=GetValue(Line,72,5,2);
        /*
        if(Tmax>Tmin) strcpy(Type,"R");
        else strcpy(Type,"T");
        */
        strcpy(Type,"T");
      }
      GetStr(Line,13,1,1,Ckt);
      if(!strcmp(Ckt," ")) strcpy(Ckt,"0");
      Eptr=ElemInList2(ACptr,ACptrp,NacEl,Type,Ckt,TRUE);
      if (Eptr!=NULL) {
        Zb=KVp*KVp/Sn;
        R=GetValue(Line,45,7,5)/Zb;
        X=GetValue(Line,52,8,5)/Zb;
        if (fabs(R)<1e-10 && fabs(X)<1e-10) {
          fCustomPrint(stderr,"Input Line-> %d\n%s",LineNum,Line);
          ErrorHalt("AC element is a short circuit. Try eliminating it.");
          G=B=0;
        } else {
          G=R/(R*R+X*X);
          B=-X/(R*R+X*X);
        }
        B2=B1=GetValue(Line,60,7,3)*Zb/1000000./2;
        if (Line[79]=='T') {
          Taps=KV/KVp;
          Tap=GetValue(Line,29,5,0)*Taps;
        }
        else Tap=Taps=1;
        Ang=0;
        Eptr->Imax=GetValue(Line,23,5,0)/(Sn*1000/(sqrt(3.0)*KVp));
        Eptr->G=G;
        Eptr->B=B;
        Eptr->B1=B1;
        Eptr->B2=B2;
        Eptr->Tap=Tap;
        Eptr->Taps=Taps;
        Eptr->Ang=Ang;
        GetStr(Line,4,1,1,Eptr->Zone);
        strcpy(Eptr->Owner,Eptr->Zone);
        NacEl++;
        Eptr->Meter=ACptr;
        strcpy(Eptr->Type,Type);
        /*
        if (!strcmp(Type,"R")) {
          Eptr->Tmax=Tmax*Taps;
          Eptr->Tmin=Tmin*Taps;
          NregV++;
          if (!strcmp(ACptrp->Type,"B")) strcpy(ACptrp->Type,"BT");
          ACptrp->Reg=AddElemToList(ACptrp->Reg,Eptr);
          Eptr->Cont=ACptrp;
        }
        */
      }
    }


  /* -------------------- DC data -------------------------------- */
    else if (!strncmp(Line,"BD ",3) || !strncmp(Line,"BZ ",3)||
             !strncmp(Line,"LD ",3)) ReadEPRIdc(Line);


                        /* FACTS */

 /* ---------------------- SVC data ------------------------ */
    else if (!strncmp(Line,"FS ",3)) ReadSVC(Line);

 /* ---------------------- TCSC data ------------------------ */
    else if (!strncmp(Line,"FC ",3)) ReadTCSC(Line);

 /* ---------------------- STATCOM data ------------------------ */
    else if (!strncmp(Line,"FT ",3)) ReadSTATCOM(Line);

                       /* END FACTS */

    else if (Line[79]=='F') break;
    else {
      fCustomPrint(stderr,"Input Line-> %d\n%s",LineNum,Line);
      fCustomPrint(stderr,"***Warning: The program will ignore this line.\n");
    }
  }
  fclose(InputDataFile);
  MaxIter=50;

  if (!ReadADDfile()) for (ACptr=dataPtr->ACbus; ACptr!=NULL; ACptr=ACptr->Next){
    Aptr=(AreaData *) AreaInList(0,ACptr->Zone,Narea);
    if (Aptr->N==0) {
      Narea++;
      Aptr->N=Narea;
      strcpy(Aptr->Zone[1],ACptr->Zone);
    }
    ACptr->Area=Aptr;
    if (strpbrk(ACptr->Type,"S")) Aptr->Slack=Aptr->BSptr=ACptr;
    else if (strpbrk(ACptr->Type,"Q") && Aptr->Slack==NULL) Aptr->Slack=Aptr->BSptr=ACptr;
  }

}


/* --------------- ReadADDfile ------------------- */
BOOLEAN ReadADDfile()
/* Read ADD COLAS file with new voltage information, areas, and
   generator and load data for voltage collapse studies. */
{
  FILE *InputFile;
  char *Name,Line[BUFLEN],Code[20],Ckt[2],Zone[5];
  ACbusData *ACptr,*ACptrp,*ACptrs;
  AreaData *Aptr;
  ElementList *ELptr;
  ElementData *Eptr;
  VALUETYPE KV,KVp,KVs,KVmax,KVmin,val,Tap,Taps,Tmax,Tmin,Q,Qmax,Qmin,Pn,Qn,Sum=0;
  BOOLEAN flagAreas=FALSE,flagPrint=TRUE,flagScards=FALSE;
  INDEX N,NJcard=0,N2SVCarea=0,KVl;



  Name=NameParameter('6');
  if (!NullName(Name) && (InputFile=OpenInput(Name))!=NULL) {
    LineNum=0;
    for (;;) {
      if (fgets(Line,BUFLEN,InputFile)==NULL) break;
      LineNum++;

      /* --------------- Comments ----------------------------- */
      if (Line[79]=='C') continue;

      /* --------------- AC buses ----------------------------- */
      else if (Line[79]=='K') {
        N=GetInt(Line,1,3);
        GetStr(Line,4,1,1,Zone);
        KVl=GetInt(Line,5,1);
        sprintf(Code,"%3d%1s%1d",N,Zone,KVl);
        ACptr=ACbusInList2(0,Code,Nac,1,FALSE);
        if (ACptr==NULL) {
          fCustomPrint(stderr,"Input Line-> %d\n%s",LineNum,Line);
          fCustomPrint(stderr,"***Warning: This bus has not been defined on the main input data file.\n");
          fCustomPrint(stderr,"            This line in the ADD file will be ignored.\n");
        } else {
          KV=GetValue(Line,6,7,0);
          if (KV==0){
            fCustomPrint(stderr,"Input Line-> %d\n%s",LineNum,Line);
            fCustomPrint(stderr,"***Warning: This bus has a zero base bus voltage.\n");
            fCustomPrint(stderr,"            This line in the ADD file will be ignored.\n");
          } else {
            sprintf(ACptr->Name,"%5s %6.0lf",Code,KV);
            KVp=ACptr->KV;
            val=KVp/KV;
            ACptr->VCont=ACptr->V=ACptr->V*val;
            ACptr->KV=KV;
            KVmin=GetValue(Line,13,7,0);
            KVmax=GetValue(Line,20,7,0);
            if (KVmax>KVmin && KVmax>=KV && KVmin<=KV) {
              ACptr->Vlmax=KVmax/KV;
              ACptr->Vlmin=KVmin/KV;
            }
            val=KV*KV/(KVp*KVp);
            ACptr->B=ACptr->B*val;
            for(ELptr=ACptr->Elem; ELptr!=NULL; ELptr=ELptr->Next) {
              Eptr=ELptr->Eptr;
              if (Eptr->To==ACptr) {
                val=KV*KV/(KVp*KVp);
                Eptr->G=Eptr->G*val;
                Eptr->B=Eptr->B*val;
                Eptr->B1=Eptr->B1*val;
                Eptr->B2=Eptr->B2*val;
                val=KVp/KV;
                Eptr->Imax=Eptr->Imax/val;
                Eptr->Tap=Eptr->Tap*val;
                Eptr->Taps=Eptr->Taps*val;
                Eptr->Tmax=Eptr->Tmax*val;
                Eptr->Tmin=Eptr->Tmin*val;
              } else {
                val=KV/KVp;
                Eptr->Tap=Eptr->Tap*val;
                Eptr->Taps=Eptr->Taps*val;
                Eptr->Tmax=Eptr->Tmax*val;
                Eptr->Tmin=Eptr->Tmin*val;
              }
            }
          }
          GetStr(Line,50,1,1,Code);