53 template <
typename Derived1,
typename Derived2>
54 dyn_mat<typename Derived1::Scalar>
56 const Eigen::MatrixBase<Derived2>& A,
const std::vector<idx>& ctrl,
57 const std::vector<idx>& subsys,
const std::vector<idx>& dims) {
58 const typename Eigen::MatrixBase<Derived1>::EvalReturnType& rstate =
65 if (!std::is_same<
typename Derived1::Scalar,
66 typename Derived2::Scalar>::value)
82 idx d = ctrl.size() > 0 ? dims[ctrl[0]] : 1;
83 for (
idx i = 1; i < ctrl.size(); ++i)
84 if (dims[ctrl[i]] != d)
96 std::vector<idx> subsys_dims(subsys.size());
97 for (
idx i = 0; i < subsys.size(); ++i)
98 subsys_dims[i] = dims[subsys[i]];
102 std::vector<idx> ctrlgate = ctrl;
103 ctrlgate.insert(std::end(ctrlgate), std::begin(subsys), std::end(subsys));
104 std::sort(std::begin(ctrlgate), std::end(ctrlgate));
113 std::vector<dyn_mat<typename Derived1::Scalar>> Ai;
114 std::vector<dyn_mat<typename Derived1::Scalar>> Aidagger;
115 for (
idx i = 0; i < std::max(d, static_cast<idx>(2)); ++i) {
116 Ai.push_back(
powm(rA, i));
120 idx D =
static_cast<idx>(rstate.rows());
122 idx ctrlsize = ctrl.size();
123 idx ctrlgatesize = ctrlgate.size();
124 idx subsyssize = subsys.size();
126 idx Dctrl =
static_cast<idx>(std::llround(std::pow(d, ctrlsize)));
127 idx DA =
static_cast<idx>(rA.rows());
135 std::vector<idx> ctrlgate_bar =
complement(ctrlgate, n);
137 idx ctrlgate_barsize = ctrlgate_bar.size();
140 for (
idx i = 0; i < ctrlgate_barsize; ++i)
141 DCTRLA_bar *= dims[ctrlgate_bar[i]];
143 for (
idx k = 0; k < n; ++k)
145 for (
idx k = 0; k < subsyssize; ++k)
146 CdimsA[k] = dims[subsys[k]];
147 for (
idx k = 0; k < ctrlsize; ++k)
149 for (
idx k = 0; k < ctrlgate_barsize; ++k)
150 CdimsCTRLA_bar[k] = dims[ctrlgate_bar[k]];
154 auto coeff_idx_ket = [&](
idx i_,
idx m_,
idx r_) noexcept
155 ->std::pair<typename Derived1::Scalar, idx> {
157 typename Derived1::Scalar coeff = 0;
166 for (
idx k = 0; k < ctrlsize; ++k) {
173 for (
idx k = 0; k < n - ctrlgatesize; ++k) {
174 Cmidx[ctrlgate_bar[k]] = CmidxCTRLA_bar[k];
179 for (
idx k = 0; k < subsyssize; ++k) {
180 Cmidx[subsys[k]] = CmidxA[k];
187 for (
idx n_ = 0; n_ < DA; ++n_) {
189 for (
idx k = 0; k < subsyssize; ++k) {
190 Cmidx[subsys[k]] = CmidxA[k];
196 return std::make_pair(coeff, indx);
204 ->std::tuple<typename Derived1::Scalar, idx, idx> {
207 typename Derived1::Scalar coeff = 0, lhs = 1, rhs = 1;
224 for (
idx k = 0; k < ctrlsize; ++k) {
225 Cmidxrow[ctrl[k]] = CmidxCTRLrow[k];
226 Cmidxcol[ctrl[k]] = CmidxCTRLcol[k];
234 for (
idx k = 0; k < n - ctrlgatesize; ++k) {
235 Cmidxrow[ctrlgate_bar[k]] = CmidxCTRLA_barrow[k];
236 Cmidxcol[ctrlgate_bar[k]] = CmidxCTRLA_barcol[k];
242 for (
idx k = 0; k < subsys.size(); ++k) {
243 Cmidxrow[subsys[k]] = CmidxArow[k];
244 Cmidxcol[subsys[k]] = CmidxAcol[k];
252 bool all_ctrl_rows_equal =
true;
253 bool all_ctrl_cols_equal =
true;
255 idx first_ctrl_row, first_ctrl_col;
257 first_ctrl_row = CmidxCTRLrow[0];
258 first_ctrl_col = CmidxCTRLcol[0];
260 first_ctrl_row = first_ctrl_col = 1;
263 for (
idx k = 1; k < ctrlsize; ++k) {
264 if (CmidxCTRLrow[k] != first_ctrl_row) {
265 all_ctrl_rows_equal =
false;
269 for (
idx k = 1; k < ctrlsize; ++k) {
270 if (CmidxCTRLcol[k] != first_ctrl_col) {
271 all_ctrl_cols_equal =
false;
277 for (
idx n1_ = 0; n1_ < DA; ++n1_) {
279 for (
idx k = 0; k < subsyssize; ++k) {
280 Cmidxrow[subsys[k]] = CmidxArow[k];
284 if (all_ctrl_rows_equal) {
285 lhs = Ai[first_ctrl_row](m1_, n1_);
287 lhs = (m1_ == n1_) ? 1 : 0;
290 for (
idx n2_ = 0; n2_ < DA; ++n2_) {
292 for (
idx k = 0; k < subsyssize; ++k) {
293 Cmidxcol[subsys[k]] = CmidxAcol[k];
296 if (all_ctrl_cols_equal) {
297 rhs = Aidagger[first_ctrl_col](n2_, m2_);
299 rhs = (n2_ == m2_) ? 1 : 0;
304 coeff += lhs * rstate(idxrowtmp, idxcoltmp) * rhs;
308 return std::make_tuple(coeff, idxrow, idxcol);
323 #pragma omp parallel for collapse(2) 324 #endif // WITH_OPENMP_ 325 for (
idx m = 0; m < DA; ++m)
326 for (
idx r = 0; r < DCTRLA_bar; ++r) {
329 result(coeff_idx_ket(1, m, r).second) =
330 coeff_idx_ket(1, m, r).first;
332 for (
idx i = 0; i < d; ++i) {
333 result(coeff_idx_ket(i, m, r).second) =
334 coeff_idx_ket(i, m, r).first;
353 #pragma omp parallel for collapse(4) 354 #endif // WITH_OPENMP_ 355 for (
idx m1 = 0; m1 < DA; ++m1)
356 for (
idx r1 = 0; r1 < DCTRLA_bar; ++r1)
357 for (
idx m2 = 0; m2 < DA; ++m2)
358 for (
idx r2 = 0; r2 < DCTRLA_bar; ++r2)
362 coeff_idx_rho(1, m1, r1, 1, m2, r2);
363 result(std::get<1>(coeff_idxes),
364 std::get<2>(coeff_idxes)) =
365 std::get<0>(coeff_idxes);
367 for (
idx i1 = 0; i1 < Dctrl; ++i1)
368 for (
idx i2 = 0; i2 < Dctrl; ++i2) {
370 coeff_idx_rho(i1, m1, r1, i2, m2, r2);
371 result(std::get<1>(coeff_idxes),
372 std::get<2>(coeff_idxes)) =
373 std::get<0>(coeff_idxes);
399 template <
typename Derived1,
typename Derived2>
402 const Eigen::MatrixBase<Derived2>& A,
const std::vector<idx>& ctrl,
403 const std::vector<idx>& subsys,
idx d = 2) {
404 const typename Eigen::MatrixBase<Derived1>::EvalReturnType& rstate =
420 std::vector<idx> dims(n, d);
422 return applyCTRL(rstate, rA, ctrl, subsys, dims);
438 template <
typename Derived1,
typename Derived2>
440 apply(
const Eigen::MatrixBase<Derived1>& state,
441 const Eigen::MatrixBase<Derived2>& A,
const std::vector<idx>& subsys,
442 const std::vector<idx>& dims) {
443 const typename Eigen::MatrixBase<Derived1>::EvalReturnType& rstate =
450 if (!std::is_same<
typename Derived1::Scalar,
451 typename Derived2::Scalar>::value)
475 std::vector<idx> subsys_dims(subsys.size());
476 for (
idx i = 0; i < subsys.size(); ++i)
477 subsys_dims[i] = dims[subsys[i]];
489 return applyCTRL(rstate, rA, {}, subsys, dims);
498 return applyCTRL(rstate, rA, {}, subsys, dims);
518 template <
typename Derived1,
typename Derived2>
520 apply(
const Eigen::MatrixBase<Derived1>& state,
521 const Eigen::MatrixBase<Derived2>& A,
const std::vector<idx>& subsys,
523 const typename Eigen::MatrixBase<Derived1>::EvalReturnType& rstate =
539 std::vector<idx> dims(n, d);
541 return apply(rstate, rA, subsys, dims);
552 template <
typename Derived>
553 cmat apply(
const Eigen::MatrixBase<Derived>& A,
const std::vector<cmat>& Ks) {
554 const cmat& rA = A.derived();
566 if (Ks[0].rows() != rA.rows())
569 if (it.rows() != Ks[0].rows() || it.cols() != Ks[0].rows())
573 cmat result = cmat::Zero(rA.rows(), rA.rows());
576 #pragma omp parallel for 577 #endif // WITH_OPENMP_ 578 for (
idx i = 0; i < Ks.size(); ++i) {
581 #endif // WITH_OPENMP_ 582 { result += Ks[i] * rA *
adjoint(Ks[i]); }
598 template <
typename Derived>
599 cmat apply(
const Eigen::MatrixBase<Derived>& A,
const std::vector<cmat>& Ks,
600 const std::vector<idx>& subsys,
const std::vector<idx>& dims) {
601 const cmat& rA = A.derived();
625 std::vector<idx> subsys_dims(subsys.size());
626 for (
idx i = 0; i < subsys.size(); ++i)
627 subsys_dims[i] = dims[subsys[i]];
637 if (it.rows() != Ks[0].rows() || it.cols() != Ks[0].rows())
641 cmat result = cmat::Zero(rA.rows(), rA.rows());
643 for (
idx i = 0; i < Ks.size(); ++i)
644 result +=
apply(rA, Ks[i], subsys, dims);
659 template <
typename Derived>
660 cmat apply(
const Eigen::MatrixBase<Derived>& A,
const std::vector<cmat>& Ks,
661 const std::vector<idx>& subsys,
idx d = 2) {
662 const cmat& rA = A.derived();
676 std::vector<idx> dims(n, d);
678 return apply(rA, Ks, subsys, dims);
702 if (it.rows() != Ks[0].rows() || it.cols() != Ks[0].rows())
706 idx D =
static_cast<idx>(Ks[0].rows());
708 cmat result(D * D, D * D);
709 cmat MN = cmat::Zero(D, D);
710 bra A = bra::Zero(D);
711 ket B = ket::Zero(D);
712 cmat EMN = cmat::Zero(D, D);
715 #pragma omp parallel for collapse(2) 716 #endif // WITH_OPENMP_ 717 for (
idx m = 0; m < D; ++m) {
718 for (
idx n = 0; n < D; ++n) {
721 #endif // WITH_OPENMP_ 725 for (
idx i = 0; i < Ks.size(); ++i)
726 EMN += Ks[i] * MN *
adjoint(Ks[i]);
729 for (
idx a = 0; a < D; ++a) {
731 for (
idx b = 0; b < D; ++b) {
734 result(a * D + b, m * D + n) =
735 static_cast<cmat>(A * EMN * B).value();
740 EMN = cmat::Zero(D, D);
773 if (it.rows() != Ks[0].rows() || it.cols() != Ks[0].rows())
777 idx D =
static_cast<idx>(Ks[0].rows());
781 cmat MES = cmat::Zero(D * D, 1);
782 for (
idx a = 0; a < D; ++a)
787 cmat result = cmat::Zero(D * D, D * D);
790 #pragma omp parallel for 791 #endif // WITH_OPENMP_ 792 for (
idx i = 0; i < Ks.size(); ++i) {
795 #endif // WITH_OPENMP_ 797 result +=
kron(cmat::Identity(D, D), Ks[i]) * Omega *
827 if (D * D != static_cast<idx>(A.rows()))
833 std::vector<cmat> result;
835 for (
idx i = 0; i < D * D; ++i) {
836 if (std::abs(ev(i)) >
eps)
837 result.push_back(std::sqrt(std::abs(ev(i))) *
860 if (D * D != static_cast<idx>(A.rows()))
864 cmat result(D * D, D * D);
867 #pragma omp parallel for collapse(4) 868 #endif // WITH_OPENMP_ 869 for (
idx a = 0; a < D; ++a)
870 for (
idx b = 0; b < D; ++b)
871 for (
idx m = 0; m < D; ++m)
872 for (
idx n = 0; n < D; ++n)
873 result(a * D + b, m * D + n) = A(m * D + a, n * D + b);
894 if (D * D != static_cast<idx>(A.rows()))
898 cmat result(D * D, D * D);
901 #pragma omp parallel for collapse(4) 902 #endif // WITH_OPENMP_ 903 for (
idx a = 0; a < D; ++a)
904 for (
idx b = 0; b < D; ++b)
905 for (
idx m = 0; m < D; ++m)
906 for (
idx n = 0; n < D; ++n)
907 result(m * D + a, n * D + b) = A(a * D + b, m * D + n);
925 template <
typename Derived>
927 const std::vector<idx>& dims) {
928 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
941 if (dims.size() != 2)
958 auto worker = [&](
idx i,
idx j) noexcept->typename Derived::Scalar {
959 typename Derived::Scalar
sum = 0;
960 for (
idx m = 0; m < DA; ++m)
961 sum += rA(m * DB + i) * std::conj(rA(m * DB + j));
967 #pragma omp parallel for collapse(2) 968 #endif // WITH_OPENMP_ 970 for (
idx j = 0; j < DB; ++j)
971 for (
idx i = 0; i < DB; ++i)
972 result(i, j) = worker(i, j);
983 auto worker = [&](
idx i,
idx j) noexcept->typename Derived::Scalar {
984 typename Derived::Scalar
sum = 0;
985 for (
idx m = 0; m < DA; ++m)
986 sum += rA(m * DB + i, m * DB + j);
992 #pragma omp parallel for collapse(2) 993 #endif // WITH_OPENMP_ 995 for (
idx j = 0; j < DB; ++j)
996 for (
idx i = 0; i < DB; ++i)
997 result(i, j) = worker(i, j);
1019 template <
typename Derived>
1035 std::vector<idx> dims(2, d);
1053 template <
typename Derived>
1055 const std::vector<idx>& dims) {
1056 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
1069 if (dims.size() != 2)
1086 auto worker = [&](
idx i,
idx j) noexcept->typename Derived::Scalar {
1087 typename Derived::Scalar
sum = 0;
1088 for (
idx m = 0; m < DB; ++m)
1089 sum += rA(i * DB + m) * std::conj(rA(j * DB + m));
1095 #pragma omp parallel for collapse(2) 1096 #endif // WITH_OPENMP_ 1098 for (
idx j = 0; j < DA; ++j)
1099 for (
idx i = 0; i < DA; ++i)
1100 result(i, j) = worker(i, j);
1112 #pragma omp parallel for collapse(2) 1113 #endif // WITH_OPENMP_ 1115 for (
idx j = 0; j < DA; ++j)
1116 for (
idx i = 0; i < DA; ++i)
1117 result(i, j) =
trace(rA.block(i * DB, j * DB, DB, DB));
1139 template <
typename Derived>
1155 std::vector<idx> dims(2, d);
1174 template <
typename Derived>
1176 const std::vector<idx>& subsys,
1177 const std::vector<idx>& dims) {
1178 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
1195 idx D =
static_cast<idx>(rA.rows());
1196 idx n = dims.size();
1197 idx n_subsys = subsys.size();
1198 idx n_subsys_bar = n - n_subsys;
1200 for (
idx i = 0; i < n_subsys; ++i)
1201 Dsubsys *= dims[subsys[i]];
1202 idx Dsubsys_bar = D / Dsubsys;
1212 std::vector<idx> subsys_bar =
complement(subsys, n);
1213 std::copy(std::begin(subsys_bar), std::end(subsys_bar),
1214 std::begin(Csubsys_bar));
1216 for (
idx i = 0; i < n; ++i) {
1219 for (
idx i = 0; i < n_subsys; ++i) {
1220 Csubsys[i] = subsys[i];
1221 Cdimssubsys[i] = dims[subsys[i]];
1223 for (
idx i = 0; i < n_subsys_bar; ++i) {
1224 Cdimssubsys_bar[i] = dims[subsys_bar[i]];
1237 if (subsys.size() == dims.size()) {
1238 result(0, 0) = (
adjoint(rA) * rA).value();
1242 if (subsys.size() == 0)
1245 auto worker = [&](
idx i) noexcept->typename Derived::Scalar {
1255 Cmidxrowsubsys_bar);
1257 for (
idx k = 0; k < n_subsys_bar; ++k) {
1258 Cmidxrow[Csubsys_bar[k]] = Cmidxrowsubsys_bar[k];
1259 Cmidxcol[Csubsys_bar[k]] = Cmidxcolsubsys_bar[k];
1261 typename Derived::Scalar sm = 0;
1262 for (
idx a = 0; a < Dsubsys; ++a) {
1266 for (
idx k = 0; k < n_subsys; ++k)
1267 Cmidxrow[Csubsys[k]] = Cmidxcol[Csubsys[k]] =
1278 for (
idx j = 0; j < Dsubsys_bar; ++j)
1282 Cmidxcolsubsys_bar);
1284 #pragma omp parallel for 1285 #endif // WITH_OPENMP_ 1286 for (
idx i = 0; i < Dsubsys_bar; ++i) {
1287 result(i, j) = worker(i);
1300 if (subsys.size() == dims.size()) {
1301 result(0, 0) = rA.trace();
1305 if (subsys.size() == 0)
1308 auto worker = [&](
idx i) noexcept->typename Derived::Scalar {
1318 Cmidxrowsubsys_bar);
1320 for (
idx k = 0; k < n_subsys_bar; ++k) {
1321 Cmidxrow[Csubsys_bar[k]] = Cmidxrowsubsys_bar[k];
1322 Cmidxcol[Csubsys_bar[k]] = Cmidxcolsubsys_bar[k];
1324 typename Derived::Scalar sm = 0;
1325 for (
idx a = 0; a < Dsubsys; ++a) {
1329 for (
idx k = 0; k < n_subsys; ++k)
1330 Cmidxrow[Csubsys[k]] = Cmidxcol[Csubsys[k]] =
1341 for (
idx j = 0; j < Dsubsys_bar; ++j)
1345 Cmidxcolsubsys_bar);
1347 #pragma omp parallel for 1348 #endif // WITH_OPENMP_ 1349 for (
idx i = 0; i < Dsubsys_bar; ++i) {
1350 result(i, j) = worker(i);
1375 template <
typename Derived>
1377 const std::vector<idx>& subsys,
1379 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
1393 std::vector<idx> dims(n, d);
1395 return ptrace(rA, subsys, dims);
1411 template <
typename Derived>
1413 ptranspose(
const Eigen::MatrixBase<Derived>& A,
const std::vector<idx>& subsys,
1414 const std::vector<idx>& dims) {
1415 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
1432 idx D =
static_cast<idx>(rA.rows());
1433 idx n = dims.size();
1434 idx n_subsys = subsys.size();
1440 for (
idx i = 0; i < n; ++i)
1442 for (
idx i = 0; i < n_subsys; ++i)
1443 Csubsys[i] = subsys[i];
1454 if (subsys.size() == dims.size())
1457 if (subsys.size() == 0)
1460 auto worker = [&](
idx i) noexcept->typename Derived::Scalar {
1465 for (
idx k = 0; k < n; ++k)
1466 midxcoltmp[k] = Cmidxcol[k];
1471 for (
idx k = 0; k < n_subsys; ++k)
1472 std::swap(midxcoltmp[Csubsys[k]], midxrow[Csubsys[k]]);
1479 for (
idx j = 0; j < D; ++j) {
1484 #pragma omp parallel for 1485 #endif // WITH_OPENMP_ 1486 for (
idx i = 0; i < D; ++i)
1487 result(i, j) = worker(i);
1499 if (subsys.size() == dims.size())
1500 return rA.transpose();
1502 if (subsys.size() == 0)
1505 auto worker = [&](
idx i) noexcept->typename Derived::Scalar {
1510 for (
idx k = 0; k < n; ++k)
1511 midxcoltmp[k] = Cmidxcol[k];
1516 for (
idx k = 0; k < n_subsys; ++k)
1517 std::swap(midxcoltmp[Csubsys[k]], midxrow[Csubsys[k]]);
1524 for (
idx j = 0; j < D; ++j) {
1529 #pragma omp parallel for 1530 #endif // WITH_OPENMP_ 1531 for (
idx i = 0; i < D; ++i)
1532 result(i, j) = worker(i);
1555 template <
typename Derived>
1557 ptranspose(
const Eigen::MatrixBase<Derived>& A,
const std::vector<idx>& subsys,
1559 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
1573 std::vector<idx> dims(n, d);
1590 template <
typename Derived>
1592 syspermute(
const Eigen::MatrixBase<Derived>& A,
const std::vector<idx>& perm,
1593 const std::vector<idx>& dims) {
1594 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
1611 if (perm.size() != dims.size())
1615 idx D =
static_cast<idx>(rA.rows());
1616 idx n = dims.size();
1631 for (
idx i = 0; i < n; ++i) {
1635 result.resize(D, 1);
1637 auto worker = [&Cdims, &Cperm, n ](
idx i) noexcept->idx {
1647 for (
idx k = 0; k < n; ++k) {
1648 permdims[k] = Cdims[Cperm[k]];
1649 midxtmp[k] = midx[Cperm[k]];
1656 #pragma omp parallel for 1657 #endif // WITH_OPENMP_ 1658 for (
idx i = 0; i < D; ++i)
1659 result(worker(i)) = rA(i);
1674 for (
idx i = 0; i < n; ++i) {
1676 Cdims[i + n] = dims[i];
1678 Cperm[i + n] = perm[i] + n;
1680 result.resize(D * D, 1);
1683 Eigen::Map<dyn_mat<typename Derived::Scalar>>(
1684 const_cast<typename Derived::Scalar*
>(rA.data()), D * D, 1);
1686 auto worker = [&Cdims, &Cperm, n ](
idx i) noexcept->idx {
1696 for (
idx k = 0; k < 2 * n; ++k) {
1697 permdims[k] = Cdims[Cperm[k]];
1698 midxtmp[k] = midx[Cperm[k]];
1705 #pragma omp parallel for 1706 #endif // WITH_OPENMP_ 1707 for (
idx i = 0; i < D * D; ++i)
1708 result(worker(i)) = rA(i);
1729 template <
typename Derived>
1731 syspermute(
const Eigen::MatrixBase<Derived>& A,
const std::vector<idx>& perm,
1733 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
1747 std::vector<idx> dims(n, d);
1763 template <
typename Derived>
1765 const std::vector<idx>& subsys,
1782 std::vector<idx> dims(n, d);
1808 idx n_subsys = subsys.size();
1812 for (
idx i = 0; i < n_subsys; ++i) {
1816 for (
idx j = 2; j <= n_subsys - i; ++j) {
1819 Rj << 1, 0, 0,
omega(std::pow(2, j));
1821 applyCTRL(result, Rj, {subsys[i + j - 1]}, {subsys[i]});
1825 for (
idx i = 0; i < n_subsys / 2; ++i) {
1827 {subsys[i], subsys[n_subsys - i - 1]});
1831 for (
idx i = 0; i < n_subsys; ++i) {
1835 for (
idx j = 2; j <= n_subsys - i; ++j) {
1837 cmat Rj = cmat::Zero(d, d);
1838 for (
idx m = 0; m < d; ++m) {
1839 Rj(m, m) = exp(2.0 *
pi * m * 1_i / std::pow(d, j));
1842 applyCTRL(result, Rj, {subsys[i + j - 1]}, {subsys[i]}, d);
1846 for (
idx i = 0; i < n_subsys / 2; ++i) {
1848 {subsys[i], subsys[n_subsys - i - 1]}, d);
1863 template <
typename Derived>
1881 std::vector<idx> dims(n, d);
1902 std::vector<idx> subsys(n);
1903 std::iota(std::begin(subsys), std::end(subsys), 0);
bool check_nonzero_size(const T &x) noexcept
Definition: util.h:123
Dimensions not equal exception.
Definition: exception.h:284
Dimension(s) mismatch matrix size exception.
Definition: exception.h:300
constexpr idx maxn
Maximum number of allowed qubits/qudits (subsystems)
Definition: constants.h:75
bool check_subsys_match_dims(const std::vector< idx > &subsys, const std::vector< idx > &dims)
Definition: util.h:210
bool check_cvector(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:117
std::vector< cmat > choi2kraus(const cmat &A)
Orthogonal Kraus operators from Choi matrix.
Definition: operations.h:818
Eigen::MatrixXd dmat
Real (double precision) dynamic Eigen matrix.
Definition: types.h:69
cplx omega(idx D)
D-th root of unity.
Definition: constants.h:97
constexpr double eps
Used to decide whether a number or expression in double precision is zero or not. ...
Definition: constants.h:68
Eigen::RowVectorXcd bra
Complex (double precision) dynamic Eigen row vector.
Definition: types.h:59
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:1592
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
idx multiidx2n(const idx *const midx, idx numdims, const idx *const dims) noexcept
Definition: util.h:71
void n2multiidx(idx n, idx numdims, const idx *const dims, idx *result) noexcept
Definition: util.h:43
Subsystems mismatch dimensions exception.
Definition: exception.h:365
Eigen::VectorXcd ket
Complex (double precision) dynamic Eigen column vector.
Definition: types.h:54
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:440
Quantum++ main namespace.
Definition: codes.h:35
bool check_perm(const std::vector< idx > &perm)
Definition: util.h:259
idx get_num_subsys(idx sz, idx d)
Definition: util.h:366
Matrix is not square nor column vector exception.
Definition: exception.h:209
dyn_col_vect< double > hevals(const Eigen::MatrixBase< Derived > &A)
Hermitian eigenvalues.
Definition: functions.h:376
Invalid dimension(s) exception.
Definition: exception.h:269
dyn_col_vect< typename Derived::Scalar > QFT(const Eigen::MatrixBase< Derived > &A, idx d=2)
Qudit quantum Fourier transform.
Definition: operations.h:1864
cmat choi2super(const cmat &A)
Converts Choi matrix to superoperator matrix.
Definition: operations.h:851
cmat super2choi(const cmat &A)
Converts superoperator matrix to Choi matrix.
Definition: operations.h:885
std::vector< T > complement(std::vector< T > subsys, idx N)
Constructs the complement of a subsystem vector.
Definition: functions.h:1733
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:486
dyn_mat< typename Derived::Scalar > applyQFT(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &subsys, idx d=2)
Applies the qudit quantum Fourier transform to the part subsys of the multi-partite state vector or d...
Definition: operations.h:1764
dyn_mat< typename Derived::Scalar > ptrace1(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Partial trace.
Definition: operations.h:926
Invalid permutation exception.
Definition: exception.h:380
dyn_mat< typename Derived::Scalar > adjoint(const Eigen::MatrixBase< Derived > &A)
Adjoint.
Definition: functions.h:83
dyn_mat< typename T::Scalar > kron(const T &head)
Kronecker product.
Definition: functions.h:868
dyn_mat< typename Derived::Scalar > ptrace2(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Partial trace.
Definition: operations.h:1054
Matrix mismatch subsystems exception.
Definition: exception.h:254
Dimension(s) mismatch column vector size exception.
Definition: exception.h:316
dyn_mat< typename Derived::Scalar > transpose(const Eigen::MatrixBase< Derived > &A)
Transpose.
Definition: functions.h:46
dyn_mat< typename Derived::Scalar > powm(const Eigen::MatrixBase< Derived > &A, idx n)
Fast matrix power based on the SQUARE-AND-MULTIPLY algorithm.
Definition: functions.h:751
Type mismatch exception.
Definition: exception.h:530
idx get_dim_subsys(idx sz, idx N)
Definition: util.h:378
cmat hevects(const Eigen::MatrixBase< Derived > &A)
Hermitian eigenvectors.
Definition: functions.h:401
cmat kraus2choi(const std::vector< cmat > &Ks)
Choi matrix.
Definition: operations.h:763
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:1175
bool check_square_mat(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:99
Derived::Scalar trace(const Eigen::MatrixBase< Derived > &A)
Trace.
Definition: functions.h:124
constexpr double pi
Definition: constants.h:80
cmat kraus2super(const std::vector< cmat > &Ks)
Superoperator matrix.
Definition: operations.h:692
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:55
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:1413
Eigen::MatrixXcd cmat
Complex (double precision) dynamic Eigen matrix.
Definition: types.h:64
bool check_dims(const std::vector< idx > &dims)
Definition: util.h:134
Permutation mismatch dimensions exception.
Definition: exception.h:396
std::size_t idx
Non-negative integer index.
Definition: types.h:39
Matrix is not square exception.
Definition: exception.h:149
Derived::Scalar sum(const Eigen::MatrixBase< Derived > &A)
Element-wise sum of A.
Definition: functions.h:201
dyn_mat< typename Derived::Scalar > reshape(const Eigen::MatrixBase< Derived > &A, idx rows, idx cols)
Reshape.
Definition: functions.h:1095
Object has zero size exception.
Definition: exception.h:134