e2nIEE / pandapower
 1 ```# Copyright (c) 1996-2015 PSERC. All rights reserved. ``` 2 ```# Use of this source code is governed by a BSD-style ``` 3 ```# license that can be found in the LICENSE file. ``` 4 5 1 ```"""Evaluates nonlinear constraints and their Jacobian for OPF. ``` 6 ```""" ``` 7 8 1 ```from numpy import zeros, ones, conj, exp, r_, Inf, arange ``` 9 10 1 ```from scipy.sparse import lil_matrix, vstack, hstack, csr_matrix as sparse ``` 11 12 1 ```from pandapower.pypower.idx_gen import GEN_BUS, PG, QG ``` 13 1 ```from pandapower.pypower.idx_brch import F_BUS, T_BUS, RATE_A ``` 14 15 1 ```from pandapower.pypower.makeSbus import makeSbus ``` 16 1 ```from pandapower.pypower.dSbus_dV import dSbus_dV ``` 17 1 ```from pandapower.pypower.dIbr_dV import dIbr_dV ``` 18 1 ```from pandapower.pypower.dSbr_dV import dSbr_dV ``` 19 1 ```from pandapower.pypower.dAbr_dV import dAbr_dV ``` 20 21 22 1 ```def opf_consfcn(x, om, Ybus, Yf, Yt, ppopt, il=None, *args): ``` 23 ``` """Evaluates nonlinear constraints and their Jacobian for OPF. ``` 24 25 ``` Constraint evaluation function for AC optimal power flow, suitable ``` 26 ``` for use with L{pips}. Computes constraint vectors and their gradients. ``` 27 28 ``` @param x: optimization vector ``` 29 ``` @param om: OPF model object ``` 30 ``` @param Ybus: bus admittance matrix ``` 31 ``` @param Yf: admittance matrix for "from" end of constrained branches ``` 32 ``` @param Yt: admittance matrix for "to" end of constrained branches ``` 33 ``` @param ppopt: PYPOWER options vector ``` 34 ``` @param il: (optional) vector of branch indices corresponding to ``` 35 ``` branches with flow limits (all others are assumed to be ``` 36 ``` unconstrained). The default is C{range(nl)} (all branches). ``` 37 ``` C{Yf} and C{Yt} contain only the rows corresponding to C{il}. ``` 38 39 ``` @return: C{h} - vector of inequality constraint values (flow limits) ``` 40 ``` limit^2 - flow^2, where the flow can be apparent power real power or ``` 41 ``` current, depending on value of C{OPF_FLOW_LIM} in C{ppopt} (only for ``` 42 ``` constrained lines). C{g} - vector of equality constraint values (power ``` 43 ``` balances). C{dh} - (optional) inequality constraint gradients, column ``` 44 ``` j is gradient of h(j). C{dg} - (optional) equality constraint gradients. ``` 45 46 ``` @see: L{opf_costfcn}, L{opf_hessfcn} ``` 47 48 ``` @author: Carlos E. Murillo-Sanchez (PSERC Cornell & Universidad ``` 49 ``` Autonoma de Manizales) ``` 50 ``` @author: Ray Zimmerman (PSERC Cornell) ``` 51 ``` """ ``` 52 ``` ##----- initialize ----- ``` 53 54 ``` ## unpack data ``` 55 1 ``` ppc = om.get_ppc() ``` 56 1 ``` baseMVA, bus, gen, branch = \ ``` 57 ``` ppc["baseMVA"], ppc["bus"], ppc["gen"], ppc["branch"] ``` 58 1 ``` vv, _, _, _ = om.get_idx() ``` 59 60 ``` ## problem dimensions ``` 61 1 ``` nb = bus.shape[0] ## number of buses ``` 62 1 ``` nl = branch.shape[0] ## number of branches ``` 63 1 ``` ng = gen.shape[0] ## number of dispatchable injections ``` 64 1 ``` nxyz = len(x) ## total number of control vars of all types ``` 65 66 ``` ## set default constrained lines ``` 67 1 ``` if il is None: ``` 68 0 ``` il = arange(nl) ## all lines have limits by default ``` 69 1 ``` nl2 = len(il) ## number of constrained lines ``` 70 71 ``` ## grab Pg & Qg ``` 72 1 ``` Pg = x[vv["i1"]["Pg"]:vv["iN"]["Pg"]] ## active generation in p.u. ``` 73 1 ``` Qg = x[vv["i1"]["Qg"]:vv["iN"]["Qg"]] ## reactive generation in p.u. ``` 74 75 ``` ## put Pg & Qg back in gen ``` 76 1 ``` gen[:, PG] = Pg * baseMVA ## active generation in MW ``` 77 1 ``` gen[:, QG] = Qg * baseMVA ## reactive generation in MVAr ``` 78 79 ``` ## rebuild Sbus ``` 80 1 ``` Sbus = makeSbus(baseMVA, bus, gen) ## net injected power in p.u. ``` 81 82 ``` ## ----- evaluate constraints ----- ``` 83 ``` ## reconstruct V ``` 84 1 ``` Va = x[vv["i1"]["Va"]:vv["iN"]["Va"]] ``` 85 1 ``` Vm = x[vv["i1"]["Vm"]:vv["iN"]["Vm"]] ``` 86 1 ``` V = Vm * exp(1j * Va) ``` 87 88 ``` ## evaluate power flow equations ``` 89 1 ``` mis = V * conj(Ybus * V) - Sbus ``` 90 91 ``` ##----- evaluate constraint function values ----- ``` 92 ``` ## first, the equality constraints (power flow) ``` 93 1 ``` g = r_[ mis.real, ## active power mismatch for all buses ``` 94 ``` mis.imag ] ## reactive power mismatch for all buses ``` 95 96 ``` ## then, the inequality constraints (branch flow limits) ``` 97 1 ``` if nl2 > 0: ``` 98 1 ``` flow_max = (branch[il, RATE_A] / baseMVA)**2 ``` 99 1 ``` flow_max[flow_max == 0] = Inf ``` 100 1 ``` if ppopt['OPF_FLOW_LIM'] == 2: ## current magnitude limit, |I| ``` 101 1 ``` If = Yf * V ``` 102 1 ``` It = Yt * V ``` 103 1 ``` h = r_[ If * conj(If) - flow_max, ## branch I limits (from bus) ``` 104 ``` It * conj(It) - flow_max ].real ## branch I limits (to bus) ``` 105 ``` else: ``` 106 ``` ## compute branch power flows ``` 107 ``` ## complex power injected at "from" bus (p.u.) ``` 108 0 ``` Sf = V[ branch[il, F_BUS].astype(int) ] * conj(Yf * V) ``` 109 ``` ## complex power injected at "to" bus (p.u.) ``` 110 0 ``` St = V[ branch[il, T_BUS].astype(int) ] * conj(Yt * V) ``` 111 0 ``` if ppopt['OPF_FLOW_LIM'] == 1: ## active power limit, P (Pan Wei) ``` 112 0 ``` h = r_[ Sf.real**2 - flow_max, ## branch P limits (from bus) ``` 113 ``` St.real**2 - flow_max ] ## branch P limits (to bus) ``` 114 ``` else: ## apparent power limit, |S| ``` 115 0 ``` h = r_[ Sf * conj(Sf) - flow_max, ## branch S limits (from bus) ``` 116 ``` St * conj(St) - flow_max ].real ## branch S limits (to bus) ``` 117 ``` else: ``` 118 1 ``` h = zeros((0,1)) ``` 119 120 ``` ##----- evaluate partials of constraints ----- ``` 121 ``` ## index ranges ``` 122 1 ``` iVa = arange(vv["i1"]["Va"], vv["iN"]["Va"]) ``` 123 1 ``` iVm = arange(vv["i1"]["Vm"], vv["iN"]["Vm"]) ``` 124 1 ``` iPg = arange(vv["i1"]["Pg"], vv["iN"]["Pg"]) ``` 125 1 ``` iQg = arange(vv["i1"]["Qg"], vv["iN"]["Qg"]) ``` 126 1 ``` iVaVmPgQg = r_[iVa, iVm, iPg, iQg].T ``` 127 128 ``` ## compute partials of injected bus powers ``` 129 1 ``` dSbus_dVm, dSbus_dVa = dSbus_dV(Ybus, V) ## w.r.t. V ``` 130 ``` ## Pbus w.r.t. Pg, Qbus w.r.t. Qg ``` 131 1 ``` neg_Cg = sparse((-ones(ng), (gen[:, GEN_BUS], range(ng))), (nb, ng)) ``` 132 133 ``` ## construct Jacobian of equality constraints (power flow) and transpose it ``` 134 1 ``` dg = lil_matrix((2 * nb, nxyz)) ``` 135 1 ``` blank = sparse((nb, ng)) ``` 136 1 ``` dg[:, iVaVmPgQg] = vstack([ ``` 137 ``` ## P mismatch w.r.t Va, Vm, Pg, Qg ``` 138 ``` hstack([dSbus_dVa.real, dSbus_dVm.real, neg_Cg, blank]), ``` 139 ``` ## Q mismatch w.r.t Va, Vm, Pg, Qg ``` 140 ``` hstack([dSbus_dVa.imag, dSbus_dVm.imag, blank, neg_Cg]) ``` 141 ``` ], "csr") ``` 142 1 ``` dg = dg.T ``` 143 144 1 ``` if nl2 > 0: ``` 145 ``` ## compute partials of Flows w.r.t. V ``` 146 1 ``` if ppopt['OPF_FLOW_LIM'] == 2: ## current ``` 147 1 ``` dFf_dVa, dFf_dVm, dFt_dVa, dFt_dVm, Ff, Ft = \ ``` 148 ``` dIbr_dV(branch[il, :], Yf, Yt, V) ``` 149 ``` else: ## power ``` 150 0 ``` dFf_dVa, dFf_dVm, dFt_dVa, dFt_dVm, Ff, Ft = \ ``` 151 ``` dSbr_dV(branch[il, :], Yf, Yt, V) ``` 152 1 ``` if ppopt['OPF_FLOW_LIM'] == 1: ## real part of flow (active power) ``` 153 0 ``` dFf_dVa = dFf_dVa.real ``` 154 0 ``` dFf_dVm = dFf_dVm.real ``` 155 0 ``` dFt_dVa = dFt_dVa.real ``` 156 0 ``` dFt_dVm = dFt_dVm.real ``` 157 0 ``` Ff = Ff.real ``` 158 0 ``` Ft = Ft.real ``` 159 160 ``` ## squared magnitude of flow (of complex power or current, or real power) ``` 161 1 ``` df_dVa, df_dVm, dt_dVa, dt_dVm = \ ``` 162 ``` dAbr_dV(dFf_dVa, dFf_dVm, dFt_dVa, dFt_dVm, Ff, Ft) ``` 163 164 ``` ## construct Jacobian of inequality constraints (branch limits) ``` 165 ``` ## and transpose it. ``` 166 1 ``` dh = lil_matrix((2 * nl2, nxyz)) ``` 167 1 ``` dh[:, r_[iVa, iVm].T] = vstack([ ``` 168 ``` hstack([df_dVa, df_dVm]), ## "from" flow limit ``` 169 ``` hstack([dt_dVa, dt_dVm]) ## "to" flow limit ``` 170 ``` ], "csr") ``` 171 1 ``` dh = dh.T ``` 172 ``` else: ``` 173 1 ``` dh = None ``` 174 175 1 ``` return h, g, dh, dg ```

Read our documentation on viewing source code .

Loading