e2nIEE / pandapower
1
# -*- coding: utf-8 -*-
2

3
# Copyright 1996-2015 PSERC. All rights reserved.
4
# Use of this source code is governed by a BSD-style
5
# license that can be found in the LICENSE file.
6

7
# Copyright (c) 2016-2021 by University of Kassel and Fraunhofer Institute for Energy Economics
8
# and Energy System Technology (IEE), Kassel. All rights reserved.
9

10

11 1
"""Constructs an OPF model object from a PYPOWER case dict.
12
"""
13

14 1
from sys import stdout, stderr
15

16 1
from numpy import array, any, delete, unique, arange, nonzero, pi, r_, ones, Inf, flatnonzero as find
17 1
from scipy.sparse import hstack, csr_matrix as sparse
18 1
from pandapower.pypower.idx_brch import RATE_A
19 1
from pandapower.pypower.idx_bus import BUS_TYPE, REF, VA, VM, PD, GS, VMAX, VMIN
20 1
from pandapower.pypower.idx_cost import MODEL, NCOST, PW_LINEAR, COST, POLYNOMIAL
21 1
from pandapower.pypower.idx_gen import GEN_BUS, VG, PG, QG, PMAX, PMIN, QMAX, QMIN
22 1
from pandapower.pypower.makeAang import makeAang
23 1
from pandapower.pypower.makeApq import makeApq
24 1
from pandapower.pypower.makeAvl import makeAvl
25 1
from pandapower.pypower.makeAy import makeAy
26 1
from pandapower.pypower.makeBdc import makeBdc
27 1
from pandapower.pypower.opf_args import opf_args
28 1
from pandapower.pypower.pqcost import pqcost
29 1
from pandapower.pypower.run_userfcn import run_userfcn
30 1
from pandapower.pypower.opf_model import opf_model
31

32

33 1
def opf_setup(ppc, ppopt):
34
    """Constructs an OPF model object from a PYPOWER case dict.
35

36
    Assumes that ppc is a PYPOWER case dict with internal indexing,
37
    all equipment in-service, etc.
38

39
    @see: L{opf}, L{ext2int}, L{opf_execute}
40

41
    @author: Ray Zimmerman (PSERC Cornell)
42
    @author: Carlos E. Murillo-Sanchez (PSERC Cornell & Universidad
43
    Autonoma de Manizales)
44
    @author: Richard Lincoln
45

46
    Modified by University of Kassel (Friederike Meier): Bugfix in line 110
47
    """
48
    ## options
49 1
    dc  = ppopt['PF_DC']        ## 1 = DC OPF, 0 = AC OPF
50 1
    alg = ppopt['OPF_ALG']
51 1
    verbose = ppopt['VERBOSE']
52

53
    ## data dimensions
54 1
    nb = ppc['bus'].shape[0]    ## number of buses
55 1
    nl = ppc['branch'].shape[0] ## number of branches
56 1
    ng = ppc['gen'].shape[0]    ## number of dispatchable injections
57 1
    if 'A' in ppc:
58 0
        nusr = ppc['A'].shape[0]    ## number of linear user constraints
59
    else:
60 1
        nusr = 0
61

62 1
    if 'N' in ppc:
63 0
        nw = ppc['N'].shape[0]      ## number of general cost vars, w
64
    else:
65 1
        nw = 0
66

67 1
    if dc:
68
        ## ignore reactive costs for DC
69 1
        ppc['gencost'], _ = pqcost(ppc['gencost'], ng)
70

71
        ## reduce A and/or N from AC dimensions to DC dimensions, if needed
72
        if nusr or nw: # pragma: no cover
73
            acc = r_[nb + arange(nb), 2 * nb + ng + arange(ng)]   ## Vm and Qg columns
74

75
            if nusr and (ppc['A'].shape[1] >= 2*nb + 2*ng):
76
                ## make sure there aren't any constraints on Vm or Qg
77
                if ppc['A'][:, acc].nnz > 0:
78
                    stderr.write('opf_setup: attempting to solve DC OPF with user constraints on Vm or Qg\n')
79

80
                # FIXME: delete sparse matrix columns
81
                bcc = delete(arange(ppc['A'].shape[1]), acc)
82
                ppc['A'] = ppc['A'].tolil()[:, bcc].tocsr()           ## delete Vm and Qg columns
83

84
            if nw and (ppc['N'].shape[1] >= 2*nb + 2*ng):
85
                ## make sure there aren't any costs on Vm or Qg
86
                if ppc['N'][:, acc].nnz > 0:
87
                    ii, _ = nonzero(ppc['N'][:, acc])
88
                    _, ii = unique(ii, return_index=True)    ## indices of w with potential non-zero cost terms from Vm or Qg
89
                    if any(ppc['Cw'][ii]) | ( ('H' in ppc) & (len(ppc['H']) > 0) &
90
                            any(any(ppc['H'][:, ii])) ):
91
                        stderr.write('opf_setup: attempting to solve DC OPF with user costs on Vm or Qg\n')
92

93
                # FIXME: delete sparse matrix columns
94
                bcc = delete(arange(ppc['N'].shape[1]), acc)
95
                ppc['N'] = ppc['N'].tolil()[:, bcc].tocsr()               ## delete Vm and Qg columns
96

97
#    ## convert single-block piecewise-linear costs into linear polynomial cost
98 1
    pwl1 = find((ppc['gencost'][:, MODEL] == PW_LINEAR) & (ppc['gencost'][:, NCOST] == 2))
99
    # p1 = array([])
100 1
    if len(pwl1) > 0:
101 1
        x0 = ppc['gencost'][pwl1, COST]
102 1
        y0 = ppc['gencost'][pwl1, COST + 1]
103 1
        x1 = ppc['gencost'][pwl1, COST + 2]
104 1
        y1 = ppc['gencost'][pwl1, COST + 3]
105 1
        m = (y1 - y0) / (x1 - x0)
106 1
        b = y0 - m * x0
107 1
        ppc['gencost'][pwl1, MODEL] = POLYNOMIAL
108 1
        ppc['gencost'][pwl1, NCOST] = 2
109 1
        ppc['gencost'][pwl1, COST:COST + 2] = r_['1',m.reshape(len(m),1), b.reshape(len(b),1)] # changed from ppc['gencost'][pwl1, COST:COST + 2] = r_[m, b] because we need to make sure, that m and b have the same shape, resulted in a value error due to shape mismatch before
110

111
    ## create (read-only) copies of individual fields for convenience
112 1
    baseMVA, bus, gen, branch, gencost, _, lbu, ubu, ppopt, \
113
            _, fparm, H, Cw, z0, zl, zu, userfcn, _ = opf_args(ppc, ppopt)
114

115
    ## warn if there is more than one reference bus
116 1
    refs = find(bus[:, BUS_TYPE] == REF)
117 1
    if len(refs) > 1 and verbose > 0:
118 0
        errstr = '\nopf_setup: Warning: Multiple reference buses.\n' + \
119
            '           For a system with islands, a reference bus in each island\n' + \
120
            '           may help convergence, but in a fully connected system such\n' + \
121
            '           a situation is probably not reasonable.\n\n'
122 0
        stdout.write(errstr)
123

124
    ## set up initial variables and bounds
125 1
    gbus = gen[:, GEN_BUS].astype(int)
126 1
    Va   = bus[:, VA] * (pi / 180.0)
127 1
    Vm   = bus[:, VM].copy()
128 1
    Vm[gbus] = gen[:, VG]   ## buses with gens, init Vm from gen data
129 1
    Pg   = gen[:, PG] / baseMVA
130 1
    Qg   = gen[:, QG] / baseMVA
131 1
    Pmin = gen[:, PMIN] / baseMVA
132 1
    Pmax = gen[:, PMAX] / baseMVA
133 1
    Qmin = gen[:, QMIN] / baseMVA
134 1
    Qmax = gen[:, QMAX] / baseMVA
135

136 1
    if dc:               ## DC model
137
        ## more problem dimensions
138 1
        nv    = 0            ## number of voltage magnitude vars
139 1
        nq    = 0            ## number of Qg vars
140 1
        q1    = array([])    ## index of 1st Qg column in Ay
141

142
        ## power mismatch constraints
143 1
        B, Bf, Pbusinj, Pfinj = makeBdc(bus, branch)
144 1
        neg_Cg = sparse((-ones(ng), (gen[:, GEN_BUS], arange(ng))), (nb, ng))   ## Pbus w.r.t. Pg
145 1
        Amis = hstack([B, neg_Cg], 'csr')
146 1
        bmis = -(bus[:, PD] + bus[:, GS]) / baseMVA - Pbusinj
147

148
        ## branch flow constraints
149 1
        il = find((branch[:, RATE_A] != 0) & (branch[:, RATE_A] < 1e10))
150 1
        nl2 = len(il)         ## number of constrained lines
151 1
        lpf = -Inf * ones(nl2)
152 1
        upf = branch[il, RATE_A] / baseMVA - Pfinj[il]
153 1
        upt = branch[il, RATE_A] / baseMVA + Pfinj[il]
154

155 1
        user_vars = ['Va', 'Pg']
156 1
        ycon_vars = ['Pg', 'y']
157
    else:                ## AC model
158
        ## more problem dimensions
159 1
        nv    = nb           ## number of voltage magnitude vars
160 1
        nq    = ng           ## number of Qg vars
161 1
        q1    = ng           ## index of 1st Qg column in Ay
162

163
        ## dispatchable load, constant power factor constraints
164 1
        Avl, lvl, uvl, _  = makeAvl(baseMVA, gen)
165

166
        ## generator PQ capability curve constraints
167 1
        Apqh, ubpqh, Apql, ubpql, Apqdata = makeApq(baseMVA, gen)
168

169 1
        user_vars = ['Va', 'Vm', 'Pg', 'Qg']
170 1
        ycon_vars = ['Pg', 'Qg', 'y']
171

172
    ## voltage angle reference constraints
173 1
    Vau = Inf * ones(nb)
174 1
    Val = -Vau
175 1
    Vau[refs] = Va[refs]
176 1
    Val[refs] = Va[refs]
177

178
    ## branch voltage angle difference limits
179 1
    Aang, lang, uang, iang  = makeAang(baseMVA, branch, nb, ppopt)
180

181
    ## basin constraints for piece-wise linear gen cost variables
182
    if alg == 545 or alg == 550:     ## SC-PDIPM or TRALM, no CCV cost vars # pragma: no cover
183
        ny = 0
184
        Ay = None
185
        by = array([])
186
    else:
187 1
        ipwl = find(gencost[:, MODEL] == PW_LINEAR)  ## piece-wise linear costs
188 1
        ny = ipwl.shape[0]   ## number of piece-wise linear cost vars
189 1
        Ay, by = makeAy(baseMVA, ng, gencost, 1, q1, 1+ng+nq)
190

191 1
    if any((gencost[:, MODEL] != POLYNOMIAL) & (gencost[:, MODEL] != PW_LINEAR)):
192 0
        stderr.write('opf_setup: some generator cost rows have invalid MODEL value\n')
193

194
    ## more problem dimensions
195 1
    nx = nb+nv + ng+nq  ## number of standard OPF control variables
196
    if nusr: # pragma: no cover
197
        nz = ppc['A'].shape[1] - nx  ## number of user z variables
198
        if nz < 0:
199
            stderr.write('opf_setup: user supplied A matrix must have at least %d columns.\n' % nx)
200
    else:
201 1
        nz = 0               ## number of user z variables
202 1
        if nw:               ## still need to check number of columns of N
203 0
            if ppc['N'].shape[1] != nx:
204 0
                stderr.write('opf_setup: user supplied N matrix must have %d columns.\n' % nx)
205

206
    ## construct OPF model object
207 1
    om = opf_model(ppc)
208 1
    if len(pwl1) > 0:
209 1
        om.userdata('pwl1', pwl1)
210

211 1
    if dc:
212 1
        om.userdata('Bf', Bf)
213 1
        om.userdata('Pfinj', Pfinj)
214 1
        om.userdata('iang', iang)
215 1
        om.add_vars('Va', nb, Va, Val, Vau)
216 1
        om.add_vars('Pg', ng, Pg, Pmin, Pmax)
217 1
        om.add_constraints('Pmis', Amis, bmis, bmis, ['Va', 'Pg']) ## nb
218 1
        om.add_constraints('Pf',  Bf[il, :], lpf, upf, ['Va'])     ## nl
219 1
        om.add_constraints('Pt', -Bf[il, :], lpf, upt, ['Va'])     ## nl
220 1
        om.add_constraints('ang', Aang, lang, uang, ['Va'])        ## nang
221
    else:
222 1
        om.userdata('Apqdata', Apqdata)
223 1
        om.userdata('iang', iang)
224 1
        om.add_vars('Va', nb, Va, Val, Vau)
225 1
        om.add_vars('Vm', nb, Vm, bus[:, VMIN], bus[:, VMAX])
226 1
        om.add_vars('Pg', ng, Pg, Pmin, Pmax)
227 1
        om.add_vars('Qg', ng, Qg, Qmin, Qmax)
228 1
        om.add_constraints('Pmis', nb, 'nonlinear')
229 1
        om.add_constraints('Qmis', nb, 'nonlinear')
230 1
        om.add_constraints('Sf', nl, 'nonlinear')
231 1
        om.add_constraints('St', nl, 'nonlinear')
232 1
        om.add_constraints('PQh', Apqh, array([]), ubpqh, ['Pg', 'Qg'])   ## npqh
233 1
        om.add_constraints('PQl', Apql, array([]), ubpql, ['Pg', 'Qg'])   ## npql
234 1
        om.add_constraints('vl',  Avl, lvl, uvl,   ['Pg', 'Qg'])   ## nvl
235 1
        om.add_constraints('ang', Aang, lang, uang, ['Va'])        ## nang
236

237
    ## y vars, constraints for piece-wise linear gen costs
238 1
    if ny > 0:
239 1
        om.add_vars('y', ny)
240 1
        om.add_constraints('ycon', Ay, array([]), by, ycon_vars)          ## ncony
241

242
    ## execute userfcn callbacks for 'formulation' stage
243 1
    run_userfcn(userfcn, 'formulation', om)
244

245 1
    return om

Read our documentation on viewing source code .

Loading