diff --git a/src/common/field.hh b/src/common/field.hh
index a0e8c07..f6ed3c2 100644
--- a/src/common/field.hh
+++ b/src/common/field.hh
@@ -1,412 +1,413 @@
/**
* file field.hh
*
* @author Till Junge <till.junge@epfl.ch>
*
* @date 07 Sep 2017
*
* @brief header-only implementation of a field for field collections
*
* @section LICENCE
*
* Copyright (C) 2017 Till Junge
*
* µSpectre is free software; you can redistribute it and/or
* modify it under the terms of the GNU General Public License as
* published by the Free Software Foundation, either version 3, or (at
* your option) any later version.
*
* µSpectre is distributed in the hope that it will be useful, but
* WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
* General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with GNU Emacs; see the file COPYING. If not, write to the
* Free Software Foundation, Inc., 59 Temple Place - Suite 330,
* Boston, MA 02111-1307, USA.
*/
#ifndef FIELD_H
#define FIELD_H
#include <string>
#include <utility>
#include <typeinfo>
#include <vector>
#include <algorithm>
#include <Eigen/Dense>
#include <cmath>
#include <memory>
#include <type_traits>
namespace muSpectre {
namespace internal {
/* ---------------------------------------------------------------------- */
template <class FieldCollection>
class FieldBase
{
protected:
//! constructor
//! unique name (whithin Collection)
//! number of components
//! collection to which this field belongs (eg, material, system)
FieldBase(std::string unique_name,
size_t nb_components,
FieldCollection & collection);
public:
using collection_t = FieldCollection;
//! Copy constructor
FieldBase(const FieldBase &other) = delete;
//! Move constructor
FieldBase(FieldBase &&other) noexcept = delete;
//! Destructor
virtual ~FieldBase() noexcept = default;
//! Copy assignment operator
FieldBase& operator=(const FieldBase &other) = delete;
//! Move assignment operator
FieldBase& operator=(FieldBase &&other) noexcept = delete;
/* ---------------------------------------------------------------------- */
//!Identifying accessors
//! return field name
inline const std::string & get_name() const;
//! return field type
//inline const Field_t & get_type() const;
//! return my collection (for iterating)
inline const FieldCollection & get_collection() const;
//! return my collection (for iterating)
inline const size_t & get_nb_components() const;
//! return type_id of stored type
virtual const std::type_info & get_stored_typeid() const = 0;
virtual size_t size() const = 0;
//! initialise field to zero (do more complicated initialisations through
//! fully typed maps)
virtual void set_zero() = 0;
//! give access to collections
friend FieldCollection;
protected:
/* ---------------------------------------------------------------------- */
//! allocate memory etc
virtual void resize(size_t size) = 0;
const std::string name;
const size_t nb_components;
const FieldCollection & collection;
private:
};
/* ---------------------------------------------------------------------- */
//! declaraton for friending
- template <class FieldCollection, typename T, Dim_t NbComponents>
+ template <class FieldCollection, typename T, Dim_t NbComponents, bool isConst>
class FieldMap;
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t NbComponents>
class TypedFieldBase: public FieldBase<FieldCollection>
{
- friend class FieldMap<FieldCollection, T, NbComponents>;
+ friend class FieldMap<FieldCollection, T, NbComponents, true>;
+ friend class FieldMap<FieldCollection, T, NbComponents, false>;
public:
constexpr static auto nb_components{NbComponents};
using Parent = FieldBase<FieldCollection>;
using Parent::collection_t;
using Scalar = T;
using Base = Parent;
//using storage_type = Eigen::Array<T, Eigen::Dynamic, NbComponents>;
using StoredType = Eigen::Array<T, NbComponents, 1>;
using StorageType = std::vector<StoredType,
Eigen::aligned_allocator<StoredType>>;
TypedFieldBase(std::string unique_name,
FieldCollection& collection);
virtual ~TypedFieldBase() = default;
//! return type_id of stored type
virtual const std::type_info & get_stored_typeid() const;
//! initialise field to zero (do more complicated initialisations through
//! fully typed maps)
inline void set_zero() override final;
inline void push_back(const StoredType & value);
template< class... Args>
inline void emplace_back(Args&&... args);
size_t size() const override final;
protected:
inline T* get_ptr_to_entry(const size_t&& index);
inline T& get_ref_to_entry(const size_t&& index);
inline virtual void resize(size_t size) override final;
StorageType array{};
};
} // internal
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t order, Dim_t dim>
class TensorField: public internal::TypedFieldBase<FieldCollection,
T,
ipow(dim,order)>
{
public:
using Parent = internal::TypedFieldBase<FieldCollection,
T,
ipow(dim,order)>;
using Base = typename Parent::Base;
using Field_p = typename FieldCollection::Field_p;
using component_type = T;
//! Copy constructor
TensorField(const TensorField &other) = delete;
//! Move constructor
TensorField(TensorField &&other) noexcept = delete;
//! Destructor
virtual ~TensorField() noexcept = default;
//! Copy assignment operator
TensorField& operator=(const TensorField &other) = delete;
//! Move assignment operator
TensorField& operator=(TensorField &&other) noexcept = delete;
//! accessors
inline Dim_t get_order() const;
inline Dim_t get_dim() const;
//! factory function
template<class FieldType, class CollectionType, typename... Args>
friend typename FieldType::Base& make_field(std::string unique_name,
CollectionType & collection,
Args&&... args);
protected:
//! constructor protected!
TensorField(std::string unique_name, FieldCollection & collection);
private:
};
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t NbRow, Dim_t NbCol=NbRow>
class MatrixField: public internal::TypedFieldBase<FieldCollection,
T,
NbRow*NbCol>
{
public:
using Parent = internal::TypedFieldBase<FieldCollection,
T,
NbRow*NbCol>;
using Base = typename Parent::Base;
using Field_p = std::unique_ptr<internal::FieldBase<FieldCollection>>;
using component_type = T;
//! Copy constructor
MatrixField(const MatrixField &other) = delete;
//! Move constructor
MatrixField(MatrixField &&other) noexcept = delete;
//! Destructor
virtual ~MatrixField() noexcept = default;
//! Copy assignment operator
MatrixField& operator=(const MatrixField &other) = delete;
//! Move assignment operator
MatrixField& operator=(MatrixField &&other) noexcept = delete;
//! accessors
inline Dim_t get_nb_row() const;
inline Dim_t get_nb_col() const;
//! factory function
template<class FieldType, class CollectionType, typename... Args>
friend typename FieldType::Base& make_field(std::string unique_name,
CollectionType & collection,
Args&&... args);
protected:
//! constructor protected!
MatrixField(std::string unique_name, FieldCollection & collection);
private:
};
/* ---------------------------------------------------------------------- */
//! convenience alias
template <class FieldCollection, typename T>
using ScalarField = TensorField<FieldCollection, T, 1, 1>;
/* ---------------------------------------------------------------------- */
// Implementations
namespace internal {
/* ---------------------------------------------------------------------- */
template <class FieldCollection>
FieldBase<FieldCollection>::FieldBase(std::string unique_name,
size_t nb_components_,
FieldCollection & collection_)
:name(unique_name), nb_components(nb_components_),
collection(collection_) {}
/* ---------------------------------------------------------------------- */
template <class FieldCollection>
inline const std::string & FieldBase<FieldCollection>::get_name() const {
return this->name;
}
/* ---------------------------------------------------------------------- */
template <class FieldCollection>
inline const FieldCollection & FieldBase<FieldCollection>::
get_collection() const {
return this->collection;
}
/* ---------------------------------------------------------------------- */
template <class FieldCollection>
inline const size_t & FieldBase<FieldCollection>::
get_nb_components() const {
return this->nb_components;
}
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t NbComponents>
TypedFieldBase<FieldCollection, T, NbComponents>::
TypedFieldBase(std::string unique_name, FieldCollection & collection)
:FieldBase<FieldCollection>(unique_name, NbComponents, collection){
static_assert
((std::is_arithmetic<T>::value ||
std::is_same<T, Complex>::value),
"Use TypedFieldBase for integer, real or complex scalars for T");
static_assert(NbComponents > 0,
"Only fields with more than 0 components");
}
/* ---------------------------------------------------------------------- */
//! return type_id of stored type
template <class FieldCollection, typename T, Dim_t NbComponents>
const std::type_info & TypedFieldBase<FieldCollection, T, NbComponents>::
get_stored_typeid() const {
return typeid(T);
}
/* ---------------------------------------------------------------------- */
template<class FieldCollection, typename T, Dim_t NbComponents>
void
TypedFieldBase<FieldCollection, T, NbComponents>::
set_zero() {
std::fill(this->array.begin(), this->array.end(), T{});
}
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t NbComponents>
size_t TypedFieldBase<FieldCollection, T, NbComponents>::
size() const {
return this->array.size();
}
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t NbComponents>
T* TypedFieldBase<FieldCollection, T, NbComponents>::
get_ptr_to_entry(const size_t&& index) {
return &this->array[std::move(index)](0, 0);
}
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t NbComponents>
T& TypedFieldBase<FieldCollection, T, NbComponents>::
get_ref_to_entry(const size_t && index) {
return this->array[std::move(index)](0, 0);
}
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t NbComponents>
void TypedFieldBase<FieldCollection, T, NbComponents>::
resize(size_t size) {
this->array.resize(size);
}
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t NbComponents>
void TypedFieldBase<FieldCollection, T, NbComponents>::
push_back(const StoredType & value) {
this->array.push_back(value);
}
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t NbComponents>
template <class... Args>
void TypedFieldBase<FieldCollection, T, NbComponents>::
emplace_back(Args&&... args) {
this->array.emplace_back(std::move(args...));
}
} // internal
/* ---------------------------------------------------------------------- */
//! Factory function, guarantees that only fields get created that are
//! properly registered and linked to a collection.
template<class FieldType, class FieldCollection, typename... Args>
typename FieldType::Base &
make_field(std::string unique_name,
FieldCollection & collection,
Args&&... args) {
auto && ptr = std::unique_ptr<FieldType>{
new FieldType(unique_name, collection, args...)};
auto& retref = *ptr;
collection.register_field(std::move(ptr));
return retref;
}
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t order, Dim_t dim>
TensorField<FieldCollection, T, order, dim>::
TensorField(std::string unique_name, FieldCollection & collection)
:Parent(unique_name, collection) {}
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t order, Dim_t dim>
Dim_t TensorField<FieldCollection, T, order, dim>::
get_order() const {
return order;
}
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t order, Dim_t dim>
Dim_t TensorField<FieldCollection, T, order, dim>::
get_dim() const {
return dim;
}
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t NbRow, Dim_t NbCol>
MatrixField<FieldCollection, T, NbRow, NbCol>::
MatrixField(std::string unique_name, FieldCollection & collection)
:Parent(unique_name, collection) {}
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t NbRow, Dim_t NbCol>
Dim_t MatrixField<FieldCollection, T, NbRow, NbCol>::
get_nb_col() const {
return NbCol;
}
/* ---------------------------------------------------------------------- */
template <class FieldCollection, typename T, Dim_t NbRow, Dim_t NbCol>
Dim_t MatrixField<FieldCollection, T, NbRow, NbCol>::
get_nb_row() const {
return NbRow;
}
} // muSpectre
#endif /* FIELD_H */
diff --git a/src/common/field_map_base.hh b/src/common/field_map_base.hh
index 9c355a3..b9d1dc2 100644
--- a/src/common/field_map_base.hh
+++ b/src/common/field_map_base.hh
@@ -1,529 +1,548 @@
/**
* file field_map.hh
*
* @author Till Junge <till.junge@epfl.ch>
*
* @date 12 Sep 2017
*
* @brief Defined a strongly defines proxy that iterates efficiently over a field
*
* @section LICENCE
*
* Copyright (C) 2017 Till Junge
*
* µSpectre is free software; you can redistribute it and/or
* modify it under the terms of the GNU General Public License as
* published by the Free Software Foundation, either version 3, or (at
* your option) any later version.
*
* µSpectre is distributed in the hope that it will be useful, but
* WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
* General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with GNU Emacs; see the file COPYING. If not, write to the
* Free Software Foundation, Inc., 59 Temple Place - Suite 330,
* Boston, MA 02111-1307, USA.
*/
#ifndef FIELD_MAP_H
#define FIELD_MAP_H
#include <array>
#include <string>
#include <memory>
#include <type_traits>
#include <unsupported/Eigen/CXX11/Tensor>
#include "common/field.hh"
#include "common/common.hh"
namespace muSpectre {
namespace internal {
//----------------------------------------------------------------------------//
//! little helper to automate creation of const maps without duplication
template<class T, bool isConst>
struct const_corrector {
using type = typename T::reference;
};
template<class T>
struct const_corrector<T, true> {
using type = typename T::const_reference;
};
template<class T, bool isConst>
using const_corrector_t = typename const_corrector<T, isConst>::type;
//----------------------------------------------------------------------------//
- template <class FieldCollection, typename T, Dim_t NbComponents>
+ template <class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
class FieldMap
{
public:
constexpr static auto nb_components{NbComponents};
using TypedField = ::muSpectre::internal::TypedFieldBase
<FieldCollection, T, NbComponents>;
using Field = typename TypedField::Parent;
using size_type = std::size_t;
//! Default constructor
FieldMap() = delete;
- FieldMap(Field & field);
+ template <bool isntConst=!ConstField>
+ FieldMap(std::enable_if_t<isntConst, Field &> field);
+ template <bool isConst=ConstField>
+ FieldMap(std::enable_if_t<isConst, const Field &> field);
template<class FC, typename T2, Dim_t NbC>
FieldMap(TypedFieldBase<FC, T2, NbC> & field);
//! Copy constructor
FieldMap(const FieldMap &other) = delete;
//! Move constructor
FieldMap(FieldMap &&other) noexcept = delete;
//! Destructor
virtual ~FieldMap() noexcept = default;
//! Copy assignment operator
FieldMap& operator=(const FieldMap &other) = delete;
//! Move assignment operator
FieldMap& operator=(FieldMap &&other) noexcept = delete;
//! give human-readable field map type
virtual std::string info_string() const = 0;
//! return field name
inline const std::string & get_name() const;
//! return my collection (for iterating)
inline const FieldCollection & get_collection() const;
//! member access needs to be implemented by inheriting classes
//inline value_type operator[](size_t index);
//inline value_type operator[](Ccoord ccord);
//! check compatibility (must be called by all inheriting classes at the
//! end of their constructors
inline void check_compatibility();
//! convenience call to collection's size method
inline size_t size() const;
//! compile-time compatibility check
template<class Field>
struct is_compatible;
- template <class FullyTypedFieldMap, bool isConst=false>
+ template <class FullyTypedFieldMap, bool ConstIter=false>
class iterator
{
+ static_assert(!((ConstIter==false) && (ConstField==true)),
+ "You can't have a non-const iterator over a const "
+ "field");
public:
using value_type =
- const_corrector_t<FullyTypedFieldMap, isConst>;
+ const_corrector_t<FullyTypedFieldMap, ConstIter>;
using pointer = typename FullyTypedFieldMap::pointer;
using difference_type = std::ptrdiff_t;
using iterator_category = std::random_access_iterator_tag;
using Ccoord = typename FieldCollection::Ccoord;
+ using reference = typename FullyTypedFieldMap::reference;
//! Default constructor
iterator() = delete;
//! constructor
inline iterator(FullyTypedFieldMap & fieldmap, bool begin=true);
//! constructor for random access
inline iterator(FullyTypedFieldMap & fieldmap, size_t index);
//! Copy constructor
iterator(const iterator &other)= default;
//! Move constructor
iterator(iterator &&other) noexcept = default;
//! Destructor
virtual ~iterator() noexcept = default;
//! Copy assignment operator
iterator& operator=(const iterator &other) = default;
//! Move assignment operator
iterator& operator=(iterator &&other) noexcept = default;
//! pre-increment
inline iterator & operator++();
//! post-increment
inline iterator operator++(int);
//! dereference
inline value_type operator*();
//! member of pointer
inline pointer operator->();
//! pre-decrement
inline iterator & operator--();
//! post-decrement
inline iterator operator--(int);
//! access subscripting
inline value_type operator[](difference_type diff);
//! equality
inline bool operator==(const iterator & other) const;
//! inequality
inline bool operator!=(const iterator & other) const;
//! div. comparisons
inline bool operator<(const iterator & other) const;
inline bool operator<=(const iterator & other) const;
inline bool operator>(const iterator & other) const;
inline bool operator>=(const iterator & other) const;
//! additions, subtractions and corresponding assignments
inline iterator operator+(difference_type diff) const;
inline iterator operator-(difference_type diff) const;
inline iterator& operator+=(difference_type diff);
inline iterator& operator-=(difference_type diff);
//! get pixel coordinates
inline Ccoord get_ccoord() const;
//! ostream operator (mainly for debug
friend std::ostream & operator<<(std::ostream & os,
const iterator& it) {
- if (isConst) {
+ if (ConstIter) {
os << "const ";
}
os << "iterator on field '"
<< it.fieldmap.get_name()
<< "', entry " << it.index;
return os;
}
protected:
const FieldCollection & collection;
FullyTypedFieldMap & fieldmap;
size_t index;
private:
};
protected:
inline T* get_ptr_to_entry(size_t index);
inline T& get_ref_to_entry(size_t index);
const FieldCollection & collection;
TypedField & field;
private:
};
} // internal
namespace internal {
/* ---------------------------------------------------------------------- */
- template<class FieldCollection, typename T, Dim_t NbComponents>
- FieldMap<FieldCollection, T, NbComponents>::
- FieldMap(Field & field)
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template <bool isntConst>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::
+ FieldMap(std::enable_if_t<isntConst, Field &> field)
:collection(field.get_collection()), field(static_cast<TypedField&>(field)) {
static_assert(NbComponents>0,
"Only fields with more than 0 components allowed");
}
/* ---------------------------------------------------------------------- */
- template <class FieldCollection, typename T, Dim_t NbComponents>
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template <bool isConst>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::
+ FieldMap(std::enable_if_t<isConst, const Field &> field)
+ :collection(field.get_collection()), field(static_cast<TypedField&>(field)) {
+ static_assert(NbComponents>0,
+ "Only fields with more than 0 components allowed");
+ }
+
+ /* ---------------------------------------------------------------------- */
+ template <class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
template <class FC, typename T2, Dim_t NbC>
- FieldMap<FieldCollection, T, NbComponents>::
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::
FieldMap(TypedFieldBase<FC, T2, NbC> & field)
:collection(field.get_collection()), field(static_cast<TypedField&>(field)) {
static_assert(std::is_same<FC, FieldCollection>::value,
"The field does not have the expected FieldCollection type");
static_assert(std::is_same<T2, T>::value,
"The field does not have the expected Scalar type");
static_assert((NbC == NbComponents),
"The field does not have the expected number of components");
}
/* ---------------------------------------------------------------------- */
- template<class FieldCollection, typename T, Dim_t NbComponents>
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
void
- FieldMap<FieldCollection, T, NbComponents>::
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::
check_compatibility() {
if (typeid(T).hash_code() !=
this->field.get_stored_typeid().hash_code()) {
throw FieldInterpretationError
("Cannot create a Map of type '" +
this->info_string() +
"' for field '" + this->field.get_name() + "' of type '" +
this->field.get_stored_typeid().name() + "'");
}
//check size compatibility
if (NbComponents != this->field.get_nb_components()) {
throw FieldInterpretationError
("Cannot create a Map of type '" +
this->info_string() +
"' for field '" + this->field.get_name() + "' with " +
std::to_string(this->field.get_nb_components()) + " components");
}
}
/* ---------------------------------------------------------------------- */
- template<class FieldCollection, typename T, Dim_t NbComponents>
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
size_t
- FieldMap<FieldCollection, T, NbComponents>::
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::
size() const {
return this->collection.size();
}
/* ---------------------------------------------------------------------- */
- template <class FieldCollection, typename T, Dim_t NbComponents>
+ template <class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
template <class Field>
- struct FieldMap<FieldCollection, T, NbComponents>::is_compatible {
+ struct FieldMap<FieldCollection, T, NbComponents, ConstField>::is_compatible {
constexpr static bool explain() {
static_assert
(std::is_same<typename Field::collection_t, FieldCollection>::value,
"The field does not have the expected FieldCollection type");
static_assert
(std::is_same<typename Field::Scalar, T>::value,
"The // FIXME: eld does not have the expected Scalar type");
static_assert((Field::nb_components == NbComponents),
"The field does not have the expected number of components");
//The static asserts wouldn't pass in the incompatible case, so this is it
return true;
}
constexpr static bool value{std::is_base_of<TypedField, Field>::value};
};
/* ---------------------------------------------------------------------- */
- template<class FieldCollection, typename T, Dim_t NbComponents>
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
const std::string &
- FieldMap<FieldCollection, T, NbComponents>::
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::
get_name() const {
return this->field.get_name();
}
/* ---------------------------------------------------------------------- */
- template<class FieldCollection, typename T, Dim_t NbComponents>
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
const FieldCollection &
- FieldMap<FieldCollection, T, NbComponents>::
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::
get_collection() const {
return this->collection;
}
/* ---------------------------------------------------------------------- */
/* ---------------------------------------------------------------------- */
// Iterator implementations
//! constructor
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator
+ <FullyTypedFieldMap, ConstIter>::
iterator(FullyTypedFieldMap & fieldmap, bool begin)
:collection(fieldmap.get_collection()), fieldmap(fieldmap),
index(begin ? 0 : fieldmap.field.size()) {}
/* ---------------------------------------------------------------------- */
//! constructor for random access
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
iterator(FullyTypedFieldMap & fieldmap, size_t index)
:collection(fieldmap.collection), fieldmap(fieldmap),
index(index) {}
/* ---------------------------------------------------------------------- */
//! pre-increment
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst> &
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter> &
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator++() {
this->index++;
return *this;
}
/* ---------------------------------------------------------------------- */
//! post-increment
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator++(int) {
iterator current = *this;
this->index++;
return current;
}
/* ---------------------------------------------------------------------- */
//! dereference
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
- typename FieldMap<FieldCollection, T, NbComponents>::template iterator<FullyTypedFieldMap, isConst>::value_type
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
+ typename FieldMap<FieldCollection, T, NbComponents, ConstField>::template iterator<FullyTypedFieldMap, ConstIter>::value_type
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator*() {
return this->fieldmap.template operator[]<value_type>(this->index);
}
/* ---------------------------------------------------------------------- */
//! member of pointer
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
typename FullyTypedFieldMap::pointer
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator->() {
return this->fieldmap.ptr_to_val_t(this->index);
}
/* ---------------------------------------------------------------------- */
//! pre-decrement
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst> &
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter> &
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator--() {
this->index--;
return *this;
}
/* ---------------------------------------------------------------------- */
//! post-decrement
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator--(int) {
iterator current = *this;
this->index--;
return current;
}
/* ---------------------------------------------------------------------- */
//! Access subscripting
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
- typename FieldMap<FieldCollection, T, NbComponents>::template iterator<FullyTypedFieldMap, isConst>::value_type
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
+ typename FieldMap<FieldCollection, T, NbComponents, ConstField>::template iterator<FullyTypedFieldMap, ConstIter>::value_type
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator[](difference_type diff) {
return this->fieldmap[this->index+diff];
}
/* ---------------------------------------------------------------------- */
//! equality
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
bool
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator==(const iterator & other) const {
return (this->index == other.index &&
&this->fieldmap == &other.fieldmap);
}
/* ---------------------------------------------------------------------- */
//! inquality
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
bool
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator!=(const iterator & other) const {
return !(*this == other);
}
/* ---------------------------------------------------------------------- */
//! div. comparisons
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
bool
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator<(const iterator & other) const {
return (this->index < other.index);
}
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
bool
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator<=(const iterator & other) const {
return (this->index <= other.index);
}
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
bool
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator>(const iterator & other) const {
return (this->index > other.index);
}
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
bool
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator>=(const iterator & other) const {
return (this->index >= other.index);
}
/* ---------------------------------------------------------------------- */
//! additions, subtractions and corresponding assignments
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator+(difference_type diff) const {
return iterator(this->fieldmap, this->index + diff);
}
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator-(difference_type diff) const {
return iterator(this->fieldmap, this->index - diff);
}
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst> &
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter> &
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator+=(difference_type diff) {
this->index += diff;
return *this;
}
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst> &
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter> &
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
operator-=(difference_type diff) {
this->index -= diff;
return *this;
}
/* ---------------------------------------------------------------------- */
//! get pixel coordinates
- template<class FieldCollection, typename T, Dim_t NbComponents>
- template<class FullyTypedFieldMap, bool isConst>
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
+ template<class FullyTypedFieldMap, bool ConstIter>
typename FieldCollection::Ccoord
- FieldMap<FieldCollection, T, NbComponents>::iterator<FullyTypedFieldMap, isConst>::
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::iterator<FullyTypedFieldMap, ConstIter>::
get_ccoord() const {
return this->collection.get_ccoord(this->index);
}
////----------------------------------------------------------------------------//
//template<class FieldCollection, typename T, Dim_t NbComponents, class FullyTypedFieldMap>
//std::ostream & operator <<
//(std::ostream &os,
-// const typename FieldMap<FieldCollection, T, NbComponents>::
-// template iterator<FullyTypedFieldMap, isConst> & it) {
+// const typename FieldMap<FieldCollection, T, NbComponents, ConstField>::
+// template iterator<FullyTypedFieldMap, ConstIter> & it) {
// os << "iterator on field '"
// << it.field.get_name()
// << "', entry " << it.index;
// return os;
//}
/* ---------------------------------------------------------------------- */
- template<class FieldCollection, typename T, Dim_t NbComponents>
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
T*
- FieldMap<FieldCollection, T, NbComponents>::get_ptr_to_entry(size_t index) {
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::get_ptr_to_entry(size_t index) {
return this->field.get_ptr_to_entry(std::move(index));
}
/* ---------------------------------------------------------------------- */
- template<class FieldCollection, typename T, Dim_t NbComponents>
+ template<class FieldCollection, typename T, Dim_t NbComponents, bool ConstField>
T&
- FieldMap<FieldCollection, T, NbComponents>::get_ref_to_entry(size_t index) {
+ FieldMap<FieldCollection, T, NbComponents, ConstField>::get_ref_to_entry(size_t index) {
return this->field.get_ref_to_entry(std::move(index));
}
} // internal
} // muSpectre
#endif /* FIELD_MAP_H */
diff --git a/src/common/field_map_matrixlike.hh b/src/common/field_map_matrixlike.hh
index 055b53f..c4f7492 100644
--- a/src/common/field_map_matrixlike.hh
+++ b/src/common/field_map_matrixlike.hh
@@ -1,227 +1,268 @@
/**
* file field_map_matrixlike.hh
*
* @author Till Junge <till.junge@epfl.ch>
*
* @date 26 Sep 2017
*
* @brief Eigen-Matrix and -Array maps over strongly typed fields
*
* @section LICENCE
*
* Copyright (C) 2017 Till Junge
*
* µSpectre is free software; you can redistribute it and/or
* modify it under the terms of the GNU General Public License as
* published by the Free Software Foundation, either version 3, or (at
* your option) any later version.
*
* µSpectre is distributed in the hope that it will be useful, but
* WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
* General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with GNU Emacs; see the file COPYING. If not, write to the
* Free Software Foundation, Inc., 59 Temple Place - Suite 330,
* Boston, MA 02111-1307, USA.
*/
#ifndef FIELD_MAP_MATRIXLIKE_H
#define FIELD_MAP_MATRIXLIKE_H
#include <memory>
#include <Eigen/Dense>
#include "common/field_map_base.hh"
#include "common/T4_map_proxy.hh"
namespace muSpectre {
namespace internal {
/* ---------------------------------------------------------------------- */
enum class Map_t{Matrix, Array, T4Matrix};
template<Map_t map_type>
struct NameWrapper {
const static std::string field_info_root;
};
template<>
const std::string NameWrapper<Map_t::Array>::field_info_root{"Array"};
template<>
const std::string NameWrapper<Map_t::Matrix>::field_info_root{"Matrix"};
template<>
const std::string NameWrapper<Map_t::T4Matrix>::field_info_root{"T4Matrix"};
/* ---------------------------------------------------------------------- */
- template <class FieldCollection, class EigenArray, Map_t map_type>
+ template <class FieldCollection, class EigenArray, Map_t map_type,
+ bool ConstField>
class MatrixLikeFieldMap: public FieldMap<FieldCollection,
typename EigenArray::Scalar,
- EigenArray::SizeAtCompileTime>
+ EigenArray::SizeAtCompileTime,
+ ConstField>
{
public:
using Parent = FieldMap<FieldCollection,
typename EigenArray::Scalar,
- EigenArray::SizeAtCompileTime>;
+ EigenArray::SizeAtCompileTime, ConstField>;
using Ccoord = Ccoord_t<FieldCollection::spatial_dim()>;
using value_type = EigenArray;
- using reference = value_type; // since it's a resource handle
using const_reference = const value_type;
+ using reference = std::conditional_t<ConstField,
+ const_reference,
+ value_type>; // since it's a resource handle
using size_type = typename Parent::size_type;
using pointer = std::unique_ptr<EigenArray>;
using TypedField = typename Parent::TypedField;
using Field = typename TypedField::Parent;
- using iterator = typename Parent::template iterator<MatrixLikeFieldMap>;
using const_iterator= typename Parent::template iterator<MatrixLikeFieldMap, true>;
+ using iterator = std::conditional_t<
+ ConstField,
+ const_iterator,
+ typename Parent::template iterator<MatrixLikeFieldMap>>;
using reverse_iterator = std::reverse_iterator<iterator>;
using const_reverse_iterator = std::reverse_iterator<const_iterator>;
friend iterator;
//! Default constructor
MatrixLikeFieldMap() = delete;
/**
* Constructor using a (non-typed) field. Compatibility is enforced at
* runtime. This should not be a performance concern, as this constructor
* will not be called in anny inner loops (if used as intended).
*/
- MatrixLikeFieldMap(Field & field);
+ template <bool isntConst=!ConstField>
+ MatrixLikeFieldMap(std::enable_if_t<isntConst, Field &> field);
+ template <bool isConst=ConstField>
+ MatrixLikeFieldMap(std::enable_if_t<isConst, const Field &> field);
/**
* Constructor using a typed field. Compatibility is enforced
* statically. It is not always possible to call this constructor, as the
* type of the field might not be known at compile time.
*/
template<class FC, typename T2, Dim_t NbC>
MatrixLikeFieldMap(TypedFieldBase<FC, T2, NbC> & field);
//! Copy constructor
MatrixLikeFieldMap(const MatrixLikeFieldMap &other) = delete;
//! Move constructor
MatrixLikeFieldMap(MatrixLikeFieldMap &&other) noexcept = delete;
//! Destructor
virtual ~MatrixLikeFieldMap() noexcept = default;
//! Copy assignment operator
MatrixLikeFieldMap& operator=(const MatrixLikeFieldMap &other) = delete;
//! Move assignment operator
MatrixLikeFieldMap& operator=(MatrixLikeFieldMap &&other) noexcept = delete;
//! give human-readable field map type
inline std::string info_string() const override final;
//! member access
template <class ref = reference>
inline ref operator[](size_type index);
inline reference operator[](const Ccoord& ccoord);
//! return an iterator to head of field for ranges
inline iterator begin(){return iterator(*this);}
inline const_iterator cbegin(){return const_iterator(*this);}
//! return an iterator to tail of field for ranges
inline iterator end(){return iterator(*this, false);};
inline const_iterator cend(){return const_iterator(*this, false);}
protected:
//! for sad, legacy iterator use
inline pointer ptr_to_val_t(size_type index);
const static std::string field_info_root;
private:
};
- template <class FieldCollection, class EigenArray, Map_t map_type>
- const std::string MatrixLikeFieldMap<FieldCollection, EigenArray, map_type>::
+ /* ---------------------------------------------------------------------- */
+ template <class FieldCollection, class EigenArray, Map_t map_type,
+ bool ConstField>
+ const std::string MatrixLikeFieldMap<FieldCollection, EigenArray, map_type, ConstField>::
field_info_root{NameWrapper<map_type>::field_info_root};
/* ---------------------------------------------------------------------- */
- template<class FieldCollection, class EigenArray, Map_t map_type>
- MatrixLikeFieldMap<FieldCollection, EigenArray, map_type>::
- MatrixLikeFieldMap(Field & field)
+ template<class FieldCollection, class EigenArray, Map_t map_type,
+ bool ConstField>
+ template <bool isntConst>
+ MatrixLikeFieldMap<FieldCollection, EigenArray, map_type, ConstField>::
+ MatrixLikeFieldMap(std::enable_if_t<isntConst, Field &> field)
+ :Parent(field) {
+ static_assert((isntConst != ConstField),
+ "let the compiler deduce isntConst, this is a SFINAE "
+ "parameter");
+ this->check_compatibility();
+ }
+
+ /* ---------------------------------------------------------------------- */
+ template<class FieldCollection, class EigenArray, Map_t map_type,
+ bool ConstField>
+ template <bool isConst>
+ MatrixLikeFieldMap<FieldCollection, EigenArray, map_type, ConstField>::
+ MatrixLikeFieldMap(std::enable_if_t<isConst, const Field &> field)
:Parent(field) {
+ static_assert((isConst == ConstField),
+ "let the compiler deduce isntConst, this is a SFINAE "
+ "parameter");
this->check_compatibility();
}
/* ---------------------------------------------------------------------- */
- template<class FieldCollection, class EigenArray, Map_t map_type>
+ template<class FieldCollection, class EigenArray, Map_t map_type,
+ bool ConstField>
template<class FC, typename T2, Dim_t NbC>
- MatrixLikeFieldMap<FieldCollection, EigenArray, map_type>::
+ MatrixLikeFieldMap<FieldCollection, EigenArray, map_type, ConstField>::
MatrixLikeFieldMap(TypedFieldBase<FC, T2, NbC> & field)
:Parent(field) {
}
/* ---------------------------------------------------------------------- */
//! human-readable field map type
- template<class FieldCollection, class EigenArray, Map_t map_type>
+ template<class FieldCollection, class EigenArray, Map_t map_type,
+ bool ConstField>
std::string
- MatrixLikeFieldMap<FieldCollection, EigenArray, map_type>::
+ MatrixLikeFieldMap<FieldCollection, EigenArray, map_type, ConstField>::
info_string() const {
std::stringstream info;
info << field_info_root << "("
<< typeid(typename EigenArray::value_type).name() << ", "
<< EigenArray::RowsAtCompileTime << "x"
<< EigenArray::ColsAtCompileTime << ")";
return info.str();
}
/* ---------------------------------------------------------------------- */
//! member access
- template<class FieldCollection, class EigenArray, Map_t map_type>
+ template<class FieldCollection, class EigenArray, Map_t map_type,
+ bool ConstField>
template <class ref>
ref
- MatrixLikeFieldMap<FieldCollection, EigenArray, map_type>::
+ MatrixLikeFieldMap<FieldCollection, EigenArray, map_type, ConstField>::
operator[](size_type index) {
return ref(this->get_ptr_to_entry(index));
}
- template<class FieldCollection, class EigenArray, Map_t map_type>
- typename MatrixLikeFieldMap<FieldCollection, EigenArray, map_type>::reference
- MatrixLikeFieldMap<FieldCollection, EigenArray, map_type>::
+ template<class FieldCollection, class EigenArray, Map_t map_type,
+ bool ConstField>
+ typename MatrixLikeFieldMap<FieldCollection, EigenArray, map_type, ConstField>::reference
+ MatrixLikeFieldMap<FieldCollection, EigenArray, map_type, ConstField>::
operator[](const Ccoord & ccoord) {
auto && index = this->collection.get_index(ccoord);
return reference(this->get_ptr_to_entry(std::move(index)));
}
/* ---------------------------------------------------------------------- */
- template<class FieldCollection, class EigenArray, Map_t map_type>
- typename MatrixLikeFieldMap<FieldCollection, EigenArray, map_type>::pointer
- MatrixLikeFieldMap<FieldCollection, EigenArray, map_type>::
+ template<class FieldCollection, class EigenArray, Map_t map_type,
+ bool ConstField>
+ typename MatrixLikeFieldMap<FieldCollection, EigenArray, map_type, ConstField>::pointer
+ MatrixLikeFieldMap<FieldCollection, EigenArray, map_type, ConstField>::
ptr_to_val_t(size_type index) {
return std::make_unique<value_type>
(this->get_ptr_to_entry(std::move(index)));
}
} // internal
/* ---------------------------------------------------------------------- */
//! short-hand for an Eigen matrix map as iterate
- template <class FieldCollection, typename T, Dim_t NbRows, Dim_t NbCols>
+ template <class FieldCollection, typename T, Dim_t NbRows, Dim_t NbCols,
+ bool ConstField=false>
using MatrixFieldMap = internal::MatrixLikeFieldMap
<FieldCollection,
- Eigen::Map<Eigen::Matrix<T, NbRows, NbCols>>,
- internal::Map_t::Matrix>;
+ std::conditional_t<ConstField,
+ Eigen::Map<const Eigen::Matrix<T, NbRows, NbCols>>,
+ Eigen::Map<Eigen::Matrix<T, NbRows, NbCols>>>,
+ internal::Map_t::Matrix, ConstField>;
/* ---------------------------------------------------------------------- */
//! short-hand for an Eigen matrix map as iterate
template <class FieldCollection, typename T, Dim_t Dim,
bool MapConst=false, bool Symmetric=false>
using T4MatrixFieldMap = internal::MatrixLikeFieldMap
<FieldCollection,
T4Map<T, Dim, MapConst, Symmetric>,
- internal::Map_t::T4Matrix>;
+ internal::Map_t::T4Matrix, MapConst>;
/* ---------------------------------------------------------------------- */
//! short-hand for an Eigen array map as iterate
- template <class FieldCollection, typename T, Dim_t NbRows, Dim_t NbCols>
+ template <class FieldCollection, typename T, Dim_t NbRows, Dim_t NbCols,
+ bool ConstField=false>
using ArrayFieldMap = internal::MatrixLikeFieldMap
<FieldCollection,
- Eigen::Map<Eigen::Array<T, NbRows, NbCols>>,
- internal::Map_t::Array>;
+ std::conditional_t<ConstField,
+ Eigen::Map<const Eigen::Array<T, NbRows, NbCols>>,
+ Eigen::Map<Eigen::Array<T, NbRows, NbCols>>>,
+ internal::Map_t::Array, ConstField>;
} // muSpectre
#endif /* FIELD_MAP_MATRIXLIKE_H */
diff --git a/src/common/field_map_scalar.hh b/src/common/field_map_scalar.hh
index 8f74a98..4ad0897 100644
--- a/src/common/field_map_scalar.hh
+++ b/src/common/field_map_scalar.hh
@@ -1,139 +1,165 @@
/**
* file field_map_scalar.hh
*
* @author Till Junge <till.junge@epfl.ch>
*
* @date 26 Sep 2017
*
* @brief maps over scalar fields
*
* @section LICENCE
*
* Copyright (C) 2017 Till Junge
*
* µSpectre is free software; you can redistribute it and/or
* modify it under the terms of the GNU General Public License as
* published by the Free Software Foundation, either version 3, or (at
* your option) any later version.
*
* µSpectre is distributed in the hope that it will be useful, but
* WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
* General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with GNU Emacs; see the file COPYING. If not, write to the
* Free Software Foundation, Inc., 59 Temple Place - Suite 330,
* Boston, MA 02111-1307, USA.
*/
#ifndef FIELD_MAP_SCALAR_H
#define FIELD_MAP_SCALAR_H
#include "common/field_map_base.hh"
namespace muSpectre {
- template <class FieldCollection, typename T>
+ template <class FieldCollection, typename T, bool ConstField=false>
class ScalarFieldMap
- : public internal::FieldMap<FieldCollection, T, 1>
+ : public internal::FieldMap<FieldCollection, T, 1, ConstField>
{
public:
- using parent = internal::FieldMap<FieldCollection, T, 1>;
+ using parent = internal::FieldMap<FieldCollection, T, 1, ConstField>;
using Ccoord = Ccoord_t<FieldCollection::spatial_dim()>;
using value_type = T;
- using reference = value_type &;
using const_reference = const value_type &;
+ using reference = std::conditional_t<ConstField,
+ const_reference,
+ value_type &>;
using size_type = typename parent::size_type;
using pointer = T*;
using Field = typename parent::Field;
- using iterator = typename parent::template iterator<ScalarFieldMap>;
using const_iterator= typename parent::template iterator<ScalarFieldMap, true>;
+ using iterator = std::conditional_t
+ <ConstField,
+ const_iterator,
+ typename parent::template iterator<ScalarFieldMap>>;
using reverse_iterator = std::reverse_iterator<iterator>;
using const_reverse_iterator = std::reverse_iterator<const_iterator>;
friend iterator;
//! Default constructor
ScalarFieldMap() = delete;
- ScalarFieldMap(Field & field);
+ template <bool isntConst=!ConstField>
+ ScalarFieldMap(std::enable_if_t<isntConst, Field &> field);
+ template <bool isConst=ConstField>
+ ScalarFieldMap(std::enable_if_t<isConst, const Field &> field);
//! Copy constructor
ScalarFieldMap(const ScalarFieldMap &other) = delete;
//! Move constructor
ScalarFieldMap(ScalarFieldMap &&other) noexcept = delete;
//! Destructor
virtual ~ScalarFieldMap() noexcept = default;
//! Copy assignment operator
ScalarFieldMap& operator=(const ScalarFieldMap &other) = delete;
//! Move assignment operator
ScalarFieldMap& operator=(ScalarFieldMap &&other) noexcept = delete;
//! give human-readable field map type
inline std::string info_string() const override final;
//! member access
template <class ref_t = reference>
inline ref_t operator[](size_type index);
inline reference operator[](const Ccoord& ccoord);
//! return an iterator to head of field for ranges
inline iterator begin(){return iterator(*this);}
inline const_iterator cbegin(){return const_iterator(*this);}
//! return an iterator to tail of field for ranges
inline iterator end(){return iterator(*this, false);}
inline const_iterator cend(){return const_iterator(*this, false);}
protected:
//! for sad, legacy iterator use
inline pointer ptr_to_val_t(size_type index);
const static std::string field_info_root;
private:
};
+
+ /* ---------------------------------------------------------------------- */
+ template <class FieldCollection, typename T, bool ConstField>
+ template <bool isntConst>
+ ScalarFieldMap<FieldCollection, T, ConstField>::
+ ScalarFieldMap(std::enable_if_t<isntConst, Field &> field)
+ :parent(field) {
+ static_assert((isntConst != ConstField),
+ "let the compiler deduce isntConst, this is a SFINAE "
+ "parameter");
+ this->check_compatibility();
+ }
+
/* ---------------------------------------------------------------------- */
- template <class FieldCollection, typename T>
- ScalarFieldMap<FieldCollection, T>::ScalarFieldMap(Field & field)
+ template <class FieldCollection, typename T, bool ConstField>
+ template <bool isConst>
+ ScalarFieldMap<FieldCollection, T, ConstField>::
+ ScalarFieldMap(std::enable_if_t<isConst, const Field &> field)
:parent(field) {
+ static_assert((isConst == ConstField),
+ "let the compiler deduce isntConst, this is a SFINAE "
+ "parameter");
this->check_compatibility();
}
/* ---------------------------------------------------------------------- */
//! human-readable field map type
- template<class FieldCollection, typename T>
+ template<class FieldCollection, typename T, bool ConstField>
std::string
- ScalarFieldMap<FieldCollection, T>::info_string() const {
+ ScalarFieldMap<FieldCollection, T, ConstField>::info_string() const {
std::stringstream info;
info << "Scalar(" << typeid(T).name() << ")";
return info.str();
}
/* ---------------------------------------------------------------------- */
- template <class FieldCollection, typename T>
- const std::string ScalarFieldMap<FieldCollection, T>::field_info_root{
+ template <class FieldCollection, typename T, bool ConstField>
+ const std::string ScalarFieldMap<FieldCollection, T, ConstField>::field_info_root{
"Scalar"};
/* ---------------------------------------------------------------------- */
//! member access
- template <class FieldCollection, typename T>
+ template <class FieldCollection, typename T, bool ConstField>
template <class ref_t>
ref_t
- ScalarFieldMap<FieldCollection, T>::operator[](size_type index) {
+ ScalarFieldMap<FieldCollection, T, ConstField>::operator[](size_type index) {
return this->get_ref_to_entry(std::move(index));
}
/* ---------------------------------------------------------------------- */
//! member access
- template <class FieldCollection, typename T>
- typename ScalarFieldMap<FieldCollection, T>::reference
- ScalarFieldMap<FieldCollection, T>::operator[](const Ccoord& ccoord) {
+ template <class FieldCollection, typename T, bool ConstField>
+ typename ScalarFieldMap<FieldCollection, T, ConstField>::reference
+ ScalarFieldMap<FieldCollection, T, ConstField>::operator[](const Ccoord& ccoord) {
auto && index = this->collection.get_index(std::move(ccoord));
return this->get_ref_to_entry(std::move(index));
}
} // muSpectre
#endif /* FIELD_MAP_SCALAR_H */
diff --git a/src/common/field_map_tensor.hh b/src/common/field_map_tensor.hh
index 6c4204c..3a9fc87 100644
--- a/src/common/field_map_tensor.hh
+++ b/src/common/field_map_tensor.hh
@@ -1,169 +1,195 @@
/**
* file field_map_tensor.hh
*
* @author Till Junge <till.junge@epfl.ch>
*
* @date 26 Sep 2017
*
* @brief Defines an Eigen-Tensor map over strongly typed fields
*
* @section LICENCE
*
* Copyright (C) 2017 Till Junge
*
* µSpectre is free software; you can redistribute it and/or
* modify it under the terms of the GNU General Public License as
* published by the Free Software Foundation, either version 3, or (at
* your option) any later version.
*
* µSpectre is distributed in the hope that it will be useful, but
* WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
* General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with GNU Emacs; see the file COPYING. If not, write to the
* Free Software Foundation, Inc., 59 Temple Place - Suite 330,
* Boston, MA 02111-1307, USA.
*/
#ifndef FIELD_MAP_TENSOR_H
#define FIELD_MAP_TENSOR_H
#include <sstream>
#include <memory>
#include "common/eigen_tools.hh"
#include "common/field_map_base.hh"
namespace muSpectre {
/* ---------------------------------------------------------------------- */
- template <class FieldCollection, typename T, Dim_t order, Dim_t dim>
+ template <class FieldCollection, typename T, Dim_t order, Dim_t dim,
+ bool ConstField=false>
class TensorFieldMap: public
internal::FieldMap<FieldCollection,
T,
- SizesByOrder<order, dim>::Sizes::total_size
- >
+ SizesByOrder<order, dim>::Sizes::total_size,
+ ConstField>
{
public:
using Parent = internal::FieldMap<FieldCollection, T,
- TensorFieldMap::nb_components>;
+ TensorFieldMap::nb_components, ConstField>;
using Ccoord = Ccoord_t<FieldCollection::spatial_dim()>;
using Sizes = typename SizesByOrder<order, dim>::Sizes;
using T_t = Eigen::TensorFixedSize<T, Sizes>;
using value_type = Eigen::TensorMap<T_t>;
- using reference = value_type; // since it's a resource handle
using const_reference = Eigen::TensorMap<const T_t>;
+ using reference = std::conditional_t<ConstField,
+ const_reference,
+ value_type>; // since it's a resource handle
using size_type = typename Parent::size_type;
using pointer = std::unique_ptr<value_type>;
using TypedField = typename Parent::TypedField;
using Field = typename TypedField::Parent;
- using iterator = typename Parent::template iterator<TensorFieldMap>;
using const_iterator = typename Parent::template iterator<TensorFieldMap, true>;
+ using iterator = std::conditional_t<
+ ConstField,
+ const_iterator,
+ typename Parent::template iterator<TensorFieldMap>>;
using reverse_iterator = std::reverse_iterator<iterator>;
using const_reverse_iterator = std::reverse_iterator<const_iterator>;
friend iterator;
//! Default constructor
TensorFieldMap() = delete;
- TensorFieldMap(Field & field);
+ template <bool isntConst=!ConstField>
+ TensorFieldMap(std::enable_if_t<isntConst, Field &> field);
+ template <bool isConst=ConstField>
+ TensorFieldMap(std::enable_if_t<isConst, const Field &> field);
//! Copy constructor
TensorFieldMap(const TensorFieldMap &other) = delete;
//! Move constructor
TensorFieldMap(TensorFieldMap &&other) noexcept = delete;
//! Destructor
virtual ~TensorFieldMap() noexcept = default;
//! Copy assignment operator
TensorFieldMap& operator=(const TensorFieldMap &other) = delete;
//! Move assignment operator
TensorFieldMap& operator=(TensorFieldMap &&other) noexcept = delete;
//! give human-readable field map type
inline std::string info_string() const override final;
//! member access
template <class ref_t = reference>
inline ref_t operator[](size_type index);
inline reference operator[](const Ccoord & ccoord);
//! return an iterator to head of field for ranges
inline iterator begin(){return iterator(*this);}
inline const_iterator cbegin(){return const_iterator(*this);}
//! return an iterator to tail of field for ranges
inline iterator end(){return iterator(*this, false);}
inline const_iterator cend(){return const_iterator(*this, false);}
protected:
//! for sad, legacy iterator use
inline pointer ptr_to_val_t(size_type index);
private:
};
/* ---------------------------------------------------------------------- */
- template<class FieldCollection, typename T, Dim_t order, Dim_t dim>
- TensorFieldMap<FieldCollection, T, order, dim>::
- TensorFieldMap(Field & field)
+ template<class FieldCollection, typename T, Dim_t order, Dim_t dim,
+ bool ConstField>
+ template <bool isntConst>
+ TensorFieldMap<FieldCollection, T, order, dim, ConstField>::
+ TensorFieldMap(std::enable_if_t<isntConst, Field &> field)
:Parent(field) {
+ static_assert((isntConst != ConstField),
+ "let the compiler deduce isntConst, this is a SFINAE "
+ "parameter");
this->check_compatibility();
}
+ /* ---------------------------------------------------------------------- */
+ template<class FieldCollection, typename T, Dim_t order, Dim_t dim,
+ bool ConstField>
+ template <bool isConst>
+ TensorFieldMap<FieldCollection, T, order, dim, ConstField>::
+ TensorFieldMap(std::enable_if_t<isConst, const Field &> field)
+ :Parent(field) {
+ static_assert((isConst == ConstField),
+ "let the compiler deduce isntConst, this is a SFINAE "
+ "parameter");
+ this->check_compatibility();
+ }
/* ---------------------------------------------------------------------- */
//! human-readable field map type
- template<class FieldCollection, typename T, Dim_t order, Dim_t dim>
+ template<class FieldCollection, typename T, Dim_t order, Dim_t dim, bool ConstField>
std::string
- TensorFieldMap<FieldCollection, T, order, dim>::info_string() const {
+ TensorFieldMap<FieldCollection, T, order, dim, ConstField>::info_string() const {
std::stringstream info;
info << "Tensor(" << typeid(T).name() << ", " << order
<< "_o, " << dim << "_d)";
return info.str();
}
/* ---------------------------------------------------------------------- */
//! member access
- template <class FieldCollection, typename T, Dim_t order, Dim_t dim>
+ template <class FieldCollection, typename T, Dim_t order, Dim_t dim, bool ConstField>
template <class ref_t>
ref_t
- TensorFieldMap<FieldCollection, T, order, dim>::operator[](size_type index) {
+ TensorFieldMap<FieldCollection, T, order, dim, ConstField>::operator[](size_type index) {
auto && lambda = [this, &index](auto&&...tens_sizes) {
return ref_t(this->get_ptr_to_entry(index), tens_sizes...);
};
return call_sizes<order, dim>(lambda);
}
- template<class FieldCollection, typename T, Dim_t order, Dim_t dim>
- typename TensorFieldMap<FieldCollection, T, order, dim>::reference
- TensorFieldMap<FieldCollection, T, order, dim>::
+ template<class FieldCollection, typename T, Dim_t order, Dim_t dim, bool ConstField>
+ typename TensorFieldMap<FieldCollection, T, order, dim, ConstField>::reference
+ TensorFieldMap<FieldCollection, T, order, dim, ConstField>::
operator[](const Ccoord & ccoord) {
auto && index = this->collection.get_index(ccoord);
auto && lambda = [this, &index](auto&&...sizes) {
return ref_t(this->get_ptr_to_entry(index), sizes...);
};
return call_sizes<order, dim>(lambda);
}
/* ---------------------------------------------------------------------- */
//! for sad, legacy iterator use. Don't use unless you have to.
- template<class FieldCollection, typename T, Dim_t order, Dim_t dim>
- typename TensorFieldMap<FieldCollection, T, order, dim>::pointer
- TensorFieldMap<FieldCollection, T, order, dim>::
+ template<class FieldCollection, typename T, Dim_t order, Dim_t dim, bool ConstField>
+ typename TensorFieldMap<FieldCollection, T, order, dim, ConstField>::pointer
+ TensorFieldMap<FieldCollection, T, order, dim, ConstField>::
ptr_to_val_t(size_type index) {
auto && lambda = [this, &index](auto&&... tens_sizes) {
return std::make_unique<value_type>(this->get_ptr_to_entry(index),
tens_sizes...);
};
return call_sizes<order, dim>(lambda);
}
} // muSpectre
#endif /* FIELD_MAP_TENSOR_H */
diff --git a/src/materials/material_hyper_elastic1.hh b/src/materials/material_hyper_elastic1.hh
index ff108d9..dd9739c 100644
--- a/src/materials/material_hyper_elastic1.hh
+++ b/src/materials/material_hyper_elastic1.hh
@@ -1,104 +1,104 @@
/**
* file material_hyper_elastic1.hh
*
* @author Till Junge <till.junge@epfl.ch>
*
* @date 13 Nov 2017
*
* @brief Implementation for hyperelastic reference material like in de Geus
* 2017. This follows the simplest and likely not most efficient
* implementation (with exception of the Python law)
*
* @section LICENCE
*
* Copyright (C) 2017 Till Junge
*
* µSpectre is free software; you can redistribute it and/or
* modify it under the terms of the GNU General Public License as
* published by the Free Software Foundation, either version 3, or (at
* your option) any later version.
*
* µSpectre is distributed in the hope that it will be useful, but
* WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
* General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with GNU Emacs; see the file COPYING. If not, write to the
* Free Software Foundation, Inc., 59 Temple Place - Suite 330,
* Boston, MA 02111-1307, USA.
*/
#include "materials/material_muSpectre_base.hh"
#ifndef MATERIAL_HYPER_ELASTIC1_H
#define MATERIAL_HYPER_ELASTIC1_H
namespace muSpectre {
//! DimS spatial dimension (dimension of problem
//! DimM material_dimension (dimension of constitutive law)
template<Dim_t DimS, Dim_t DimM>
class MaterialHyperElastic1: public MaterialMuSpectre<MaterialHyperElastic1<DimS, DimM>, DimS, DimM>
{
public:
using Parent = MaterialMuSpectre<MaterialHyperElastic1, DimS, DimM>;
using NeedTangent = typename Parent::NeedTangent;
using GFieldCollection_t = typename Parent::GFieldCollection_t;
// declare what type of strain measure your law takes as input
constexpr static auto strain_measure{StrainMeasure::GreenLagrange};
// declare what type of stress measure your law yields as output
constexpr static auto stress_measure{StressMeasure::PK2};
// declare whether the derivative of stress with respect to strain is uniform
constexpr static bool uniform_stiffness = true;
// declare the type in which you wish to receive your strain measure
- using StrainMap_t = MatrixFieldMap<GFieldCollection_t, Real, DimM, DimM>;
- using StressMap_t = StrainMap_t;
+ using StrainMap_t = MatrixFieldMap<GFieldCollection_t, Real, DimM, DimM, true>;
+ using StressMap_t = MatrixFieldMap<GFieldCollection_t, Real, DimM, DimM>;
using TangentMap_t = T4MatrixFieldMap<GFieldCollection_t, Real, DimM>;
using Strain_t = typename StrainMap_t::const_reference;
using Stress_t = typename StressMap_t::reference;
using Tangent_t = typename TangentMap_t::reference::ConstType;
using Stiffness_t = Eigen::TensorFixedSize
<Real, Eigen::Sizes<DimM, DimM, DimM, DimM>>;
using InternalVariables = typename Parent::DefaultInternalVariables;
//! Default constructor
MaterialHyperElastic1() = delete;
//! Copy constructor
MaterialHyperElastic1(const MaterialHyperElastic1 &other) = delete;
//! Construct by name, Young's modulus and Poisson's ratio
MaterialHyperElastic1(std::string name, Real young, Real poisson);
//! Move constructor
MaterialHyperElastic1(MaterialHyperElastic1 &&other) noexcept = delete;
//! Destructor
virtual ~MaterialHyperElastic1() noexcept = default;
//! Copy assignment operator
MaterialHyperElastic1& operator=(const MaterialHyperElastic1 &other) = delete;
//! Move assignment operator
MaterialHyperElastic1& operator=(MaterialHyperElastic1 &&other) noexcept = delete;
template <class s_t>
decltype(auto) evaluate_stress(s_t && E);
template <class s_t>
decltype(auto) evaluate_stress_tangent(s_t && E);
const Tangent_t & get_stiffness() const;
protected:
const Real young, poisson, lambda, mu;
const Stiffness_t C;
private:
};
} // muSpectre
#endif /* MATERIAL_HYPER_ELASTIC1_H */
diff --git a/src/materials/material_muSpectre_base.hh b/src/materials/material_muSpectre_base.hh
index 72cdeb5..d1ad361 100644
--- a/src/materials/material_muSpectre_base.hh
+++ b/src/materials/material_muSpectre_base.hh
@@ -1,503 +1,508 @@
/**
* file material_muSpectre_base.hh
*
* @author Till Junge <till.junge@epfl.ch>
*
* @date 25 Oct 2017
*
* @brief Base class for materials written for µSpectre specifically. These
* can take full advantage of the configuration-change utilities of
* µSpectre. The user can inherit from them to define new constitutive
* laws and is merely required to provide the methods for computing the
* second Piola-Kirchhoff stress and Tangent. This class uses the
* "curiously recurring template parameter" to avoid virtual calls.
*
* @section LICENCE
*
* Copyright (C) 2017 Till Junge
*
* µSpectre is free software; you can redistribute it and/or
* modify it under the terms of the GNU General Public License as
* published by the Free Software Foundation, either version 3, or (at
* your option) any later version.
*
* µSpectre is distributed in the hope that it will be useful, but
* WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
* General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with GNU Emacs; see the file COPYING. If not, write to the
* Free Software Foundation, Inc., 59 Temple Place - Suite 330,
* Boston, MA 02111-1307, USA.
*/
#include <tuple>
#include <type_traits>
#include <iterator>
#include <stdexcept>
#include <boost/range/combine.hpp>
#include "materials/material_base.hh"
#include "materials/materials_toolbox.hh"
#include "common/field_collection.hh"
#include "common/field.hh"
#ifndef MATERIAL_MUSPECTRE_BASE_H
#define MATERIAL_MUSPECTRE_BASE_H
namespace muSpectre {
//! 'Material' is a CRTP
template<class Material, Dim_t DimS, Dim_t DimM>
class MaterialMuSpectre: public MaterialBase<DimS, DimM>
{
public:
enum class NeedTangent {yes, no};
using Parent = MaterialBase<DimS, DimM>;
using GFieldCollection_t = typename Parent::GFieldCollection_t;
using MFieldCollection_t = typename Parent::MFieldCollection_t;
using StressField_t = typename Parent::StressField_t;
using StrainField_t = typename Parent::StrainField_t;
using TangentField_t = typename Parent::TangentField_t;
using DefaultInternalVariables = std::tuple<>;
//! Default constructor
MaterialMuSpectre() = delete;
//! Construct by name
MaterialMuSpectre(std::string name);
//! Copy constructor
MaterialMuSpectre(const MaterialMuSpectre &other) = delete;
//! Move constructor
MaterialMuSpectre(MaterialMuSpectre &&other) noexcept = delete;
//! Destructor
virtual ~MaterialMuSpectre() noexcept = default;
//! Copy assignment operator
MaterialMuSpectre& operator=(const MaterialMuSpectre &other) = delete;
//! Move assignment operator
MaterialMuSpectre& operator=(MaterialMuSpectre &&other) noexcept = delete;
//! computes stress
virtual void compute_stresses(const StrainField_t & F,
StressField_t & P,
Formulation form) override final;
//! computes stress and tangent modulus
virtual void compute_stresses_tangent(const StrainField_t & F,
StressField_t & P,
TangentField_t & K,
Formulation form) override final;
protected:
//! computes stress with the formulation available at compile time
template <Formulation Form>
inline void compute_stresses_worker(const StrainField_t & F,
StressField_t & P);
//! computes stress with the formulation available at compile time
template <Formulation Form>
inline void compute_stresses_worker(const StrainField_t & F,
StressField_t & P,
TangentField_t & K);
//! this iterable class is a default for simple laws that just take a strain
//! the iterable is just a templated wrapper to provide a range to iterate over
//! that does or does not include tangent moduli
//template<NeedTangent need_tgt = NeedTangent::no>
//class iterable_proxy; TODO: clean this
/**
* get a combined range that can be used to iterate over all fields
*/
inline decltype(auto) get_zipped_fields(const StrainField_t & F,
StressField_t & P);
inline decltype(auto) get_zipped_fields(const StrainField_t & F,
StressField_t & P,
TangentField_t & K);
template <class... Args>
inline decltype(auto) get_zipped_fields_worker(Args && ... args);
/**
* inheriting classes with internal variables need to overload this function
*/
decltype(auto) get_internals() {
return typename Material::InternalVariables{};}
private:
};
/* ---------------------------------------------------------------------- */
template <class Material, Dim_t DimS, Dim_t DimM>
MaterialMuSpectre<Material, DimS, DimM>::
MaterialMuSpectre(std::string name):Parent(name) {
using stress_compatible = typename Material::StressMap_t::
template is_compatible<StressField_t>;
using strain_compatible = typename Material::StrainMap_t::
template is_compatible<StrainField_t>;
using tangent_compatible = typename Material::TangentMap_t::
template is_compatible<TangentField_t>;
static_assert((stress_compatible::value &&
stress_compatible::explain()),
"The material's declared stress map is not compatible "
"with the stress field. More info in previously shown "
"assert.");
static_assert((strain_compatible::value &&
strain_compatible::explain()),
"The material's declared strain map is not compatible "
"with the strain field. More info in previously shown "
"assert.");
static_assert((tangent_compatible::value &&
tangent_compatible::explain()),
"The material's declared tangent map is not compatible "
"with the tangent field. More info in previously shown "
"assert.");
}
/* ---------------------------------------------------------------------- */
template <class Material, Dim_t DimS, Dim_t DimM>
void MaterialMuSpectre<Material, DimS, DimM>::
compute_stresses(const StrainField_t &F, StressField_t &P,
muSpectre::Formulation form) {
switch (form) {
case Formulation::finite_strain: {
this->template compute_stresses_worker<Formulation::finite_strain>(F, P);
break;
}
case Formulation::small_strain: {
this->template compute_stresses_worker<Formulation::small_strain>(F, P);
break;
}
default:
throw std::runtime_error("Unknown formulation");
break;
}
}
/* ---------------------------------------------------------------------- */
template <class Material, Dim_t DimS, Dim_t DimM>
void MaterialMuSpectre<Material, DimS, DimM>::
compute_stresses_tangent(const StrainField_t & F, StressField_t & P,
TangentField_t & K,
muSpectre::Formulation form) {
switch (form) {
case Formulation::finite_strain: {
this->template compute_stresses_worker<Formulation::finite_strain>(F, P, K);
break;
}
case Formulation::small_strain: {
this->template compute_stresses_worker<Formulation::small_strain>(F, P, K);
break;
}
default:
throw std::runtime_error("Unknown formulation");
break;
}
}
/* ---------------------------------------------------------------------- */
template <class Material, Dim_t DimS, Dim_t DimM>
template <Formulation Form>
void MaterialMuSpectre<Material, DimS, DimM>::
compute_stresses_worker(const StrainField_t & F,
StressField_t & P,
TangentField_t & K){
//! in small strain problems, we always want Cauchy stress as a function
//! of the infinitesimal strain tensor
// using IterableSmallStrain_t = iterable_proxy<NeedTangent::yes>;
// using IterableFiniteStrain_t = iterable_proxy<NeedTangent::yes>;
// using iterable =
// std::conditional_t
// <(Form==Formulation::small_strain),
// IterableSmallStrain_t,
// IterableFiniteStrain_t>;
// TODO: CLEAN THIS
/* These lambdas are executed for every integration point.
F contains the transformation gradient for finite strain calculations and
the infinitesimal strain tensor in small strain problems
The internal_variables tuple contains whatever internal variables
Material declared (e.g., eigenstrain, strain rate, etc.)
*/
auto constitutive_law_small_strain = [this]
(auto && F, auto && ... internal_variables) {
constexpr StrainMeasure stored_strain_m{get_stored_strain_type(Form)};
constexpr StrainMeasure expected_strain_m{
get_formulation_strain_type(Form, Material::strain_measure)};
auto && strain = MatTB::convert_strain<stored_strain_m, expected_strain_m>(F);
// return value contains a tuple of rvalue_refs to both stress and tangent moduli
auto && stress_tgt = static_cast<Material&>(*this).evaluate_stress_tangent(std::move(strain), internal_variables...);
return std::move(stress_tgt);
};
auto constitutive_law_finite_strain = [this]
(auto && F, auto && ... internal_variables) {
constexpr StrainMeasure stored_strain_m{get_stored_strain_type(Form)};
constexpr StrainMeasure expected_strain_m{
get_formulation_strain_type(Form, Material::strain_measure)};
auto && strain = MatTB::convert_strain<stored_strain_m, expected_strain_m>(F);
// return value contains a tuple of rvalue_refs to both stress and tangent moduli
auto && stress_tgt = static_cast<Material&>(*this).evaluate_stress_tangent(std::move(strain), internal_variables...);
auto && stress = std::get<0>(stress_tgt);
auto && tangent = std::get<1>(stress_tgt);
return MatTB::PK1_stress<Material::stress_measure, Material::strain_measure>
(F, stress, tangent);
};
- MatrixFieldMap<GFieldCollection_t, Real, DimM, DimM> Pmap(P);
- const MatrixFieldMap<GFieldCollection_t, Real, DimM, DimM> Fmap(F);
auto it{this->get_zipped_fields(F, P, K)};
for (auto && tuples: it) {
// the iterator yields a pair of tuples. this first tuple contains
// references to stress and stiffness in the global arrays, the second
// contains references to the deformation gradient and internal variables
// (some of them are const).
auto && stress_tgt = std::get<0>(tuples);
auto && inputs = std::get<1>(tuples);
switch (Form) {
case Formulation::small_strain: {
stress_tgt = std::apply(constitutive_law_small_strain, std::move(inputs));
break;
}
case Formulation::finite_strain: {
stress_tgt = std::apply(constitutive_law_finite_strain, std::move(inputs));
break;
}
}
}
}
/* ---------------------------------------------------------------------- */
template <class Material, Dim_t DimS, Dim_t DimM>
decltype(auto)
MaterialMuSpectre<Material, DimS, DimM>::
get_zipped_fields(const StrainField_t & F, StressField_t & P) {
- return this->get_zipped_fields_worker(F, P);
+ typename Material::StrainMap_t Fmap(F);//TODO: think about whether F and P should be std::forward<>()'ed
+ typename Material::StressMap_t Pmap(P);
+
+ return this->get_zipped_fields_worker(std::move(Fmap),
+ std::move(Pmap));
}
/* ---------------------------------------------------------------------- */
template <class Material, Dim_t DimS, Dim_t DimM>
decltype(auto)
MaterialMuSpectre<Material, DimS, DimM>::
get_zipped_fields(const StrainField_t & F, StressField_t & P,
TangentField_t & K) {
- auto test = boost::begin(F);
- std::cout << test ;
+ typename Material::StrainMap_t Fmap(F);
+ typename Material::StressMap_t Pmap(P);
+ typename Material::TangentMap_t Kmap(K);
- return this->get_zipped_fields_worker(F, P, K);
+ return this->get_zipped_fields_worker(std::move(Fmap),
+ std::move(Pmap),
+ std::move(Kmap));
}
/* ---------------------------------------------------------------------- */
template <class Material, Dim_t DimS, Dim_t DimM>
template <class... Args>
decltype(auto)
MaterialMuSpectre<Material, DimS, DimM>::
get_zipped_fields_worker(Args && ... args) {
return std::apply
([this, &args...] (auto & ... internal_fields) {
return boost::combine(args..., internal_fields...);},
this->get_internals());
}
//TODO: CLEAN BELOW
///* ---------------------------------------------------------------------- */
////! this iterator class is a default for simple laws that just take a strain
//template <class Material, Dim_t DimS, Dim_t DimM>
//template <typename MaterialMuSpectre<Material, DimS, DimM>::NeedTangent NeedTgt>
//class MaterialMuSpectre<Material, DimS, DimM>::iterable_proxy {
//public:
// //! Default constructor
// iterable_proxy() = delete;
// using NeedTangent =
// typename MaterialMuSpectre<Material, DimS, DimM>::NeedTangent;
// /** Iterator uses the material's internal variables field
// collection to iterate selectively over the global fields
// (such as the transformation gradient F and first
// Piola-Kirchhoff stress P.
// **/
// template<bool DoNeedTgt=(NeedTgt == NeedTangent::yes)>
// iterable_proxy(const MaterialMuSpectre & mat,
// const StrainField_t & F,
// StressField_t & P,
// std::enable_if_t<DoNeedTgt, TangentField_t> & K);
//
// template<bool DontNeedTgt=(NeedTgt == NeedTangent::no)>
// iterable_proxy(const MaterialMuSpectre & mat,
// const StrainField_t & F,
// std::enable_if_t<DontNeedTgt, StressField_t> & P);
//
// using Strain_t = typename Material::Strain_t;
// using Stress_t = typename Material::Stress_t;
// using Tangent_t = typename Material::Tangent_t;
//
// //! Copy constructor
// iterable_proxy(const iterable_proxy &other) = default;
//
// //! Move constructor
// iterable_proxy(iterable_proxy &&other) noexcept = default;
//
// //! Destructor
// virtual ~iterable_proxy() noexcept = default;
//
// //! Copy assignment operator
// iterable_proxy& operator=(const iterable_proxy &other) = default;
//
// //! Move assignment operator
// iterable_proxy& operator=(iterable_proxy &&other) noexcept = default;
//
// class iterator
// {
// public:
// using value_type =
// std::conditional_t
// <(NeedTgt == NeedTangent::yes),
// std::tuple<std::tuple<Stress_t, Tangent_t>, std::tuple<Strain_t>>,
// std::tuple<std::tuple<Stress_t>, std::tuple<Strain_t>>>;
// using iterator_category = std::forward_iterator_tag;
//
// //! Default constructor
// iterator() = delete;
//
// /** Iterator uses the material's internal variables field
// collection to iterate selectively over the global fields
// (such as the transformation gradient F and first
// Piola-Kirchhoff stress P.
// **/
// iterator(const iterable_proxy & it, bool begin = true)
// : it{it}, index{begin ? 0:it.mat.internal_fields.size()}{}
//
//
// //! Copy constructor
// iterator(const iterator &other) = default;
//
// //! Move constructor
// iterator(iterator &&other) noexcept = default;
//
// //! Destructor
// virtual ~iterator() noexcept = default;
//
// //! Copy assignment operator
// iterator& operator=(const iterator &other) = default;
//
// //! Move assignment operator
// iterator& operator=(iterator &&other) noexcept = default;
//
// //! pre-increment
// inline iterator & operator++();
// //! dereference
// inline value_type operator*();
// //! inequality
// inline bool operator!=(const iterator & other) const;
//
//
// protected:
// const iterable_proxy & it;
// size_t index;
// private:
// };
//
// iterator begin() {return iterator(*this);}
// iterator end() {return iterator(*this, false);}
//
//protected:
// using StressTup = std::conditional_t
// <(NeedTgt == NeedTangent::yes),
// std::tuple<StressField_t&, TangentField_t&>,
// std::tuple<StressField_t&>>;
// const MaterialMuSpectre & material;
// const StrainField_t & strain_field;
// auto ZippedFields;
// std::tupleStrainField_t
// StressTup stress_tup;
//private:
//};
//
///* ---------------------------------------------------------------------- */
//template <class Material, Dim_t DimS, Dim_t DimM>
//template <typename MaterialMuSpectre<Material, DimS, DimM>::NeedTangent Tangent>
//template <bool DoNeedTgt>
//MaterialMuSpectre<Material, DimS, DimM>::iterable_proxy<Tangent>::
//iterable_proxy(const MaterialMuSpectre<Material, DimS, DimM> & mat,
// const StrainField_t & F,
// StressField_t & P,
// std::enable_if_t<DoNeedTgt, TangentField_t> & K)
// : material{mat}, strain_field{F}, stress_tup(P, K),
// ZippedFields{std::apply(boost::combine,
// std::tuple_cat(std::tuple<F>,
// stress_tup,
// get_internals()))}
//{
// static_assert((Tangent == NeedTangent::yes),
// "Explicitly choosing the template parameter 'DoNeedTgt' "
// "breaks the SFINAE intention of thi code. Let it be deduced "
// "by the compiler.");
//}
//
///* ---------------------------------------------------------------------- */
//template <class Material, Dim_t DimS, Dim_t DimM>
//template <typename MaterialMuSpectre<Material, DimS, DimM>::NeedTangent Tangent>
//template <bool DontNeedTgt>
//MaterialMuSpectre<Material, DimS, DimM>::iterable_proxy<Tangent>::
//iterable_proxy(const MaterialMuSpectre<Material, DimS, DimM> & mat,
// const StrainField_t & F,
// std::enable_if_t<DontNeedTgt, StressField_t> & P)
// : material{mat}, strain_field{F}, stress_tup(P),
// ZippedFields{std::apply(boost::combine,
// std::tuple_cat(std::tuple<F>,
// stress_tup,
// get_internals()))} {
// static_assert((Tangent == NeedTangent::no),
// "Explicitly choosing the template parameter 'DontNeedTgt' "
// "breaks the SFINAE intention of thi code. Let it be deduced "
// "by the compiler.");
//}
//
///* ---------------------------------------------------------------------- */
////! pre-increment
//template<class Material, Dim_t DimS, Dim_t DimM>
//template<typename MaterialMuSpectre<Material, DimS, DimM>::NeedTangent Tangent>
//
//typename MaterialMuSpectre<Material, DimS, DimM>::
//template iterable_proxy<Tangent>::iterator &
//
//MaterialMuSpectre<Material, DimS, DimM>::
//template iterable_proxy<Tangent>::iterator::
//operator++() {
// ++this->index;
// return *this;
//}
//
///* ---------------------------------------------------------------------- */
////! dereference
//template<class Material, Dim_t DimS, Dim_t DimM>
//template<typename MaterialMuSpectre<Material, DimS, DimM>::NeedTangent Tangent>
//
//typename MaterialMuSpectre<Material, DimS, DimM>::
//template iterable_proxy<Tangent>::iterator::value_type
//
//MaterialMuSpectre<Material, DimS, DimM>::
//template iterable_proxy<Tangent>::iterator::
//operator*() {
// static_assert(Tangent, "Template specialisation failed");
// return std::make_tuple
// (this->material.internal_fields);
//}
} // muSpectre
#endif /* MATERIAL_MUSPECTRE_BASE_H */