Coverage for picos/expressions/exp_quantkeydist.py: 88.28%

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

4# This file is part of PICOS.

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

19"""Implements :class:`QuantumKeyDistribution`."""

21# TODO: Common base class for QuantumEntropy and NegativeQuantumEntropy.

23import math

24import operator

25from collections import namedtuple

26from functools import reduce

28import cvxopt

29import numpy

31from .. import glyphs

32from ..apidoc import api_end, api_start

33from ..caching import cached_unary_operator

34from ..constraints import (

35 QuantKeyDistributionConstraint,

36 ComplexQuantKeyDistributionConstraint,

37)

38from .data import convert_and_refine_arguments, convert_operands, cvx2np

39from .exp_affine import AffineExpression, ComplexAffineExpression

40from .expression import Expression, refine_operands, validate_prediction

42_API_START = api_start(globals())

43# -------------------------------

46class QuantumKeyDistribution(Expression):

47 r"""Slice of quantum relative entropy used to compute quantum key rates.

49 :Definition:

51 Let :math:`X` be an :math:`n \times n`-dimensional symmetric or hermitian

52 matrix. Let :math:`\mathcal{Z}` be the pinching map which maps off-diagonal

53 blocks of a given block structure to zero, i.e., for a bipartite state

55 .. math::

57 \mathcal{Z}(X) = \sum_{i} Z_i X Z_i^\dagger,

59 where :math:`Z_i= | i \rangle \langle i | \otimes \mathbb{I}` if we block-

60 diagonalize over the first subsystem, and :math:`Z_i= \mathbb{I} \otimes

61 | i \rangle \langle i |` if we block-diagonalize over the second subsystem.

62 We also generalize this definition to multipartite systems where we block-

63 diagonalize over any number of subsystems.

65 1. In general, this is the expression

67 .. math::

69 -S(\mathcal{G}(X)) + S(\mathcal{Z}(\mathcal{G}(X))),

71 where :math:`S(X)=-\operatorname{Tr}(X \log(X))` is the quantum entropy,

72 :math:`mathcal{G}` is a positive linear map given by Kraus operators

74 .. math::

76 \mathcal{G}(X) = \sum_{i} K_i X K_i^\dagger.

78 2. If ``K_list=None``, then :math:`\mathcal{G}` is assumed to be the

79 identity map, and then this expression is simplified to

81 .. math::

83 -S(X) + S(\mathcal{Z}(X)).

85 .. warning::

87 When you pose an upper bound on this expression, then PICOS enforces

88 :math:`X \succeq 0` through an auxiliary constraint during solution

89 search.

90 """

92 # --------------------------------------------------------------------------

93 # Initialization and factory methods.

94 # --------------------------------------------------------------------------

96 @convert_and_refine_arguments("X")

97 def __init__(self, X, subsystems=0, dimensions=2, K_list=None):

98 r"""Construct an :class:`QuantumKeyDistribution`.

100 :param X: The affine expression :math:`X`.

101 :type X: ~picos.expressions.AffineExpression

102

103 :param subsystems: A collection of or a single subystem number, indexed

104 from zero, corresponding to subsystems that will be block-

105 diagonalized over. The value :math:`-1` refers to the last

106 subsystem.

107 :type subsystems: int or tuple or list

108

109 :param dimensions: Either an integer :math:`d` so that the subsystems

110 are assumed to be all of shape :math:`d \times d`, or a sequence of

111 subsystem shapes where an integer :math:`d` within the sequence is

112 read as :math:`d \times d`. In any case, the elementwise product

113 over all subsystem shapes must equal the expression's shape.

114 :type dimensions: int or tuple or list

115

116 :param K_list: A list of Kraus operators representing the linear map

117 :math:`\mathcal{G}`. If ``K_list=None``, then :math:`\mathcal{G}`

118 is defined as the identity map.

119 :type K_list: None or list(numpy.ndarray)

120 """

121 # Check that X is an affine Hermitian expression

122 if not isinstance(X, ComplexAffineExpression):

123 raise TypeError(

124 "Can only take the matrix logarithm of a real "

125 "or complex affine expression, not of {}.".format(

126 type(X).__name__

127 )

128 )

129 if not X.hermitian:

130 raise TypeError(

131 "Can only take the matrix logarithm of a symmetric "

132 "or Hermitian expression, not of {}.".format(type(X).__name__)

133 )

134

135 self._X = X

136

137 self._iscomplex = (not isinstance(X, AffineExpression)) or (

138 K_list is not None

139 and any([numpy.iscomplexobj(Ki) for Ki in K_list])

140 )

141

142 # Check that subsystems and dimension are compatible with X.

143 if isinstance(dimensions, int):

144 dimensions = self._square_equal_subsystem_dims(dimensions)

145 else:

146 dimensions = [

147 (d, d) if isinstance(d, int) else d for d in dimensions

148 ]

149

150 if (

151 reduce(lambda x, y: (x[0] * y[0], x[1] * y[1]), dimensions)

152 != X.shape

153 ):

154 raise TypeError("Subsystem dimensions do not match expression.")

155

156 if isinstance(subsystems, int):

157 subsystems = (subsystems,)

158

159 numSys = len(dimensions)

160 subsystems = set(numSys - 1 if sys == -1 else sys for sys in subsystems)

161

162 for sys in subsystems:

163 if not isinstance(sys, int):

164 raise IndexError(

165 "Subsystem indices must be integer, not {}.".format(

166 type(sys).__name__

167 )

168 )

169 elif sys < 0:

170 raise IndexError("Subsystem indices must be nonnegative.")

171 elif sys >= numSys:

172 raise IndexError(

173 "Subsystem index {} out of range for {} "

174 "systems total.".format(sys, numSys)

175 )

176 elif dimensions[sys][0] != dimensions[sys][1]:

177 raise TypeError(

178 "Subsystem index {} refers to a non-square subsystem that "

179 "cannot be traced over.".format(sys)

180 )

181

182 self._subsystems = subsystems

183 self._dimensions = dimensions

184 self._K_list = K_list

185

186 self._build_Z_list()

187

188 typeStr = "Quantum Key Distribution"

189 if K_list is None:

190 zxStr = "Z(" + X.string + ")"

191 symbStr = glyphs.qre(X.string, zxStr)

192 else:

193 gxStr = "G(" + X.string + ")"

194 gzxStr = "Z(G(" + X.string + "))"

195 symbStr = glyphs.qre(gxStr, gzxStr)

196

197 Expression.__init__(self, typeStr, symbStr)

198

199 def _square_equal_subsystem_dims(self, diagLen):

200 m, n = self._X.shape

201 k = math.log(m, diagLen)

202

203 if m != n or int(k) != k:

204 raise TypeError(

205 "The expression has shape {} so it cannot be "

206 "decomposed into subsystems of shape {}.".format(

207 glyphs.shape(self._X.shape), glyphs.shape((diagLen,) * 2)

208 )

209 )

210

211 return ((diagLen,) * 2,) * int(k)

212

213 def _build_Z_list(self):

214 r = numpy.meshgrid(

215 *[range(self.dimensions[k][0]) for k in list(self.subsystems)]

216 )

217 r = list(numpy.array(r).reshape(len(self.subsystems), -1).T)

218 self._Z_list = []

219 for i in range(len(r)):

220 Z_i = numpy.array([1])

221 counter = 0

222 for k, dimk in enumerate(self.dimensions):

223 if k in self.subsystems:

224 Z_ik = numpy.zeros(dimk[0])

225 Z_ik[r[i][counter]] = 1

226 Z_i = numpy.kron(Z_i, Z_ik)

227 counter += 1

228 else:

229 Z_i = numpy.kron(Z_i, numpy.ones(dimk[0]))

230 self._Z_list += [numpy.diag(Z_i)]

231 return self._Z_list

232

233 # --------------------------------------------------------------------------

234 # Abstract method implementations and method overridings, except _predict.

235 # --------------------------------------------------------------------------

236

237 def _get_refined(self):

238 if self._X.constant:

239 return AffineExpression.from_constant(self.value, 1, self._symbStr)

240 else:

241 return self

242

243 Subtype = namedtuple("Subtype", ("argdim", "iscomplex"))

244

245 def _get_subtype(self):

246 return self.Subtype(len(self._X), self._iscomplex)

247

248 def _get_value(self):

249 X = cvx2np(self._X._get_value())

250 if self.K_list is not None:

251 X = sum([K @ X @ K.conj().T for K in self.K_list])

252 eigX = numpy.linalg.eigvalsh(X)

253 eigX = eigX[eigX > 1e-12]

254

255 ZX = sum([Z @ X @ Z.conj().T for Z in self.Z_list])

256 eigZX = numpy.linalg.eigvalsh(ZX)

257 eigZX = eigZX[eigZX > 1e-12]

258

259 s = numpy.dot(eigX, numpy.log(eigX))

260 s -= numpy.dot(eigZX, numpy.log(eigZX))

261

262 return cvxopt.matrix(s)

263

264 @cached_unary_operator

265 def _get_mutables(self):

266 return self._X._get_mutables()

267

268 def _is_convex(self):

269 return True

270

271 def _is_concave(self):

272 return False

273

274 def _replace_mutables(self, mapping):

275 return self.__class__(self._X._replace_mutables(mapping))

276

277 def _freeze_mutables(self, freeze):

278 return self.__class__(self._X._freeze_mutables(freeze))

279

280 # --------------------------------------------------------------------------

281 # Methods and properties that return expressions.

282 # --------------------------------------------------------------------------

283

284 @property

285 def X(self):

286 """The expression :math:`X`."""

287 return self._X

288

289 # --------------------------------------------------------------------------

290 # Methods and properties that describe the expression.

291 # --------------------------------------------------------------------------

292

293 @property

294 def subsystems(self):

295 """The subsystems being block-diagonalized of :math:`X`."""

296 return self._subsystems

297

298 @property

299 def dimensions(self):

300 """The dimensions of the subsystems of :math:`X`."""

301 return self._dimensions

302

303 @property

304 def K_list(self):

305 r"""The Kraus operators :math:`K_i` of :math:`\mathcal{G}`."""

306 return self._K_list

307

308 @property

309 def Z_list(self):

310 r"""The Kraus operators :math:`Z_i` of :math:`\mathcal{Z}`."""

311 return self._Z_list

312

313 @property

314 def n(self):

315 """Length of :attr:`X`."""

316 return self._X.shape[0]

317

318 @property

319 def iscomplex(self):

320 """Whether :attr:`X` is a complex expression or not."""

321 return self._iscomplex

322

323 # --------------------------------------------------------------------------

324 # Constraint-creating operators, and _predict.

325 # --------------------------------------------------------------------------

326

327 @classmethod

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

329 assert isinstance(subtype, cls.Subtype)

330

331 if relation == operator.__le__:

332 if (

333 issubclass(other.clstype, AffineExpression)

334 and other.subtype.dim == 1

335 ):

336 if subtype.iscomplex:

337 return ComplexQuantKeyDistributionConstraint.make_type(

338 argdim=subtype.argdim

339 )

340 else:

341 return QuantKeyDistributionConstraint.make_type(

342 argdim=subtype.argdim

343 )

344 return NotImplemented

345

346 @convert_operands(scalarRHS=True)

347 @validate_prediction

348 @refine_operands()

349 def __le__(self, other):

350 if isinstance(other, AffineExpression):

351 if self.iscomplex:

352 return ComplexQuantKeyDistributionConstraint(self, other)

353 else:

354 return QuantKeyDistributionConstraint(self, other)

355 else:

356 return NotImplemented

357

358

359# --------------------------------------

360__all__ = api_end(_API_START, globals())