Quantum++  v1.0.0-beta2
C++11 quantum computing library
entanglement.h
Go to the documentation of this file.
1 /*
2  * Quantum++
3  *
4  * Copyright (c) 2013 - 2017 Vlad Gheorghiu (vgheorgh@gmail.com)
5  *
6  * This file is part of Quantum++.
7  *
8  * Quantum++ is free software: you can redistribute it and/or modify
9  * it under the terms of the GNU General Public License as published by
10  * the Free Software Foundation, either version 3 of the License, or
11  * (at your option) any later version.
12  *
13  * Quantum++ is distributed in the hope that it will be useful,
14  * but WITHOUT ANY WARRANTY; without even the implied warranty of
15  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  * GNU General Public License for more details.
17  *
18  * You should have received a copy of the GNU General Public License
19  * along with Quantum++. If not, see <http://www.gnu.org/licenses/>.
20  */
21 
27 #ifndef ENTANGLEMENT_H_
28 #define ENTANGLEMENT_H_
29 
30 // entanglement
31 
32 namespace qpp
33 {
34 
45 template<typename Derived>
46 dyn_col_vect<double> schmidtcoeffs(const Eigen::MatrixBase<Derived>& A,
47  const std::vector<idx>& dims)
48 {
49  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
50 
51  // EXCEPTION CHECKS
52 
53  // check zero-size
55  throw Exception("qpp::schmidtcoeffs()", Exception::Type::ZERO_SIZE);
56  // check bi-partite
57  if (dims.size() != 2)
58  throw Exception("qpp::schmidtcoeffs()", Exception::Type::NOT_BIPARTITE);
59  // check column vector
60  if (!internal::check_cvector(rA))
61  throw Exception("qpp::schmidtcoeffs()",
63  // check matching dimensions
64  if (!internal::check_dims_match_mat(dims, rA))
65  throw Exception("qpp::schmidtcoeffs()",
67  // END EXCEPTION CHECKS
68 
69  return svals(transpose(reshape(rA, dims[1], dims[0])));
70 }
71 
82 template<typename Derived>
83 dyn_col_vect<double> schmidtcoeffs(const Eigen::MatrixBase<Derived>& A,
84  idx d = 2)
85 {
86  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
87 
88  // EXCEPTION CHECKS
89 
90  // check zero size
92  throw Exception("qpp::schmidtcoeffs()", Exception::Type::ZERO_SIZE);
93 
94  // check valid dims
95  if (d == 0)
96  throw Exception("qpp::schmidtcoeffs()",
98  // END EXCEPTION CHECKS
99 
100  idx N = internal::get_num_subsys(static_cast<idx>(rA.rows()), d);
101  std::vector<idx> dims(N, d); // local dimensions vector
102 
103  return schmidtcoeffs(A, dims);
104 }
105 
114 template<typename Derived>
115 cmat schmidtA(const Eigen::MatrixBase<Derived>& A,
116  const std::vector<idx>& dims)
117 {
118  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
119 
120  // EXCEPTION CHECKS
121 
122  // check zero-size
124  throw Exception("qpp::schmidtU()", Exception::Type::ZERO_SIZE);
125  // check bi-partite
126  if (dims.size() != 2)
127  throw Exception("qpp::schmidtU()", Exception::Type::NOT_BIPARTITE);
128  // check column vector
129  if (!internal::check_cvector(rA))
130  throw Exception("qpp::schmidtU()",
132  // check matching dimensions
133  if (!internal::check_dims_match_mat(dims, rA))
134  throw Exception("qpp::schmidtU()",
136  // END EXCEPTION CHECKS
137 
138  return svdU(transpose(reshape(rA, dims[1], dims[0])));
139 }
140 
149 template<typename Derived>
150 cmat schmidtA(const Eigen::MatrixBase<Derived>& A, idx d = 2)
151 {
152  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
153 
154  // EXCEPTION CHECKS
155 
156  // check zero size
158  throw Exception("qpp::schmidtA()", Exception::Type::ZERO_SIZE);
159 
160  // check valid dims
161  if (d == 0)
162  throw Exception("qpp::schmidtA()",
164  // END EXCEPTION CHECKS
165 
166  idx N = internal::get_num_subsys(static_cast<idx>(rA.rows()), d);
167  std::vector<idx> dims(N, d); // local dimensions vector
168 
169  return schmidtA(A, dims);
170 }
171 
180 template<typename Derived>
181 cmat schmidtB(const Eigen::MatrixBase<Derived>& A,
182  const std::vector<idx>& dims)
183 {
184  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
185 
186  // EXCEPTION CHECKS
187 
188  // check zero-size
190  throw Exception("qpp::schmidtV()", Exception::Type::ZERO_SIZE);
191  // check bi-partite
192  if (dims.size() != 2)
193  throw Exception("qpp::schmidtV()", Exception::Type::NOT_BIPARTITE);
194  // check column vector
195  if (!internal::check_cvector(rA))
196  throw Exception("qpp::schmidtV()",
198  // check matching dimensions
199  if (!internal::check_dims_match_mat(dims, rA))
200  throw Exception("qpp::schmidtV()",
202  // END EXCEPTION CHECKS
203 
204  // by default returns V^*, we need V, i.e. the complex conjugate,
205  // i.e. adjoint(transpose(V))
206  return svdV(transpose(reshape(conjugate(rA), dims[1], dims[0])));
207 }
208 
217 template<typename Derived>
218 cmat schmidtB(const Eigen::MatrixBase<Derived>& A, idx d = 2)
219 {
220  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
221 
222  // EXCEPTION CHECKS
223 
224  // check zero size
226  throw Exception("qpp::schmidtB()", Exception::Type::ZERO_SIZE);
227 
228  // check valid dims
229  if (d == 0)
230  throw Exception("qpp::schmidtB()",
232  // END EXCEPTION CHECKS
233 
234  idx N = internal::get_num_subsys(static_cast<idx>(rA.rows()), d);
235  std::vector<idx> dims(N, d); // local dimensions vector
236 
237  return schmidtB(A, dims);
238 }
239 
251 template<typename Derived>
252 std::vector<double> schmidtprobs(const Eigen::MatrixBase<Derived>& A,
253  const std::vector<idx>& dims)
254 {
255  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
256 
257  // EXCEPTION CHECKS
258 
259  // check zero-size
261  throw Exception("qpp::schmidtprobs()", Exception::Type::ZERO_SIZE);
262  // check bi-partite
263  if (dims.size() != 2)
264  throw Exception("qpp::schmidtprobs()", Exception::Type::NOT_BIPARTITE);
265  // check column vector
266  if (!internal::check_cvector(rA))
267  throw Exception("qpp::schmidtprobs()",
269  // check matching dimensions
270  if (!internal::check_dims_match_mat(dims, rA))
271  throw Exception("qpp::schmidtprobs()",
273  // END EXCEPTION CHECKS
274 
275  std::vector<double> result;
276  dyn_col_vect<double> scf = schmidtcoeffs(rA, dims);
277  for (idx i = 0; i < static_cast<idx>(scf.rows()); ++i)
278  result.push_back(std::pow(scf(i), 2));
279 
280  return result;
281 }
282 
294 template<typename Derived>
295 std::vector<double> schmidtprobs(const Eigen::MatrixBase<Derived>& A, idx d = 2)
296 {
297  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
298 
299  // EXCEPTION CHECKS
300 
301  // check zero size
303  throw Exception("qpp::schmidtprobs()", Exception::Type::ZERO_SIZE);
304 
305  // check valid dims
306  if (d == 0)
307  throw Exception("qpp::schmidtprobs()",
309  // END EXCEPTION CHECKS
310 
311  idx N = internal::get_num_subsys(static_cast<idx>(rA.rows()), d);
312  std::vector<idx> dims(N, d); // local dimensions vector
313 
314  return schmidtprobs(A, dims);
315 }
316 
328 template<typename Derived>
329 double entanglement(const Eigen::MatrixBase<Derived>& A,
330  const std::vector<idx>& dims)
331 {
332  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
333 
334  // EXCEPTION CHECKS
335 
336  // check zero-size
338  throw Exception("qpp::entanglement()", Exception::Type::ZERO_SIZE);
339  // check bi-partite
340  if (dims.size() != 2)
341  throw Exception("qpp::entanglement()", Exception::Type::NOT_BIPARTITE);
342  // check column vector
343  if (!internal::check_cvector(rA))
344  throw Exception("qpp::entanglement()",
346  // check matching dimensions
347  if (!internal::check_dims_match_mat(dims, rA))
348  throw Exception("qpp::entanglement()",
350  // END EXCEPTION CHECKS
351 
352  return entropy(schmidtprobs(rA, dims));
353 }
354 
366 template<typename Derived>
367 double entanglement(const Eigen::MatrixBase<Derived>& A, idx d = 2)
368 {
369  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
370 
371  // EXCEPTION CHECKS
372 
373  // check zero size
375  throw Exception("qpp::entanglement()", Exception::Type::ZERO_SIZE);
376 
377  // check valid dims
378  if (d == 0)
379  throw Exception("qpp::entanglement()",
381  // END EXCEPTION CHECKS
382 
383  idx N = internal::get_num_subsys(static_cast<idx>(rA.rows()), d);
384  std::vector<idx> dims(N, d); // local dimensions vector
385 
386  return entanglement(A, dims);
387 }
388 
400 // the G-concurrence
401 template<typename Derived>
402 double gconcurrence(const Eigen::MatrixBase<Derived>& A)
403 {
404  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
405 
406  // EXCEPTION CHECKS
407 
408  // check zero-size
410  throw Exception("qpp::gconcurrence()", Exception::Type::ZERO_SIZE);
411  // check column vector
412  if (!internal::check_cvector(rA))
413  throw Exception("qpp::gconcurrence()",
415 
416  idx d = internal::get_dim_subsys(static_cast<idx>(rA.rows()), 2);
417 
418  // check equal local dimensions
419  if (d * d != static_cast<idx>(rA.rows()))
420  throw Exception("qpp::gconcurrence()",
422  // END EXCEPTION CHECKS
423 
424  // we compute exp(logdet()) to avoid underflow
425  return d * std::abs(std::exp(2. / d * logdet(reshape(rA, d, d))));
426 }
427 
435 template<typename Derived>
436 double negativity(const Eigen::MatrixBase<Derived>& A,
437  const std::vector<idx>& dims)
438 {
439  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
440 
441  // EXCEPTION CHECKS
442 
443  // check zero-size
445  throw Exception("qpp::negativity()", Exception::Type::ZERO_SIZE);
446  // check bi-partite
447  if (dims.size() != 2)
448  throw Exception("qpp::negativity()", Exception::Type::NOT_BIPARTITE);
449  // check square matrix vector
451  throw Exception("qpp::negativity()",
453  // check matching dimensions
454  if (!internal::check_dims_match_mat(dims, rA))
455  throw Exception("qpp::negativity()",
457  // END EXCEPTION CHECKS
458 
459  return (schatten(ptranspose(rA, {0}, dims), 1) - 1.) / 2.;
460 }
461 
469 template<typename Derived>
470 double negativity(const Eigen::MatrixBase<Derived>& A, idx d = 2)
471 {
472  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
473 
474  // EXCEPTION CHECKS
475 
476  // check zero size
478  throw Exception("qpp::negativity()", Exception::Type::ZERO_SIZE);
479 
480  // check valid dims
481  if (d == 0)
482  throw Exception("qpp::negativity()",
484  // END EXCEPTION CHECKS
485 
486  idx N = internal::get_num_subsys(static_cast<idx>(rA.rows()), d);
487  std::vector<idx> dims(N, d); // local dimensions vector
488 
489  return negativity(A, dims);
490 }
491 
499 template<typename Derived>
500 double lognegativity(const Eigen::MatrixBase<Derived>& A,
501  const std::vector<idx>& dims)
502 {
503  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
504 
505  // EXCEPTION CHECKS
506 
507  // check zero-size
509  throw Exception("qpp::lognegativity()", Exception::Type::ZERO_SIZE);
510  // check bi-partite
511  if (dims.size() != 2)
512  throw Exception("qpp::lognegativity()",
514  // check square matrix vector
516  throw Exception("qpp::lognegativity()",
518  // check matching dimensions
519  if (!internal::check_dims_match_mat(dims, rA))
520  throw Exception("qpp::lognegativity()",
522  // END EXCEPTION CHECKS
523 
524  return std::log2(2 * negativity(rA, dims) + 1);
525 }
526 
534 template<typename Derived>
535 double lognegativity(const Eigen::MatrixBase<Derived>& A, idx d = 2)
536 {
537  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
538 
539  // EXCEPTION CHECKS
540 
541  // check zero size
543  throw Exception("qpp::lognegativity()", Exception::Type::ZERO_SIZE);
544 
545  // check valid dims
546  if (d == 0)
547  throw Exception("qpp::lognegativity()",
549  // END EXCEPTION CHECKS
550 
551  idx N = internal::get_num_subsys(static_cast<idx>(rA.rows()), d);
552  std::vector<idx> dims(N, d); // local dimensions vector
553 
554  return lognegativity(A, dims);
555 }
556 
563 template<typename Derived>
564 double concurrence(const Eigen::MatrixBase<Derived>& A)
565 {
566  const dyn_mat<typename Derived::Scalar>& rA = A.derived();
567 
568  // EXCEPTION CHECKS
569 
570  // check zero-size
572  throw Exception("qpp::concurrence()", Exception::Type::ZERO_SIZE);
573  // check square matrix vector
575  throw Exception("qpp::concurrence()",
577  // check that the state is a 2-qubit state
578  if (rA.rows() != 4)
579  throw Exception("qpp::concurrence()",
581  // END EXCEPTION CHECKS
582 
583  cmat sigmaY = Gates::get_instance().Y;
584  dyn_col_vect<double> lambdas =
585  evals(rA * kron(sigmaY, sigmaY) * conjugate(rA)
586  * kron(sigmaY, sigmaY)).real();
587 
588  std::vector<double> lambdas_sorted(lambdas.data(),
589  lambdas.data() + lambdas.size());
590 
591  std::sort(std::begin(lambdas_sorted), std::end(lambdas_sorted),
592  std::greater<double>());
593  std::transform(std::begin(lambdas_sorted), std::end(lambdas_sorted),
594  std::begin(lambdas_sorted), [](double elem)
595  {
596  return std::sqrt(std::abs(elem));
597  }); // chop tiny negatives
598 
599  return std::max(0.,
600  lambdas_sorted[0] - lambdas_sorted[1] - lambdas_sorted[2]
601  - lambdas_sorted[3]);
602 }
603 
604 } /* namespace qpp */
605 
606 #endif /* ENTANGLEMENT_H_ */
bool check_nonzero_size(const T &x) noexcept
Definition: util.h:127
Derived::Scalar logdet(const Eigen::MatrixBase< Derived > &A)
Logarithm of the determinant.
Definition: functions.h:173
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:252
bool check_cvector(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:120
double gconcurrence(const Eigen::MatrixBase< Derived > &A)
G-concurrence of the bi-partite pure state A.
Definition: entanglement.h:402
double lognegativity(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Logarithmic negativity of the bi-partite mixed state A.
Definition: entanglement.h:500
Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic > dyn_mat
Dynamic Eigen matrix over the field specified by Scalar.
Definition: types.h:78
bool check_dims_match_mat(const std::vector< idx > &dims, const Eigen::MatrixBase< Derived > &A)
Definition: util.h:156
Quantum++ main namespace.
Definition: codes.h:30
idx get_num_subsys(idx sz, idx d)
Definition: util.h:370
double schatten(const Eigen::MatrixBase< Derived > &A, double p)
Schatten matrix norm.
Definition: functions.h:826
cmat svdU(const Eigen::MatrixBase< Derived > &A)
Left singular vectors.
Definition: functions.h:494
cmat schmidtA(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Schmidt basis on Alice side.
Definition: entanglement.h:115
dyn_col_vect< double > svals(const Eigen::MatrixBase< Derived > &A)
Singular values.
Definition: functions.h:470
dyn_mat< typename T::Scalar > kron(const T &head)
Kronecker product.
Definition: functions.h:895
double entropy(const Eigen::MatrixBase< Derived > &A)
von-Neumann entropy of the density matrix A
Definition: entropies.h:43
Generates custom exceptions, used when validating function parameters.
Definition: exception.h:39
dyn_mat< typename Derived::Scalar > transpose(const Eigen::MatrixBase< Derived > &A)
Transpose.
Definition: functions.h:44
dyn_col_vect< cplx > evals(const Eigen::MatrixBase< Derived > &A)
Eigenvalues.
Definition: functions.h:304
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:64
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:181
cmat svdV(const Eigen::MatrixBase< Derived > &A)
Right singular vectors.
Definition: functions.h:519
double entanglement(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Entanglement of the bi-partite pure state A.
Definition: entanglement.h:329
Eigen::Matrix< Scalar, Eigen::Dynamic, 1 > dyn_col_vect
Dynamic Eigen column vector over the field specified by Scalar.
Definition: types.h:90
static const Gates & get_instance() noexcept(std::is_nothrow_constructible< const Gates >::value)
Definition: singleton.h:90
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:1540
Eigen::MatrixXcd cmat
Complex (double precision) dynamic Eigen matrix.
Definition: types.h:61
std::size_t idx
Non-negative integer index.
Definition: types.h:36
double concurrence(const Eigen::MatrixBase< Derived > &A)
Wootters concurrence of the bi-partite qubit mixed state A.
Definition: entanglement.h:564
cmat Y
Pauli Sigma-Y gate.
Definition: gates.h:46
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:46
dyn_mat< typename Derived::Scalar > reshape(const Eigen::MatrixBase< Derived > &A, idx rows, idx cols)
Reshape.
Definition: functions.h:1135
double negativity(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Negativity of the bi-partite mixed state A.
Definition: entanglement.h:436