Quantum++  v1.0-rc4
A modern C++11 quantum computing library
entanglement.h
Go to the documentation of this file.
1 /*
2  * This file is part of Quantum++.
3  *
4  * MIT License
5  *
6  * Copyright (c) 2013 - 2018 Vlad Gheorghiu (vgheorgh@gmail.com)
7  *
8  * Permission is hereby granted, free of charge, to any person obtaining a copy
9  * of this software and associated documentation files (the "Software"), to deal
10  * in the Software without restriction, including without limitation the rights
11  * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
12  * copies of the Software, and to permit persons to whom the Software is
13  * furnished to do so, subject to the following conditions:
14  *
15  * The above copyright notice and this permission notice shall be included in
16  * all copies or substantial portions of the Software.
17  *
18  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
19  * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
20  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
21  * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
22  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
23  * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
24  * SOFTWARE.
25  */
26 
32 #ifndef ENTANGLEMENT_H_
33 #define ENTANGLEMENT_H_
34 
35 namespace qpp {
47 template <typename Derived>
48 dyn_col_vect<double> schmidtcoeffs(const Eigen::MatrixBase<Derived>& A,
49  const std::vector<idx>& dims) {
50  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
51 
52  // EXCEPTION CHECKS
53 
54  // check zero-size
56  throw exception::ZeroSize("qpp::schmidtcoeffs()");
57  // check bi-partite
58  if (dims.size() != 2)
59  throw exception::NotBipartite("qpp::schmidtcoeffs()");
60  // check column vector
61  if (!internal::check_cvector(rA))
62  throw exception::MatrixNotCvector("qpp::schmidtcoeffs()");
63  // check matching dimensions
64  if (!internal::check_dims_match_cvect(dims, rA))
65  throw exception::DimsMismatchCvector("qpp::schmidtcoeffs()");
66  // END EXCEPTION CHECKS
67 
68  return svals(transpose(reshape(rA, dims[1], dims[0])));
69 }
70 
82 template <typename Derived>
83 dyn_col_vect<double> schmidtcoeffs(const Eigen::MatrixBase<Derived>& A,
84  idx d = 2) {
85  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
86 
87  // EXCEPTION CHECKS
88 
89  // check zero size
91  throw exception::ZeroSize("qpp::schmidtcoeffs()");
92 
93  // check valid dims
94  if (d < 2)
95  throw exception::DimsInvalid("qpp::schmidtcoeffs()");
96  // END EXCEPTION CHECKS
97 
98  idx N = internal::get_num_subsys(static_cast<idx>(rA.rows()), d);
99  std::vector<idx> dims(N, d); // local dimensions vector
100 
101  return schmidtcoeffs(A, dims);
102 }
103 
112 template <typename Derived>
113 cmat schmidtA(const Eigen::MatrixBase<Derived>& A,
114  const std::vector<idx>& dims) {
115  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
116 
117  // EXCEPTION CHECKS
118 
119  // check zero-size
121  throw exception::ZeroSize("qpp::schmidtA()");
122  // check bi-partite
123  if (dims.size() != 2)
124  throw exception::NotBipartite("qpp::schmidtA()");
125  // check column vector
126  if (!internal::check_cvector(rA))
127  throw exception::MatrixNotCvector("qpp::schmidtA()");
128  // check matching dimensions
129  if (!internal::check_dims_match_cvect(dims, rA))
130  throw exception::DimsMismatchCvector("qpp::schmidtA()");
131  // END EXCEPTION CHECKS
132 
133  return svdU(transpose(reshape(rA, dims[1], dims[0])));
134 }
135 
144 template <typename Derived>
145 cmat schmidtA(const Eigen::MatrixBase<Derived>& A, idx d = 2) {
146  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
147 
148  // EXCEPTION CHECKS
149 
150  // check zero size
152  throw exception::ZeroSize("qpp::schmidtA()");
153 
154  // check valid dims
155  if (d < 2)
156  throw exception::DimsInvalid("qpp::schmidtA()");
157  // END EXCEPTION CHECKS
158 
159  idx N = internal::get_num_subsys(static_cast<idx>(rA.rows()), d);
160  std::vector<idx> dims(N, d); // local dimensions vector
161 
162  return schmidtA(A, dims);
163 }
164 
173 template <typename Derived>
174 cmat schmidtB(const Eigen::MatrixBase<Derived>& A,
175  const std::vector<idx>& dims) {
176  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
177 
178  // EXCEPTION CHECKS
179 
180  // check zero-size
182  throw exception::ZeroSize("qpp::schmidtB()");
183  // check bi-partite
184  if (dims.size() != 2)
185  throw exception::NotBipartite("qpp::schmidtB()");
186  // check column vector
187  if (!internal::check_cvector(rA))
188  throw exception::MatrixNotCvector("qpp::schmidtB()");
189  // check matching dimensions
190  if (!internal::check_dims_match_cvect(dims, rA))
191  throw exception::DimsMismatchCvector("qpp::schmidtB()");
192  // END EXCEPTION CHECKS
193 
194  // by default returns U_B^*, we need U_B, i.e. the complex conjugate,
195  // i.e. adjoint(transpose(U_B))
196  return svdV(transpose(reshape(conjugate(rA), dims[1], dims[0])));
197 }
198 
207 template <typename Derived>
208 cmat schmidtB(const Eigen::MatrixBase<Derived>& A, idx d = 2) {
209  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
210 
211  // EXCEPTION CHECKS
212 
213  // check zero size
215  throw exception::ZeroSize("qpp::schmidtB()");
216 
217  // check valid dims
218  if (d < 2)
219  throw exception::DimsInvalid("qpp::schmidtB()");
220  // END EXCEPTION CHECKS
221 
222  idx N = internal::get_num_subsys(static_cast<idx>(rA.rows()), d);
223  std::vector<idx> dims(N, d); // local dimensions vector
224 
225  return schmidtB(A, dims);
226 }
227 
240 template <typename Derived>
241 std::vector<double> schmidtprobs(const Eigen::MatrixBase<Derived>& A,
242  const std::vector<idx>& dims) {
243  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
244 
245  // EXCEPTION CHECKS
246 
247  // check zero-size
249  throw exception::ZeroSize("qpp::schmidtprobs()");
250  // check bi-partite
251  if (dims.size() != 2)
252  throw exception::NotBipartite("qpp::schmidtprobs()");
253  // check column vector
254  if (!internal::check_cvector(rA))
255  throw exception::MatrixNotCvector("qpp::schmidtprobs()");
256  // check matching dimensions
257  if (!internal::check_dims_match_cvect(dims, rA))
258  throw exception::DimsMismatchCvector("qpp::schmidtprobs()");
259  // END EXCEPTION CHECKS
260 
261  std::vector<double> result;
262  dyn_col_vect<double> scf = schmidtcoeffs(rA, dims);
263  for (idx i = 0; i < static_cast<idx>(scf.rows()); ++i)
264  result.push_back(std::pow(scf(i), 2));
265 
266  return result;
267 }
268 
281 template <typename Derived>
282 std::vector<double> schmidtprobs(const Eigen::MatrixBase<Derived>& A,
283  idx d = 2) {
284  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
285 
286  // EXCEPTION CHECKS
287 
288  // check zero size
290  throw exception::ZeroSize("qpp::schmidtprobs()");
291 
292  // check valid dims
293  if (d < 2)
294  throw exception::DimsInvalid("qpp::schmidtprobs()");
295  // END EXCEPTION CHECKS
296 
297  idx N = internal::get_num_subsys(static_cast<idx>(rA.rows()), d);
298  std::vector<idx> dims(N, d); // local dimensions vector
299 
300  return schmidtprobs(A, dims);
301 }
302 
314 template <typename Derived>
315 double entanglement(const Eigen::MatrixBase<Derived>& A,
316  const std::vector<idx>& dims) {
317  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
318 
319  // EXCEPTION CHECKS
320 
321  // check zero-size
323  throw exception::ZeroSize("qpp::entanglement()");
324  // check bi-partite
325  if (dims.size() != 2)
326  throw exception::NotBipartite("qpp::entanglement()");
327  // check column vector
328  if (!internal::check_cvector(rA))
329  throw exception::MatrixNotCvector("qpp::entanglement()");
330  // check matching dimensions
331  if (!internal::check_dims_match_cvect(dims, rA))
332  throw exception::DimsMismatchCvector("qpp::entanglement()");
333  // END EXCEPTION CHECKS
334 
335  return entropy(schmidtprobs(rA, dims));
336 }
337 
349 template <typename Derived>
350 double entanglement(const Eigen::MatrixBase<Derived>& A, idx d = 2) {
351  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
352 
353  // EXCEPTION CHECKS
354 
355  // check zero size
357  throw exception::ZeroSize("qpp::entanglement()");
358 
359  // check valid dims
360  if (d < 2)
361  throw exception::DimsInvalid("qpp::entanglement()");
362  // END EXCEPTION CHECKS
363 
364  idx N = internal::get_num_subsys(static_cast<idx>(rA.rows()), d);
365  std::vector<idx> dims(N, d); // local dimensions vector
366 
367  return entanglement(A, dims);
368 }
369 
381 // the G-concurrence
382 template <typename Derived>
383 double gconcurrence(const Eigen::MatrixBase<Derived>& A) {
384  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
385 
386  // EXCEPTION CHECKS
387 
388  // check zero-size
390  throw exception::ZeroSize("qpp::gconcurrence()");
391  // check column vector
392  if (!internal::check_cvector(rA))
393  throw exception::MatrixNotCvector("qpp::gconcurrence()");
394 
395  idx d = internal::get_dim_subsys(static_cast<idx>(rA.rows()), 2);
396 
397  // check equal local dimensions
398  if (d * d != static_cast<idx>(rA.rows()))
399  throw exception::DimsNotEqual("qpp::gconcurrence()");
400  // END EXCEPTION CHECKS
401 
402  // we compute exp(logdet()) to avoid underflow
403  return d * std::abs(std::exp(2. / d * logdet(reshape(rA, d, d))));
404 }
405 
413 template <typename Derived>
414 double negativity(const Eigen::MatrixBase<Derived>& A,
415  const std::vector<idx>& dims) {
416  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
417 
418  // EXCEPTION CHECKS
419 
420  // check zero-size
422  throw exception::ZeroSize("qpp::negativity()");
423  // check bi-partite
424  if (dims.size() != 2)
425  throw exception::NotBipartite("qpp::negativity()");
426  // check square matrix vector
428  throw exception::MatrixNotSquare("qpp::negativity()");
429  // check matching dimensions
430  if (!internal::check_dims_match_mat(dims, rA))
431  throw exception::DimsMismatchMatrix("qpp::negativity()");
432  // END EXCEPTION CHECKS
433 
434  return (schatten(ptranspose(rA, {0}, dims), 1) - 1.) / 2.;
435 }
436 
444 template <typename Derived>
445 double negativity(const Eigen::MatrixBase<Derived>& A, idx d = 2) {
446  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
447 
448  // EXCEPTION CHECKS
449 
450  // check zero size
452  throw exception::ZeroSize("qpp::negativity()");
453 
454  // check valid dims
455  if (d < 2)
456  throw exception::DimsInvalid("qpp::negativity()");
457  // END EXCEPTION CHECKS
458 
459  idx N = internal::get_num_subsys(static_cast<idx>(rA.rows()), d);
460  std::vector<idx> dims(N, d); // local dimensions vector
461 
462  return negativity(A, dims);
463 }
464 
472 template <typename Derived>
473 double lognegativity(const Eigen::MatrixBase<Derived>& A,
474  const std::vector<idx>& dims) {
475  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
476 
477  // EXCEPTION CHECKS
478 
479  // check zero-size
481  throw exception::ZeroSize("qpp::lognegativity()");
482  // check bi-partite
483  if (dims.size() != 2)
484  throw exception::NotBipartite("qpp::lognegativity()");
485  // check square matrix vector
487  throw exception::MatrixNotSquare("qpp::lognegativity()");
488  // check matching dimensions
489  if (!internal::check_dims_match_mat(dims, rA))
490  throw exception::DimsMismatchMatrix("qpp::lognegativity()");
491  // END EXCEPTION CHECKS
492 
493  return std::log2(2 * negativity(rA, dims) + 1);
494 }
495 
503 template <typename Derived>
504 double lognegativity(const Eigen::MatrixBase<Derived>& A, idx d = 2) {
505  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
506 
507  // EXCEPTION CHECKS
508 
509  // check zero size
511  throw exception::ZeroSize("qpp::lognegativity()");
512 
513  // check valid dims
514  if (d < 2)
515  throw exception::DimsInvalid("qpp::lognegativity()");
516  // END EXCEPTION CHECKS
517 
518  idx N = internal::get_num_subsys(static_cast<idx>(rA.rows()), d);
519  std::vector<idx> dims(N, d); // local dimensions vector
520 
521  return lognegativity(A, dims);
522 }
523 
530 template <typename Derived>
531 double concurrence(const Eigen::MatrixBase<Derived>& A) {
532  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
533 
534  // EXCEPTION CHECKS
535 
536  // check zero-size
538  throw exception::ZeroSize("qpp::concurrence()");
539  // check square matrix vector
541  throw exception::MatrixNotSquare("qpp::concurrence()");
542  // check that the state is a 2-qubit state
543  if (rA.rows() != 4)
544  throw exception::NotQubitSubsys("qpp::concurrence()");
545  // END EXCEPTION CHECKS
546 
547  cmat sigmaY = Gates::get_instance().Y;
548  dyn_col_vect<double> lambdas =
549  evals(rA * kron(sigmaY, sigmaY) * conjugate(rA) * kron(sigmaY, sigmaY))
550  .real();
551 
552  std::vector<double> lambdas_sorted(lambdas.data(),
553  lambdas.data() + lambdas.size());
554 
555  std::sort(std::begin(lambdas_sorted), std::end(lambdas_sorted),
556  std::greater<double>());
557  std::transform(std::begin(lambdas_sorted), std::end(lambdas_sorted),
558  std::begin(lambdas_sorted), [](double elem) {
559  return std::sqrt(std::abs(elem));
560  }); // chop tiny negatives
561 
562  return std::max(0., lambdas_sorted[0] - lambdas_sorted[1] -
563  lambdas_sorted[2] - lambdas_sorted[3]);
564 }
565 
566 } /* namespace qpp */
567 
568 #endif /* ENTANGLEMENT_H_ */
bool check_nonzero_size(const T &x) noexcept
Definition: util.h:123
Dimensions not equal exception.
Definition: exception.h:282
Dimension(s) mismatch matrix size exception.
Definition: exception.h:298
Derived::Scalar logdet(const Eigen::MatrixBase< Derived > &A)
Logarithm of the determinant.
Definition: functions.h:168
std::vector< double > schmidtprobs(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Schmidt probabilities of the bi-partite pure state A.
Definition: entanglement.h:241
bool check_cvector(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:117
double gconcurrence(const Eigen::MatrixBase< Derived > &A)
G-concurrence of the bi-partite pure state A.
Definition: entanglement.h:383
double lognegativity(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Logarithmic negativity of the bi-partite mixed state A.
Definition: entanglement.h:473
Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic > dyn_mat
Dynamic Eigen matrix over the field specified by Scalar.
Definition: types.h:81
bool check_dims_match_mat(const std::vector< idx > &dims, const Eigen::MatrixBase< Derived > &A)
Definition: util.h:150
Matrix is not a column vector exception.
Definition: exception.h:162
Quantum++ main namespace.
Definition: codes.h:35
idx get_num_subsys(idx sz, idx d)
Definition: util.h:366
Invalid dimension(s) exception.
Definition: exception.h:267
Subsystems are not qubits exception.
Definition: exception.h:469
double schatten(const Eigen::MatrixBase< Derived > &A, double p)
Schatten matrix norm.
Definition: functions.h:796
bool check_dims_match_cvect(const std::vector< idx > &dims, const Eigen::MatrixBase< Derived > &A)
Definition: util.h:165
Not bi-partite exception.
Definition: exception.h:484
cmat svdU(const Eigen::MatrixBase< Derived > &A)
Left singular vectors.
Definition: functions.h:478
cmat schmidtA(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Schmidt basis on Alice side.
Definition: entanglement.h:113
dyn_col_vect< double > svals(const Eigen::MatrixBase< Derived > &A)
Singular values.
Definition: functions.h:455
dyn_mat< typename T::Scalar > kron(const T &head)
Kronecker product.
Definition: functions.h:868
double entropy(const Eigen::MatrixBase< Derived > &A)
von-Neumann entropy of the density matrix A
Definition: entropies.h:43
Dimension(s) mismatch column vector size exception.
Definition: exception.h:314
dyn_mat< typename Derived::Scalar > transpose(const Eigen::MatrixBase< Derived > &A)
Transpose.
Definition: functions.h:46
dyn_col_vect< cplx > evals(const Eigen::MatrixBase< Derived > &A)
Eigenvalues.
Definition: functions.h:294
idx get_dim_subsys(idx sz, idx N)
Definition: util.h:378
dyn_mat< typename Derived::Scalar > conjugate(const Eigen::MatrixBase< Derived > &A)
Complex conjugate.
Definition: functions.h:65
bool check_square_mat(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:99
cmat schmidtB(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Schmidt basis on Bob side.
Definition: entanglement.h:174
cmat svdV(const Eigen::MatrixBase< Derived > &A)
Right singular vectors.
Definition: functions.h:502
double entanglement(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Entanglement of the bi-partite pure state A.
Definition: entanglement.h:315
Eigen::Matrix< Scalar, Eigen::Dynamic, 1 > dyn_col_vect
Dynamic Eigen column vector over the field specified by Scalar.
Definition: types.h:93
static const Gates & get_instance() noexcept(std::is_nothrow_constructible< const Gates >::value)
Definition: singleton.h:92
dyn_mat< typename Derived::Scalar > ptranspose(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &subsys, const std::vector< idx > &dims)
Partial transpose.
Definition: operations.h:1414
Eigen::MatrixXcd cmat
Complex (double precision) dynamic Eigen matrix.
Definition: types.h:64
std::size_t idx
Non-negative integer index.
Definition: types.h:39
Matrix is not square exception.
Definition: exception.h:147
double concurrence(const Eigen::MatrixBase< Derived > &A)
Wootters concurrence of the bi-partite qubit mixed state A.
Definition: entanglement.h:531
cmat Y
Pauli Sigma-Y gate.
Definition: gates.h:49
dyn_col_vect< double > schmidtcoeffs(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Schmidt coefficients of the bi-partite pure state A.
Definition: entanglement.h:48
dyn_mat< typename Derived::Scalar > reshape(const Eigen::MatrixBase< Derived > &A, idx rows, idx cols)
Reshape.
Definition: functions.h:1095
double negativity(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Negativity of the bi-partite mixed state A.
Definition: entanglement.h:414
Object has zero size exception.
Definition: exception.h:132