r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1# ------------------------------------------------------------------------------

2# Copyright (C) 2012-2017, 2020 Guillaume Sagnol

3# Copyright (C) 2018-2019 Maximilian Stahlberg

4#

5# This file is part of PICOS.

6#

7# PICOS is free software: you can redistribute it and/or modify it under the

9# Foundation, either version 3 of the License, or (at your option) any later

10# version.

11#

12# PICOS is distributed in the hope that it will be useful, but WITHOUT ANY

13# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR

14# A PARTICULAR PURPOSE. See the GNU General Public License for more details.

15#

16# You should have received a copy of the GNU General Public License along with

17# this program. If not, see <http://www.gnu.org/licenses/>.

18# ------------------------------------------------------------------------------

20"""Implementation of matrix norm constraints."""

22import operator

23from collections import namedtuple

25import cvxopt as cvx

27from .. import glyphs

28from ..apidoc import api_end, api_start

29from .constraint import Constraint, ConstraintConversion

31_API_START = api_start(globals())

32# -------------------------------

35class MatrixNormConstraint(Constraint):

36 """Upper bound on a matrix :math:`(p,q)`-norm."""

38 class VectorNormConversion(ConstraintConversion):

39 """Upper bound on a :math:`(p,q)`-norm constraint conversion."""

41 @classmethod

42 def predict(cls, subtype, options):

43 """Implement :meth:`~.constraint.ConstraintConversion.predict`."""

44 from ..expressions import AffineExpression, RealVariable, Norm

46 shape, pNum, pDen, qNum, qDen = subtype

47 rows, cols = shape

49 # HACK: Whether the bound is constant is irrelevant.

50 bound = AffineExpression.make_type((1, 1), None, None)

51 pNorm = Norm.make_type((rows, 1), pNum, pDen, pNum, pDen)

52 qNorm = Norm.make_type((cols, 1), qNum, qDen, qNum, qDen)

54 yield ("var", RealVariable.make_var_type(dim=cols, bnd=0), 1)

55 yield ("con", pNorm.predict(operator.__le__, bound), cols)

56 yield ("con", qNorm.predict(operator.__le__, bound), 1)

58 @classmethod

59 def convert(cls, con, options):

60 """Implement :meth:`~.constraint.ConstraintConversion.convert`."""

61 from ..expressions import no_refinement, Norm, RealVariable

62 from ..modeling import Problem

64 norm = con.norm

65 p, q = norm.p, norm.q

66 cols = norm.x.size[1]

68 P = Problem()

70 u = RealVariable("__u", cols)

72 # Bound the p-norm of every column.

73 with no_refinement():

74 for j in range(cols):

75 P.add_constraint(Norm(norm.x[:, j], p) <= u[j])

77 # Bound the q-norm of the column norms.

80 return P

82 def __init__(self, norm, upperBound):

83 """Construct a :class:`MatrixNormConstraint`.

85 :param ~picos.expressions.Norm norm:

86 The norm.

87 :param ~picos.expressions.AffineExpression upperBound:

88 The scalar upper bound.

89 """

90 from ..expressions import AffineExpression, Norm

92 assert isinstance(norm, Norm)

93 assert isinstance(upperBound, AffineExpression)

94 assert len(upperBound) == 1

96 assert norm.qnum is not None or norm.pnum != norm.qnum, \

97 "Won't create a (p,q)-norm constraint for a p-norm."

99 self.norm = norm

100 self.upperBound = upperBound

102 super(MatrixNormConstraint, self).__init__(norm._typeStr)

104 Subtype = namedtuple("Subtype", ("shape", "pNum", "pDen", "qNum", "qDen"))

106 @classmethod

107 def _cost(cls, subtype):

108 return subtype.shape[0] * subtype.shape[1] + 1

110 def _subtype(self):

111 return self.Subtype(self.norm.x.size, self.norm.pnum, self.norm.pden,

112 self.norm.qnum, self.norm.qden)

114 def _expression_names(self):

115 yield "norm"

116 yield "upperBound"

118 def _str(self):

119 return glyphs.le(self.norm.string, self.upperBound.string)

121 def _get_slack(self):

122 return self.upperBound.safe_value - self.norm.safe_value

125class SpectralNormConstraint(Constraint):

126 """Spectral norm of a matrix."""

128 class Conversion(ConstraintConversion):

129 """Spectral norm constraint conversion."""

131 @classmethod

132 def predict(cls, subtype, options):

133 """Implement :meth:`~.constraint.ConstraintConversion.predict`."""

134 from . import (ComplexLMIConstraint, LMIConstraint)

135 m, n = subtype.shape

137 if subtype.hermitian:

138 if subtype.complex:

139 yield ("con", ComplexLMIConstraint.make_type(diag=n), 2)

140 else:

141 yield ("con", LMIConstraint.make_type(diag=n), 2)

142 else:

143 if subtype.complex:

144 yield ("con", ComplexLMIConstraint.make_type(diag=n+m), 1)

145 else:

146 yield ("con", LMIConstraint.make_type(diag=n+m), 1)

148 @classmethod

149 def convert(cls, con, options):

150 """Implement :meth:`~.constraint.ConstraintConversion.convert`."""

151 from ..modeling import Problem

152 from ..expressions import block, Constant

154 x = con.norm.x

155 m, n = x.shape

156 t = con.upperBound

158 In = Constant('I', cvx.spdiag([1.] * n))

159 P = Problem()

160 if x.hermitian:

162 P.add_constraint(x >> -t * In)

163 else:

164 Im = Constant('I', cvx.spdiag([1.] * m))

165 P.add_constraint(block([[t*Im, x], [x.H, t*In]]) >> 0)

166 return P

168 def __init__(self, norm, upperBound):

169 """Construct a :class:`SpectralNormConstraint`.

171 :param ~picos.expressions.SpectralNorm norm:

172 Constrained spectral norm

173 :param ~picos.expressions.AffineExpression upperBound:

174 Upper bound on the expression.

175 """

176 from ..expressions import AffineExpression, SpectralNorm

178 assert isinstance(norm, SpectralNorm)

179 assert isinstance(upperBound, AffineExpression)

180 assert len(upperBound) == 1

182 self.norm = norm

183 self.upperBound = upperBound

185 super(SpectralNormConstraint, self).__init__(norm._typeStr)

187 Subtype = namedtuple("Subtype", ("shape", "complex", "hermitian"))

189 def _subtype(self):

190 x = self.norm.x

191 return self.Subtype(x.shape, x.complex, x.hermitian)

193 @classmethod

194 def _cost(cls, subtype):

195 m = subtype.shape[0]

196 n = subtype.shape[1]

198 if subtype.hermitian:

199 if subtype.complex:

200 return 2 * n**2 + 1

201 else:

202 return 2 * n*(n+1) // 2 + 1

203 else:

204 if subtype.complex:

205 return (n + m) ** 2 + 1

206 else:

207 return (n + m) * (n + m + 1) // 2 + 1

209 def _expression_names(self):

210 yield "norm"

211 yield "upperBound"

213 def _str(self):

214 return glyphs.le(self.norm.string, self.upperBound.string)

216 def _get_slack(self):

217 return self.upperBound.safe_value - self.norm.safe_value

220class NuclearNormConstraint(Constraint):

221 """Nuclear norm of a matrix."""

223 class Conversion(ConstraintConversion):

224 """Nuclear norm constraint conversion."""

226 @classmethod

227 def predict(cls, subtype, options):

228 """Implement :meth:`~.constraint.ConstraintConversion.predict`."""

229 from ..expressions import SymmetricVariable, HermitianVariable

230 from . import (ComplexLMIConstraint, LMIConstraint,

231 AffineConstraint, ComplexAffineConstraint)

233 m, n = subtype.shape

235 if subtype.complex:

236 MatrixVariable = HermitianVariable

237 dm, dn = m**2, n**2

238 else:

239 MatrixVariable = SymmetricVariable

240 dm = (m * (m + 1)) // 2

241 dn = (n * (n + 1)) // 2

243 yield ("var", MatrixVariable.make_var_type(dim=dm, bnd=0), 1)

244 yield ("var", MatrixVariable.make_var_type(dim=dn, bnd=0), 1)

246 if subtype.hermitian:

247 l = (n * (n + 1)) // 2 # Number of lower triangular elements.

248 if subtype.complex:

249 yield ("con", ComplexAffineConstraint.make_type(dim=l), 1)

250 yield ("con", ComplexLMIConstraint.make_type(diag=n), 2)

251 else:

252 yield ("con", AffineConstraint.make_type(dim=l, eq=True), 1)

253 yield ("con", LMIConstraint.make_type(diag=n), 2)

254 else:

255 if subtype.complex:

256 yield ("con", ComplexLMIConstraint.make_type(diag=n+m), 1)

257 else:

258 yield ("con", LMIConstraint.make_type(diag=n+m), 1)

260 yield ("con", AffineConstraint.make_type(dim=1, eq=False), 1)

262 @classmethod

263 def convert(cls, con, options):

264 """Implement :meth:`~.constraint.ConstraintConversion.convert`."""

265 from ..modeling import Problem

266 from ..expressions import SymmetricVariable, HermitianVariable

267 from ..expressions.algebra import block, trace

269 x = con.norm.x

270 m, n = x.shape

271 t = con.upperBound

273 if con.norm._complex:

274 Y = HermitianVariable("__Y", (m, m))

275 Z = HermitianVariable("__Z", (n, n))

276 else:

277 Y = SymmetricVariable("__Y", (m, m))

278 Z = SymmetricVariable("__Z", (n, n))

280 P = Problem()

281 if x.hermitian:

282 P.add_constraint(x.trilvec == (Y - Z).trilvec)

285 P.add_constraint(trace(Y) + trace(Z) <= t)

286 else:

287 P.add_constraint(block([[Y, x], [x.H, Z]]) >> 0)

288 P.add_constraint(trace(Y) + trace(Z) <= 2 * t)

289 return P

291 def __init__(self, norm, upperBound):

292 """Construct a :class:`NuclearNormConstraint`.

294 :param ~picos.expressions.NuclearNorm norm:

295 Constrained nuclear norm

296 :param ~picos.expressions.AffineExpression upperBound:

297 Upper bound on the expression.

298 """

299 from ..expressions import AffineExpression, NuclearNorm

301 assert isinstance(norm, NuclearNorm)

302 assert isinstance(upperBound, AffineExpression)

303 assert len(upperBound) == 1

305 self.norm = norm

306 self.upperBound = upperBound

308 super(NuclearNormConstraint, self).__init__(norm._typeStr)

310 Subtype = namedtuple("Subtype", ("shape", "complex", "hermitian"))

312 def _subtype(self):

313 x = self.norm.x

314 return self.Subtype(x.shape, x.complex, x.hermitian)

316 @classmethod

317 def _cost(cls, subtype):

318 m, n = subtype.shape

320 if subtype.hermitian:

321 if subtype.complex:

322 return 2 * n ** 2 + 1

323 else:

324 return 2 * n * (n+1) // 2 + 1

325 else:

326 if subtype.complex:

327 return (n + m) ** 2 + 1

328 else:

329 return (n + m) * (n + m + 1) // 2 + 1

331 def _expression_names(self):

332 yield "norm"

333 yield "upperBound"

335 def _str(self):

336 return glyphs.le(self.norm.string, self.upperBound.string)

338 def _get_slack(self):

339 return self.upperBound.safe_value - self.norm.safe_value

342# --------------------------------------

343__all__ = api_end(_API_START, globals())