XACC
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 - 2020 Vlad Gheorghiu.
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
55  if (!internal::check_nonzero_size(rA))
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
90  if (!internal::check_nonzero_size(A))
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
120  if (!internal::check_nonzero_size(rA))
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
151  if (!internal::check_nonzero_size(A))
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
181  if (!internal::check_nonzero_size(rA))
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
214  if (!internal::check_nonzero_size(A))
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
248  if (!internal::check_nonzero_size(rA))
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.emplace_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
289  if (!internal::check_nonzero_size(A))
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
322  if (!internal::check_nonzero_size(rA))
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
356  if (!internal::check_nonzero_size(A))
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
389  if (!internal::check_nonzero_size(rA))
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
421  if (!internal::check_nonzero_size(rA))
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
427  if (!internal::check_square_mat(rA))
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
451  if (!internal::check_nonzero_size(A))
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
480  if (!internal::check_nonzero_size(rA))
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
486  if (!internal::check_square_mat(rA))
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
510  if (!internal::check_nonzero_size(A))
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
537  if (!internal::check_nonzero_size(rA))
538  throw exception::ZeroSize("qpp::concurrence()");
539  // check square matrix vector
540  if (!internal::check_square_mat(rA))
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_ */
Quantum++ main namespace.
Definition: circuits.h:35