27 #ifndef INTERNAL_UTIL_H_ 28 #define INTERNAL_UTIL_H_ 48 for (
idx i = 0; i < numdims; ++i)
54 for (
idx i = 0; i < numdims; ++i)
56 result[numdims - i - 1] = n % (dims[numdims - i - 1]);
57 n /= (dims[numdims - i - 1]);
63 #if (__GNUC__ && !__clang__) 64 #pragma GCC diagnostic push 65 #pragma GCC diagnostic ignored "-Warray-bounds" 66 #pragma GCC diagnostic ignored "-Wmaybe-uninitialized" 85 part_prod[numdims - 1] = 1;
86 for (
idx i = 1; i < numdims; ++i)
88 part_prod[numdims - i - 1] =
89 part_prod[numdims - i] * dims[numdims - i];
90 result += midx[numdims - i - 1] * part_prod[numdims - i - 1];
93 return result + midx[numdims - 1];
95 #if (__GNUC__ && !__clang__) 96 #pragma GCC diagnostic pop 100 template<
typename Derived>
103 return A.rows() == A.cols();
107 template<
typename Derived>
110 return A.rows() == 1 || A.cols() == 1;
114 template<
typename Derived>
117 return A.rows() == 1;
121 template<
typename Derived>
124 return A.cols() == 1;
131 return x.size() != 0;
135 template<
typename T1,
typename T2>
138 return lhs.size() == rhs.size();
144 if (dims.size() == 0)
147 return std::find_if(std::begin(dims), std::end(dims),
148 [dims](
idx i) ->
bool 150 if (i == 0)
return true;
152 }) == std::end(dims);
157 template<
typename Derived>
159 const Eigen::MatrixBase<Derived>& A)
163 assert(dims.size() > 0);
164 assert(A.rows() == A.cols());
166 idx proddim = std::accumulate(std::begin(dims), std::end(dims),
167 static_cast<idx>(1), std::multiplies<idx>());
169 return proddim ==
static_cast<idx>(A.rows());
173 template<
typename Derived>
175 const Eigen::MatrixBase<Derived>& A)
179 assert(dims.size() > 0);
180 assert(A.rows() > 0);
181 assert(A.cols() == 1);
183 idx proddim = std::accumulate(std::begin(dims), std::end(dims),
184 static_cast<idx>(1), std::multiplies<idx>());
186 return proddim ==
static_cast<idx>(A.rows());
190 template<
typename Derived>
192 const Eigen::MatrixBase<Derived>& A)
196 assert(dims.size() > 0);
197 assert(A.cols() > 0);
198 assert(A.rows() == 1);
200 idx proddim = std::accumulate(std::begin(dims), std::end(dims),
201 static_cast<idx>(1), std::multiplies<idx>());;
203 return proddim ==
static_cast<idx>(A.cols());
211 assert(dims.size() > 0);
222 const std::vector<idx>& dims)
227 if (subsys.size() > dims.size())
231 std::vector<idx> subsyssort = subsys;
232 std::sort(std::begin(subsyssort), std::end(subsyssort));
235 if (std::unique(std::begin(subsyssort), std::end(subsyssort))
236 != std::end(subsyssort))
240 return std::find_if(std::begin(subsyssort), std::end(subsyssort),
241 [dims](
idx i) ->
bool 243 return i > dims.size() - 1;
244 }) == std::end(subsyssort);
248 template<
typename Derived>
251 return A.rows() == 2 && A.cols() == 2;
255 template<
typename Derived>
258 return A.rows() == 2 && A.cols() == 1;
262 template<
typename Derived>
265 return A.rows() == 1 && A.cols() == 2;
269 template<
typename Derived>
272 return (A.rows() == 1 && A.cols() == 2) ||
273 (A.rows() == 2 && A.cols() == 1);
280 if (perm.size() == 0)
283 std::vector<idx> ordered(perm.size());
284 std::iota(std::begin(ordered), std::end(ordered), 0);
286 return std::is_permutation(std::begin(ordered), std::end(ordered),
292 template<
typename Derived1,
typename Derived2>
294 const Eigen::MatrixBase<Derived2>& B)
302 if (!std::is_same<
typename Derived1::Scalar,
303 typename Derived2::Scalar>::value)
315 idx Acols =
static_cast<idx>(rA.cols());
316 idx Arows =
static_cast<idx>(rA.rows());
317 idx Bcols =
static_cast<idx>(rB.cols());
318 idx Brows =
static_cast<idx>(rB.rows());
321 result.resize(Arows * Brows, Acols * Bcols);
324 #pragma omp parallel for collapse(2) 325 #endif // WITH_OPENMP_ 326 for (
idx j = 0; j < Acols; ++j)
327 for (
idx i = 0; i < Arows; ++i)
328 result.block(i * Brows, j * Bcols, Brows, Bcols) = rA(i, j) * rB;
335 template<
typename Derived1,
typename Derived2>
337 const Eigen::MatrixBase<Derived1>& A,
338 const Eigen::MatrixBase<Derived2>& B)
346 if (!std::is_same<
typename Derived1::Scalar,
347 typename Derived2::Scalar>::value)
359 idx Acols =
static_cast<idx>(rA.cols());
360 idx Arows =
static_cast<idx>(rA.rows());
361 idx Bcols =
static_cast<idx>(rB.cols());
362 idx Brows =
static_cast<idx>(rB.rows());
368 result.block(0, 0, Arows, Acols) = rA;
369 result.block(Arows, Acols, Brows, Bcols) = rB;
382 template<
typename T,
typename First,
typename ... Args>
385 v.emplace_back(std::forward<First>(first));
398 return static_cast<idx>(std::llround(std::log2(sz) / std::log2(d)));
412 return static_cast<idx>(std::llround(std::sqrt(sz)));
414 return static_cast<idx>(std::llround(std::pow(sz, 1. / N)));
426 std::ostringstream ostr;
429 std::vector<std::string> vstr;
432 for (
idx i = 0; i < static_cast<idx>(A.rows()); ++i)
435 j < static_cast<idx>(A.cols()); ++j)
439 ostr.str(std::string {});
442 double re =
static_cast<cplx>(A(i, j)).real();
443 double im =
static_cast<cplx>(A(i, j)).imag();
445 if (std::abs(re) <
chop && std::abs(im) <
chop)
450 vstr.push_back(ostr.str());
451 }
else if (std::abs(re) <
chop)
454 vstr.push_back(ostr.str() +
"i");
455 }
else if (std::abs(im) <
chop)
458 vstr.push_back(ostr.str() +
" ");
464 strA += (im > 0 ?
" + " :
" - ");
466 ostr.str(std::string());
467 ostr << std::abs(im);
470 vstr.push_back(strA);
476 std::vector<idx> maxlengthcols(A.cols(), 0);
478 for (
idx i = 0; i < static_cast<idx>(A.rows());
481 j < static_cast<idx>(A.cols()); ++j)
482 if (vstr[i * A.cols() + j].size() > maxlengthcols[j])
483 maxlengthcols[j] = vstr[i * A.cols() + j].size();
486 for (
idx i = 0; i < static_cast<idx>(A.rows()); ++i)
488 os << std::setw(static_cast<int>(maxlengthcols[0])) << std::right
489 << vstr[i * A.cols()];
492 j < static_cast<idx>(A.cols()); ++j)
493 os << std::setw(static_cast<int>(maxlengthcols[j] + 2))
494 << std::right << vstr[i * A.cols() + j];
496 if (i < static_cast<idx>(A.rows()) - 1)
bool check_nonzero_size(const T &x) noexcept
Definition: util.h:129
constexpr double chop
Used in qpp::disp() for setting to zero numbers that have their absolute value smaller than qpp::chop...
Definition: constants.h:56
constexpr idx maxn
Maximum number of allowed qubits/qudits (subsystems)
Definition: constants.h:71
bool check_subsys_match_dims(const std::vector< idx > &subsys, const std::vector< idx > &dims)
Definition: util.h:221
dyn_mat< typename Derived1::Scalar > kron2(const Eigen::MatrixBase< Derived1 > &A, const Eigen::MatrixBase< Derived2 > &B)
Definition: util.h:293
bool check_cvector(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:122
std::ostream & display_impl_(const T &A, std::ostream &os, double chop=qpp::chop) const
Definition: util.h:422
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
idx multiidx2n(const idx *const midx, idx numdims, const idx *const dims) noexcept
Definition: util.h:70
void n2multiidx(idx n, idx numdims, const idx *const dims, idx *result) noexcept
Definition: util.h:40
Quantum++ main namespace.
Definition: codes.h:30
bool check_perm(const std::vector< idx > &perm)
Definition: util.h:278
idx get_num_subsys(idx sz, idx d)
Definition: util.h:391
bool check_rvector(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:115
bool check_eq_dims(const std::vector< idx > &dims, idx dim) noexcept
Definition: util.h:207
bool check_vector(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:108
bool check_dims_match_cvect(const std::vector< idx > &dims, const Eigen::MatrixBase< Derived > &A)
Definition: util.h:174
bool check_dims_match_rvect(const std::vector< idx > &dims, const Eigen::MatrixBase< Derived > &A)
Definition: util.h:191
bool check_matching_sizes(const T1 &lhs, const T2 &rhs) noexcept
Definition: util.h:136
std::complex< double > cplx
Complex number in double precision.
Definition: types.h:45
void variadic_vector_emplace(std::vector< T > &)
Definition: util.h:377
Type mismatch exception.
Definition: exception.h:585
idx get_dim_subsys(idx sz, idx N)
Definition: util.h:404
bool check_square_mat(const Eigen::MatrixBase< Derived > &A)
Definition: util.h:101
dyn_mat< typename Derived1::Scalar > dirsum2(const Eigen::MatrixBase< Derived1 > &A, const Eigen::MatrixBase< Derived2 > &B)
Definition: util.h:336
bool check_qubit_rvector(const Eigen::MatrixBase< Derived > &A) noexcept
Definition: util.h:263
bool check_dims(const std::vector< idx > &dims)
Definition: util.h:142
bool check_qubit_vector(const Eigen::MatrixBase< Derived > &A) noexcept
Definition: util.h:270
std::size_t idx
Non-negative integer index.
Definition: types.h:35
bool check_qubit_cvector(const Eigen::MatrixBase< Derived > &A) noexcept
Definition: util.h:256
Object has zero size exception.
Definition: exception.h:134
bool check_qubit_matrix(const Eigen::MatrixBase< Derived > &A) noexcept
Definition: util.h:249