42 template<
typename Derived>
43 double entropy(
const Eigen::MatrixBase<Derived>& A)
57 for (
idx i = 0; i < static_cast<idx>(ev.rows()); ++i)
59 result -= ev(i) * std::log2(ev(i));
70 double entropy(
const std::vector<double>& prob)
77 for (
idx i = 0; i < prob.size(); ++i)
78 if (std::abs(prob[i]) != 0)
79 result -= std::abs(prob[i]) * std::log2(std::abs(prob[i]));
96 template<
typename Derived>
97 double renyi(
const Eigen::MatrixBase<Derived>& A,
double alpha)
113 return std::log2(static_cast<double>( rA.rows()));
119 return -std::log2(
svals(rA)[0]);
123 for (
idx i = 0; i < static_cast<idx>(sv.rows()); ++i)
125 result += std::pow(sv(i), alpha);
127 return std::log2(result) / (1 - alpha);
142 double renyi(
const std::vector<double>& prob,
double alpha)
152 return std::log2(static_cast<double>( prob.size()));
160 for (
idx i = 0; i < prob.size(); ++i)
161 if (std::abs(prob[i]) > max)
162 max = std::abs(prob[i]);
164 return -std::log2(max);
168 for (
idx i = 0; i < prob.size(); ++i)
169 if (std::abs(prob[i]) != 0)
170 result += std::pow(std::abs(prob[i]), alpha);
172 return std::log2(result) / (1 - alpha);
187 template<
typename Derived>
188 double tsallis(
const Eigen::MatrixBase<Derived>& A,
double q)
204 return entropy(rA) * std::log(2.);
208 for (
idx i = 0; i < static_cast<idx>(ev.rows()); ++i)
210 result += std::pow(ev(i), q);
212 return (result - 1) / (1 - q);
227 double tsallis(
const std::vector<double>& prob,
double q)
237 return entropy(prob) * std::log(2.);
240 for (
idx i = 0; i < prob.size(); ++i)
241 if (std::abs(prob[i]) != 0)
242 result += std::pow(std::abs(prob[i]), q);
244 return (result - 1) / (1 - q);
256 template<
typename Derived>
258 const std::vector<idx>& subsysA,
259 const std::vector<idx>& subsysB,
260 const std::vector<idx>& dims)
291 std::vector<idx> full_system(dims.size());
292 std::iota(std::begin(full_system), std::end(full_system), 0);
295 std::vector<idx> subsysAsorted{subsysA};
296 std::vector<idx> subsysBsorted{subsysB};
299 std::sort(std::begin(subsysAsorted), std::end(subsysAsorted));
300 std::sort(std::begin(subsysBsorted), std::end(subsysBsorted));
303 std::vector<idx> subsysAB;
304 std::set_union(std::begin(subsysAsorted), std::end(subsysAsorted),
305 std::begin(subsysBsorted), std::end(subsysBsorted),
306 std::back_inserter(subsysAB));
309 std::vector<idx> subsysAbar =
complement(subsysA, dims.size());
310 std::vector<idx> subsysBbar =
complement(subsysB, dims.size());;
311 std::vector<idx> subsysABbar =
complement(subsysAB, dims.size());
329 template<
typename Derived>
331 const std::vector<idx>& subsysA,
332 const std::vector<idx>& subsysB,
344 static_cast<idx>(std::llround(std::log2(rA.rows()) /
346 std::vector<idx> dims(n, d);
double qmutualinfo(const Eigen::MatrixBase< Derived > &A, const std::vector< idx > &subsysA, const std::vector< idx > &subsysB, const std::vector< idx > &dims)
Quantum mutual information between 2 subsystems of a composite system.
Definition: entropies.h:257
double tsallis(const Eigen::MatrixBase< Derived > &A, double q)
Tsallis- entropy of the density matrix A, for .
Definition: entropies.h:188
bool _check_subsys_match_dims(const std::vector< idx > &subsys, const std::vector< idx > &dims)
Definition: util.h:189
bool _check_dims_match_mat(const std::vector< idx > &dims, const Eigen::MatrixBase< Derived > &A)
Definition: util.h:141
Eigen::MatrixXd dmat
Real (double precision) dynamic Eigen matrix.
Definition: types.h:61
Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic > dyn_mat
Dynamic Eigen matrix over the field specified by Scalar.
Definition: types.h:73
Quantum++ main namespace.
Definition: codes.h:30
bool _check_square_mat(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:83
std::vector< T > complement(std::vector< T > subsys, idx N)
Constructs the complement of a subsystem vector.
Definition: functions.h:1652
dyn_col_vect< double > svals(const Eigen::MatrixBase< Derived > &A)
Singular values.
Definition: functions.h:426
double entropy(const Eigen::MatrixBase< Derived > &A)
von-Neumann entropy of the density matrix A
Definition: entropies.h:43
Generates custom exceptions, used when validating function parameters.
Definition: exception.h:39
double renyi(const Eigen::MatrixBase< Derived > &A, double alpha)
Renyi- entropy of the density matrix A, for .
Definition: entropies.h:97
bool _check_nonzero_size(const T &x) noexcept
Definition: util.h:119
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:1079
constexpr double infty
Used to denote infinity in double precision.
Definition: constants.h:96
Eigen::MatrixXcd cmat
Complex (double precision) dynamic Eigen matrix.
Definition: types.h:56
std::size_t idx
Non-negative integer index.
Definition: types.h:36
bool _check_dims(const std::vector< idx > &dims)
Definition: util.h:125