How to use the cvxopt.spdiag function in cvxopt
To help you get started, we’ve selected a few cvxopt examples, based on popular ways it is used in public projects.
Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
gsagnol / picos / mytests / complexSDP.py View on Github 
re2 = pic.Problem()
XX = re2.add_variable('XX',(6,6),'symmetric')
PP = pic.tools._cplx_mat_to_real_mat(P)
QQ = pic.tools._cplx_mat_to_real_mat(Q)
I = cvx.spdiag([1]*3)
II = pic.tools._cplx_mat_to_real_mat(I)
re2.add_constraint(PP|XX>2)
re2.add_constraint(QQ|XX>2)
re2.add_constraint(XX>>0)
re2.set_objective('min',II|XX)
re3 = pic.Problem()
X = re3.add_variable('X',(3,3),'symmetric')
Y = re3.add_variable('Y',(3,3))
I = cvx.spdiag([1]*3)
re3.add_constraint((P.real()|X)-(P.imag()|Y)>1)
re3.add_constraint((Q.real()|X)-(Q.imag()|Y)>1)
re3.add_constraint(((X & -Y)// (Y & X))>>0)
re3.set_objective('min',I|X)
c2 = pic.Problem()
Z = c2.add_variable('Z',(3,3),'hermitian')
u = c2.add_variable('u',2)
#c2.set_objective('min','I'|Z+3*u[0]+2*u[1])
c2.set_objective('min',Q*Q.H|Z)
c2.add_constraint(P|Z>1)
#c2.add_constraint(Q|Z>1)
#c2.add_constraint(P*Z+Z*P.H+ u[0]*Q*Q.H + u[1]*P*P.H >> P+P.H+Q+Q.H)
c2.add_constraint(R*Z+Z*R.H>>0)
c2.add_constraint('I'|Z==1)def d2AIbr_dV2(dIbr_dVa, dIbr_dVm, Ibr, Ybr, V, lam):
""" Computes 2nd derivatives of |complex current|**2 w.r.t. V.
"""
diaglam = spdiag(lam)
diagIbr_conj = spdiag(conj(Ibr))
Iaa, Iav, Iva, Ivv = d2Ibr_dV2(Ybr, V, diagIbr_conj * lam)
Haa = 2 * ( Iaa + dIbr_dVa.T * diaglam * conj(dIbr_dVa) ).real()
Hva = 2 * ( Iva + dIbr_dVm.T * diaglam * conj(dIbr_dVa) ).real()
Hav = 2 * ( Iav + dIbr_dVa.T * diaglam * conj(dIbr_dVm) ).real()
Hvv = 2 * ( Ivv + dIbr_dVm.T * diaglam * conj(dIbr_dVm) ).real()
return Haa, Hav, Hva, Hvvdef d2ASbr_dV2(dSbr_dVa, dSbr_dVm, Sbr, Cbr, Ybr, V, lam):
""" Computes 2nd derivatives of |complex power flow|**2 w.r.t. V.
"""
diaglam = spdiag(lam)
diagSbr_conj = spdiag(conj(Sbr))
Saa, Sav, Sva, Svv = d2Sbr_dV2(Cbr, Ybr, V, diagSbr_conj * lam)
Haa = 2 * ( Saa + dSbr_dVa.T * diaglam * conj(dSbr_dVa) ).real()
Hva = 2 * ( Sva + dSbr_dVm.T * diaglam * conj(dSbr_dVa) ).real()
Hav = 2 * ( Sav + dSbr_dVa.T * diaglam * conj(dSbr_dVm) ).real()
Hvv = 2 * ( Svv + dSbr_dVm.T * diaglam * conj(dSbr_dVm) ).real()
return Haa, Hav, Hva, Hvvdef __rmul__(self,fact):
selfcopy=self.copy()
if isinstance(fact,AffinExp):
if fact.isconstant():
fac,facString=cvx.sparse(fact.eval()),fact.string
else:
raise Exception('not implemented')
else:
fac,facString=_retrieve_matrix(fact,self.size[0])
if fac.size==(1,1) and selfcopy.size[0]<>1:
fac=fac[0]*cvx.spdiag([1.]*selfcopy.size[0])
if self.size==(1,1) and fac.size[1]<>1:
oldstring=selfcopy.string
selfcopy=selfcopy.diag(fac.size[1])
selfcopy.string=oldstring
if selfcopy.size[0]<>fac.size[1]:
raise Exception('incompatible dimensions')
bfac=_blocdiag(fac,selfcopy.size[1])
for k in selfcopy.factors:
newfac=bfac*selfcopy.factors[k]
selfcopy.factors[k]=newfac
if selfcopy.constant is None:
newfac=None
else:
newfac=bfac*selfcopy.constant
selfcopy.constant=newfac
selfcopy.size=(fac.size[0],selfcopy.size[1])def F(x = None, z = None):
if x is None: return 0, matrix(0.0, (3,1))
if max(abs(x)) >= 1.0: return None
u = 1 - x**2
val = -sum(log(u))
Df = div(2*x, u).T
if z is None: return val, Df
H = spdiag(2 * z[0] * div(1 + u**2, u**2))
return val, Df, Hdef F(x=None, z=None):
if x is None: return 0, cvxopt.matrix(1.0, (n,1))
if min(x) <= 0.0: return None
f = -sum(cvxopt.log(x))
Df = -(x**-1).T
if z is None: return f, Df
H = cvxopt.spdiag(z[0] * x**-2)
return f, Df, H
sol = cvxopt.solvers.cp(F, A=A, b=b)branch power flow w.r.t voltage.
"""
nl = len(branches)
nb = len(V)
f = matrix([l.from_bus._i for l in branches])
t = matrix([l.to_bus._i for l in branches])
# Compute currents.
If = Yf * V
It = Yt * V
Vnorm = div(V, abs(V))
diagVf = spdiag(V[f])
diagIf = spdiag(If)
diagVt = spdiag(V[t])
diagIt = spdiag(It)
diagV = spdiag(V)
diagVnorm = spdiag(Vnorm)
ibr = range(nl)
size = (nl, nb)
# Partial derivative of S w.r.t voltage phase angle.
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) * \def implicit_step(self):
"""
Integrate one step using trapezoidal method. Sets convergence and niter flags.
Returns
-------
None
"""
config = self.config
system = self.system
dae = self.system.dae
# constant short names
In = spdiag([1] * dae.n)
h = self.h
while self.err > config.tol and self.niter < config.maxit:
if self.t - self.t_jac >= 5:
dae.rebuild = True
self.t_jac = self.t
elif self.niter > 4:
dae.rebuild = True
elif dae.factorize:
dae.rebuild = True
elif self.config.honest is True:
dae.rebuild = True
# rebuild Jacobian
if dae.rebuild:
dae.factorize = Truedef dIbr_dV(Yf, Yt, V):
""" Computes partial derivatives of branch currents w.r.t. voltage.
Ray Zimmerman, "dIbr_dV.m", MATPOWER, version 4.0b1,
PSERC (Cornell), http://www.pserc.cornell.edu/matpower/
"""
# nb = len(V)
Vnorm = div(V, abs(V))
diagV = spdiag(V)
diagVnorm = spdiag(Vnorm)
dIf_dVa = Yf * 1j * diagV
dIf_dVm = Yf * diagVnorm
dIt_dVa = Yt * 1j * diagV
dIt_dVm = Yt * diagVnorm
# Compute currents.
If = Yf * V
It = Yt * V
return dIf_dVa, dIf_dVm, dIt_dVa, dIt_dVm, If, Itdef dSbus_dV(Y, V):
""" Computes the partial derivative of power injection w.r.t. voltage.
References:
Ray Zimmerman, "dSbus_dV.m", MATPOWER, version 3.2,
PSERC (Cornell), http://www.pserc.cornell.edu/matpower/
"""
I = Y * V
diagV = spdiag(V)
diagIbus = spdiag(I)
diagVnorm = spdiag(div(V, abs(V))) # Element-wise division.
dS_dVm = diagV * conj(Y * diagVnorm) + conj(diagIbus) * diagVnorm
dS_dVa = 1j * diagV * conj(diagIbus - Y * diagV)
return dS_dVm, dS_dVa
