44 template <
typename Derived>
45 dyn_mat<typename Derived::Scalar>
53 return rA.transpose();
63 template <
typename Derived>
72 return rA.conjugate();
82 template <
typename Derived>
103 template <
typename Derived>
123 template <
typename Derived>
124 typename Derived::Scalar
trace(
const Eigen::MatrixBase<Derived>& A) {
144 template <
typename Derived>
145 typename Derived::Scalar
det(
const Eigen::MatrixBase<Derived>& A) {
155 return rA.determinant();
167 template <
typename Derived>
168 typename Derived::Scalar
logdet(
const Eigen::MatrixBase<Derived>& A) {
182 Eigen::PartialPivLU<dyn_mat<typename Derived::Scalar>> lu(rA);
184 lu.matrixLU().template triangularView<Eigen::Upper>();
185 typename Derived::Scalar result = std::log(U(0, 0));
187 for (
idx i = 1; i < static_cast<idx>(rA.rows()); ++i)
188 result += std::log(U(i, i));
200 template <
typename Derived>
201 typename Derived::Scalar
sum(
const Eigen::MatrixBase<Derived>& A) {
221 template <
typename Derived>
222 typename Derived::Scalar
prod(
const Eigen::MatrixBase<Derived>& A) {
241 template <
typename Derived>
242 double norm(
const Eigen::MatrixBase<Derived>& A) {
253 return (rA.template cast<cplx>()).norm();
264 template <
typename Derived>
265 std::pair<dyn_col_vect<cplx>,
cmat>
267 eig(
const Eigen::MatrixBase<Derived>& A) {
281 Eigen::ComplexEigenSolver<cmat> es(rA.template cast<cplx>());
283 return std::make_pair(es.eigenvalues(), es.eigenvectors());
293 template <
typename Derived>
308 return eig(rA).first;
318 template <
typename Derived>
333 Eigen::ComplexEigenSolver<cmat> es(rA.template cast<cplx>());
335 return eig(rA).second;
346 template <
typename Derived>
347 std::pair<dyn_col_vect<double>,
cmat>
349 heig(
const Eigen::MatrixBase<Derived>& A) {
363 Eigen::SelfAdjointEigenSolver<cmat> es(rA.template cast<cplx>());
365 return std::make_pair(es.eigenvalues(), es.eigenvectors());
375 template <
typename Derived>
390 return heig(rA).first;
400 template <
typename Derived>
415 return heig(rA).second;
427 template <
typename Derived>
428 std::tuple<cmat, dyn_col_vect<double>,
cmat>
430 svd(
const Eigen::MatrixBase<Derived>& A) {
440 Eigen::JacobiSVD<dyn_mat<typename Derived::Scalar>> sv(
441 rA, Eigen::DecompositionOptions::ComputeFullU |
442 Eigen::DecompositionOptions::ComputeFullV);
444 return std::make_tuple(sv.matrixU(), sv.singularValues(), sv.matrixV());
454 template <
typename Derived>
465 Eigen::JacobiSVD<dyn_mat<typename Derived::Scalar>> sv(rA);
467 return sv.singularValues();
477 template <
typename Derived>
488 Eigen::JacobiSVD<dyn_mat<typename Derived::Scalar>> sv(
489 rA, Eigen::DecompositionOptions::ComputeFullU);
501 template <
typename Derived>
512 Eigen::JacobiSVD<dyn_mat<typename Derived::Scalar>> sv(
513 rA, Eigen::DecompositionOptions::ComputeFullV);
527 template <
typename Derived>
542 Eigen::ComplexEigenSolver<cmat> es(rA.template cast<cplx>());
545 for (
idx i = 0; i < static_cast<idx>(evals.rows()); ++i)
548 cmat evalsdiag = evals.asDiagonal();
550 return evects * evalsdiag * evects.inverse();
559 template <
typename Derived>
574 return funm(rA, &std::sqrt);
583 template <
typename Derived>
607 template <
typename Derived>
622 return funm(rA, &std::exp);
631 template <
typename Derived>
646 return funm(rA, &std::log);
655 template <
typename Derived>
670 return funm(rA, &std::sin);
679 template <
typename Derived>
694 return funm(rA, &std::cos);
708 template <
typename Derived>
724 if (real(z) == 0 && imag(z) == 0)
725 return cmat::Identity(rA.rows(), rA.rows());
727 Eigen::ComplexEigenSolver<cmat> es(rA.template cast<cplx>());
730 for (
idx i = 0; i < static_cast<idx>(evals.rows()); ++i)
733 cmat evalsdiag = evals.asDiagonal();
735 return evects * evalsdiag * evects.inverse();
750 template <
typename Derived>
778 for (; n > 0; n /= 2) {
780 result = (result * cA).eval();
781 cA = (cA * cA).eval();
795 template <
typename Derived>
796 double schatten(
const Eigen::MatrixBase<Derived>& A,
double p) {
813 for (
idx i = 0; i < static_cast<idx>(sv.rows()); ++i)
814 result += std::pow(sv[i], p);
816 return std::pow(result, 1. / p);
829 template <
typename OutputScalar,
typename Derived>
831 cwise(
const Eigen::MatrixBase<Derived>& A,
832 OutputScalar (*f)(
const typename Derived::Scalar&)) {
845 #pragma omp parallel for collapse(2) 846 #endif // WITH_OPENMP_ 848 for (
idx j = 0; j < static_cast<idx>(rA.cols()); ++j)
849 for (
idx i = 0; i < static_cast<idx>(rA.rows()); ++i)
850 result(i, j) = (*f)(rA(i, j));
867 template <
typename T>
882 template <
typename T,
typename... Args>
896 template <
typename Derived>
909 for (
idx i = 1; i < As.size(); ++i) {
910 result =
kron(result, As[i]);
928 template <
typename Derived>
930 kron(
const std::initializer_list<Derived>& As) {
931 return kron(std::vector<Derived>(As));
943 template <
typename Derived>
959 std::vector<dyn_mat<typename Derived::Scalar>> As(n, rA);
976 template <
typename T>
991 template <
typename T,
typename... Args>
1005 template <
typename Derived>
1012 for (
auto&& it : As)
1017 idx total_rows = 0, total_cols = 0;
1018 for (
idx i = 0; i < As.size(); ++i) {
1019 total_rows +=
static_cast<idx>(As[i].rows());
1020 total_cols +=
static_cast<idx>(As[i].cols());
1025 idx cur_row = 0, cur_col = 0;
1026 for (
idx i = 0; i < As.size(); ++i) {
1027 result.block(cur_row, cur_col, As[i].rows(), As[i].cols()) = As[i];
1028 cur_row +=
static_cast<idx>(As[i].rows());
1029 cur_col +=
static_cast<idx>(As[i].cols());
1047 template <
typename Derived>
1049 dirsum(
const std::initializer_list<Derived>& As) {
1050 return dirsum(std::vector<Derived>(As));
1062 template <
typename Derived>
1078 std::vector<dyn_mat<typename Derived::Scalar>> As(n, rA);
1094 template <
typename Derived>
1099 idx Arows =
static_cast<idx>(rA.rows());
1100 idx Acols =
static_cast<idx>(rA.cols());
1108 if (Arows * Acols != rows * cols)
1112 return Eigen::Map<dyn_mat<typename Derived::Scalar>>(
1113 const_cast<typename Derived::Scalar*
>(rA.data()), rows, cols);
1128 template <
typename Derived1,
typename Derived2>
1130 const Eigen::MatrixBase<Derived2>& B) {
1137 if (!std::is_same<
typename Derived1::Scalar,
1138 typename Derived2::Scalar>::value)
1150 if (rA.rows() != rB.rows())
1154 return rA * rB - rB * rA;
1169 template <
typename Derived1,
typename Derived2>
1172 const Eigen::MatrixBase<Derived2>& B) {
1179 if (!std::is_same<
typename Derived1::Scalar,
1180 typename Derived2::Scalar>::value)
1192 if (rA.rows() != rB.rows())
1196 return rA * rB + rB * rA;
1209 template <
typename Derived>
1224 double normA =
norm(rA);
1226 return rA *
adjoint(rA) / (normA * normA);
1238 template <
typename Derived>
1246 for (
auto&& it : As)
1255 for (
auto&& it : As)
1256 if (it.rows() != As[0].rows() || it.cols() != 1)
1266 std::vector<dyn_mat<typename Derived::Scalar>> outvecs;
1269 for (pos = 0; pos < As.size(); ++pos) {
1272 outvecs.push_back(As[pos]);
1278 for (
idx i = pos + 1; i < As.size(); ++i) {
1279 cut -=
prj(outvecs[i - 1 - pos]);
1281 outvecs.push_back(vi);
1287 for (
auto&& it : outvecs) {
1288 double normA =
norm(it);
1291 result.col(cnt) = it / normA;
1296 return result.block(0, 0, As[0].rows(), cnt);
1307 template <
typename Derived>
1309 grams(
const std::initializer_list<Derived>& As) {
1310 return grams(std::vector<Derived>(As));
1320 template <
typename Derived>
1330 std::vector<dyn_mat<typename Derived::Scalar>> input;
1332 for (
idx i = 0; i < static_cast<idx>(rA.cols()); ++i)
1333 input.push_back(rA.col(i));
1335 return grams<dyn_mat<typename Derived::Scalar>>(input);
1354 if (n >= std::accumulate(std::begin(dims), std::end(dims),
1355 static_cast<idx>(1), std::multiplies<idx>()))
1363 return std::vector<idx>(result, result + dims.size());
1377 const std::vector<idx>& dims) {
1383 for (
idx i = 0; i < dims.size(); ++i)
1384 if (midx[i] >= dims[i])
1405 inline ket mket(
const std::vector<idx>& mask,
const std::vector<idx>& dims) {
1406 idx N = mask.size();
1408 idx D = std::accumulate(std::begin(dims), std::end(dims),
1409 static_cast<idx>(1), std::multiplies<idx>());
1420 if (mask.size() != dims.size())
1423 for (
idx i = 0; i < N; ++i)
1424 if (mask[i] >= dims[i])
1428 ket result = ket::Zero(D);
1449 idx N = mask.size();
1450 idx D =
static_cast<idx>(std::llround(std::pow(d, N)));
1463 for (
idx i = 0; i < N; ++i)
1468 ket result = ket::Zero(D);
1469 std::vector<idx> dims(N, d);
1491 inline cmat mprj(
const std::vector<idx>& mask,
const std::vector<idx>& dims) {
1492 idx N = mask.size();
1494 idx D = std::accumulate(std::begin(dims), std::end(dims),
1495 static_cast<idx>(1), std::multiplies<idx>());
1506 if (mask.size() != dims.size())
1509 for (
idx i = 0; i < N; ++i)
1510 if (mask[i] >= dims[i])
1514 cmat result = cmat::Zero(D, D);
1516 result(pos, pos) = 1;
1537 idx N = mask.size();
1538 idx D =
static_cast<idx>(std::llround(std::pow(d, N)));
1551 for (
idx i = 0; i < N; ++i)
1556 cmat result = cmat::Zero(D, D);
1557 std::vector<idx> dims(N, d);
1559 result(pos, pos) = 1;
1572 template <
typename InputIterator>
1573 std::vector<double>
abssq(InputIterator first, InputIterator last) {
1574 std::vector<double> weights(std::distance(first, last));
1575 std::transform(first, last, std::begin(weights),
1587 template <
typename Container>
1598 return abssq(std::begin(c), std::end(c));
1607 template <
typename Derived>
1608 std::vector<double>
abssq(
const Eigen::MatrixBase<Derived>& A) {
1618 return abssq(rA.data(), rA.data() + rA.size());
1629 template <
typename InputIterator>
1630 typename std::iterator_traits<InputIterator>::value_type
1631 sum(InputIterator first, InputIterator last) {
1632 using value_type =
typename std::iterator_traits<InputIterator>::value_type;
1634 return std::accumulate(first, last, static_cast<value_type>(0));
1644 template <
typename Container>
1645 typename Container::value_type
1648 return sum(std::begin(c), std::end(c));
1659 template <
typename InputIterator>
1660 typename std::iterator_traits<InputIterator>::value_type
1661 prod(InputIterator first, InputIterator last) {
1662 using value_type =
typename std::iterator_traits<InputIterator>::value_type;
1664 return std::accumulate(first, last, static_cast<value_type>(1),
1665 std::multiplies<value_type>());
1675 template <
typename Container>
1676 typename Container::value_type
1679 return prod(std::begin(c), std::end(c));
1694 template <
typename Derived>
1714 for (
idx k = 0; k < static_cast<idx>(rA.rows()); ++k) {
1715 if (std::abs(tmp_evals(k)) >
eps) {
1716 result = tmp_evects.col(k);
1732 template <
typename T>
1736 if (N < subsys.size())
1740 std::vector<T> all(N);
1741 std::vector<T> subsys_bar(N - subsys.size());
1743 std::iota(std::begin(all), std::end(all), 0);
1744 std::sort(std::begin(subsys), std::end(subsys));
1745 std::set_difference(std::begin(all), std::end(all), std::begin(subsys),
1746 std::end(subsys), std::begin(subsys_bar));
1761 template <
typename Derived>
1762 std::vector<double>
rho2bloch(
const Eigen::MatrixBase<Derived>& A) {
1772 std::vector<double> result(3);
1773 cmat X(2, 2), Y(2, 2), Z(2, 2);
1775 Y << 0, -1_i, 1_i, 0;
1777 result[0] = std::real(
trace(rA * X));
1778 result[1] = std::real(
trace(rA * Y));
1779 result[2] = std::real(
trace(rA * Z));
1798 "r is not a 3-dimensional vector!");
1801 cmat X(2, 2), Y(2, 2), Z(2, 2), Id2(2, 2);
1803 Y << 0, -1_i, 1_i, 0;
1807 return (Id2 + r[0] * X + r[1] * Y + r[2] * Z) / 2.;
1810 inline namespace literals {
1821 template <
char... Bits>
1823 constexpr
idx n =
sizeof...(Bits);
1824 constexpr
char bits[n + 1] = {Bits...,
'\0'};
1825 qpp::ket q = qpp::ket::Zero(std::pow(2, n));
1831 for (
idx i = 0; i < n; ++i) {
1832 if (bits[i] !=
'0' && bits[i] !=
'1')
1837 pos = std::stoi(bits,
nullptr, 2);
1852 template <
char... Bits>
1854 constexpr
idx n =
sizeof...(Bits);
1855 constexpr
char bits[n + 1] = {Bits...,
'\0'};
1856 qpp::bra q = qpp::ket::Zero(std::pow(2, n));
1862 for (
idx i = 0; i < n; ++i) {
1863 if (bits[i] !=
'0' && bits[i] !=
'1')
1868 pos = std::stoi(bits,
nullptr, 2);
1885 template <
char... Bits>
1887 constexpr
idx n =
sizeof...(Bits);
1888 constexpr
char bits[n + 1] = {Bits...,
'\0'};
1893 for (
idx i = 0; i < n; ++i) {
1894 if (bits[i] !=
'0' && bits[i] !=
'1')
1899 return kron(
operator""_ket<Bits...>(),
operator""_bra<Bits...>());
bool check_nonzero_size(const T &x) noexcept
Definition: util.h:123
std::pair< dyn_col_vect< double >, cmat > heig(const Eigen::MatrixBase< Derived > &A)
Full eigen decomposition of Hermitian expression.
Definition: functions.h:349
Dimensions not equal exception.
Definition: exception.h:284
Dimension(s) mismatch matrix size exception.
Definition: exception.h:300
Derived::Scalar logdet(const Eigen::MatrixBase< Derived > &A)
Logarithm of the determinant.
Definition: functions.h:168
constexpr idx maxn
Maximum number of allowed qubits/qudits (subsystems)
Definition: constants.h:75
cmat absm(const Eigen::MatrixBase< Derived > &A)
Matrix absolute value.
Definition: functions.h:584
dyn_mat< typename Derived1::Scalar > kron2(const Eigen::MatrixBase< Derived1 > &A, const Eigen::MatrixBase< Derived2 > &B)
Definition: util.h:273
bool check_cvector(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:117
Custom exception.
Definition: exception.h:571
cmat expm(const Eigen::MatrixBase< Derived > &A)
Matrix exponential.
Definition: functions.h:608
Derived::Scalar prod(const Eigen::MatrixBase< Derived > &A)
Element-wise product of A.
Definition: functions.h:222
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 T::Scalar > dirsum(const T &head)
Direct sum.
Definition: functions.h:977
Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic > dyn_mat
Dynamic Eigen matrix over the field specified by Scalar.
Definition: types.h:81
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
std::vector< double > rho2bloch(const Eigen::MatrixBase< Derived > &A)
Computes the 3-dimensional real Bloch vector corresponding to the qubit density matrix A...
Definition: functions.h:1762
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 > anticomm(const Eigen::MatrixBase< Derived1 > &A, const Eigen::MatrixBase< Derived2 > &B)
Anti-commutator.
Definition: functions.h:1171
Matrix is not a column vector exception.
Definition: exception.h:164
Quantum++ main namespace.
Definition: codes.h:35
std::pair< dyn_col_vect< cplx >, cmat > eig(const Eigen::MatrixBase< Derived > &A)
Full eigen decomposition.
Definition: functions.h:267
double norm(const Eigen::MatrixBase< Derived > &A)
Frobenius norm.
Definition: functions.h:242
cmat spectralpowm(const Eigen::MatrixBase< Derived > &A, const cplx z)
Matrix power.
Definition: functions.h:709
dyn_col_vect< double > hevals(const Eigen::MatrixBase< Derived > &A)
Hermitian eigenvalues.
Definition: functions.h:376
dyn_mat< OutputScalar > cwise(const Eigen::MatrixBase< Derived > &A, OutputScalar(*f)(const typename Derived::Scalar &))
Functor.
Definition: functions.h:831
Invalid dimension(s) exception.
Definition: exception.h:269
Checks whether T is compatible with an STL-like iterable container.
Definition: traits.h:71
dyn_col_vect< typename Derived::Scalar > rho2pure(const Eigen::MatrixBase< Derived > &A)
Finds the pure state representation of a matrix proportional to a projector onto a pure state...
Definition: functions.h:1696
dyn_mat< typename Derived::Scalar > dirsumpow(const Eigen::MatrixBase< Derived > &A, idx n)
Direct sum power.
Definition: functions.h:1063
std::vector< T > complement(std::vector< T > subsys, idx N)
Constructs the complement of a subsystem vector.
Definition: functions.h:1733
cmat evects(const Eigen::MatrixBase< Derived > &A)
Eigenvectors.
Definition: functions.h:319
cmat funm(const Eigen::MatrixBase< Derived > &A, cplx(*f)(const cplx &))
Functional calculus f(A)
Definition: functions.h:528
double schatten(const Eigen::MatrixBase< Derived > &A, double p)
Schatten matrix norm.
Definition: functions.h:796
cmat svdU(const Eigen::MatrixBase< Derived > &A)
Left singular vectors.
Definition: functions.h:478
cmat bloch2rho(const std::vector< double > &r)
Computes the density matrix corresponding to the 3-dimensional real Bloch vector r.
Definition: functions.h:1792
std::vector< idx > n2multiidx(idx n, const std::vector< idx > &dims)
Non-negative integer index to multi-index.
Definition: functions.h:1348
dyn_mat< typename Derived::Scalar > adjoint(const Eigen::MatrixBase< Derived > &A)
Adjoint.
Definition: functions.h:83
idx multiidx2n(const std::vector< idx > &midx, const std::vector< idx > &dims)
Multi-index to non-negative integer index.
Definition: functions.h:1376
dyn_col_vect< double > svals(const Eigen::MatrixBase< Derived > &A)
Singular values.
Definition: functions.h:455
dyn_mat< typename T::Scalar > kron(const T &head)
Kronecker product.
Definition: functions.h:868
std::complex< double > cplx
Complex number in double precision.
Definition: types.h:49
cmat mprj(const std::vector< idx > &mask, const std::vector< idx > &dims)
Projector onto multi-partite qudit ket.
Definition: functions.h:1491
Derived::Scalar det(const Eigen::MatrixBase< Derived > &A)
Determinant.
Definition: functions.h:145
dyn_mat< typename Derived::Scalar > transpose(const Eigen::MatrixBase< Derived > &A)
Transpose.
Definition: functions.h:46
dyn_mat< typename Derived::Scalar > kronpow(const Eigen::MatrixBase< Derived > &A, idx n)
Kronecker power.
Definition: functions.h:944
dyn_mat< typename Derived::Scalar > prj(const Eigen::MatrixBase< Derived > &A)
Projector.
Definition: functions.h:1210
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
dyn_col_vect< cplx > evals(const Eigen::MatrixBase< Derived > &A)
Eigenvalues.
Definition: functions.h:294
cmat hevects(const Eigen::MatrixBase< Derived > &A)
Hermitian eigenvectors.
Definition: functions.h:401
dyn_mat< typename Derived::Scalar > conjugate(const Eigen::MatrixBase< Derived > &A)
Complex conjugate.
Definition: functions.h:65
constexpr double infty
Used to denote infinity in double precision.
Definition: constants.h:89
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
cmat svdV(const Eigen::MatrixBase< Derived > &A)
Right singular vectors.
Definition: functions.h:502
dyn_mat< typename Derived1::Scalar > dirsum2(const Eigen::MatrixBase< Derived1 > &A, const Eigen::MatrixBase< Derived2 > &B)
Definition: util.h:317
cmat cosm(const Eigen::MatrixBase< Derived > &A)
Matrix cos.
Definition: functions.h:680
Parameter out of range exception.
Definition: exception.h:515
Eigen::Matrix< Scalar, Eigen::Dynamic, 1 > dyn_col_vect
Dynamic Eigen column vector over the field specified by Scalar.
Definition: types.h:93
cmat sinm(const Eigen::MatrixBase< Derived > &A)
Matrix sin.
Definition: functions.h:656
std::vector< double > abssq(InputIterator first, InputIterator last)
Computes the absolute values squared of an STL-like range of complex numbers.
Definition: functions.h:1573
Matrix is not 2 x 2 exception.
Definition: exception.h:411
dyn_mat< typename Derived::Scalar > grams(const std::vector< Derived > &As)
Gram-Schmidt orthogonalization.
Definition: functions.h:1239
dyn_mat< typename Derived1::Scalar > comm(const Eigen::MatrixBase< Derived1 > &A, const Eigen::MatrixBase< Derived2 > &B)
Commutator.
Definition: functions.h:1129
Eigen::MatrixXcd cmat
Complex (double precision) dynamic Eigen matrix.
Definition: types.h:64
ket mket(const std::vector< idx > &mask, const std::vector< idx > &dims)
Multi-partite qudit ket.
Definition: functions.h:1405
dyn_mat< typename Derived::Scalar > inverse(const Eigen::MatrixBase< Derived > &A)
Inverse.
Definition: functions.h:104
cmat logm(const Eigen::MatrixBase< Derived > &A)
Matrix logarithm.
Definition: functions.h:632
bool check_dims(const std::vector< idx > &dims)
Definition: util.h:134
cmat sqrtm(const Eigen::MatrixBase< Derived > &A)
Matrix square root.
Definition: functions.h:560
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
std::tuple< cmat, dyn_col_vect< double >, cmat > svd(const Eigen::MatrixBase< Derived > &A)
Full singular value decomposition.
Definition: functions.h:430
Object has zero size exception.
Definition: exception.h:134
bool check_qubit_matrix(const Eigen::MatrixBase< Derived > &A) noexcept
Definition: util.h:236