Source Code for Module pypower.opf_setup

  1  # Copyright (C) 1996-2011 Power System Engineering Research Center (PSERC) 
  2  # Copyright (C) 2011 Richard Lincoln 
  3  # 
  4  # PYPOWER is free software: you can redistribute it and/or modify 
  5  # it under the terms of the GNU General Public License as published 
  6  # by the Free Software Foundation, either version 3 of the License, 
  7  # or (at your option) any later version. 
  8  # 
  9  # PYPOWER is distributed in the hope that it will be useful, 
 10  # but WITHOUT ANY WARRANTY; without even the implied warranty of 
 11  # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 
 12  # GNU General Public License for more details. 
 13  # 
 14  # You should have received a copy of the GNU General Public License 
 15  # along with PYPOWER. If not, see <http://www.gnu.org/licenses/>. 
 16   
 17  """Constructs an OPF model object from a PYPOWER case dict. 
 18  """ 
 19   
 20  from sys import stdout, stderr 
 21   
 22  from numpy import array, any, delete, unique, arange, nonzero, pi, \ 
 23      r_, ones, Inf 
 24  from numpy import flatnonzero as find 
 25   
 26  from scipy.sparse import hstack, csr_matrix as sparse 
 27   
 28  from pypower.pqcost import pqcost 
 29  from pypower.opf_args import opf_args 
 30  from pypower.makeBdc import makeBdc 
 31  from pypower.makeAvl import makeAvl 
 32  from pypower.makeApq import makeApq 
 33  from pypower.makeAang import makeAang 
 34  from pypower.makeAy import makeAy 
 35  from pypower.opf_model import opf_model 
 36  from pypower.run_userfcn import run_userfcn 
 37   
 38  from pypower.idx_cost import MODEL, NCOST, PW_LINEAR, COST, POLYNOMIAL 
 39  from pypower.idx_bus import BUS_TYPE, REF, VA, VM, PD, GS, VMAX, VMIN 
 40  from pypower.idx_gen import GEN_BUS, VG, PG, QG, PMAX, PMIN, QMAX, QMIN 
 41  from pypower.idx_brch import RATE_A 
 42   
 43   
 44 -def opf_setup(ppc, ppopt): 
 45      """Constructs an OPF model object from a PYPOWER case dict. 
 46   
 47      Assumes that ppc is a PYPOWER case dict with internal indexing, 
 48      all equipment in-service, etc. 
 49   
 50      @see: L{opf}, L{ext2int}, L{opf_execute} 
 51   
 52      @author: Ray Zimmerman (PSERC Cornell) 
 53      @author: Carlos E. Murillo-Sanchez (PSERC Cornell & Universidad 
 54      Autonoma de Manizales) 
 55      @author: Richard Lincoln 
 56      """ 
 57      ## options 
 58      dc  = ppopt['PF_DC']        ## 1 = DC OPF, 0 = AC OPF 
 59      alg = ppopt['OPF_ALG'] 
 60      verbose = ppopt['VERBOSE'] 
 61   
 62      ## data dimensions 
 63      nb = ppc['bus'].shape[0]    ## number of buses 
 64      nl = ppc['branch'].shape[0] ## number of branches 
 65      ng = ppc['gen'].shape[0]    ## number of dispatchable injections 
 66      if 'A' in ppc: 
 67          nusr = ppc['A'].shape[0]    ## number of linear user constraints 
 68      else: 
 69          nusr = 0 
 70   
 71      if 'N' in ppc: 
 72          nw = ppc['N'].shape[0]      ## number of general cost vars, w 
 73      else: 
 74          nw = 0 
 75   
 76      if dc: 
 77          ## ignore reactive costs for DC 
 78          ppc['gencost'], _ = pqcost(ppc['gencost'], ng) 
 79   
 80          ## reduce A and/or N from AC dimensions to DC dimensions, if needed 
 81          if nusr or nw: 
 82              acc = r_[nb + arange(nb), 2 * nb + ng + arange(ng)]   ## Vm and Qg columns 
 83   
 84              if nusr and (ppc['A'].shape[1] >= 2*nb + 2*ng): 
 85                  ## make sure there aren't any constraints on Vm or Qg 
 86                  if ppc['A'][:, acc].nnz > 0: 
 87                      stderr.write('opf_setup: attempting to solve DC OPF with user constraints on Vm or Qg\n') 
 88   
 89                  # FIXME: delete sparse matrix columns 
 90                  bcc = delete(arange(ppc['A'].shape[1]), acc) 
 91                  ppc['A'] = ppc['A'].tolil()[:, bcc].tocsr()           ## delete Vm and Qg columns 
 92   
 93              if nw and (ppc['N'].shape[1] >= 2*nb + 2*ng): 
 94                  ## make sure there aren't any costs on Vm or Qg 
 95                  if ppc['N'][:, acc].nnz > 0: 
 96                      ii, _ = nonzero(ppc['N'][:, acc]) 
 97                      _, ii = unique(ii, return_index=True)    ## indices of w with potential non-zero cost terms from Vm or Qg 
 98                      if any(ppc['Cw'][ii]) | ( ('H' in ppc) & (len(ppc['H']) > 0) & 
 99                              any(any(ppc['H'][:, ii])) ): 
100                          stderr.write('opf_setup: attempting to solve DC OPF with user costs on Vm or Qg\n') 
101   
102                  # FIXME: delete sparse matrix columns 
103                  bcc = delete(arange(ppc['N'].shape[1]), acc) 
104                  ppc['N'] = ppc['N'].tolil()[:, bcc].tocsr()               ## delete Vm and Qg columns 
105   
106      ## convert single-block piecewise-linear costs into linear polynomial cost 
107      pwl1 = find((ppc['gencost'][:, MODEL] == PW_LINEAR) & (ppc['gencost'][:, NCOST] == 2)) 
108      # p1 = array([]) 
109      if len(pwl1) > 0: 
110          x0 = ppc['gencost'][pwl1, COST] 
111          y0 = ppc['gencost'][pwl1, COST + 1] 
112          x1 = ppc['gencost'][pwl1, COST + 2] 
113          y1 = ppc['gencost'][pwl1, COST + 3] 
114          m = (y1 - y0) / (x1 - x0) 
115          b = y0 - m * x0 
116          ppc['gencost'][pwl1, MODEL] = POLYNOMIAL 
117          ppc['gencost'][pwl1, NCOST] = 2 
118          ppc['gencost'][pwl1, COST:COST + 2] = r_[m, b] 
119   
120      ## create (read-only) copies of individual fields for convenience 
121      baseMVA, bus, gen, branch, gencost, _, lbu, ubu, ppopt, \ 
122              _, fparm, H, Cw, z0, zl, zu, userfcn, _ = opf_args(ppc, ppopt) 
123   
124      ## warn if there is more than one reference bus 
125      refs = find(bus[:, BUS_TYPE] == REF) 
126      if len(refs) > 1 and verbose > 0: 
127          errstr = '\nopf_setup: Warning: Multiple reference buses.\n' + \ 
128              '           For a system with islands, a reference bus in each island\n' + \ 
129              '           may help convergence, but in a fully connected system such\n' + \ 
130              '           a situation is probably not reasonable.\n\n' 
131          stdout.write(errstr) 
132   
133      ## set up initial variables and bounds 
134      gbus = gen[:, GEN_BUS].astype(int) 
135      Va   = bus[:, VA] * (pi / 180.0) 
136      Vm   = bus[:, VM].copy() 
137      Vm[gbus] = gen[:, VG]   ## buses with gens, init Vm from gen data 
138      Pg   = gen[:, PG] / baseMVA 
139      Qg   = gen[:, QG] / baseMVA 
140      Pmin = gen[:, PMIN] / baseMVA 
141      Pmax = gen[:, PMAX] / baseMVA 
142      Qmin = gen[:, QMIN] / baseMVA 
143      Qmax = gen[:, QMAX] / baseMVA 
144   
145      if dc:               ## DC model 
146          ## more problem dimensions 
147          nv    = 0            ## number of voltage magnitude vars 
148          nq    = 0            ## number of Qg vars 
149          q1    = array([])    ## index of 1st Qg column in Ay 
150   
151          ## power mismatch constraints 
152          B, Bf, Pbusinj, Pfinj = makeBdc(baseMVA, bus, branch) 
153          neg_Cg = sparse((-ones(ng), (gen[:, GEN_BUS], arange(ng))), (nb, ng))   ## Pbus w.r.t. Pg 
154          Amis = hstack([B, neg_Cg], 'csr') 
155          bmis = -(bus[:, PD] + bus[:, GS]) / baseMVA - Pbusinj 
156   
157          ## branch flow constraints 
158          il = find((branch[:, RATE_A] != 0) & (branch[:, RATE_A] < 1e10)) 
159          nl2 = len(il)         ## number of constrained lines 
160          lpf = -Inf * ones(nl2) 
161          upf = branch[il, RATE_A] / baseMVA - Pfinj[il] 
162          upt = branch[il, RATE_A] / baseMVA + Pfinj[il] 
163   
164          user_vars = ['Va', 'Pg'] 
165          ycon_vars = ['Pg', 'y'] 
166      else:                ## AC model 
167          ## more problem dimensions 
168          nv    = nb           ## number of voltage magnitude vars 
169          nq    = ng           ## number of Qg vars 
170          q1    = ng           ## index of 1st Qg column in Ay 
171   
172          ## dispatchable load, constant power factor constraints 
173          Avl, lvl, uvl, _  = makeAvl(baseMVA, gen) 
174   
175          ## generator PQ capability curve constraints 
176          Apqh, ubpqh, Apql, ubpql, Apqdata = makeApq(baseMVA, gen) 
177   
178          user_vars = ['Va', 'Vm', 'Pg', 'Qg'] 
179          ycon_vars = ['Pg', 'Qg', 'y'] 
180   
181      ## voltage angle reference constraints 
182      Vau = Inf * ones(nb) 
183      Val = -Vau 
184      Vau[refs] = Va[refs] 
185      Val[refs] = Va[refs] 
186   
187      ## branch voltage angle difference limits 
188      Aang, lang, uang, iang  = makeAang(baseMVA, branch, nb, ppopt) 
189   
190      ## basin constraints for piece-wise linear gen cost variables 
191      if alg == 545 or alg == 550:     ## SC-PDIPM or TRALM, no CCV cost vars 
192          ny = 0 
193          Ay = None 
194          by = array([]) 
195      else: 
196          ipwl = find(gencost[:, MODEL] == PW_LINEAR)  ## piece-wise linear costs 
197          ny = ipwl.shape[0]   ## number of piece-wise linear cost vars 
198          Ay, by = makeAy(baseMVA, ng, gencost, 1, q1, 1+ng+nq) 
199   
200      if any((gencost[:, MODEL] != POLYNOMIAL) & (gencost[:, MODEL] != PW_LINEAR)): 
201          stderr.write('opf_setup: some generator cost rows have invalid MODEL value\n') 
202   
203      ## more problem dimensions 
204      nx = nb+nv + ng+nq;  ## number of standard OPF control variables 
205      if nusr: 
206          nz = ppc['A'].shape[1] - nx  ## number of user z variables 
207          if nz < 0: 
208              stderr.write('opf_setup: user supplied A matrix must have at least %d columns.\n' % nx) 
209      else: 
210          nz = 0               ## number of user z variables 
211          if nw:               ## still need to check number of columns of N 
212              if ppc['N'].shape[1] != nx: 
213                  stderr.write('opf_setup: user supplied N matrix must have %d columns.\n' % nx) 
214   
215      ## construct OPF model object 
216      om = opf_model(ppc) 
217      if len(pwl1) > 0: 
218          om.userdata('pwl1', pwl1) 
219   
220      if dc: 
221          om.userdata('Bf', Bf) 
222          om.userdata('Pfinj', Pfinj) 
223          om.userdata('iang', iang) 
224          om.add_vars('Va', nb, Va, Val, Vau) 
225          om.add_vars('Pg', ng, Pg, Pmin, Pmax) 
226          om.add_constraints('Pmis', Amis, bmis, bmis, ['Va', 'Pg']) ## nb 
227          om.add_constraints('Pf',  Bf[il, :], lpf, upf, ['Va'])     ## nl 
228          om.add_constraints('Pt', -Bf[il, :], lpf, upt, ['Va'])     ## nl 
229          om.add_constraints('ang', Aang, lang, uang, ['Va'])        ## nang 
230      else: 
231          om.userdata('Apqdata', Apqdata) 
232          om.userdata('iang', iang) 
233          om.add_vars('Va', nb, Va, Val, Vau) 
234          om.add_vars('Vm', nb, Vm, bus[:, VMIN], bus[:, VMAX]) 
235          om.add_vars('Pg', ng, Pg, Pmin, Pmax) 
236          om.add_vars('Qg', ng, Qg, Qmin, Qmax) 
237          om.add_constraints('Pmis', nb, 'nonlinear') 
238          om.add_constraints('Qmis', nb, 'nonlinear') 
239          om.add_constraints('Sf', nl, 'nonlinear') 
240          om.add_constraints('St', nl, 'nonlinear') 
241          om.add_constraints('PQh', Apqh, array([]), ubpqh, ['Pg', 'Qg'])   ## npqh 
242          om.add_constraints('PQl', Apql, array([]), ubpql, ['Pg', 'Qg'])   ## npql 
243          om.add_constraints('vl',  Avl, lvl, uvl,   ['Pg', 'Qg'])   ## nvl 
244          om.add_constraints('ang', Aang, lang, uang, ['Va'])        ## nang 
245   
246      ## y vars, constraints for piece-wise linear gen costs 
247      if ny > 0: 
248          om.add_vars('y', ny) 
249          om.add_constraints('ycon', Ay, array([]), by, ycon_vars)          ## ncony 
250   
251      ## add user vars, constraints and costs (as specified via A, ..., N, ...) 
252      if nz > 0: 
253          om.add_vars('z', nz, z0, zl, zu) 
254          user_vars.append('z') 
255   
256      if nusr: 
257          om.add_constraints('usr', ppc['A'], lbu, ubu, user_vars)      ## nusr 
258   
259      if nw: 
260          user_cost = {} 
261          user_cost['N'] = ppc['N'] 
262          user_cost['Cw'] = Cw 
263          if len(fparm) > 0: 
264              user_cost['dd'] = fparm[:, 0] 
265              user_cost['rh'] = fparm[:, 1] 
266              user_cost['kk'] = fparm[:, 2] 
267              user_cost['mm'] = fparm[:, 3] 
268   
269  #        if len(H) > 0: 
270          user_cost['H'] = H 
271   
272          om.add_costs('usr', user_cost, user_vars) 
273   
274      ## execute userfcn callbacks for 'formulation' stage 
275      run_userfcn(userfcn, 'formulation', om) 
276   
277      return om 
278