// This file is part of Eigen, a lightweight C++ template library
// for linear algebra.
//
// Copyright (C) 2013 Christian Seiler <christian@iwakd.de>
//
// This Source Code Form is subject to the terms of the Mozilla
// Public License v. 2.0. If a copy of the MPL was not distributed
// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
#ifndef EIGEN_CXX11_TENSORSYMMETRY_SYMMETRY_H
#define EIGEN_CXX11_TENSORSYMMETRY_SYMMETRY_H
namespace Eigen {
enum {
NegationFlag = 0x01,
ConjugationFlag = 0x02
};
enum {
GlobalRealFlag = 0x01,
GlobalImagFlag = 0x02,
GlobalZeroFlag = 0x03
};
namespace internal {
template<std::size_t NumIndices, typename... Sym> struct tensor_symmetry_pre_analysis;
template<std::size_t NumIndices, typename... Sym> struct tensor_static_symgroup;
template<bool instantiate, std::size_t NumIndices, typename... Sym> struct tensor_static_symgroup_if;
template<typename Tensor_> struct tensor_symmetry_calculate_flags;
template<typename Tensor_> struct tensor_symmetry_assign_value;
template<typename... Sym> struct tensor_symmetry_num_indices;
} // end namespace internal
template<int One_, int Two_>
struct Symmetry
{
static_assert(One_ != Two_, "Symmetries must cover distinct indices.");
constexpr static int One = One_;
constexpr static int Two = Two_;
constexpr static int Flags = 0;
};
template<int One_, int Two_>
struct AntiSymmetry
{
static_assert(One_ != Two_, "Symmetries must cover distinct indices.");
constexpr static int One = One_;
constexpr static int Two = Two_;
constexpr static int Flags = NegationFlag;
};
template<int One_, int Two_>
struct Hermiticity
{
static_assert(One_ != Two_, "Symmetries must cover distinct indices.");
constexpr static int One = One_;
constexpr static int Two = Two_;
constexpr static int Flags = ConjugationFlag;
};
template<int One_, int Two_>
struct AntiHermiticity
{
static_assert(One_ != Two_, "Symmetries must cover distinct indices.");
constexpr static int One = One_;
constexpr static int Two = Two_;
constexpr static int Flags = ConjugationFlag | NegationFlag;
};
/** \class DynamicSGroup
* \ingroup TensorSymmetry_Module
*
* \brief Dynamic symmetry group
*
* The %DynamicSGroup class represents a symmetry group that need not be known at
* compile time. It is useful if one wants to support arbitrary run-time defineable
* symmetries for tensors, but it is also instantiated if a symmetry group is defined
* at compile time that would be either too large for the compiler to reasonably
* generate (using templates to calculate this at compile time is very inefficient)
* or that the compiler could generate the group but that it wouldn't make sense to
* unroll the loop for setting coefficients anymore.
*/
class DynamicSGroup;
/** \internal
*
* \class DynamicSGroupFromTemplateArgs
* \ingroup TensorSymmetry_Module
*
* \brief Dynamic symmetry group, initialized from template arguments
*
* This class is a child class of DynamicSGroup. It uses the template arguments
* specified to initialize itself.
*/
template<typename... Gen>
class DynamicSGroupFromTemplateArgs;
/** \class StaticSGroup
* \ingroup TensorSymmetry_Module
*
* \brief Static symmetry group
*
* This class represents a symmetry group that is known and resolved completely
* at compile time. Ideally, no run-time penalty is incurred compared to the
* manual unrolling of the symmetry.
*
* <b><i>CAUTION:</i></b>
*
* Do not use this class directly for large symmetry groups. The compiler
* may run into a limit, or segfault or in the very least will take a very,
* very, very long time to compile the code. Use the SGroup class instead
* if you want a static group. That class contains logic that will
* automatically select the DynamicSGroup class instead if the symmetry
* group becomes too large. (In that case, unrolling may not even be
* beneficial.)
*/
template<typename... Gen>
class StaticSGroup;
/** \class SGroup
* \ingroup TensorSymmetry_Module
*
* \brief Symmetry group, initialized from template arguments
*
* This class represents a symmetry group whose generators are already
* known at compile time. It may or may not be resolved at compile time,
* depending on the estimated size of the group.
*
* \sa StaticSGroup
* \sa DynamicSGroup
*/
template<typename... Gen>
class SGroup : public internal::tensor_symmetry_pre_analysis<internal::tensor_symmetry_num_indices<Gen...>::value, Gen...>::root_type
{
public:
constexpr static std::size_t NumIndices = internal::tensor_symmetry_num_indices<Gen...>::value;
typedef typename internal::tensor_symmetry_pre_analysis<NumIndices, Gen...>::root_type Base;
// make standard constructors + assignment operators public
inline SGroup() : Base() { }
inline SGroup(const SGroup<Gen...>& other) : Base(other) { }
inline SGroup(SGroup<Gen...>&& other) : Base(other) { }
inline SGroup<Gen...>& operator=(const SGroup<Gen...>& other) { Base::operator=(other); return *this; }
inline SGroup<Gen...>& operator=(SGroup<Gen...>&& other) { Base::operator=(other); return *this; }
// all else is defined in the base class
};
namespace internal {
template<typename... Sym> struct tensor_symmetry_num_indices
{
constexpr static std::size_t value = 1;
};
template<int One_, int Two_, typename... Sym> struct tensor_symmetry_num_indices<Symmetry<One_, Two_>, Sym...>
{
private:
constexpr static std::size_t One = static_cast<std::size_t>(One_);
constexpr static std::size_t Two = static_cast<std::size_t>(Two_);
constexpr static std::size_t Three = tensor_symmetry_num_indices<Sym...>::value;
// don't use std::max, since it's not constexpr until C++14...
constexpr static std::size_t maxOneTwoPlusOne = ((One > Two) ? One : Two) + 1;
public:
constexpr static std::size_t value = (maxOneTwoPlusOne > Three) ? maxOneTwoPlusOne : Three;
};
template<int One_, int Two_, typename... Sym> struct tensor_symmetry_num_indices<AntiSymmetry<One_, Two_>, Sym...>
: public tensor_symmetry_num_indices<Symmetry<One_, Two_>, Sym...> {};
template<int One_, int Two_, typename... Sym> struct tensor_symmetry_num_indices<Hermiticity<One_, Two_>, Sym...>
: public tensor_symmetry_num_indices<Symmetry<One_, Two_>, Sym...> {};
template<int One_, int Two_, typename... Sym> struct tensor_symmetry_num_indices<AntiHermiticity<One_, Two_>, Sym...>
: public tensor_symmetry_num_indices<Symmetry<One_, Two_>, Sym...> {};
/** \internal
*
* \class tensor_symmetry_pre_analysis
* \ingroup TensorSymmetry_Module
*
* \brief Pre-select whether to use a static or dynamic symmetry group
*
* When a symmetry group could in principle be determined at compile time,
* this template implements the logic whether to actually do that or whether
* to rather defer that to runtime.
*
* The logic is as follows:
* <dl>
* <dt><b>No generators (trivial symmetry):</b></dt>
* <dd>Use a trivial static group. Ideally, this has no performance impact
* compared to not using symmetry at all. In practice, this might not
* be the case.</dd>
* <dt><b>More than 4 generators:</b></dt>
* <dd>Calculate the group at run time, it is likely far too large for the
* compiler to be able to properly generate it in a realistic time.</dd>
* <dt><b>Up to and including 4 generators:</b></dt>
* <dd>Actually enumerate all group elements, but then check how many there
* are. If there are more than 16, it is unlikely that unrolling the
* loop (as is done in the static compile-time case) is sensible, so
* use a dynamic group instead. If there are at most 16 elements, actually
* use that static group. Note that the largest group with 4 generators
* still compiles with reasonable resources.</dd>
* </dl>
*
* Note: Example compile time performance with g++-4.6 on an Intenl Core i5-3470
* with 16 GiB RAM (all generators non-redundant and the subgroups don't
* factorize):
*
* # Generators -O0 -ggdb -O2
* -------------------------------------------------------------------
* 1 0.5 s / 250 MiB 0.45s / 230 MiB
* 2 0.5 s / 260 MiB 0.5 s / 250 MiB
* 3 0.65s / 310 MiB 0.62s / 310 MiB
* 4 2.2 s / 860 MiB 1.7 s / 770 MiB
* 5 130 s / 13000 MiB 120 s / 11000 MiB
*
* It is clear that everything is still very efficient up to 4 generators, then
* the memory and CPU requirements become unreasonable. Thus we only instantiate
* the template group theory logic if the number of generators supplied is 4 or
* lower, otherwise this will be forced to be done during runtime, where the
* algorithm is reasonably fast.
*/
template<std::size_t NumIndices>
struct tensor_symmetry_pre_analysis<NumIndices>
{
typedef StaticSGroup<> root_type;
};
template<std::size_t NumIndices, typename Gen_, typename... Gens_>
struct tensor_symmetry_pre_analysis<NumIndices, Gen_, Gens_...>
{
constexpr static std::size_t max_static_generators = 4;
constexpr static std::size_t max_static_elements = 16;
typedef tensor_static_symgroup_if<(sizeof...(Gens_) + 1 <= max_static_generators), NumIndices, Gen_, Gens_...> helper;
constexpr static std::size_t possible_size = helper::size;
typedef typename conditional<
possible_size == 0 || possible_size >= max_static_elements,
DynamicSGroupFromTemplateArgs<Gen_, Gens_...>,
typename helper::type
>::type root_type;
};
template<bool instantiate, std::size_t NumIndices, typename... Gens>
struct tensor_static_symgroup_if
{
constexpr static std::size_t size = 0;
typedef void type;
};
template<std::size_t NumIndices, typename... Gens>
struct tensor_static_symgroup_if<true, NumIndices, Gens...> : tensor_static_symgroup<NumIndices, Gens...> {};
template<typename Tensor_>
struct tensor_symmetry_assign_value
{
typedef typename Tensor_::Index Index;
typedef typename Tensor_::Scalar Scalar;
constexpr static std::size_t NumIndices = Tensor_::NumIndices;
static inline int run(const std::array<Index, NumIndices>& transformed_indices, int transformation_flags, int dummy, Tensor_& tensor, const Scalar& value_)
{
Scalar value(value_);
if (transformation_flags & ConjugationFlag)
value = numext::conj(value);
if (transformation_flags & NegationFlag)
value = -value;
tensor.coeffRef(transformed_indices) = value;
return dummy;
}
};
template<typename Tensor_>
struct tensor_symmetry_calculate_flags
{
typedef typename Tensor_::Index Index;
constexpr static std::size_t NumIndices = Tensor_::NumIndices;
static inline int run(const std::array<Index, NumIndices>& transformed_indices, int transform_flags, int current_flags, const std::array<Index, NumIndices>& orig_indices)
{
if (transformed_indices == orig_indices) {
if (transform_flags & (ConjugationFlag | NegationFlag))
return current_flags | GlobalImagFlag; // anti-hermitian diagonal
else if (transform_flags & ConjugationFlag)
return current_flags | GlobalRealFlag; // hermitian diagonal
else if (transform_flags & NegationFlag)
return current_flags | GlobalZeroFlag; // anti-symmetric diagonal
}
return current_flags;
}
};
template<typename Tensor_, typename Symmetry_, int Flags = 0>
class tensor_symmetry_value_setter
{
public:
typedef typename Tensor_::Index Index;
typedef typename Tensor_::Scalar Scalar;
constexpr static std::size_t NumIndices = Tensor_::NumIndices;
inline tensor_symmetry_value_setter(Tensor_& tensor, Symmetry_ const& symmetry, std::array<Index, NumIndices> const& indices)
: m_tensor(tensor), m_symmetry(symmetry), m_indices(indices) { }
inline tensor_symmetry_value_setter<Tensor_, Symmetry_, Flags>& operator=(Scalar const& value)
{
doAssign(value);
return *this;
}
private:
Tensor_& m_tensor;
Symmetry_ m_symmetry;
std::array<Index, NumIndices> m_indices;
inline void doAssign(Scalar const& value)
{
#ifdef EIGEN_TENSOR_SYMMETRY_CHECK_VALUES
int value_flags = m_symmetry.template apply<internal::tensor_symmetry_calculate_flags<Tensor_>, int>(m_indices, m_symmetry.globalFlags(), m_indices);
if (value_flags & GlobalRealFlag)
eigen_assert(numext::imag(value) == 0);
if (value_flags & GlobalImagFlag)
eigen_assert(numext::real(value) == 0);
#endif
m_symmetry.template apply<internal::tensor_symmetry_assign_value<Tensor_>, int>(m_indices, 0, m_tensor, value);
}
};
} // end namespace internal
} // end namespace Eigen
#endif // EIGEN_CXX11_TENSORSYMMETRY_SYMMETRY_H
/*
* kate: space-indent on; indent-width 2; mixedindent off; indent-mode cstyle;
*/