31 #if (__GNUC__ && !__clang__)
32 #pragma GCC diagnostic ignored "-Wunused-but-set-variable"
53 template<
typename Derived1,
typename Derived2>
55 const Eigen::MatrixBase<Derived1>& state,
56 const Eigen::MatrixBase<Derived2>& A,
57 const std::vector<idx>& ctrl,
58 const std::vector<idx>& subsys,
59 const std::vector<idx>& dims)
67 if (!std::is_same<
typename Derived1::Scalar,
68 typename Derived2::Scalar>::value)
85 idx d = ctrl.size() > 0 ? dims[ctrl[0]] : 1;
86 for (
idx i = 1; i < ctrl.size(); ++i)
87 if (dims[ctrl[i]] != d)
101 std::vector<idx> subsys_dims(subsys.size());
102 for (
idx i = 0; i < subsys.size(); ++i)
103 subsys_dims[i] = dims[subsys[i]];
108 std::vector<idx> ctrlgate = ctrl;
109 ctrlgate.insert(std::end(ctrlgate), std::begin(subsys), std::end(subsys));
110 std::sort(std::begin(ctrlgate), std::end(ctrlgate));
120 std::vector<dyn_mat<typename Derived1::Scalar>> Ai;
121 std::vector<dyn_mat<typename Derived1::Scalar>> Aidagger;
122 for (
idx i = 0; i < std::max(d, static_cast<idx>(2)); ++i)
124 Ai.push_back(
powm(rA, i));
128 idx D =
static_cast<idx>(rstate.rows());
130 idx ctrlsize = ctrl.size();
131 idx ctrlgatesize = ctrlgate.size();
132 idx subsyssize = subsys.size();
134 idx Dctrl =
static_cast<idx>(std::llround(std::pow(d, ctrlsize)));
135 idx DA =
static_cast<idx>(rA.rows());
143 std::vector<idx> ctrlgatebar =
complement(ctrlgate, n);
145 idx ctrlgatebarsize = ctrlgatebar.size();
148 for (
idx i = 0; i < ctrlgatebarsize; ++i)
149 DCTRLAbar *= dims[ctrlgatebar[i]];
151 for (
idx k = 0; k < n; ++k)
153 for (
idx k = 0; k < subsyssize; ++k)
154 CdimsA[k] = dims[subsys[k]];
155 for (
idx k = 0; k < ctrlsize; ++k)
157 for (
idx k = 0; k < ctrlgatebarsize; ++k)
158 CdimsCTRLAbar[k] = dims[ctrlgatebar[k]];
162 auto coeff_idx_ket = [=](
idx _i,
idx _m,
idx _r) noexcept
163 -> std::pair<typename Derived1::Scalar, idx>
166 typename Derived1::Scalar coeff = 0;
175 for (
idx k = 0; k < ctrlsize; ++k)
182 CdimsCTRLAbar, CmidxCTRLAbar);
183 for (
idx k = 0; k < n - ctrlgatesize; ++k)
185 Cmidx[ctrlgatebar[k]] = CmidxCTRLAbar[k];
190 for (
idx k = 0; k < subsyssize; ++k)
192 Cmidx[subsys[k]] = CmidxA[k];
199 for (
idx _n = 0; _n < DA; ++_n)
202 for (
idx k = 0; k < subsyssize; ++k)
204 Cmidx[subsys[k]] = CmidxA[k];
206 coeff += Ai[_i](_m, _n) *
210 return std::make_pair(coeff, indx);
216 auto coeff_idx_rho = [&](
idx _i1,
idx _m1,
218 -> std::tuple<typename Derived1::Scalar, idx, idx>
222 typename Derived1::Scalar coeff = 0, lhs = 1, rhs = 1;
237 CdimsCTRL, CmidxCTRLrow);
239 CdimsCTRL, CmidxCTRLcol);
241 for (
idx k = 0; k < ctrlsize; ++k)
243 Cmidxrow[ctrl[k]] = CmidxCTRLrow[k];
244 Cmidxcol[ctrl[k]] = CmidxCTRLcol[k];
249 CdimsCTRLAbar, CmidxCTRLAbarrow);
251 CdimsCTRLAbar, CmidxCTRLAbarcol);
252 for (
idx k = 0; k < n - ctrlgatesize; ++k)
254 Cmidxrow[ctrlgatebar[k]] = CmidxCTRLAbarrow[k];
255 Cmidxcol[ctrlgatebar[k]] = CmidxCTRLAbarcol[k];
261 for (
idx k = 0; k < subsys.size(); ++k)
263 Cmidxrow[subsys[k]] = CmidxArow[k];
264 Cmidxcol[subsys[k]] = CmidxAcol[k];
272 bool all_ctrl_rows_equal =
true;
273 bool all_ctrl_cols_equal =
true;
275 idx first_ctrl_row, first_ctrl_col;
278 first_ctrl_row = CmidxCTRLrow[0];
279 first_ctrl_col = CmidxCTRLcol[0];
283 first_ctrl_row = first_ctrl_col = 1;
286 for (
idx k = 1; k < ctrlsize; ++k)
288 if (CmidxCTRLrow[k] != first_ctrl_row)
290 all_ctrl_rows_equal =
false;
294 for (
idx k = 1; k < ctrlsize; ++k)
296 if (CmidxCTRLcol[k] != first_ctrl_col)
298 all_ctrl_cols_equal =
false;
304 for (
idx _n1 = 0; _n1 < DA; ++_n1)
307 for (
idx k = 0; k < subsyssize; ++k)
309 Cmidxrow[subsys[k]] = CmidxArow[k];
313 if (all_ctrl_rows_equal)
315 lhs = Ai[first_ctrl_row](_m1, _n1);
319 lhs = (_m1 == _n1) ? 1 : 0;
322 for (
idx _n2 = 0; _n2 < DA; ++_n2)
325 for (
idx k = 0; k < subsyssize; ++k)
327 Cmidxcol[subsys[k]] = CmidxAcol[k];
330 if (all_ctrl_cols_equal)
332 rhs = Aidagger[first_ctrl_col](_n2, _m2);
336 rhs = (_n2 == _m2) ? 1 : 0;
341 coeff += lhs * rstate(idxrowtmp, idxcoltmp) * rhs;
345 return std::make_tuple(coeff, idxrow, idxcol);
360 #pragma omp parallel for collapse(2)
361 for (
idx m = 0; m < DA; ++m)
362 for (
idx r = 0; r < DCTRLAbar; ++r)
366 result(coeff_idx_ket(1, m, r).second) =
367 coeff_idx_ket(1, m, r).first;
370 for (
idx i = 0; i < d; ++i)
372 result(coeff_idx_ket(i, m, r).second) =
373 coeff_idx_ket(i, m, r).first;
392 #pragma omp parallel for collapse(4)
393 for (
idx m1 = 0; m1 < DA; ++m1)
394 for (
idx r1 = 0; r1 < DCTRLAbar; ++r1)
395 for (
idx m2 = 0; m2 < DA; ++m2)
396 for (
idx r2 = 0; r2 < DCTRLAbar; ++r2)
399 auto coeff_idxes = coeff_idx_rho(1, m1, r1,
401 result(std::get<1>(coeff_idxes),
402 std::get<2>(coeff_idxes)) =
403 std::get<0>(coeff_idxes);
407 for (
idx i1 = 0; i1 < Dctrl; ++i1)
408 for (
idx i2 = 0; i2 < Dctrl; ++i2)
410 auto coeff_idxes = coeff_idx_rho(
413 result(std::get<1>(coeff_idxes),
414 std::get<2>(coeff_idxes)) =
415 std::get<0>(coeff_idxes);
442 template<
typename Derived1,
typename Derived2>
444 const Eigen::MatrixBase<Derived1>& state,
445 const Eigen::MatrixBase<Derived2>& A,
446 const std::vector<idx>& ctrl,
447 const std::vector<idx>& subsys,
465 static_cast<idx>(std::llround(std::log2(rstate.rows()) /
467 std::vector<idx> dims(n, d);
469 return applyCTRL(rstate, rA, ctrl, subsys, dims);
485 template<
typename Derived1,
typename Derived2>
487 const Eigen::MatrixBase<Derived1>& state,
488 const Eigen::MatrixBase<Derived2>& A,
489 const std::vector<idx>& subsys,
490 const std::vector<idx>& dims)
498 if (!std::is_same<
typename Derived1::Scalar,
499 typename Derived2::Scalar>::value)
523 std::vector<idx> subsys_dims(subsys.size());
524 for (
idx i = 0; i < subsys.size(); ++i)
525 subsys_dims[i] = dims[subsys[i]];
539 return applyCTRL(rstate, rA, {}, subsys, dims);
550 return applyCTRL(rstate, rA, {}, subsys, dims);
571 template<
typename Derived1,
typename Derived2>
573 const Eigen::MatrixBase<Derived1>& state,
574 const Eigen::MatrixBase<Derived2>& A,
575 const std::vector<idx>& subsys,
593 static_cast<idx>(std::llround(std::log2(rstate.rows()) /
595 std::vector<idx> dims(n, d);
597 return apply(rstate, rA, subsys, dims);
608 template<
typename Derived>
610 const std::vector<cmat>& Ks)
612 const cmat& rrho = rho;
624 if (Ks[0].rows() != rrho.rows())
628 if (it.rows() != Ks[0].rows() || it.cols() != Ks[0].rows())
632 cmat result = cmat::Zero(rrho.rows(), rrho.rows());
634 #pragma omp parallel for
635 for (
idx i = 0; i < Ks.size(); ++i)
639 result += Ks[i] * rrho *
adjoint(Ks[i]);
656 template<
typename Derived>
658 const std::vector<cmat>& Ks,
659 const std::vector<idx>& subsys,
660 const std::vector<idx>& dims)
662 const cmat& rrho = rho;
688 std::vector<idx> subsys_dims(subsys.size());
689 for (
idx i = 0; i < subsys.size(); ++i)
690 subsys_dims[i] = dims[subsys[i]];
701 if (it.rows() != Ks[0].rows() || it.cols() != Ks[0].rows())
705 cmat result = cmat::Zero(rrho.rows(), rrho.rows());
707 for (
idx i = 0; i < Ks.size(); ++i)
708 result +=
apply(rrho, Ks[i], subsys, dims);
723 template<
typename Derived>
725 const std::vector<cmat>& Ks,
726 const std::vector<idx>& subsys,
729 const cmat& rrho = rho;
743 static_cast<idx>(std::llround(std::log2(rrho.rows()) /
745 std::vector<idx> dims(n, d);
747 return apply(rrho, Ks, subsys, dims);
773 if (it.rows() != Ks[0].rows() || it.cols() != Ks[0].rows())
778 idx D =
static_cast<idx>(Ks[0].rows());
780 cmat result(D * D, D * D);
781 cmat MN = cmat::Zero(D, D);
782 bra A = bra::Zero(D);
783 ket B = ket::Zero(D);
784 cmat EMN = cmat::Zero(D, D);
786 #pragma omp parallel for collapse(2)
787 for (
idx m = 0; m < D; ++m)
789 for (
idx n = 0; n < D; ++n)
795 for (
idx i = 0; i < Ks.size(); ++i)
796 EMN += Ks[i] * MN *
adjoint(Ks[i]);
799 for (
idx a = 0; a < D; ++a)
802 for (
idx b = 0; b < D; ++b)
806 result(a * D + b, m * D + n) =
807 static_cast<cmat>(A * EMN * B).value();
812 EMN = cmat::Zero(D, D);
847 if (it.rows() != Ks[0].rows() || it.cols() != Ks[0].rows())
852 idx D =
static_cast<idx>(Ks[0].rows());
856 cmat MES = cmat::Zero(D * D, 1);
857 for (
idx a = 0; a < D; ++a)
862 cmat result = cmat::Zero(D * D, D * D);
864 #pragma omp parallel for
865 for (
idx i = 0; i < Ks.size(); ++i)
869 result +=
kron(cmat::Identity(D, D), Ks[i]) * Omega
899 idx D =
static_cast<idx>(std::llround(
900 std::sqrt(static_cast<double>(A.rows()))));
901 if (D * D != static_cast<idx>(A.rows()))
907 std::vector<cmat> result;
909 for (
idx i = 0; i < D * D; ++i)
911 if (std::abs(ev(i)) >
eps)
913 std::sqrt(std::abs(ev(i))) *
reshape(evec.col(i), D, D));
935 idx D =
static_cast<idx>(std::llround(
936 std::sqrt(static_cast<double>(A.rows()))));
937 if (D * D != static_cast<idx>(A.rows()))
941 cmat result(D * D, D * D);
943 #pragma omp parallel for collapse(4)
944 for (
idx a = 0; a < D; ++a)
945 for (
idx b = 0; b < D; ++b)
946 for (
idx m = 0; m < D; ++m)
947 for (
idx n = 0; n < D; ++n)
948 result(a * D + b, m * D + n) = A(m * D + a, n * D + b);
969 idx D =
static_cast<idx>(std::llround(
970 std::sqrt(static_cast<double>(A.rows()))));
971 if (D * D != static_cast<idx>(A.rows()))
975 cmat result(D * D, D * D);
977 #pragma omp parallel for collapse(4)
978 for (
idx a = 0; a < D; ++a)
979 for (
idx b = 0; b < D; ++b)
980 for (
idx m = 0; m < D; ++m)
981 for (
idx n = 0; n < D; ++n)
982 result(m * D + a, n * D + b) = A(a * D + b, m * D + n);
1001 template<
typename Derived>
1003 const std::vector<idx>& dims)
1018 if (dims.size() != 2)
1036 auto worker = [=](
idx i,
idx j) noexcept
1037 ->
typename Derived::Scalar
1039 typename Derived::Scalar
sum = 0;
1040 for (
idx m = 0; m < DA; ++m)
1041 sum += rA(m * DB + i) * std::conj(rA(m * DB + j));
1046 #pragma omp parallel for collapse(2)
1047 for (
idx j = 0; j < DB; ++j)
1048 for (
idx i = 0; i < DB; ++i)
1049 result(i, j) = worker(i, j);
1061 auto worker = [=](
idx i,
idx j) noexcept
1062 ->
typename Derived::Scalar
1064 typename Derived::Scalar
sum = 0;
1065 for (
idx m = 0; m < DA; ++m)
1066 sum += rA(m * DB + i, m * DB + j);
1071 #pragma omp parallel for collapse(2)
1072 for (
idx j = 0; j < DB; ++j)
1073 for (
idx i = 0; i < DB; ++i)
1074 result(i, j) = worker(i, j);
1098 template<
typename Derived>
1100 const std::vector<idx>& dims)
1115 if (dims.size() != 2)
1133 auto worker = [=](
idx i,
idx j) noexcept
1134 ->
typename Derived::Scalar
1136 typename Derived::Scalar
sum = 0;
1137 for (
idx m = 0; m < DB; ++m)
1138 sum += rA(i * DB + m) * std::conj(rA(j * DB + m));
1143 #pragma omp parallel for collapse(2)
1144 for (
idx j = 0; j < DA; ++j)
1145 for (
idx i = 0; i < DA; ++i)
1146 result(i, j) = worker(i, j);
1158 #pragma omp parallel for collapse(2)
1159 for (
idx j = 0; j < DA; ++j)
1160 for (
idx i = 0; i < DA; ++i)
1161 result(i, j) =
trace(rA.block(i * DB, j * DB, DB, DB));
1185 template<
typename Derived>
1187 const std::vector<idx>& subsys,
1188 const std::vector<idx>& dims)
1208 idx D =
static_cast<idx>(rA.rows());
1209 idx n = dims.size();
1210 idx nsubsys = subsys.size();
1211 idx nsubsysbar = n - nsubsys;
1213 for (
idx i = 0; i < nsubsys; ++i)
1214 Dsubsys *= dims[subsys[i]];
1215 idx Dsubsysbar = D / Dsubsys;
1225 std::vector<idx> subsys_bar =
complement(subsys, n);
1226 std::copy(std::begin(subsys_bar), std::end(subsys_bar),
1227 std::begin(Csubsysbar));
1229 for (
idx i = 0; i < n; ++i)
1233 for (
idx i = 0; i < nsubsys; ++i)
1235 Csubsys[i] = subsys[i];
1236 Cdimssubsys[i] = dims[subsys[i]];
1238 for (
idx i = 0; i < nsubsysbar; ++i)
1240 Cdimssubsysbar[i] = dims[subsys_bar[i]];
1254 if (subsys.size() == dims.size())
1256 result(0, 0) = (
adjoint(rA) * rA).value();
1260 if (subsys.size() == 0)
1263 auto worker = [=, &Cmidxcolsubsysbar](
idx i) noexcept
1264 ->
typename Derived::Scalar
1275 Cdimssubsysbar, Cmidxrowsubsysbar);
1277 for (
idx k = 0; k < nsubsysbar; ++k)
1279 Cmidxrow[Csubsysbar[k]] = Cmidxrowsubsysbar[k];
1280 Cmidxcol[Csubsysbar[k]] = Cmidxcolsubsysbar[k];
1282 typename Derived::Scalar sm = 0;
1283 for (
idx a = 0; a < Dsubsys; ++a)
1288 for (
idx k = 0; k < nsubsys; ++k)
1289 Cmidxrow[Csubsys[k]] = Cmidxcol[Csubsys[k]]
1301 for (
idx j = 0; j < Dsubsysbar; ++j)
1305 Cdimssubsysbar, Cmidxcolsubsysbar);
1306 #pragma omp parallel for
1307 for (
idx i = 0; i < Dsubsysbar; ++i)
1309 result(i, j) = worker(i);
1323 if (subsys.size() == dims.size())
1325 result(0, 0) = rA.trace();
1329 if (subsys.size() == 0)
1332 auto worker = [=, &Cmidxcolsubsysbar](
idx i) noexcept
1333 ->
typename Derived::Scalar
1344 Cdimssubsysbar, Cmidxrowsubsysbar);
1346 for (
idx k = 0; k < nsubsysbar; ++k)
1348 Cmidxrow[Csubsysbar[k]] = Cmidxrowsubsysbar[k];
1349 Cmidxcol[Csubsysbar[k]] = Cmidxcolsubsysbar[k];
1351 typename Derived::Scalar sm = 0;
1352 for (
idx a = 0; a < Dsubsys; ++a)
1357 for (
idx k = 0; k < nsubsys; ++k)
1358 Cmidxrow[Csubsys[k]] = Cmidxcol[Csubsys[k]]
1369 for (
idx j = 0; j < Dsubsysbar; ++j)
1373 Cdimssubsysbar, Cmidxcolsubsysbar);
1374 #pragma omp parallel for
1375 for (
idx i = 0; i < Dsubsysbar; ++i)
1377 result(i, j) = worker(i);
1403 template<
typename Derived>
1405 const std::vector<idx>& subsys,
1422 static_cast<idx>(std::llround(std::log2(rA.rows()) /
1424 std::vector<idx> dims(n, d);
1426 return ptrace(rA, subsys, dims);
1442 template<
typename Derived>
1444 const Eigen::MatrixBase<Derived>& A,
1445 const std::vector<idx>& subsys,
1446 const std::vector<idx>& dims)
1466 idx D =
static_cast<idx>(rA.rows());
1467 idx numdims = dims.size();
1468 idx numsubsys = subsys.size();
1474 for (
idx i = 0; i < numdims; ++i)
1476 for (
idx i = 0; i < numsubsys; ++i)
1477 Csubsys[i] = subsys[i];
1489 if (subsys.size() == dims.size())
1492 if (subsys.size() == 0)
1495 auto worker = [=, &Cmidxcol](
idx i) noexcept
1496 ->
typename Derived::Scalar
1502 for (
idx k = 0; k < numdims; ++k)
1503 midxcoltmp[k] = Cmidxcol[k];
1508 for (
idx k = 0; k < numsubsys; ++k)
1509 std::swap(midxcoltmp[Csubsys[k]], midxrow[Csubsys[k]]);
1517 for (
idx j = 0; j < D; ++j)
1521 #pragma omp parallel for
1522 for (
idx i = 0; i < D; ++i)
1523 result(i, j) = worker(i);
1536 if (subsys.size() == dims.size())
1537 return rA.transpose();
1539 if (subsys.size() == 0)
1542 auto worker = [=, &Cmidxcol](
idx i) noexcept
1543 ->
typename Derived::Scalar
1549 for (
idx k = 0; k < numdims; ++k)
1550 midxcoltmp[k] = Cmidxcol[k];
1555 for (
idx k = 0; k < numsubsys; ++k)
1556 std::swap(midxcoltmp[Csubsys[k]], midxrow[Csubsys[k]]);
1563 for (
idx j = 0; j < D; ++j)
1567 #pragma omp parallel for
1568 for (
idx i = 0; i < D; ++i)
1569 result(i, j) = worker(i);
1593 template<
typename Derived>
1595 const Eigen::MatrixBase<Derived>& A,
1596 const std::vector<idx>& subsys,
1613 static_cast<idx>(std::llround(std::log2(rA.rows()) /
1615 std::vector<idx> dims(n, d);
1632 template<
typename Derived>
1634 const Eigen::MatrixBase<Derived>& A,
1635 const std::vector<idx>& perm,
1636 const std::vector<idx>& dims)
1655 if (perm.size() != dims.size())
1660 idx D =
static_cast<idx>(rA.rows());
1661 idx numdims = dims.size();
1677 for (
idx i = 0; i < numdims; ++i)
1682 result.resize(D, 1);
1684 auto worker = [&Cdims, &Cperm, numdims](
idx i) noexcept
1696 for (
idx k = 0; k < numdims; ++k)
1698 permdims[k] = Cdims[Cperm[k]];
1699 midxtmp[k] = midx[Cperm[k]];
1705 #pragma omp parallel for
1706 for (
idx i = 0; i < D; ++i)
1707 result(worker(i)) = rA(i);
1723 for (
idx i = 0; i < numdims; ++i)
1726 Cdims[i + numdims] = dims[i];
1728 Cperm[i + numdims] = perm[i] + numdims;
1730 result.resize(D * D, 1);
1733 Eigen::Map<dyn_mat<typename Derived::Scalar>>(
1734 const_cast<typename Derived::Scalar*
>(rA.data()), D * D,
1737 auto worker = [&Cdims, &Cperm, numdims](
idx i) noexcept
1749 for (
idx k = 0; k < 2 * numdims; ++k)
1751 permdims[k] = Cdims[Cperm[k]];
1752 midxtmp[k] = midx[Cperm[k]];
1758 #pragma omp parallel for
1759 for (
idx i = 0; i < D * D; ++i)
1760 result(worker(i)) = rA(i);
1782 template<
typename Derived>
1784 const Eigen::MatrixBase<Derived>& A,
1785 const std::vector<idx>& perm,
1802 static_cast<idx>(std::llround(std::log2(rA.rows()) /
1804 std::vector<idx> dims(n, d);
bool _check_cvector(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:104
constexpr idx maxn
Maximum number of allowed qu(d)its (subsystems)
Definition: constants.h:80
std::vector< cmat > choi2kraus(const cmat &A)
Orthogonal Kraus operators from Choi matrix.
Definition: operations.h:890
bool _check_subsys_match_dims(const std::vector< idx > &subsys, const std::vector< idx > &dims)
Definition: util.h:190
bool _check_dims_match_mat(const std::vector< idx > &dims, const Eigen::MatrixBase< Derived > &A)
Definition: util.h:142
Eigen::MatrixXd dmat
Real (double precision) dynamic Eigen matrix.
Definition: types.h:71
bool _check_dims_match_cvect(const std::vector< idx > &dims, const Eigen::MatrixBase< Derived > &V)
Definition: util.h:155
constexpr double eps
Used to decide whether a number or expression in double precision is zero or not. ...
Definition: constants.h:73
Eigen::RowVectorXcd bra
Complex (double precision) dynamic Eigen row vector.
Definition: types.h:61
dyn_mat< typename Derived::Scalar > syspermute(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &perm, const std::vector< idx > &dims)
Subsystem permutation.
Definition: operations.h:1633
Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic > dyn_mat
Dynamic Eigen matrix over the field specified by Scalar.
Definition: types.h:83
Eigen::VectorXcd ket
Complex (double precision) dynamic Eigen column vector.
Definition: types.h:56
dyn_mat< typename Derived1::Scalar > apply(const Eigen::MatrixBase< Derived1 > &state, const Eigen::MatrixBase< Derived2 > &A, const std::vector< idx > &subsys, const std::vector< idx > &dims)
Applies the gate A to the part subsys of the multi-partite state vector or density matrix state...
Definition: operations.h:486
Quantum++ main namespace.
Definition: codes.h:30
bool _check_square_mat(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:77
dyn_col_vect< double > hevals(const Eigen::MatrixBase< Derived > &A)
Hermitian eigenvalues.
Definition: functions.h:391
cmat choi2super(const cmat &A)
Converts Choi matrix to superoperator matrix.
Definition: operations.h:926
cmat super2choi(const cmat &A)
Converts superoperator matrix to Choi matrix.
Definition: operations.h:960
std::vector< T > complement(std::vector< T > subsys, idx N)
Constructs the complement of a subsystem vector.
Definition: functions.h:1806
dyn_mat< typename Derived::Scalar > ptrace1(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Partial trace.
Definition: operations.h:1002
dyn_mat< typename Derived::Scalar > adjoint(const Eigen::MatrixBase< Derived > &A)
Adjoint.
Definition: functions.h:84
dyn_mat< typename T::Scalar > kron(const T &head)
Kronecker product.
Definition: functions.h:890
dyn_mat< typename Derived::Scalar > ptrace2(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Partial trace.
Definition: operations.h:1099
bool _check_perm(const std::vector< idx > &perm)
Definition: util.h:257
idx _multiidx2n(const idx *midx, idx numdims, const idx *dims) noexcept
Definition: util.h:54
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_mat< typename Derived::Scalar > powm(const Eigen::MatrixBase< Derived > &A, idx n)
Matrix power.
Definition: functions.h:783
cmat hevects(const Eigen::MatrixBase< Derived > &A)
Hermitian eigenvectors.
Definition: functions.h:417
bool _check_nonzero_size(const T &x) noexcept
Definition: util.h:113
cmat kraus2choi(const std::vector< cmat > &Ks)
Choi matrix.
Definition: operations.h:835
dyn_mat< typename Derived::Scalar > ptrace(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &subsys, const std::vector< idx > &dims)
Partial trace.
Definition: operations.h:1186
Derived::Scalar trace(const Eigen::MatrixBase< Derived > &A)
Trace.
Definition: functions.h:127
cmat kraus2super(const std::vector< cmat > &Ks)
Superoperator matrix.
Definition: operations.h:761
dyn_mat< typename Derived1::Scalar > applyCTRL(const Eigen::MatrixBase< Derived1 > &state, const Eigen::MatrixBase< Derived2 > &A, const std::vector< idx > &ctrl, const std::vector< idx > &subsys, const std::vector< idx > &dims)
Applies the controlled-gate A to the part subsys of the multi-partite state vector or density matrix ...
Definition: operations.h:54
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:1443
Eigen::MatrixXcd cmat
Complex (double precision) dynamic Eigen matrix.
Definition: types.h:66
void _n2multiidx(idx n, idx numdims, const idx *dims, idx *result) noexcept
Definition: util.h:41
std::size_t idx
Non-negative integer index.
Definition: types.h:36
bool _check_dims(const std::vector< idx > &dims)
Definition: util.h:126
Derived::Scalar sum(const Eigen::MatrixBase< Derived > &A)
Element-wise sum of A.
Definition: functions.h:209
dyn_mat< typename Derived::Scalar > reshape(const Eigen::MatrixBase< Derived > &A, idx rows, idx cols)
Reshape.
Definition: functions.h:1130