xref: /petsc/src/snes/tutorials/network/power/PFReadData.c (revision fffbe07892e3829f28d4ec4d935109d0474c54e5)
1 #include "petscmat.h"
2 #include "power.h"
3 #include <string.h>
4 #include <ctype.h>
5 
6 PetscErrorCode PFReadMatPowerData(PFDATA *pf,char *filename)
7 {
8   FILE           *fp;
9   PetscErrorCode ierr;
10   VERTEX_Power   Bus;
11   LOAD           Load;
12   GEN            Gen;
13   EDGE_Power     Branch;
14   PetscInt       line_counter=0;
15   PetscInt       bus_start_line=-1,bus_end_line=-1; /* xx_end_line points to the next line after the record ends */
16   PetscInt       gen_start_line=-1,gen_end_line=-1;
17   PetscInt       br_start_line=-1,br_end_line=-1;
18   char           line[MAXLINE];
19   PetscInt       loadi=0,geni=0,bri=0,busi=0,i,j;
20   int            extbusnum,bustype_i;
21   double         Pd,Qd;
22   PetscInt       maxbusnum=-1,intbusnum,*busext2intmap,genj,loadj;
23   GEN            newgen;
24   LOAD           newload;
25 
26   PetscFunctionBegin;
27   fp = fopen(filename,"r");
28   /* Check for valid file */
29   PetscCheckFalse(!fp,PETSC_COMM_SELF,PETSC_ERR_FILE_OPEN,"Can't open Matpower data file %s",filename);
30   pf->nload=0;
31   while (fgets(line,MAXLINE,fp)) {
32     if (strstr(line,"mpc.bus = ["))    bus_start_line = line_counter+1; /* Bus data starts from next line */
33     if (strstr(line,"mpc.gen") && gen_start_line == -1)    gen_start_line = line_counter+1; /* Generator data starts from next line */
34     if (strstr(line,"mpc.branch")) br_start_line = line_counter+1; /* Branch data starts from next line */
35     if (strstr(line,"];")) {
36       if (bus_start_line != -1 && bus_end_line == -1) bus_end_line = line_counter;
37       if (gen_start_line != -1 && gen_end_line == -1) gen_end_line = line_counter;
38       if (br_start_line  != -1 && br_end_line == -1) br_end_line = line_counter;
39     }
40 
41     /* Count the number of pq loads */
42     if (bus_start_line != -1 && line_counter >= bus_start_line && bus_end_line == -1) {
43       sscanf(line,"%d %d %lf %lf",&extbusnum,&bustype_i,&Pd,&Qd);
44       if (!((Pd == 0.0) && (Qd == 0.0))) pf->nload++;
45       if (extbusnum > maxbusnum) maxbusnum = extbusnum;
46     }
47     line_counter++;
48   }
49   fclose(fp);
50 
51   pf->nbus    = bus_end_line - bus_start_line;
52   pf->ngen    = gen_end_line - gen_start_line;
53   pf->nbranch = br_end_line  - br_start_line;
54 
55   ierr = PetscCalloc1(pf->nbus,&pf->bus);CHKERRQ(ierr);
56   ierr = PetscCalloc1(pf->ngen,&pf->gen);CHKERRQ(ierr);
57   ierr = PetscCalloc1(pf->nload,&pf->load);CHKERRQ(ierr);
58   ierr = PetscCalloc1(pf->nbranch,&pf->branch);CHKERRQ(ierr);
59   Bus = pf->bus; Gen = pf->gen; Load = pf->load; Branch = pf->branch;
60 
61   /* Setting pf->sbase to 100 */
62   pf->sbase = 100.0;
63 
64   ierr = PetscMalloc1(maxbusnum+1,&busext2intmap);CHKERRQ(ierr);
65   for (i=0; i < maxbusnum+1; i++) busext2intmap[i] = -1;
66 
67   fp = fopen(filename,"r");
68   /* Reading data */
69   for (i=0;i<line_counter;i++) {
70     PetscCheckFalse(!fgets(line,MAXLINE,fp),PETSC_COMM_SELF,PETSC_ERR_SUP,"File is incorrectly formatted");
71 
72     if ((i >= bus_start_line) && (i < bus_end_line)) {
73       double gl,bl,vm,va,basekV;
74       int    bus_i,ide,area;
75       /* Bus data */
76       sscanf(line,"%d %d %lf %lf %lf %lf %d %lf %lf %lf",&bus_i,&ide,&Pd,&Qd,&gl,&bl,&area,&vm,&va,&basekV);
77       Bus[busi].bus_i = bus_i; Bus[busi].ide = ide; Bus[busi].area = area;
78       Bus[busi].gl = gl; Bus[busi].bl = bl;
79       Bus[busi].vm = vm; Bus[busi].va = va; Bus[busi].basekV = basekV;
80       Bus[busi].internal_i = busi;
81       busext2intmap[Bus[busi].bus_i] = busi;
82 
83       if (!((Pd == 0.0) && (Qd == 0.0))) {
84         Load[loadi].bus_i = Bus[busi].bus_i;
85         Load[loadi].status = 1;
86         Load[loadi].pl = Pd;
87         Load[loadi].ql = Qd;
88         Load[loadi].area = Bus[busi].area;
89         Load[loadi].internal_i = busi;
90         Bus[busi].lidx[Bus[busi].nload++] = loadi;
91         PetscCheckFalse(Bus[busi].nload > NLOAD_AT_BUS_MAX,PETSC_COMM_SELF,PETSC_ERR_SUP,"Exceeded maximum number of loads allowed at bus");
92         loadi++;
93       }
94       busi++;
95     }
96 
97     /* Read generator data */
98     if (i >= gen_start_line && i < gen_end_line) {
99       double pg,qg,qt,qb,vs,mbase,pt,pb;
100       int    bus_i,status;
101       sscanf(line,"%d %lf %lf %lf %lf %lf %lf %d %lf %lf",&bus_i,&pg,&qg,&qt,&qb,&vs,&mbase,&status,&pt,&pb);
102       Gen[geni].bus_i = bus_i; Gen[geni].status = status;
103       Gen[geni].pg = pg; Gen[geni].qg = qg; Gen[geni].qt = qt; Gen[geni].qb = qb;
104       Gen[geni].vs = vs; Gen[geni].mbase = mbase; Gen[geni].pt = pt; Gen[geni].pb = pb;
105 
106       intbusnum = busext2intmap[Gen[geni].bus_i];
107       Gen[geni].internal_i = intbusnum;
108       Bus[intbusnum].gidx[Bus[intbusnum].ngen++] = geni;
109 
110       Bus[intbusnum].vm = Gen[geni].vs;
111 
112       PetscCheckFalse(Bus[intbusnum].ngen > NGEN_AT_BUS_MAX,PETSC_COMM_SELF,PETSC_ERR_SUP,"Exceeded maximum number of generators allowed at bus");
113       geni++;
114     }
115 
116     if (i >= br_start_line && i < br_end_line) {
117       PetscScalar R,X,Bc,B,G,Zm,tap,shift,tap2,tapr,tapi;
118       double      r,x,b,rateA,rateB,rateC,tapratio,phaseshift;
119       int         fbus,tbus,status;
120       sscanf(line,"%d %d %lf %lf %lf %lf %lf %lf %lf %lf %d",&fbus,&tbus,&r,&x,&b,&rateA,&rateB,&rateC,&tapratio,&phaseshift,&status);
121       Branch[bri].fbus = fbus; Branch[bri].tbus = tbus; Branch[bri].status = status;
122       Branch[bri].r = r; Branch[bri].x = x; Branch[bri].b = b;
123       Branch[bri].rateA = rateA; Branch[bri].rateB = rateB; Branch[bri].rateC = rateC;
124       Branch[bri].tapratio = tapratio; Branch[bri].phaseshift = phaseshift;
125 
126       if (Branch[bri].tapratio == 0.0) Branch[bri].tapratio = 1.0;
127       Branch[bri].phaseshift *= PETSC_PI/180.0;
128 
129       intbusnum = busext2intmap[Branch[bri].fbus];
130       Branch[bri].internal_i = intbusnum;
131 
132       intbusnum = busext2intmap[Branch[bri].tbus];
133       Branch[bri].internal_j = intbusnum;
134 
135       /* Compute self and transfer admittances */
136       R = Branch[bri].r;
137       X = Branch[bri].x;
138       Bc = Branch[bri].b;
139 
140       Zm = R*R + X*X;
141       G  = R/Zm;
142       B  = -X/Zm;
143 
144       tap = Branch[bri].tapratio;
145       shift = Branch[bri].phaseshift;
146       tap2 = tap*tap;
147       tapr = tap*PetscCosScalar(shift);
148       tapi = tap*PetscSinScalar(shift);
149 
150       Branch[bri].yff[0] = G/tap2;
151       Branch[bri].yff[1] = (B+Bc/2.0)/tap2;
152 
153       Branch[bri].yft[0] = -(G*tapr - B*tapi)/tap2;
154       Branch[bri].yft[1] = -(B*tapr + G*tapi)/tap2;
155 
156       Branch[bri].ytf[0] = -(G*tapr + B*tapi)/tap2;
157       Branch[bri].ytf[1] = -(B*tapr - G*tapi)/tap2;
158 
159       Branch[bri].ytt[0] = G;
160       Branch[bri].ytt[1] = B+Bc/2.0;
161 
162       bri++;
163     }
164   }
165   fclose(fp);
166 
167   /* Reorder the generator data structure according to bus numbers */
168   genj=0; loadj=0;
169   ierr = PetscMalloc1(pf->ngen,&newgen);CHKERRQ(ierr);
170   ierr = PetscMalloc1(pf->nload,&newload);CHKERRQ(ierr);
171   for (i = 0; i < pf->nbus; i++) {
172     for (j = 0; j < pf->bus[i].ngen; j++) {
173       ierr = PetscMemcpy(&newgen[genj++],&pf->gen[pf->bus[i].gidx[j]],sizeof(struct _p_GEN));CHKERRQ(ierr);
174     }
175     for (j = 0; j < pf->bus[i].nload; j++) {
176       ierr = PetscMemcpy(&newload[loadj++],&pf->load[pf->bus[i].lidx[j]],sizeof(struct _p_LOAD));CHKERRQ(ierr);
177     }
178   }
179   ierr = PetscFree(pf->gen);CHKERRQ(ierr);
180   ierr = PetscFree(pf->load);CHKERRQ(ierr);
181   pf->gen = newgen;
182   pf->load = newload;
183 
184   ierr = PetscFree(busext2intmap);CHKERRQ(ierr);
185   PetscFunctionReturn(0);
186 }
187