Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
1 - dae.x[self.omega_m],
mul(dae.x[self.ird], self.x1) + mul(dae.y[self.isd], self.xmu))
dae.g[self.vsd] = -dae.y[self.vsd] - mul(dae.y[self.v],
sin(dae.y[self.a]))
dae.g[self.vsq] = -dae.y[self.vsq] + mul(dae.y[self.v],
cos(dae.y[self.a]))
dae.g[self.vref] = self.vref0 - dae.y[self.vref]
dae.g[self.pwa] = mmax(mmin(2 * dae.x[self.omega_m] - 1, 1),
0) - dae.y[self.pwa]
dae.hard_limit(self.pwa, 0, 1)
dae.g[self.pw] = -dae.y[self.pw] + mul(
0.5, dae.y[self.cp], self.ngen, pi, self.rho, (self.R)**2,
(self.Vwn)**3, div(1, self.mva_mega), (dae.x[self.vw])**3)
dae.g[self.cp] = -dae.y[self.cp] + mul(
-1.1 + mul(25.52, div(1, dae.y[self.ilamb])) + mul(
-0.088, dae.x[self.theta_p]),
exp(mul(-12.5, div(1, dae.y[self.ilamb]))))
dae.g[self.lamb] = -dae.y[self.lamb] + mul(
4, self.R, self.fn, self.ngb, dae.x[self.omega_m], pi,
div(1, self.Vwn), div(1, self.npole), div(1, dae.x[self.vw]))
dae.g[self.ilamb] = div(
1,
div(1, dae.y[self.lamb] + mul(0.08, dae.x[self.theta_p])) + mul(
-0.035, div(1, 1 +
(dae.x[self.theta_p])**3))) - dae.y[self.ilamb]
dae.g += spmatrix(
mul(
self.u0, -mul(dae.x[self.ird], dae.y[self.vrd]) - mul(
dae.x[self.irq], dae.y[self.vrq]) - mul(
dae.y[self.isd], dae.y[self.vsd]) - mul(
dae.y[self.isq], dae.y[self.vsq])), self.a,
def jac0(self, dae):
dae.add_jac(Gy0, mul(self.Ic5, self.u0), self.In, self.v)
dae.add_jac(Gy0, -self.u0, self.In, self.In)
dae.add_jac(Gy0, mul(self.Ic4, self.toSg, self.u0), self.In, self.pm)
dae.add_jac(Gy0, mul(self.Ic3, self.toSg, self.u0), self.In, self.p)
dae.add_jac(Gy0, -self.u0, self.v1, self.v1)
dae.add_jac(Gy0, mul(self.T12, self.u0), self.v2, self.v1)
dae.add_jac(Gy0, -self.u0, self.v2, self.v2)
dae.add_jac(Gy0, -self.u0, self.v3, self.v3)
dae.add_jac(Gy0, mul(self.T34, self.u0), self.v3, self.v2)
dae.add_jac(Gy0, -self.u0, self.v4, self.v4)
dae.add_jac(Gy0, mul(self.KsT56, self.u0), self.v4, self.v3)
dae.add_jac(Gy0, self.u0, self.vss, self.v4)
dae.add_jac(Gy0, -self.u0, self.vss, self.vss)
dae.add_jac(Gy0, -self.u0, self.vst, self.vst)
dae.add_jac(Gy0, self.u0, self.vst, self.vss)
dae.add_jac(Gy0, self.u0, self.vf, self.vst)
dae.add_jac(Gx0, mul(self.Ic2, self.u0), self.In, self.w)
dae.add_jac(Gx0, mul(self.Ic1, self.u0), self.In, self.omega)
dae.add_jac(Gx0, mul(self.A5, self.u0), self.v1, self.q1)
dae.add_jac(Gx0, self.u0, self.v1, self.q0)
def gcall(self, dae):
if sum(self.u) == 0:
return
Vm = polar(dae.y[self.v], dae.y[self.a])
Vsh = polar(dae.y[self.vsh], dae.y[self.ash])
Ish = mul(self.Ysh, Vsh - Vm)
IshC = conj(Ish)
Ssh = mul(Vm, IshC)
# check the Vsh and Ish limits during PF iterations
vupper = list(abs(Vsh) - self.vshmax)
vlower = list(abs(Vsh) - self.vshmin)
iupper = list(abs(IshC) - self.Ishmax)
# check for Vsh and Ish limit violations
if self.system.pflow.niter >= self.system.pflow.config.ipv2pq:
for i in range(self.n):
if self.u[i] and (vupper[i] > 0 or vlower[i] < 0 or iupper[i] > 0):
if i not in self.vio.keys():
self.vio[i] = list()
if vupper[i] > 0:
if 'vmax' not in self.vio[i]:
self.vio[i].append('vmax')
dae.g += spmatrix(div(mul(self.u, dae.y[self.pdc]), dae.y[self.v1] - dae.y[self.v2]), self.v1,
[0] * self.n, (dae.m, 1), 'd')
dae.g -= spmatrix(div(mul(self.u, dae.y[self.pdc]), dae.y[self.v1] - dae.y[self.v2]), self.v2,
[0] * self.n, (dae.m, 1), 'd') # negative current injection
dae.g[self.ash] = mul(self.u, Ssh.real() - dae.y[self.psh]) # (2)
dae.g[self.vsh] = mul(self.u, Ssh.imag() - dae.y[self.qsh]) # (3)
# PQ, PV or V control
# (12), (15)
dae.g[self.psh] = mul(self.u, dae.y[self.psh] + mul(self.R, dae.y[self.v1] - self.vdcref)
- self.pref0, self.PQ + self.PV) + \
mul(self.u, (dae.y[self.v1] - dae.y[self.v2]) - self.vdcref0, self.vV + self.vQ, )
dae.g[self.qsh] = mul(self.u, dae.y[self.qsh] - self.qref0, self.PQ + self.vQ) + \
mul(self.u, dae.y[self.v] - self.vref0, self.PV + self.vV) # (13), (16)
for comp, var in self.vio.items():
for count, item in enumerate(var):
if item == 'vmax':
yidx = self.vsh[comp]
ylim = self.vshmax[comp]
elif item == 'vmin':
yidx = self.vsh[comp]
ylim = self.vshmin[comp]
elif item == 'Imax':
yidx = self.Ish[comp]
ylim = self.Ishmax[comp]
else:
raise NameError('Unknown limit variable name <{0}>.'.format(item))
if count == 0:
def fcall(self, dae):
dae.f[self.q0] = mul(dae.x[self.q1], self.u0)
dae.f[self.q1] = mul(dae.x[self.q2], self.u0)
dae.f[self.q2] = mul(dae.x[self.q3], self.u0)
dae.f[self.q3] = mul(
self.u0, div(1, self.H4),
dae.y[self.In] - dae.x[self.q0] - mul(self.H1, dae.x[self.q1]) -
mul(self.H2, dae.x[self.q2]) - mul(self.H3, dae.x[self.q3]))
dae.f[self.x1] = mul(
self.u0, div(1, self.T2),
-dae.x[self.x1] + mul(dae.y[self.v1], 1 - self.T12))
dae.f[self.x2] = mul(
self.u0, div(1, self.T4),
-dae.x[self.x2] + mul(dae.y[self.v2], 1 - self.T34))
dae.f[self.x3] = mul(self.u0, div(1, self.T6),
-dae.x[self.x3] + mul(self.KsT56, dae.y[self.v3]))
def build_gy(self, dae):
"""Build line Jacobian matrix"""
if not self.n:
idx = range(dae.m)
dae.set_jac(Gy, 1e-6, idx, idx)
return
Vn = polar(1.0, dae.y[self.a])
Vc = mul(dae.y[self.v], Vn)
Ic = self.Y * Vc
diagVn = spdiag(Vn)
diagVc = spdiag(Vc)
diagIc = spdiag(Ic)
dS = self.Y * diagVn
dS = diagVc * conj(dS)
dS += conj(diagIc) * diagVn
dR = diagIc
dR -= self.Y * diagVc
dR = diagVc.H.T * dR
self.gy_store = sparse([[dR.imag(), dR.real()], [dS.real(),
dS.imag()]])
def gycall(self, dae):
dV2 = mul(self.u, 2 * dae.y[self.v])
dPdV = mul(dV2, matrix(self.g, (self.n, 1), 'd'))
dQdV = -mul(dV2, matrix(self.b, (self.n, 1), 'd'))
dae.add_jac(Gy, dPdV, self.a, self.v)
dae.add_jac(Gy, dQdV, self.v, self.v)
dSf_dVa = 1j * (conj(diagIf) *
spmatrix(V[f], ibr, f, size) - diagVf * conj(Yf * diagV))
dSt_dVa = 1j * (conj(diagIt) *
spmatrix(V[t], ibr, t, size) - diagVt * conj(Yt * diagV))
# Partial derivative of S w.r.t. voltage amplitude.
dSf_dVm = diagVf * conj(Yf * diagVnorm) + conj(diagIf) * \
spmatrix(Vnorm[f], ibr, f, size)
dSt_dVm = diagVt * conj(Yt * diagVnorm) + conj(diagIt) * \
spmatrix(Vnorm[t], ibr, t, size)
# Compute power flow vectors.
Sf = mul(V[f], conj(If))
St = mul(V[t], conj(It))
return dSf_dVa, dSf_dVm, dSt_dVa, dSt_dVm, Sf, St
dae.y[self.In1] = mul(
self.u0,
mul(self.Ic11, -1 + dae.x[self.omega]) + mul(
self.Ic12, -1 + dae.x[self.w]) + mul(self.Ic15, dae.y[self.v])
+ mul(self.Ic13, dae.y[self.p], self.toSg) + mul(
self.Ic14, dae.y[self.pm], self.toSg))
dae.y[self.In2] = mul(
self.u0,
mul(self.Ic21, -1 + dae.x[self.omega]) + mul(
self.Ic22, -1 + dae.x[self.w]) + mul(self.Ic25, dae.y[self.v])
+ mul(self.Ic23, dae.y[self.p], self.toSg) + mul(
self.Ic24, dae.y[self.pm], self.toSg))
dae.x[self.x1] = mul(self.K1, dae.y[self.In1], self.u0)
dae.x[self.x2] = mul(self.K2, dae.y[self.In2], self.u0)
dae.y[self.In] = mul(self.u0, dae.y[self.In1] + dae.y[self.In2])
dae.x[self.u3] = mul(dae.y[self.In], self.T34, self.u0)
Fy, -mul(self.Ki, div(1, self.Teq), div(1, dae.x[self.omega_m]),
div(1, self.psip - mul(dae.y[self.isd], self.xd))),
self.isq, self.dwdt)
dae.add_jac(
Fy,
mul(
div(1, self.Teq), div(1, dae.x[self.omega_m]),
div(1, self.psip - mul(dae.y[self.isd], self.xd))), self.isq,
self.pwa)
dae.add_jac(
Fy, -mul(self.Kdc, div(1, self.Teq), div(1, dae.x[self.omega_m]),
div(1, self.psip - mul(dae.y[self.isd], self.xd))),
self.isq, self.v1)
dae.add_jac(
Fy,
mul(
self.xd, div(1, self.Teq), div(1, dae.x[self.omega_m]),
(self.psip - mul(dae.y[self.isd], self.xd))**-2,
dae.y[self.pwa] - mul(self.Kcoi, dae.y[self.dwdt_coi]) - mul(
self.Kdc, dae.y[self.v1] - dae.y[self.v2]) - mul(
self.Ki, dae.y[self.dwdt])), self.isq, self.isd)
dae.add_jac(
Fy, -mul(self.Kcoi, div(1, self.Teq), div(1, dae.x[self.omega_m]),
div(1, self.psip - mul(dae.y[self.isd], self.xd))),
self.isq, self.dwdt_coi)