Coverage for picos/constraints/con_matnorm.py: 98.45%

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 

8# terms of the GNU General Public License as published by the Free Software 

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# ------------------------------------------------------------------------------ 

19 

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

21 

22import operator 

23from collections import namedtuple 

24 

25import cvxopt as cvx 

26 

27from .. import glyphs 

28from ..apidoc import api_end, api_start 

29from .constraint import Constraint, ConstraintConversion 

30 

31_API_START = api_start(globals()) 

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

33 

34 

35class MatrixNormConstraint(Constraint): 

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

37 

38 class VectorNormConversion(ConstraintConversion): 

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

40 

41 @classmethod 

42 def predict(cls, subtype, options): 

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

44 from ..expressions import AffineExpression, RealVariable, Norm 

45 

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

47 rows, cols = shape 

48 

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) 

53 

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) 

57 

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 

63 

64 norm = con.norm 

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

66 cols = norm.x.size[1] 

67 

68 P = Problem() 

69 

70 u = RealVariable("__u", cols) 

71 

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]) 

76 

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

78 P.add_constraint(Norm(u, q) <= con.upperBound) 

79 

80 return P 

81 

82 def __init__(self, norm, upperBound): 

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

84 

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 

91 

92 assert isinstance(norm, Norm) 

93 assert isinstance(upperBound, AffineExpression) 

94 assert len(upperBound) == 1 

95 

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." 

98 

99 self.norm = norm 

100 self.upperBound = upperBound 

101 

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

103 

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

105 

106 @classmethod 

107 def _cost(cls, subtype): 

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

109 

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) 

113 

114 def _expression_names(self): 

115 yield "norm" 

116 yield "upperBound" 

117 

118 def _str(self): 

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

120 

121 def _get_slack(self): 

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

123 

124 

125class SpectralNormConstraint(Constraint): 

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

127 

128 class Conversion(ConstraintConversion): 

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

130 

131 @classmethod 

132 def predict(cls, subtype, options): 

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

134 from . import (ComplexLMIConstraint, LMIConstraint) 

135 m, n = subtype.shape 

136 

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) 

147 

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 

153 

154 x = con.norm.x 

155 m, n = x.shape 

156 t = con.upperBound 

157 

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

159 P = Problem() 

160 if x.hermitian: 

161 P.add_constraint(x << t*In) 

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 

167 

168 def __init__(self, norm, upperBound): 

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

170 

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 

177 

178 assert isinstance(norm, SpectralNorm) 

179 assert isinstance(upperBound, AffineExpression) 

180 assert len(upperBound) == 1 

181 

182 self.norm = norm 

183 self.upperBound = upperBound 

184 

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

186 

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

188 

189 def _subtype(self): 

190 x = self.norm.x 

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

192 

193 @classmethod 

194 def _cost(cls, subtype): 

195 m = subtype.shape[0] 

196 n = subtype.shape[1] 

197 

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 

208 

209 def _expression_names(self): 

210 yield "norm" 

211 yield "upperBound" 

212 

213 def _str(self): 

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

215 

216 def _get_slack(self): 

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

218 

219 

220class NuclearNormConstraint(Constraint): 

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

222 

223 class Conversion(ConstraintConversion): 

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

225 

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) 

232 

233 m, n = subtype.shape 

234 

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 

242 

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

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

245 

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) 

259 

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

261 

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 

268 

269 x = con.norm.x 

270 m, n = x.shape 

271 t = con.upperBound 

272 

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)) 

279 

280 P = Problem() 

281 if x.hermitian: 

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

283 P.add_constraint(Y >> 0) 

284 P.add_constraint(Z >> 0) 

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 

290 

291 def __init__(self, norm, upperBound): 

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

293 

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 

300 

301 assert isinstance(norm, NuclearNorm) 

302 assert isinstance(upperBound, AffineExpression) 

303 assert len(upperBound) == 1 

304 

305 self.norm = norm 

306 self.upperBound = upperBound 

307 

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

309 

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

311 

312 def _subtype(self): 

313 x = self.norm.x 

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

315 

316 @classmethod 

317 def _cost(cls, subtype): 

318 m, n = subtype.shape 

319 

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 

330 

331 def _expression_names(self): 

332 yield "norm" 

333 yield "upperBound" 

334 

335 def _str(self): 

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

337 

338 def _get_slack(self): 

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

340 

341 

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

343__all__ = api_end(_API_START, globals())