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));
117 Aidagger.push_back(powm(adjoint(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);
499 return applyCTRL(rstate, rA, {}, subsys, dims);
519 template <
typename Derived1,
typename Derived2>
521 apply(
const Eigen::MatrixBase<Derived1>& state,
522 const Eigen::MatrixBase<Derived2>& A,
const std::vector<idx>& subsys,
524 const typename Eigen::MatrixBase<Derived1>::EvalReturnType& rstate =
540 std::vector<idx> dims(N, d);
542 return apply(rstate, rA, subsys, dims);
553 template <
typename Derived>
554 cmat apply(
const Eigen::MatrixBase<Derived>& A,
const std::vector<cmat>& Ks) {
555 const cmat& rA = A.derived();
567 if (Ks[0].rows() != rA.rows())
570 if (it.rows() != Ks[0].rows() || it.cols() != Ks[0].rows())
574 cmat result = cmat::Zero(rA.rows(), rA.rows());
577 #pragma omp parallel for 578 #endif // WITH_OPENMP_ 579 for (
idx i = 0; i < Ks.size(); ++i) {
582 #endif // WITH_OPENMP_ 583 { result += Ks[i] * rA * adjoint(Ks[i]); }
599 template <
typename Derived>
600 cmat apply(
const Eigen::MatrixBase<Derived>& A,
const std::vector<cmat>& Ks,
601 const std::vector<idx>& subsys,
const std::vector<idx>& dims) {
602 const cmat& rA = A.derived();
626 std::vector<idx> subsys_dims(subsys.size());
627 for (
idx i = 0; i < subsys.size(); ++i)
628 subsys_dims[i] = dims[subsys[i]];
638 if (it.rows() != Ks[0].rows() || it.cols() != Ks[0].rows())
642 cmat result = cmat::Zero(rA.rows(), rA.rows());
644 for (
idx i = 0; i < Ks.size(); ++i)
645 result +=
apply(rA, Ks[i], subsys, dims);
660 template <
typename Derived>
661 cmat apply(
const Eigen::MatrixBase<Derived>& A,
const std::vector<cmat>& Ks,
662 const std::vector<idx>& subsys,
idx d = 2) {
663 const cmat& rA = A.derived();
677 std::vector<idx> dims(N, d);
679 return apply(rA, Ks, subsys, dims);
703 if (it.rows() != Ks[0].rows() || it.cols() != Ks[0].rows())
707 idx D =
static_cast<idx>(Ks[0].rows());
709 cmat result(D * D, D * D);
710 cmat MN = cmat::Zero(D, D);
711 bra A = bra::Zero(D);
712 ket B = ket::Zero(D);
713 cmat EMN = cmat::Zero(D, D);
716 #pragma omp parallel for collapse(2) 717 #endif // WITH_OPENMP_ 718 for (
idx m = 0; m < D; ++m) {
719 for (
idx n = 0; n < D; ++n) {
722 #endif // WITH_OPENMP_ 726 for (
idx i = 0; i < Ks.size(); ++i)
727 EMN += Ks[i] * MN * adjoint(Ks[i]);
730 for (
idx a = 0; a < D; ++a) {
732 for (
idx b = 0; b < D; ++b) {
735 result(a * D + b, m * D + n) =
736 static_cast<cmat>(A * EMN * B).value();
741 EMN = cmat::Zero(D, D);
774 if (it.rows() != Ks[0].rows() || it.cols() != Ks[0].rows())
778 idx D =
static_cast<idx>(Ks[0].rows());
782 cmat MES = cmat::Zero(D * D, 1);
783 for (
idx a = 0; a < D; ++a)
786 cmat Omega = MES * adjoint(MES);
788 cmat result = cmat::Zero(D * D, D * D);
791 #pragma omp parallel for 792 #endif // WITH_OPENMP_ 793 for (
idx i = 0; i < Ks.size(); ++i) {
796 #endif // WITH_OPENMP_ 798 result += kron(cmat::Identity(D, D), Ks[i]) * Omega *
799 adjoint(kron(cmat::Identity(D, D), Ks[i]));
828 if (D * D != static_cast<idx>(A.rows()))
833 cmat evec = hevects(A);
834 std::vector<cmat> result;
836 for (
idx i = 0; i < D * D; ++i) {
837 if (std::abs(ev(i)) >
eps)
838 result.push_back(std::sqrt(std::abs(ev(i))) *
839 reshape(evec.col(i), D, D));
861 if (D * D != static_cast<idx>(A.rows()))
865 cmat result(D * D, D * D);
868 #pragma omp parallel for collapse(4) 869 #endif // WITH_OPENMP_ 870 for (
idx a = 0; a < D; ++a)
871 for (
idx b = 0; b < D; ++b)
872 for (
idx m = 0; m < D; ++m)
873 for (
idx n = 0; n < D; ++n)
874 result(a * D + b, m * D + n) = A(m * D + a, n * D + b);
895 if (D * D != static_cast<idx>(A.rows()))
899 cmat result(D * D, D * D);
902 #pragma omp parallel for collapse(4) 903 #endif // WITH_OPENMP_ 904 for (
idx a = 0; a < D; ++a)
905 for (
idx b = 0; b < D; ++b)
906 for (
idx m = 0; m < D; ++m)
907 for (
idx n = 0; n < D; ++n)
908 result(m * D + a, n * D + b) = A(a * D + b, m * D + n);
926 template <
typename Derived>
928 const std::vector<idx>& dims) {
929 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
942 if (dims.size() != 2)
959 auto worker = [&](
idx i,
idx j) noexcept->typename Derived::Scalar {
960 typename Derived::Scalar sum = 0;
961 for (
idx m = 0; m < DA; ++m)
962 sum += rA(m * DB + i) * std::conj(rA(m * DB + j));
968 #pragma omp parallel for collapse(2) 969 #endif // WITH_OPENMP_ 971 for (
idx j = 0; j < DB; ++j)
972 for (
idx i = 0; i < DB; ++i)
973 result(i, j) = worker(i, j);
984 auto worker = [&](
idx i,
idx j) noexcept->typename Derived::Scalar {
985 typename Derived::Scalar sum = 0;
986 for (
idx m = 0; m < DA; ++m)
987 sum += rA(m * DB + i, m * DB + j);
993 #pragma omp parallel for collapse(2) 994 #endif // WITH_OPENMP_ 996 for (
idx j = 0; j < DB; ++j)
997 for (
idx i = 0; i < DB; ++i)
998 result(i, j) = worker(i, j);
1020 template <
typename Derived>
1036 std::vector<idx> dims(2, d);
1054 template <
typename Derived>
1056 const std::vector<idx>& dims) {
1057 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
1070 if (dims.size() != 2)
1087 auto worker = [&](
idx i,
idx j) noexcept->typename Derived::Scalar {
1088 typename Derived::Scalar sum = 0;
1089 for (
idx m = 0; m < DB; ++m)
1090 sum += rA(i * DB + m) * std::conj(rA(j * DB + m));
1096 #pragma omp parallel for collapse(2) 1097 #endif // WITH_OPENMP_ 1099 for (
idx j = 0; j < DA; ++j)
1100 for (
idx i = 0; i < DA; ++i)
1101 result(i, j) = worker(i, j);
1113 #pragma omp parallel for collapse(2) 1114 #endif // WITH_OPENMP_ 1116 for (
idx j = 0; j < DA; ++j)
1117 for (
idx i = 0; i < DA; ++i)
1118 result(i, j) = trace(rA.block(i * DB, j * DB, DB, DB));
1140 template <
typename Derived>
1156 std::vector<idx> dims(2, d);
1175 template <
typename Derived>
1177 const std::vector<idx>& subsys,
1178 const std::vector<idx>& dims) {
1179 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
1196 idx D =
static_cast<idx>(rA.rows());
1197 idx N = dims.size();
1198 idx Nsubsys = subsys.size();
1199 idx Nsubsys_bar = N - Nsubsys;
1201 for (
idx i = 0; i < Nsubsys; ++i)
1202 Dsubsys *= dims[subsys[i]];
1203 idx Dsubsys_bar = D / Dsubsys;
1213 std::vector<idx> subsys_bar = complement(subsys, N);
1214 std::copy(std::begin(subsys_bar), std::end(subsys_bar),
1215 std::begin(Csubsys_bar));
1217 for (
idx i = 0; i < N; ++i) {
1220 for (
idx i = 0; i < Nsubsys; ++i) {
1221 Csubsys[i] = subsys[i];
1222 Cdimssubsys[i] = dims[subsys[i]];
1224 for (
idx i = 0; i < Nsubsys_bar; ++i) {
1225 Cdimssubsys_bar[i] = dims[subsys_bar[i]];
1238 if (subsys.size() == dims.size()) {
1239 result(0, 0) = (adjoint(rA) * rA).value();
1243 if (subsys.size() == 0)
1244 return rA * adjoint(rA);
1246 auto worker = [&](
idx i) noexcept->typename Derived::Scalar {
1256 Cmidxrowsubsys_bar);
1258 for (
idx k = 0; k < Nsubsys_bar; ++k) {
1259 Cmidxrow[Csubsys_bar[k]] = Cmidxrowsubsys_bar[k];
1260 Cmidxcol[Csubsys_bar[k]] = Cmidxcolsubsys_bar[k];
1262 typename Derived::Scalar sm = 0;
1263 for (
idx a = 0; a < Dsubsys; ++a) {
1267 for (
idx k = 0; k < Nsubsys; ++k)
1268 Cmidxrow[Csubsys[k]] = Cmidxcol[Csubsys[k]] =
1279 for (
idx j = 0; j < Dsubsys_bar; ++j)
1283 Cmidxcolsubsys_bar);
1285 #pragma omp parallel for 1286 #endif // WITH_OPENMP_ 1287 for (
idx i = 0; i < Dsubsys_bar; ++i) {
1288 result(i, j) = worker(i);
1301 if (subsys.size() == dims.size()) {
1302 result(0, 0) = rA.trace();
1306 if (subsys.size() == 0)
1309 auto worker = [&](
idx i) noexcept->typename Derived::Scalar {
1319 Cmidxrowsubsys_bar);
1321 for (
idx k = 0; k < Nsubsys_bar; ++k) {
1322 Cmidxrow[Csubsys_bar[k]] = Cmidxrowsubsys_bar[k];
1323 Cmidxcol[Csubsys_bar[k]] = Cmidxcolsubsys_bar[k];
1325 typename Derived::Scalar sm = 0;
1326 for (
idx a = 0; a < Dsubsys; ++a) {
1330 for (
idx k = 0; k < Nsubsys; ++k)
1331 Cmidxrow[Csubsys[k]] = Cmidxcol[Csubsys[k]] =
1342 for (
idx j = 0; j < Dsubsys_bar; ++j)
1346 Cmidxcolsubsys_bar);
1348 #pragma omp parallel for 1349 #endif // WITH_OPENMP_ 1350 for (
idx i = 0; i < Dsubsys_bar; ++i) {
1351 result(i, j) = worker(i);
1376 template <
typename Derived>
1378 const std::vector<idx>& subsys,
1380 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
1394 std::vector<idx> dims(N, d);
1396 return ptrace(rA, subsys, dims);
1412 template <
typename Derived>
1414 ptranspose(
const Eigen::MatrixBase<Derived>& A,
const std::vector<idx>& subsys,
1415 const std::vector<idx>& dims) {
1416 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
1433 idx D =
static_cast<idx>(rA.rows());
1434 idx N = dims.size();
1435 idx Nsubsys = subsys.size();
1441 for (
idx i = 0; i < N; ++i)
1443 for (
idx i = 0; i < Nsubsys; ++i)
1444 Csubsys[i] = subsys[i];
1455 if (subsys.size() == dims.size())
1456 return (rA * adjoint(rA)).transpose();
1458 if (subsys.size() == 0)
1459 return rA * adjoint(rA);
1461 auto worker = [&](
idx i) noexcept->typename Derived::Scalar {
1466 for (
idx k = 0; k < N; ++k)
1467 midxcoltmp[k] = Cmidxcol[k];
1472 for (
idx k = 0; k < Nsubsys; ++k)
1473 std::swap(midxcoltmp[Csubsys[k]], midxrow[Csubsys[k]]);
1480 for (
idx j = 0; j < D; ++j) {
1485 #pragma omp parallel for 1486 #endif // WITH_OPENMP_ 1487 for (
idx i = 0; i < D; ++i)
1488 result(i, j) = worker(i);
1500 if (subsys.size() == dims.size())
1501 return rA.transpose();
1503 if (subsys.size() == 0)
1506 auto worker = [&](
idx i) noexcept->typename Derived::Scalar {
1511 for (
idx k = 0; k < N; ++k)
1512 midxcoltmp[k] = Cmidxcol[k];
1517 for (
idx k = 0; k < Nsubsys; ++k)
1518 std::swap(midxcoltmp[Csubsys[k]], midxrow[Csubsys[k]]);
1525 for (
idx j = 0; j < D; ++j) {
1530 #pragma omp parallel for 1531 #endif // WITH_OPENMP_ 1532 for (
idx i = 0; i < D; ++i)
1533 result(i, j) = worker(i);
1556 template <
typename Derived>
1558 ptranspose(
const Eigen::MatrixBase<Derived>& A,
const std::vector<idx>& subsys,
1560 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
1574 std::vector<idx> dims(N, d);
1591 template <
typename Derived>
1593 syspermute(
const Eigen::MatrixBase<Derived>& A,
const std::vector<idx>& perm,
1594 const std::vector<idx>& dims) {
1595 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
1612 if (perm.size() != dims.size())
1616 idx D =
static_cast<idx>(rA.rows());
1617 idx N = dims.size();
1632 for (
idx i = 0; i < N; ++i) {
1636 result.resize(D, 1);
1638 auto worker = [&Cdims, &Cperm, N ](
idx i) noexcept->idx {
1648 for (
idx k = 0; k < N; ++k) {
1649 permdims[k] = Cdims[Cperm[k]];
1650 midxtmp[k] = midx[Cperm[k]];
1657 #pragma omp parallel for 1658 #endif // WITH_OPENMP_ 1659 for (
idx i = 0; i < D; ++i)
1660 result(worker(i)) = rA(i);
1675 for (
idx i = 0; i < N; ++i) {
1677 Cdims[i + N] = dims[i];
1679 Cperm[i + N] = perm[i] + N;
1681 result.resize(D * D, 1);
1684 Eigen::Map<dyn_mat<typename Derived::Scalar>>(
1685 const_cast<typename Derived::Scalar*
>(rA.data()), D * D, 1);
1687 auto worker = [&Cdims, &Cperm, N ](
idx i) noexcept->idx {
1697 for (
idx k = 0; k < 2 * N; ++k) {
1698 permdims[k] = Cdims[Cperm[k]];
1699 midxtmp[k] = midx[Cperm[k]];
1706 #pragma omp parallel for 1707 #endif // WITH_OPENMP_ 1708 for (
idx i = 0; i < D * D; ++i)
1709 result(worker(i)) = rA(i);
1711 return reshape(result, D, D);
1730 template <
typename Derived>
1732 syspermute(
const Eigen::MatrixBase<Derived>& A,
const std::vector<idx>& perm,
1734 const typename Eigen::MatrixBase<Derived>::EvalReturnType& rA = A.derived();
1748 std::vector<idx> dims(N, d);
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
constexpr idx maxn
Maximum number of allowed qubits/qudits (subsystems)
Definition: constants.h:73
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:819
Eigen::MatrixXd dmat
Real (double precision) dynamic Eigen matrix.
Definition: types.h:69
constexpr double eps
Used to decide whether a number or expression in double precision is zero or not. ...
Definition: constants.h:66
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:1593
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:363
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:207
Invalid dimension(s) exception.
Definition: exception.h:267
cmat choi2super(const cmat &A)
Converts Choi matrix to superoperator matrix.
Definition: operations.h:852
cmat super2choi(const cmat &A)
Converts superoperator matrix to Choi matrix.
Definition: operations.h:886
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
dyn_mat< typename Derived::Scalar > ptrace1(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Partial trace.
Definition: operations.h:927
Invalid permutation exception.
Definition: exception.h:378
dyn_mat< typename Derived::Scalar > ptrace2(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &dims)
Partial trace.
Definition: operations.h:1055
Matrix mismatch subsystems exception.
Definition: exception.h:252
Dimension(s) mismatch column vector size exception.
Definition: exception.h:314
Type mismatch exception.
Definition: exception.h:528
idx get_dim_subsys(idx sz, idx N)
Definition: util.h:378
cmat kraus2choi(const std::vector< cmat > &Ks)
Choi matrix.
Definition: operations.h:764
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:1176
bool check_square_mat(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:99
cmat kraus2super(const std::vector< cmat > &Ks)
Superoperator matrix.
Definition: operations.h:693
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
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
bool check_dims(const std::vector< idx > &dims)
Definition: util.h:134
Permutation mismatch dimensions exception.
Definition: exception.h:394
std::size_t idx
Non-negative integer index.
Definition: types.h:39
Matrix is not square exception.
Definition: exception.h:147
Object has zero size exception.
Definition: exception.h:132