update.c

电力系统分析计算 学习调试程序 UNIX / LINUX / CYGWIN 系统使用

/* Update AC/DC variables. */

#include <stdlib.h>
//#ifndef WINDOWS
//#include <stdio.h>
//#else
#include "pfwstdio.h"
//#endif
#include <math.h>
#include <string.h>
#include "constant.h"
#include "param.h"
#include "sparse.h"
#include "pflow.h"

#ifdef ANSIPROTO
void UpdateACvar(VALUETYPE cons,INDEX j,BOOLEAN Limits,BOOLEAN Recover);
void UpdateDCvar(VALUETYPE cons,INDEX j,BOOLEAN Limits);
#else
void UpdateACvar();
void UpdateDCvar();
#endif

/* ------- Global Variables ------ */
extern Data *dataPtr;
extern SparseMatrix *Jac;
extern INDEX MaxIter,Nac,NacEl,NregPQ,NregV,Ndc,Nslack,Nvolt,Narea,NacVar,Bl,NZbuses,NXbuses;
extern INDEX *ACvar;
extern VALUETYPE *dx,*dF,tol,Tol,Sn,lambda,*x0;
extern VALUETYPE K1,K2,MaxdFi,alpha;
extern IntegerVector *NewRow,*OldRow,*NewCol,*OldCol,*RowPartition,*ColPartition;
extern IntegerVector *RowPer,*ColPer;
extern BOOLEAN Acont,PQcont,QRcont,Rcont,Xcont,
               PQlim,Tlim,Qlim,Vlim,Elim,Ilim,Zlim,Xlim,
               flagH,flagPoC,flagL,flagR,flagBS,flagPgMax,flagSmax;


/* -----------------UpdateACvar ---------------------------- */
#ifdef ANSIPROTO
void UpdateACvar(VALUETYPE cons,INDEX j,BOOLEAN Limits,BOOLEAN Recover)
#else
void UpdateACvar(cons,j,Limits,Recover)
VALUETYPE cons;
INDEX j;
BOOLEAN Limits,Recover;
#endif
{
  ACbusData *ACptr,*ACptrp,*BEptr;
  AClist *ALptr;
  ElementList *ELptr;
  ElementData *Eptr;
  VALUETYPE Pg,DPg,Qm,Pmax,PgMax;
  INDEX k;

  for(ALptr=dataPtr->KGbus;ALptr!=NULL;ALptr=ALptr->Next) {
    ACptr=ALptr->AC;
    k=ACvar[ACptr->N];
    if (strpbrk(ACptr->Type,"S")) {
      ACptr->Kg=ACptr->Kg+cons*dx[k];
      BEptr=ACptr;
    }
    else if (Acont) ACptr->Kg=ACptr->Kg+cons*dx[k+2];
  }
  for (ACptr=dataPtr->ACbus; ACptr!=NULL; ACptr=ACptr->Next){
    if (ACptr->Area!=NULL) BEptr=ACptr->Area->Slack;
    DPg=ACptr->DPg;
    if (flagBS) {
      if (ACptr!=BEptr) Pg=ACptr->Pg+lambda*DPg;
      else              Pg=ACptr->Pg+BEptr->Kg;
    }
    else  Pg=ACptr->Pg+BEptr->Kg*DPg;
    Pmax=ACptr->Smax*ACptr->Smax-ACptr->Qg*ACptr->Qg;
    if (Pmax>0) Pmax=sqrt(Pmax);
    else        Pmax=99999999.;
    if (!flagSmax && Pmax<ACptr->PgMax) PgMax=Pmax;
    else if (!flagPgMax)                PgMax=ACptr->PgMax;
    else                                PgMax=99999999.;
    if (Pg>PgMax) {
      Pg=PgMax;
      DPg=0;
      if (ACptr->flagPgMax!=1) ACptr->flagPgMax++;
    } else ACptr->flagPgMax=0;
    ACptr->PG=Pg;
    ACptr->DPG=DPg;
    ACptr->Pmax=PgMax;
    Qm=ACptr->Smax*ACptr->Smax-Pg*Pg;
    if (!flagSmax && Qm>0) {
      Qm=sqrt(Qm);
      if (ACptr->Qmax<Qm) ACptr->Max=ACptr->Qmax;
      else ACptr->Max=Qm;
      if (ACptr->Qmin>-Qm) ACptr->Min=ACptr->Qmin;
      else ACptr->Min=-Qm;
    } else {
      ACptr->Max=ACptr->Qmax;
      ACptr->Min=ACptr->Qmin;
    }
    if (!strcmp(ACptr->Type,"B") || !strcmp(ACptr->Type,"BA")) {
      k=ACvar[ACptr->N];
      ACptr->Ang=ACptr->Ang+cons*dx[k];
      if (j==0) ACptr->val=ACptr->V;
      ACptr->val=ACptr->val+cons*dx[k+1];
      if (ACptr->val<=0) ACptr->V=0.00001;
      else ACptr->V=ACptr->val;
    }
    else if (strpbrk(ACptr->Type,"V")) {
      k=ACvar[ACptr->N];
      ACptr->Ang=ACptr->Ang+cons*dx[k];
      if (!strpbrk(ACptr->cont,"V")) {
        if (j==0) ACptr->val=ACptr->V;
        ACptr->val=ACptr->val+cons*dx[k+1];
        if (Vlim && ACptr->val>ACptr->Vmax && Limits)      ACptr->V=ACptr->Vmax;
        else if (Vlim && ACptr->val<ACptr->Vmin && Limits) ACptr->V=ACptr->Vmin;
        else if (ACptr->val<=0)                            ACptr->V=0.00001;
        else                                               ACptr->V=ACptr->val;
      } else {
        if (j==0) ACptr->val=ACptr->Qg;
        ACptr->val=ACptr->val+cons*dx[k+1];
        if (Recover && ACptr->V>=ACptr->Vmax && ACptr->val>ACptr->VCont && Limits)      ACptr->Qg=ACptr->VCont;
        else if (Recover && ACptr->V<=ACptr->Vmin && ACptr->val<ACptr->VCont && Limits) ACptr->Qg=ACptr->VCont;
        else                                                                            ACptr->Qg=ACptr->val;
      }
    }
    else if (strpbrk(ACptr->Type,"L")){
      k=ACvar[ACptr->N];
      ACptr->Ang=ACptr->Ang+cons*dx[k];
      lambda=lambda+cons*dx[k+1];
      if (flagH) {
        if (j==0) ACptr->val=ACptr->V;
        ACptr->val=ACptr->val+cons*dx[Jac->n1];
        if (ACptr->val<=0) ACptr->V=0.00001;
        else ACptr->V=ACptr->val;
      }
    }
    else if (QRcont && strpbrk(ACptr->Type,"G")){
      k=ACvar[ACptr->N];
      if (!strpbrk(ACptr->Type,"S")) ACptr->Ang=ACptr->Ang+cons*dx[k];
      if (j==0) ACptr->val=ACptr->V;
      ACptr->val=ACptr->val+cons*dx[k+1];
      if (ACptr->val<=0) ACptr->V=0.00001;
      else ACptr->V=ACptr->val;
    }
    else if (strpbrk(ACptr->Type,"Q,S") || (!QRcont && strpbrk(ACptr->Type,"G"))){
      k=ACvar[ACptr->N];
      if (!strpbrk(ACptr->Type,"S")) ACptr->Ang=ACptr->Ang+cons*dx[k];
      if (strpbrk(ACptr->cont,"V")) {
        if (j==0) ACptr->val=ACptr->Qg;
        ACptr->val=ACptr->val+cons*dx[k+1];
        if (Qlim && ACptr->val>ACptr->Max && Limits)      ACptr->Qg=ACptr->Max;
        else if (Qlim && ACptr->val<ACptr->Min && Limits) ACptr->Qg=ACptr->Min;
        else                                              ACptr->Qg=ACptr->val;
      }
      else {
        if (j==0) ACptr->val=ACptr->V;
        ACptr->val=ACptr->val+cons*dx[k+1];
        if (ACptr->val<=0) ACptr->V=0.00001;
        else if (strpbrk(ACptr->cont,"Q")) {
          if (Recover && ACptr->Qg>=ACptr->Max && ACptr->val>ACptr->VCont && Limits)      ACptr->V=ACptr->VCont;
          else if (Recover && ACptr->Qg<=ACptr->Min && ACptr->val<ACptr->VCont && Limits) ACptr->V=ACptr->VCont;
          else                                                                            ACptr->V=ACptr->val;
        }
        else if (strpbrk(ACptr->cont,"E")) {
          if (Recover && ACptr->Gen->Eq>=ACptr->Gen->EqMax && ACptr->val>ACptr->VCont && Limits)      ACptr->V=ACptr->VCont;
          else if (Recover && ACptr->Gen->Eq<=ACptr->Gen->EqMin && ACptr->val<ACptr->VCont && Limits) ACptr->V=ACptr->VCont;
          else                                                                                        ACptr->V=ACptr->val;
        }
        else if (strpbrk(ACptr->cont,"I")) {
          if (Recover && ACptr->Gen->Ia>=ACptr->Gen->IaMax && ACptr->val>ACptr->VCont && Limits) ACptr->V=ACptr->VCont;
          else                                                                                   ACptr->V=ACptr->val;
        }
        else ACptr->V=ACptr->val;
      }
    }
    else if (strpbrk(ACptr->Type,"X")) {
      k=ACvar[ACptr->N];
      ACptr->Ang=ACptr->Ang+cons*dx[k];
      if (j==0) ACptr->val=ACptr->V;
      ACptr->val=ACptr->val+cons*dx[k+1];
      if (ACptr->val<=0)                                   ACptr->V=0.00001;
      else if (strpbrk(ACptr->cont,"X")) {
        if (Xlim && ACptr->val>ACptr->Vmax && Limits)      ACptr->V=ACptr->Vmax;
        else if (Xlim && ACptr->val<ACptr->Vmin && Limits) ACptr->V=ACptr->Vmin;
        else                                               ACptr->V=ACptr->val;
      }
      else if (strpbrk(ACptr->cont,"M")) {
        if (Recover && ACptr->val<ACptr->Vmax && Limits)   ACptr->V=ACptr->Vmax;
        else                                               ACptr->V=ACptr->val;
      }
      else {
        if (Recover && ACptr->val>ACptr->Vmin && Limits)   ACptr->V=ACptr->Vmin;
        else                                               ACptr->V=ACptr->val;
      }
    }
    else if (strpbrk(ACptr->Type,"Z")){
      k=ACvar[ACptr->N];
      ACptr->Ang=ACptr->Ang+cons*dx[k];
      if (strpbrk(ACptr->cont,"V")) {
        if (j==0) ACptr->val=ACptr->Qg;
        ACptr->val=ACptr->val+cons*dx[k+1];
        if (Zlim && ACptr->val>ACptr->Qmax && Limits)      ACptr->Qg=ACptr->Qmax;
        else if (Zlim && ACptr->val<ACptr->Qmin && Limits) ACptr->Qg=ACptr->Qmin;
        else                                               ACptr->Qg=ACptr->val;
      }
      else {
        if (j==0) ACptr->val=ACptr->V;
        ACptr->val=ACptr->val+cons*dx[k+1];
        ACptr->Qg=ACptr->val*ACptr->val*ACptr->Bz;
        if (ACptr->val<=0)                                                               ACptr->V=0.00001;
        else if (Recover && ACptr->Qg>=ACptr->Qmax && ACptr->val>ACptr->VCont && Limits) ACptr->V=ACptr->VCont;
        else if (Recover && ACptr->Qg<=ACptr->Qmin && ACptr->val<ACptr->VCont && Limits) ACptr->V=ACptr->VCont;
        else                                                                             ACptr->V=ACptr->val;
        ACptr->Qg=ACptr->V*ACptr->V*ACptr->Bz;
      }
    }
    else if (strpbrk(ACptr->Type,"C")) {
      k=ACvar[ACptr->N];
      ACptr->Ang=ACptr->Ang+cons*dx[k];
      if (!QRcont) {
        if (j==0) ACptr->val=ACptr->V;
        ACptr->val=ACptr->val+cons*dx[k+1];
        if (ACptr->val<=0) ACptr->V=0.00001;
        else ACptr->V=ACptr->val;
      }
      else if (ACptr->Kbg<1) {