Coverage for picos/expressions/cone_psd.py: 83.75%

80 statements  

« prev     ^ index     » next       coverage.py v6.5.0, created at 2023-03-26 07:46 +0000

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

2# Copyright (C) 2020 Maximilian Stahlberg 

3# 

4# This file is part of PICOS. 

5# 

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

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

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

9# version. 

10# 

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

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

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

14# 

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

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

17# ------------------------------------------------------------------------------ 

18 

19"""Implements the positive semidefinite cone.""" 

20 

21import operator 

22from collections import namedtuple 

23 

24from .. import glyphs 

25from ..apidoc import api_end, api_start 

26from ..constraints.uncertain import ScenarioUncertainConicConstraint 

27from .cone import Cone 

28from .exp_affine import AffineExpression, ComplexAffineExpression 

29from .expression import ExpressionType 

30from .uncertain import ScenarioPerturbationSet, UncertainAffineExpression 

31 

32_API_START = api_start(globals()) 

33# ------------------------------- 

34 

35 

36class PositiveSemidefiniteCone(Cone): 

37 r"""The positive semidefinite cone. 

38 

39 Unlike other :class:`cones <.cone.Cone>` which are defined only on 

40 :math:`\mathbb{K}^n`, this cone accepts both symmetric and hermitian 

41 matrices as well as their 

42 :attr:`special vectorization <.exp_biaffine.BiaffineExpression.svec>` 

43 as members. 

44 

45 :Example: 

46 

47 >>> from picos import Constant, PositiveSemidefiniteCone 

48 >>> R = Constant("R", range(16), (4, 4)) 

49 >>> S = R + R.T 

50 >>> S.shape 

51 (4, 4) 

52 >>> S.svec.shape 

53 (10, 1) 

54 >>> S.svec << PositiveSemidefiniteCone() # Constrain the matrix via svec(). 

55 <4×4 LMI Constraint: R + Rᵀ ≽ 0> 

56 >>> C = S << PositiveSemidefiniteCone(); C # Constrain the matrix directly. 

57 <4×4 LMI Constraint: R + Rᵀ ≽ 0> 

58 >>> C.conic_membership_form[0] # The conic form still refers to svec(). 

59 <10×1 Real Constant: svec(R + Rᵀ)> 

60 >>> C.conic_membership_form[1] 

61 <10-dim. Positive Semidefinite Cone: {svec(A) : xᵀ·A·x ≥ 0 ∀ x}> 

62 """ 

63 

64 def __init__(self, dim=None): 

65 r"""Construct a :class:`PositiveSemidefiniteCone`. 

66 

67 If a fixed dimensionality is given, this must be the dimensiona of the 

68 special vectorization. For a :math:`n \times n` matrix, this is 

69 :math:`\frac{n(n + 1)}{2}`. 

70 """ 

71 Cone.__init__(self, dim, "Positive Semidefinite Cone", 

72 glyphs.set(glyphs.sep(glyphs.svec("A"), glyphs.forall(glyphs.ge( 

73 glyphs.mul(glyphs.mul(glyphs.transp("x"), "A"), "x"), 

74 glyphs.scalar(0)), "x")))) 

75 

76 def _get_mutables(self): 

77 return frozenset() 

78 

79 def _replace_mutables(self): 

80 return self 

81 

82 Subtype = namedtuple("Subtype", ("dim",)) 

83 

84 def _get_subtype(self): 

85 return self.Subtype(self.dim) 

86 

87 @classmethod 

88 def _predict(cls, subtype, relation, other): 

89 assert isinstance(subtype, cls.Subtype) 

90 

91 if relation == operator.__rshift__: 

92 if issubclass(other.clstype, ( 

93 ComplexAffineExpression, UncertainAffineExpression)): 

94 m, n = other.subtype.shape 

95 

96 if m == n: 

97 svec_length = int(0.5*n*(n + 1)) 

98 

99 if subtype.dim and subtype.dim != svec_length: 

100 return NotImplemented 

101 

102 if issubclass(other.clstype, ComplexAffineExpression): 

103 # Other is already a square matrix. 

104 matrix = other 

105 else: 

106 # Predict the vector svec(other). 

107 vector = other.clstype.make_type( 

108 shape=(svec_length, 1), 

109 universe_type=other.subtype.universe_type) 

110 elif 1 in other.subtype.shape: 

111 if subtype.dim and subtype.dim != other.subtype.dim: 

112 return NotImplemented 

113 

114 if issubclass(other.clstype, ComplexAffineExpression): 

115 # Predict the square matrix desvec(other). 

116 n = 0.5*((8*other.subtype.dim + 1)**0.5 - 1) 

117 if int(n) != n: 

118 return NotImplemented 

119 n = int(n) 

120 matrix = other.clstype.make_type( 

121 shape=(n, n), constant=other.subtype.constant, 

122 nonneg=other.subtype.nonneg) 

123 else: 

124 # Other is already a vector. 

125 vector = other 

126 else: 

127 return NotImplemented 

128 

129 if issubclass(other.clstype, ComplexAffineExpression): 

130 zero = AffineExpression.make_type( 

131 shape=matrix.subtype.shape, constant=True, nonneg=True) 

132 

133 return matrix.clstype._predict( 

134 matrix.subtype, operator.__rshift__, zero) 

135 elif issubclass( 

136 other.subtype.universe_type.clstype, 

137 ScenarioPerturbationSet): 

138 dim = vector.subtype.dim 

139 count = other.subtype.universe_type.subtype.scenario_count 

140 cone = ExpressionType(cls, subtype) 

141 

142 return ScenarioUncertainConicConstraint.make_type( 

143 dim=dim, scenario_count=count, cone_type=cone) 

144 else: 

145 return NotImplemented 

146 

147 return NotImplemented 

148 

149 def _rshift_implementation(self, element): 

150 if isinstance(element, ( 

151 ComplexAffineExpression, UncertainAffineExpression)): 

152 if element.square: 

153 # Mimic _check_dimension. 

154 n = element.shape[0] 

155 d = int(0.5*n*(n + 1)) 

156 if self.dim and self.dim != d: 

157 raise TypeError( 

158 "The shape {} of {} implies a {}-dimensional " 

159 "svec-representation which does not match the fixed " 

160 "dimensionality {} of the cone {}.".format( 

161 glyphs.shape(element.shape), element.string, d, 

162 self.dim, self.string)) 

163 

164 if isinstance(element, ComplexAffineExpression): 

165 return element >> 0 

166 elif isinstance(element.universe, ScenarioPerturbationSet): 

167 return ScenarioUncertainConicConstraint(element.svec, self) 

168 else: 

169 raise TypeError("LMIs with uncertainty parameterized " 

170 "through a {} are not supported.".format( 

171 element.universe.__class__.__name__)) 

172 elif 1 in element.shape: 

173 self._check_dimension(element) 

174 

175 if isinstance(element, ComplexAffineExpression): 

176 return element.desvec >> 0 

177 elif isinstance(element.universe, ScenarioPerturbationSet): 

178 return ScenarioUncertainConicConstraint(element, self) 

179 else: 

180 raise TypeError("LMIs with uncertainty parameterized " 

181 "through a {} are not supported.".format( 

182 element.universe.__class__.__name__)) 

183 else: 

184 raise TypeError("The {} expression {} is neither square nor a " 

185 "vector so it cannot be constrained to be in the positive " 

186 "semidefinite cone.".format(element.shape, element.string)) 

187 

188 return NotImplemented 

189 

190 @property 

191 def dual_cone(self): 

192 """Implement :attr:`.cone.Cone.dual_cone`.""" 

193 return self 

194 

195 

196# -------------------------------------- 

197__all__ = api_end(_API_START, globals())