38 template<
typename Derived>
39 double entropy(
const Eigen::MatrixBase<Derived>& A)
56 for (
idx i = 0; i < static_cast<idx>(ev.rows()); ++i)
58 result -= ev(i) * std::log2(ev(i));
69 inline double entropy(
const std::vector<double>& prob)
79 for (
idx i = 0; i < prob.size(); ++i)
80 if (std::abs(prob[i]) != 0)
81 result -= std::abs(prob[i]) * std::log2(std::abs(prob[i]));
98 template<
typename Derived>
99 double renyi(
const Eigen::MatrixBase<Derived>& A,
double alpha)
118 return std::log2(rA.rows());
124 return -std::log2(
svals(rA)[0]);
128 for (
idx i = 0; i < static_cast<idx>(sv.rows()); ++i)
130 result += std::pow(sv(i), alpha);
132 return std::log2(result) / (1 - alpha);
147 inline double renyi(
const std::vector<double>& prob,
double alpha)
159 return std::log2(prob.size());
168 for (
idx i = 0; i < prob.size(); ++i)
169 if (std::abs(prob[i]) > max)
170 max = std::abs(prob[i]);
172 return -std::log2(max);
176 for (
idx i = 0; i < prob.size(); ++i)
177 if (std::abs(prob[i]) != 0)
178 result += std::pow(std::abs(prob[i]), alpha);
180 return std::log2(result) / (1 - alpha);
195 template<
typename Derived>
196 double tsallis(
const Eigen::MatrixBase<Derived>& A,
double q)
215 return entropy(rA) * std::log(2);
219 for (
idx i = 0; i < static_cast<idx>(ev.rows()); ++i)
221 result += std::pow(ev(i), q);
223 return (result - 1) / (1 - q);
238 inline double tsallis(
const std::vector<double>& prob,
double q)
251 return entropy(prob) * std::log(2);
254 for (
idx i = 0; i < prob.size(); ++i)
255 if (std::abs(prob[i]) != 0)
256 result += std::pow(std::abs(prob[i]), q);
258 return (result - 1) / (1 - q);
270 template<
typename Derived>
272 const std::vector<idx>& subsysA,
273 const std::vector<idx>& subsysB,
274 const std::vector<idx>& dims)
303 std::vector<idx> full_system(dims.size());
304 std::iota(std::begin(full_system), std::end(full_system), 0);
307 std::vector<idx> subsysAsorted{subsysA};
308 std::vector<idx> subsysBsorted{subsysB};
311 std::sort(std::begin(subsysAsorted), std::end(subsysAsorted));
312 std::sort(std::begin(subsysBsorted), std::end(subsysBsorted));
315 std::vector<idx> subsysAB;
316 std::set_union(std::begin(subsysAsorted), std::end(subsysAsorted),
317 std::begin(subsysBsorted), std::end(subsysBsorted),
318 std::back_inserter(subsysAB));
321 std::vector<idx> subsysA_bar =
complement(subsysA, dims.size());
322 std::vector<idx> subsysB_bar =
complement(subsysB, dims.size());;
323 std::vector<idx> subsysAB_bar =
complement(subsysAB, dims.size());
341 template<
typename Derived>
343 const std::vector<idx>& subsysA,
344 const std::vector<idx>& subsysB,
361 std::vector<idx> dims(N, d);
bool check_nonzero_size(const T &x) noexcept
Definition: util.h:129
Dimension(s) mismatch matrix size exception.
Definition: exception.h:322
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:271
bool check_subsys_match_dims(const std::vector< idx > &subsys, const std::vector< idx > &dims)
Definition: util.h:221
double tsallis(const Eigen::MatrixBase< Derived > &A, double q)
Tsallis- entropy of the density matrix A, for .
Definition: entropies.h:196
Eigen::MatrixXd dmat
Real (double precision) dynamic Eigen matrix.
Definition: types.h:65
Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic > dyn_mat
Dynamic Eigen matrix over the field specified by Scalar.
Definition: types.h:77
bool check_dims_match_mat(const std::vector< idx > &dims, const Eigen::MatrixBase< Derived > &A)
Definition: util.h:158
Subsystems mismatch dimensions exception.
Definition: exception.h:395
Quantum++ main namespace.
Definition: codes.h:30
idx get_num_subsys(idx sz, idx d)
Definition: util.h:391
Invalid dimension(s) exception.
Definition: exception.h:287
std::vector< T > complement(std::vector< T > subsys, idx N)
Constructs the complement of a subsystem vector.
Definition: functions.h:1811
dyn_col_vect< double > svals(const Eigen::MatrixBase< Derived > &A)
Singular values.
Definition: functions.h:470
double entropy(const Eigen::MatrixBase< Derived > &A)
von-Neumann entropy of the density matrix A
Definition: entropies.h:39
double renyi(const Eigen::MatrixBase< Derived > &A, double alpha)
Renyi- entropy of the density matrix A, for .
Definition: entropies.h:99
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:1241
constexpr double infty
Used to denote infinity in double precision.
Definition: constants.h:85
bool check_square_mat(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:101
Parameter out of range exception.
Definition: exception.h:567
Eigen::MatrixXcd cmat
Complex (double precision) dynamic Eigen matrix.
Definition: types.h:60
bool check_dims(const std::vector< idx > &dims)
Definition: util.h:142
std::size_t idx
Non-negative integer index.
Definition: types.h:35
Matrix is not square exception.
Definition: exception.h:151
Object has zero size exception.
Definition: exception.h:134