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