Reference¶
-
template<class T>
class ArangeContainer¶ - #include <iterators.hh>
helper class to generate range iterators
Public Types
-
using iterator = iterators::ArangeIterator<T>¶
undocumented
Public Functions
-
using iterator = iterators::ArangeIterator<T>¶
-
template<class T>
class ArangeIterator¶ - #include <iterators.hh>
emulates python’s range iterator
Public Types
-
using iterator_category = std::input_iterator_tag¶
undocumented
Public Functions
-
constexpr ArangeIterator(const ArangeIterator&) = default¶
undocumented
-
inline constexpr ArangeIterator &operator++()¶
undocumented
-
inline constexpr bool operator==(const ArangeIterator &other) const¶
undocumented
-
inline constexpr bool operator!=(const ArangeIterator &other) const¶
undocumented
-
using iterator_category = std::input_iterator_tag¶
-
template<Dim_t order, typename Fun_t, Dim_t dim, Dim_t... args>
struct CallSizesHelper¶ - #include <eigen_tools.hh>
Call a passed lambda with the unpacked sizes as arguments.
-
template<typename Fun_t, Dim_t dim, Dim_t... args>
struct CallSizesHelper<0, Fun_t, dim, args...>¶ - #include <eigen_tools.hh>
Call a passed lambda with the unpacked sizes as arguments.
-
class Cell¶
- #include <cell.hh>
Base class for the representation of a homogenisatonion problem in µSpectre. The
muSpectre::Cell
holds the global strain, stress and (optionally) tangent moduli fields of the problem, maintains the list of materials present, as well as the projection operator.Subclassed by muSpectre::CellSplit
Public Types
-
using Material_ptr = std::unique_ptr<MaterialBase>¶
materials handled through
std::unique_ptr
s
-
using Material_sptr = std::shared_ptr<MaterialBase>¶
-
using Projection_ptr = std::unique_ptr<ProjectionBase>¶
projections handled through
std::unique_ptr
s
-
using Adaptor = CellAdaptor<Cell>¶
adaptor to represent the cell as an Eigen sparse matrix
Public Functions
-
Cell() = delete¶
Deleted default constructor.
-
explicit Cell(Projection_ptr projection, SplitCell is_cell_split = SplitCell::no)¶
Constructor from a projection operator.
-
virtual ~Cell() = default¶
Destructor.
-
bool is_initialised() const¶
for handling double initialisations right
-
Dim_t get_nb_dof() const¶
returns the number of degrees of freedom in the cell
-
size_t get_nb_pixels() const¶
number of pixels on this processor
-
const muFFT::Communicator &get_communicator() const¶
return the communicator object
-
const Formulation &get_formulation() const¶
formulation is hard set by the choice of the projection class
-
Dim_t get_material_dim() const¶
returns the material dimension of the problem
-
void set_uniform_strain(const Eigen::Ref<const Matrix_t>&)¶
set uniform strain (typically used to initialise problems
-
virtual MaterialBase &add_material(Material_ptr mat)¶
add a new material to the cell
-
void complete_material_assignment_simple(MaterialBase &material)¶
By taking a material as input this function assigns all the untouched(not-assigned) pixels to that material
-
void make_pixels_precipitate_for_laminate_material(const std::vector<DynRcoord_t> &precipitate_vertices, MaterialBase &mat_laminate, MaterialBase &mat_precipitate_cell, Material_sptr mat_precipitate, Material_sptr mat_matrix)¶
Given the vertices of polygonal/Polyhedral precipitate, this function assign pixels 1. inside precipitate->mat_precipitate_cell, material at the interface of precipitae-> to mat_precipitate & mat_matrix according to the intersection of pixels with the precipitate
-
template<Dim_t Dim>
void make_pixels_precipitate_for_laminate_material_helper(const std::vector<DynRcoord_t> &precipitate_vertices, MaterialBase &mat_laminate, MaterialBase &mat_precipitate_cell, Material_sptr mat_precipitate, Material_sptr mat_matrix)¶
-
void save_history_variables()¶
freezes all the history variables of the materials
-
std::array<Dim_t, 2> get_strain_shape() const¶
returns the number of rows and cols for the strain matrix type (for full storage, the strain is stored in material_dim × material_dim matrices, but in symmetric storage, it is a column vector)
-
Dim_t get_strain_size() const¶
returns the number of components for the strain matrix type (for full storage, the strain is stored in material_dim × material_dim matrices, but in symmetric storage, it is a column vector)
-
const Dim_t &get_spatial_dim() const¶
return the spatial dimension of the discretisation grid
-
const Dim_t &get_nb_quad() const¶
return the number of quadrature points stored per pixel
-
virtual void check_material_coverage() const¶
makes sure every pixel has been assigned to exactly one material
-
void initialise(muFFT::FFT_PlanFlags flags = muFFT::FFT_PlanFlags::estimate)¶
initialise the projection, the materials and the global fields
-
const muGrid::CcoordOps::DynamicPixels &get_pixels() const¶
return a const reference to the grids pixels iterator
-
muGrid::FieldCollection::IndexIterable get_quad_pt_indices() const¶
return an iterable proxy to this cell’s field collection, iterable by quadrature point
-
muGrid::FieldCollection::PixelIndexIterable get_pixel_indices() const¶
return an iterable proxy to this cell’s field collection, iterable by pixel
-
const muGrid::RealField &get_tangent(bool do_create = false)¶
return a const reference to the cell’s field of tangent moduli
-
virtual const muGrid::RealField &evaluate_stress()¶
evaluates and returns the stress for the currently set strain
-
Eigen_cmap evaluate_stress_eigen()¶
evaluates and returns the stress for the currently set strain
-
virtual std::tuple<const muGrid::RealField&, const muGrid::RealField&> evaluate_stress_tangent()¶
evaluates and returns the stress and tangent moduli for the currently set strain
-
std::tuple<const Eigen_cmap, const Eigen_cmap> evaluate_stress_tangent_eigen()¶
evaluates and returns the stress and tangent moduli for the currently set strain
-
muGrid::RealField &globalise_real_internal_field(const std::string &unique_name)¶
collect the real-valued fields of name
unique_name
of each material in the cell and write their values into a global field of same type and name
-
muGrid::IntField &globalise_int_internal_field(const std::string &unique_name)¶
collect the integer-valued fields of name
unique_name
of each material in the cell and write their values into a global field of same type and name
-
muGrid::UintField &globalise_uint_internal_field(const std::string &unique_name)¶
collect the unsigned integer-valued fields of name
unique_name
of each material in the cell and write their values into a global field of same type and name
-
muGrid::ComplexField &globalise_complex_internal_field(const std::string &unique_name)¶
collect the complex-valued fields of name
unique_name
of each material in the cell and write their values into a global field of same type and name
-
muGrid::GlobalFieldCollection &get_fields()¶
return a reference to the cell’s global fields
-
void apply_projection(muGrid::TypedFieldBase<Real> &field)¶
apply the cell’s projection operator to field
field
(i.e., return G:f)
-
void evaluate_projected_directional_stiffness(const muGrid::TypedFieldBase<Real> &delta_strain, muGrid::TypedFieldBase<Real> &del_stress)¶
evaluates the directional and projected stiffness (this corresponds to G:K:δF (note the negative sign in de Geus 2017, http://dx.doi.org/10.1016/j.cma.2016.12.032).
-
void add_projected_directional_stiffness(EigenCVec_t delta_strain, const Real &alpha, EigenVec_t del_stress)¶
evaluates the directional and projected stiffness (this corresponds to G:K:δF (note the negative sign in de Geus 2017, http://dx.doi.org/10.1016/j.cma.2016.12.032). and then adds it do the values already in del_stress, scaled by alpha (i.e., del_stress += alpha*Q:K:δStrain. This function should not be used directly, as it does absolutely no input checking. Rather, it is meant to be called by the scaleAndAddTo function in the CellAdaptor
-
const ProjectionBase &get_projection() const¶
return a const ref to the projection implementation
-
bool is_point_inside(const DynRcoord_t &point) const¶
check if the pixel is inside of the cell
-
bool is_pixel_inside(const DynCcoord_t &pixel) const¶
check if the point is inside of the cell
Protected Functions
-
template<typename T>
muGrid::TypedField<T> &globalise_internal_field(const std::string &unique_name)¶ helper function for the globalise_<T>_internal_field() functions
Protected Attributes
-
bool initialised = {false}¶
to handle double initialisations right
-
std::vector<Material_ptr> materials = {}¶
container of the materials present in the cell
-
Projection_ptr projection¶
handle for the projection operator
-
std::unique_ptr<muGrid::GlobalFieldCollection> fields¶
handle for the global fields associated with this cell
Protected Static Functions
-
template<Dim_t DimM>
static void apply_directional_stiffness(const muGrid::TypedFieldBase<Real> &delta_strain, const muGrid::TypedFieldBase<Real> &tangent, muGrid::TypedFieldBase<Real> &delta_stress)¶ statically dimensioned worker for evaluating the tangent operator
-
template<Dim_t DimM>
static void add_projected_directional_stiffness_helper(const muGrid::TypedFieldBase<Real> &delta_strain, const muGrid::TypedFieldBase<Real> &tangent, const Real &alpha, muGrid::TypedFieldBase<Real> &delta_stress)¶ statically dimensioned worker for evaluating the incremental tangent operator
-
using Material_ptr = std::unique_ptr<MaterialBase>¶
-
template<class Cell>
class CellAdaptor : public Eigen::EigenBase<CellAdaptor<Cell>>¶ - #include <cell.hh>
Cell adaptors implement the matrix-vector multiplication and allow the system to be used like a sparse matrix in conjugate-gradient-type solvers
lightweight resource handle wrapping a
muSpectre::Cell
or a subclass thereof intoEigen::EigenBase
, so it can be interpreted as a sparse matrix by Eigen solversPublic Types
-
enum [anonymous]¶
Values:
-
enumerator ColsAtCompileTime¶
-
enumerator MaxColsAtCompileTime¶
-
enumerator RowsAtCompileTime¶
-
enumerator MaxRowsAtCompileTime¶
-
enumerator IsRowMajor¶
-
enumerator ColsAtCompileTime¶
-
using Scalar = double¶
sparse matrix traits
-
using RealScalar = double¶
sparse matrix traits
-
using StorageIndex = int¶
sparse matrix traits
Public Functions
-
enum [anonymous]¶
-
class CellSplit : public muSpectre::Cell¶
- #include <cell_split.hh>
DimS spatial dimension (dimension of problem DimM material_dimension (dimension of constitutive law)
Public Types
-
using Projection_ptr = std::unique_ptr<ProjectionBase>¶
projections handled through
std::unique_ptr
s
Public Functions
-
CellSplit() = delete¶
Default constructor.
-
explicit CellSplit(Projection_ptr projection)¶
constructor using sizes and resolution
-
virtual ~CellSplit() = default¶
Destructor.
-
virtual MaterialBase &add_material(Material_ptr mat) final¶
add a new material to the cell
-
void complete_material_assignment(MaterialBase &material)¶
completes the assignmnet of material with a specific material so all the under-assigned pixels would be assigned to a material.
-
std::vector<Real> get_assigned_ratios()¶
-
void make_automatic_precipitate_split_pixels(const std::vector<DynRcoord_t> &preciptiate_vertices, MaterialBase &material)¶
-
std::vector<Real> get_unassigned_ratios_incomplete_pixels() const¶
-
std::vector<int> get_index_incomplete_pixels() const¶
-
std::vector<DynCcoord_t> get_unassigned_pixels()¶
-
IncompletePixels make_incomplete_pixels()¶
-
virtual void check_material_coverage() const final¶
makes sure every pixel has been assigned to materials whose ratios add up to 1.0
Protected Functions
-
void set_p_k_zero()¶
Friends
- friend class Cell
-
using Projection_ptr = std::unique_ptr<ProjectionBase>¶
-
class Communicator¶
- #include <communicator.hh>
stub communicator object that doesn’t communicate anything
Public Functions
-
inline Communicator()¶
-
inline ~Communicator()¶
-
inline int rank() const¶
get rank of present process
-
inline int size() const¶
get total number of processes
-
template<typename T>
inline Matrix_t<T> sum_mat(const Eigen::Ref<Matrix_t<T>> &arg) const¶ sum reduction on EigenMatrix types
-
template<typename T>
inline Matrix_t<T> gather(const Eigen::Ref<Matrix_t<T>> &arg) const¶ gather on EigenMatrix types
Public Static Functions
-
static inline bool has_mpi()¶
find whether the underlying communicator is mpi
-
inline Communicator()¶
-
class ConvergenceError : public muSpectre::SolverError¶
-
template<ElasticModulus Out, ElasticModulus In1, ElasticModulus In2>
struct Converter¶ - #include <materials_toolbox.hh>
Base template for elastic modulus conversion.
Public Static Functions
-
static inline constexpr Real compute(const Real&, const Real&)¶
wrapped function (raison d’être)
-
static inline constexpr Real compute(const Real&, const Real&)¶
-
template<>
struct Converter<ElasticModulus::Bulk, ElasticModulus::lambda, ElasticModulus::Shear>¶ - #include <materials_toolbox.hh>
Specialisation K(λ, µ)
Public Static Functions
-
static inline constexpr Real compute(const Real &lambda, const Real &G)¶
wrapped function (raison d’être)
-
static inline constexpr Real compute(const Real &lambda, const Real &G)¶
-
template<>
struct Converter<ElasticModulus::Bulk, ElasticModulus::Young, ElasticModulus::Poisson>¶ - #include <materials_toolbox.hh>
Specialisation K(E, ν)
Public Static Functions
-
static inline constexpr Real compute(const Real &E, const Real &nu)¶
wrapped function (raison d’être)
-
static inline constexpr Real compute(const Real &E, const Real &nu)¶
-
template<>
struct Converter<ElasticModulus::lambda, ElasticModulus::Bulk, ElasticModulus::Shear>¶ - #include <materials_toolbox.hh>
Specialisation λ(K, µ)
Public Static Functions
-
static inline constexpr Real compute(const Real &K, const Real &mu)¶
wrapped function (raison d’être)
-
static inline constexpr Real compute(const Real &K, const Real &mu)¶
-
template<>
struct Converter<ElasticModulus::lambda, ElasticModulus::Young, ElasticModulus::Poisson>¶ - #include <materials_toolbox.hh>
Specialisation λ(E, ν)
Public Static Functions
-
static inline constexpr Real compute(const Real &E, const Real &nu)¶
wrapped function (raison d’être)
-
static inline constexpr Real compute(const Real &E, const Real &nu)¶
-
template<>
struct Converter<ElasticModulus::Poisson, ElasticModulus::Bulk, ElasticModulus::Shear>¶ - #include <materials_toolbox.hh>
Specialisation ν(K, µ)
Public Static Functions
-
static inline constexpr Real compute(const Real &K, const Real &G)¶
wrapped function (raison d’être)
-
static inline constexpr Real compute(const Real &K, const Real &G)¶
-
template<>
struct Converter<ElasticModulus::Shear, ElasticModulus::Young, ElasticModulus::Poisson>¶ - #include <materials_toolbox.hh>
Specialisation μ(E, ν)
Public Static Functions
-
static inline constexpr Real compute(const Real &E, const Real &nu)¶
wrapped function (raison d’être)
-
static inline constexpr Real compute(const Real &E, const Real &nu)¶
-
template<>
struct Converter<ElasticModulus::Young, ElasticModulus::Bulk, ElasticModulus::Shear>¶ - #include <materials_toolbox.hh>
Specialisation E(K, µ)
Public Static Functions
-
static inline constexpr Real compute(const Real &K, const Real &G)¶
wrapped function (raison d’être)
-
static inline constexpr Real compute(const Real &K, const Real &G)¶
-
template<>
struct Converter<ElasticModulus::Young, ElasticModulus::lambda, ElasticModulus::Shear>¶ - #include <materials_toolbox.hh>
Specialisation E(λ, µ)
Public Static Functions
-
static inline constexpr Real compute(const Real &lambda, const Real &G)¶
wrapped function (raison d’être)
-
static inline constexpr Real compute(const Real &lambda, const Real &G)¶
-
template<ElasticModulus Out, ElasticModulus In>
struct Converter<Out, In, Out>¶ - #include <materials_toolbox.hh>
Spectialisation for when the output is the second input
Public Static Functions
-
static inline constexpr Real compute(const Real&, const Real &B)¶
wrapped function (raison d’être)
-
static inline constexpr Real compute(const Real&, const Real &B)¶
-
template<ElasticModulus Out, ElasticModulus In>
struct Converter<Out, Out, In>¶ - #include <materials_toolbox.hh>
Spectialisation for when the output is the first input
Public Static Functions
-
static inline constexpr Real compute(const Real &A, const Real&)¶
wrapped function (raison d’être)
-
static inline constexpr Real compute(const Real &A, const Real&)¶
-
template<StrainMeasure In, StrainMeasure Out = In>
struct ConvertStrain¶ - #include <materials_toolbox.hh>
Structure for functions returning one strain measure as a function of another
-
template<>
struct ConvertStrain<StrainMeasure::Gradient, StrainMeasure::GreenLagrange>¶ - #include <materials_toolbox.hh>
Specialisation for getting Green-Lagrange strain from the transformation gradient E = ¹/₂ (C - I) = ¹/₂ (Fᵀ·F - I)
-
template<>
struct ConvertStrain<StrainMeasure::Gradient, StrainMeasure::LCauchyGreen>¶ - #include <materials_toolbox.hh>
Specialisation for getting Left Cauchy-Green strain from the transformation gradient B = F·Fᵀ = V²
-
template<>
struct ConvertStrain<StrainMeasure::Gradient, StrainMeasure::Log>¶ - #include <materials_toolbox.hh>
Specialisation for getting logarithmic (Hencky) strain from the transformation gradient E₀ = ¹/₂ ln C = ¹/₂ ln (Fᵀ·F)
-
template<>
struct ConvertStrain<StrainMeasure::Gradient, StrainMeasure::RCauchyGreen>¶ - #include <materials_toolbox.hh>
Specialisation for getting Right Cauchy-Green strain from the transformation gradient C = Fᵀ·F = U²
-
template<Dim_t DimS>
class Correction¶
-
template<>
class Correction<2>¶
-
template<>
class Correction<3>¶
-
template<Dim_t Dim>
struct DefaultOrder¶ - #include <geometry.hh>
convenience structure providing the default order of rotations around (in order) the z, x, and y axis
Public Static Attributes
-
static constexpr RotationOrder value = {RotationOrder::ZXYTaitBryan}¶
holds the value of the rotation order
-
static constexpr RotationOrder value = {RotationOrder::ZXYTaitBryan}¶
-
template<>
struct DefaultOrder<twoD>¶ - #include <geometry.hh>
specialisation for two-dimensional problems
Public Static Attributes
-
static constexpr RotationOrder value = {RotationOrder::Z}¶
holds the value of the rotation order
-
static constexpr RotationOrder value = {RotationOrder::Z}¶
-
class DerivativeBase¶
- #include <derivative.hh>
Representation of a derivative
Subclassed by muFFT::DiscreteDerivative, muFFT::FourierDerivative
Public Functions
-
DerivativeBase() = delete¶
Deleted default constructor.
-
explicit DerivativeBase(Dim_t spatial_dimension)¶
constructor with spatial dimension
-
DerivativeBase(const DerivativeBase &other) = default¶
Copy constructor.
-
DerivativeBase(DerivativeBase &&other) = default¶
Move constructor.
-
virtual ~DerivativeBase() = default¶
Destructor.
-
DerivativeBase &operator=(const DerivativeBase &other) = delete¶
Copy assignment operator.
-
DerivativeBase &operator=(DerivativeBase &&other) = delete¶
Move assignment operator.
Protected Attributes
-
Dim_t spatial_dimension¶
spatial dimension of the problem
-
DerivativeBase() = delete¶
-
class DerivativeError : public runtime_error¶
- #include <derivative.hh>
base class for projection related exceptions
-
template<class Derived>
struct DimCounter¶
-
template<class Derived>
struct DimCounter<Eigen::MatrixBase<Derived>>¶ - #include <T4_map_proxy.hh>
Convenience structure to determine the spatial dimension of a tensor represented by a fixed-size
Eigen::Matrix
. used to derive spatial dimension from input arguments of template functions thus avoiding the need for redundant explicit specification.
-
class DiscreteDerivative : public muFFT::DerivativeBase¶
- #include <derivative.hh>
Representation of a finite-differences stencil
Public Types
-
using Parent = DerivativeBase¶
base class
Public Functions
-
DiscreteDerivative() = delete¶
Default constructor.
-
DiscreteDerivative(DynCcoord_t nb_pts, DynCcoord_t lbounds, const std::vector<Real> &stencil)¶
Constructor with raw stencil information
- Parameters
nb_pts – stencil size
lbounds – relative starting point of stencil
stencil – stencil coefficients
-
DiscreteDerivative(DynCcoord_t nb_pts, DynCcoord_t lbounds, const Eigen::ArrayXd &stencil)¶
Constructor with raw stencil information.
-
DiscreteDerivative(const DiscreteDerivative &other) = default¶
Copy constructor.
-
DiscreteDerivative(DiscreteDerivative &&other) = default¶
Move constructor.
-
virtual ~DiscreteDerivative() = default¶
Destructor.
-
DiscreteDerivative &operator=(const DiscreteDerivative &other) = delete¶
Copy assignment operator.
-
DiscreteDerivative &operator=(DiscreteDerivative &&other) = delete¶
Move assignment operator.
-
inline Real operator()(const DynCcoord_t &dcoord) const¶
Return stencil value.
-
inline const DynCcoord_t &get_nb_pts() const¶
Return number of grid points in stencil.
-
inline const DynCcoord_t &get_lbounds() const¶
Return lower stencil bound.
-
inline virtual Complex fourier(const Vector &phase) const¶
Any translationally invariant linear combination of grid values (as expressed through a “stencil”) becomes a multiplication with a number in Fourier space. This method returns the Fourier representation of this stencil.
-
DiscreteDerivative rollaxes(int distance = 1) const¶
Return a new stencil rolled axes. Given a stencil on a three-dimensional grid with axes (x, y, z), the stencil that has been “rolled” by distance one has axes (z, x, y). This is a simple implementation of a rotation operation. For example, given a stencil that described the derivative in the x-direction, rollaxes(1) gives the derivative in the y-direction and rollaxes(2) gives the derivative in the z-direction.
Protected Attributes
-
const DynCcoord_t nb_pts¶
Number of stencil points.
-
const DynCcoord_t lbounds¶
Lower bound of the finite-differences stencil.
-
using Parent = DerivativeBase¶
-
template<Dim_t Dim>
struct Dotter<Dim, fourthOrder, fourthOrder>¶ - #include <tensor_algebra.hh>
Double contraction between two fourth-rank tensors A and B returns a fourth-rank tensor Cᵢⱼₖₗ = Aᵢⱼₐₑ·Bₐₑₖₗ
-
template<Dim_t Dim>
struct Dotter<Dim, fourthOrder, secondOrder>¶ - #include <tensor_algebra.hh>
Tensor-product between a fourth-rank tensor A and a second-rank tensor B. Returns a fourth-rank Cᵢⱼₖₗ = Aᵢⱼₖₐ·Bₐₗ
-
template<Dim_t Dim>
struct Dotter<Dim, secondOrder, fourthOrder>¶ - #include <tensor_algebra.hh>
Tensor-product between a second-rank tensor A and a fourth-rank tensor B. Returns a fourth-rank Cᵢⱼₖₗ = Aᵢₐ·Bₐⱼₖₗ
-
template<Dim_t Dim>
struct Dotter<Dim, secondOrder, secondOrder>¶ - #include <tensor_algebra.hh>
Double contraction between two second-rank tensors A and B returns a scalar c = AᵢⱼBᵢⱼ
-
class DynamicPixels¶
- #include <ccoord_operations.hh>
Iteration over square (or cubic) discretisation grids. Duplicates capabilities of
muGrid::Ccoordops::Pixels
without needing to be templated with the spatial dimension. Iteration is slower, though.Subclassed by muGrid::CcoordOps::Pixels< Dim >
Public Functions
-
DynamicPixels()¶
-
explicit DynamicPixels(const DynCcoord_t &nb_grid_pts, const DynCcoord_t &locations = DynCcoord_t{})¶
Constructor with default strides (column-major pixel storage order)
-
DynamicPixels(const DynCcoord_t &nb_grid_pts, const DynCcoord_t &locations, const DynCcoord_t &strides)¶
Constructor with custom strides (any, including partially transposed pixel storage order)
-
template<size_t Dim>
explicit DynamicPixels(const Ccoord_t<Dim> &nb_grid_pts, const Ccoord_t<Dim> &locations = Ccoord_t<Dim>{})¶ Constructor with default strides from statically sized coords.
-
template<size_t Dim>
DynamicPixels(const Ccoord_t<Dim> &nb_grid_pts, const Ccoord_t<Dim> &locations, const Ccoord_t<Dim> &strides)¶ Constructor with custom strides from statically sized coords.
-
DynamicPixels(const DynamicPixels &other) = default¶
Copy constructor.
-
DynamicPixels(DynamicPixels &&other) = default¶
Move constructor.
-
virtual ~DynamicPixels() = default¶
Destructor.
-
DynamicPixels &operator=(const DynamicPixels &other) = default¶
Copy assignment operator.
-
DynamicPixels &operator=(DynamicPixels &&other) = default¶
Move assignment operator.
-
inline Dim_t get_index(const DynCcoord_t &ccoord) const¶
evaluate and return the linear index corresponding to dynamic
ccoord
-
template<size_t Dim>
inline Dim_t get_index(const Ccoord_t<Dim> &ccoord) const¶ evaluate and return the linear index corresponding to
ccoord
-
template<size_t Dim>
const Pixels<Dim> &get_dimensioned_pixels() const¶ return a reference to the Pixels object cast into a statically dimensioned grid. the statically dimensioned version duplicates
muGrid::Ccoordops::DynamicPixels
’s capabilities, but iterates much more efficiently.
-
size_t size() const¶
stl conformance
-
inline const DynCcoord_t &get_nb_grid_pts() const¶
return the resolution of the discretisation grid in each spatial dim
-
inline const DynCcoord_t &get_locations() const¶
return the ccoordinates of the bottom, left, (front) pixel/voxel of this processors partition of the discretisation grid. For sequential calculations, this is alvays the origin
-
inline const DynCcoord_t &get_strides() const¶
return the strides used for iterating over the pixels
-
Enumerator enumerate() const¶
iterates in tuples of pixel index ond coordinate. Useful in parallel problems, where simple enumeration of the pixels would be incorrect
Protected Attributes
-
DynCcoord_t nb_grid_pts¶
nb_grid_pts of this domain
-
DynCcoord_t locations¶
locations of this domain
-
DynCcoord_t strides¶
strides of memory layout
-
DynamicPixels()¶
-
template<size_t MaxDim, typename T = Dim_t>
class DynCcoord¶ - #include <grid_common.hh>
Class to represent integer (cell-) coordinates or real-valued coordinates. This class can dynamically accept any spatial-dimension between 1 and MaxDim, and DynCcoord references can be cast to
muGrid::Ccoord_t &
ormuGrid::Rcoord_t &
references. These are used when templating with the spatial dimension of the problem is undesireable/impossible.Public Types
Public Functions
-
inline DynCcoord()¶
default constructor
-
inline DynCcoord(std::initializer_list<T> init_list)¶
constructor from an initialiser list for compound initialisation.
- Parameters
init_list – The length of the initialiser list becomes the spatial dimension of the coordinate, therefore the list must have a length between 1 and MaxDim
-
inline explicit DynCcoord(Dim_t dim)¶
Constructor only setting the dimension. WARNING: This constructor needs regular (round) braces ‘()’, using curly braces ‘{}’ results in the initialiser list constructor being called and creating a DynCcoord with spatial dimension 1
- Parameters
dim – spatial dimension. Needs to be between 1 and MaxDim
-
template<size_t Dim>
inline explicit DynCcoord(const std::array<T, Dim> &ccoord)¶ Constructor from a statically sized coord.
-
~DynCcoord() = default¶
nonvirtual Destructor
-
template<size_t Dim2>
inline bool operator==(const std::array<T, Dim2> &other) const¶ comparison operator
-
template<typename T2>
inline DynCcoord<MaxDim, decltype(T{} / T2{})> operator/(const DynCcoord<MaxDim, T2> &other) const¶ element-wise division
-
template<Dim_t Dim>
inline std::array<T, Dim> &get()¶ cast to a reference to a statically sized array
-
template<Dim_t Dim>
inline const std::array<T, Dim> &get() const¶ cast to a const reference to a statically sized array
-
inline const_iterator begin() const¶
const iterator to the first entry for iterating over only the valid entries
-
inline const_iterator end() const¶
const iterator past the dim-th entry for iterating over only the valid entries
Protected Attributes
Private Functions
-
inline DynCcoord()¶
-
template<typename T, class EigenPlain>
struct EigenMap¶ - #include <field_map_static.hh>
Internal struct for handling the matrix-shaped iterates of
muGrid::FieldMap
Public Types
-
using PlainType = EigenPlain¶
Eigen type of the iterate.
-
using value_type = std::conditional_t<MutIter == Mapping::Const, Eigen::Map<const PlainType>, Eigen::Map<PlainType>>¶
stl (const-correct)
-
using ref_type = value_type<MutIter>¶
stl (const-correct)
-
using Return_t = value_type<MutIter>¶
for direct access through operator[]
-
using storage_type = value_type<MutIter>¶
stored type (cannot always be same as ref_type)
Public Static Functions
-
static inline constexpr bool IsValidStaticMapType()¶
check at compile time whether the type is meant to be a map with statically sized iterates.
-
static inline constexpr bool IsScalarMapType()¶
check at compiler time whether this map is scalar
-
template<Mapping MutIter>
static inline constexpr value_type<MutIter> &provide_ref(storage_type<MutIter> &storage)¶ return the return_type version of the iterate from storage_type
-
template<Mapping MutIter>
static inline constexpr const value_type<MutIter> &provide_const_ref(const storage_type<MutIter> &storage)¶ return the const return_type version of the iterate from storage_type
-
template<Mapping MutIter>
static inline constexpr value_type<MutIter> *provide_ptr(storage_type<MutIter> &storage)¶ return a pointer to the iterate from storage_type
-
template<Mapping MutIter>
static inline constexpr Return_t<MutIter> from_data_ptr(std::conditional_t<MutIter == Mapping::Const, const T*, T*> data)¶ return a return_type version of the iterate from its pointer
-
template<Mapping MutIter>
static inline constexpr storage_type<MutIter> to_storage(value_type<MutIter> &&value)¶ return a storage_type version of the iterate from its value
-
static inline constexpr Dim_t stride()¶
return the nb of components of the iterate (known at compile time)
-
static inline std::string shape()¶
return the iterate’s shape as text, mostly for error messages
-
using PlainType = EigenPlain¶
-
class Enumerator¶
- #include <ccoord_operations.hh>
enumerator class for
muSpectre::DynamicPixels
Public Functions
-
Enumerator() = delete¶
Default constructor.
-
explicit Enumerator(const DynamicPixels &pixels)¶
Constructor.
-
Enumerator(const Enumerator &other) = default¶
Copy constructor.
-
Enumerator(Enumerator &&other) = default¶
Move constructor.
-
virtual ~Enumerator() = default¶
Destructor.
-
Enumerator &operator=(const Enumerator &other) = delete¶
Copy assignment operator.
-
Enumerator &operator=(Enumerator &&other) = delete¶
Move assignment operator.
-
size_t size() const¶
stl conformance
Protected Attributes
-
const DynamicPixels &pixels¶
-
Enumerator() = delete¶
-
template<Dim_t dim>
class FFT_freqs¶ - #include <fft_utils.hh>
simple class encapsulating the creation, and retrieval of wave vectors
Public Types
Public Functions
-
FFT_freqs() = delete¶
Default constructor.
-
inline FFT_freqs(Ccoord_t<dim> nb_grid_pts, std::array<Real, dim> lengths)¶
constructor with domain length
-
virtual ~FFT_freqs() = default¶
Destructor.
-
inline Vector get_xi(const Ccoord_t<dim> ccoord) const¶
get unnormalised wave vector (in sampling units)
-
inline VectorComplex get_complex_xi(const Ccoord_t<dim> ccoord) const¶
get unnormalised complex wave vector (in sampling units)
-
inline Dim_t get_nb_grid_pts(Dim_t i) const¶
-
FFT_freqs() = delete¶
-
class FFTEngineBase¶
- #include <fft_engine_base.hh>
Virtual base class for FFT engines. To be implemented by all FFT_engine implementations.
Subclassed by muFFT::FFTWEngine, muFFT::FFTWMPIEngine
Public Types
-
using GFieldCollection_t = muGrid::GlobalFieldCollection¶
global FieldCollection
-
using Pixels = typename GFieldCollection_t::DynamicPixels¶
pixel iterator
-
using Field_t = muGrid::TypedFieldBase<Real>¶
Field type on which to apply the projection. This is a TypedFieldBase because it need to be able to hold either TypedField or a WrappedField.
-
using Workspace_t = muGrid::ComplexField¶
Field type holding a Fourier-space representation of a real-valued second-order tensor field
-
using iterator = typename GFieldCollection_t::DynamicPixels::iterator¶
iterator over Fourier-space discretisation point
Public Functions
-
FFTEngineBase() = delete¶
Default constructor.
-
FFTEngineBase(DynCcoord_t nb_grid_pts, Dim_t nb_dof_per_pixel, Communicator comm = Communicator())¶
Constructor with the domain’s number of grid points in each direciton, the number of components to transform, and the communicator
-
FFTEngineBase(const FFTEngineBase &other) = delete¶
Copy constructor.
-
FFTEngineBase(FFTEngineBase &&other) = delete¶
Move constructor.
-
virtual ~FFTEngineBase() = default¶
Destructor.
-
FFTEngineBase &operator=(const FFTEngineBase &other) = delete¶
Copy assignment operator.
-
FFTEngineBase &operator=(FFTEngineBase &&other) = delete¶
Move assignment operator.
-
virtual void initialise(FFT_PlanFlags)¶
compute the plan, etc
-
virtual Workspace_t &fft(Field_t&) = 0¶
forward transform (dummy for interface)
-
inline virtual bool is_active() const¶
return whether this engine is active
-
const Pixels &get_pixels() const¶
iterators over only those pixels that exist in frequency space (i.e. about half of all pixels, see rfft)
-
size_t size() const¶
nb of pixels (mostly for debugging)
-
size_t fourier_size() const¶
nb of pixels in Fourier space
-
size_t workspace_size() const¶
nb of pixels in the work space (may contain a padding region)
-
inline const Communicator &get_communicator() const¶
return the communicator object
-
inline const DynCcoord_t &get_nb_subdomain_grid_pts() const¶
returns the process-local number of grid points in each direction of the cell
-
inline const DynCcoord_t &get_nb_domain_grid_pts() const¶
returns the process-local number of grid points in each direction of the cell
-
inline const DynCcoord_t &get_subdomain_locations() const¶
returns the process-local locations of the cell
-
inline const DynCcoord_t &get_nb_fourier_grid_pts() const¶
returns the process-local number of grid points in each direction of the cell in Fourier space
-
inline const DynCcoord_t &get_fourier_locations() const¶
returns the process-local locations of the cell in Fourier space
-
inline GFieldCollection_t &get_field_collection()¶
only required for testing and debugging
-
inline Workspace_t &get_work_space()¶
only required for testing and debugging
-
inline Real normalisation() const¶
factor by which to multiply projection before inverse transform (this is typically 1/nb_pixels for so-called unnormalized transforms (see, e.g. http://www.fftw.org/fftw3_doc/Multi_002dDimensional-DFTs-of-Real-Data.html#Multi_002dDimensional-DFTs-of-Real-Data or https://docs.scipy.org/doc/numpy-1.13.0/reference/routines.fft.html . Rather than scaling the inverse transform (which would cost one more loop), FFT engines provide this value so it can be used in the projection operator (where no additional loop is required)
-
const Dim_t &get_nb_dof_per_pixel() const¶
return the number of components per pixel
-
const Dim_t &get_dim() const¶
return the number of spatial dimensions
-
const Dim_t &get_nb_quad() const¶
returns the number of quadrature points
-
inline bool is_initialised() const¶
has this engine been initialised?
Protected Attributes
-
Dim_t spatial_dimension¶
spatial dimension of the grid
-
Communicator comm¶
Field collection in which to store fields associated with Fourier-space pointscommunicator
-
GFieldCollection_t work_space_container¶
Field collection to store the fft workspace.
-
DynCcoord_t nb_subdomain_grid_pts¶
nb_grid_pts of the process-local (subdomain) portion of the cell
-
DynCcoord_t subdomain_locations¶
location of the process-local (subdomain) portion of the cell
-
DynCcoord_t nb_fourier_grid_pts¶
nb_grid_pts of the process-local (subdomain) portion of the Fourier transformed data
-
DynCcoord_t fourier_locations¶
location of the process-local (subdomain) portion of the Fourier transformed data
-
const DynCcoord_t nb_domain_grid_pts¶
nb_grid_pts of the full domain of the cell
-
Workspace_t &work¶
field to store the Fourier transform of P
-
const Real norm_factor¶
normalisation coefficient of fourier transform
-
Dim_t nb_dof_per_pixel¶
number of degrees of freedom per pixel. Corresponds to the number of quadrature points per pixel multiplied by the number of components per quadrature point
-
bool initialised = {false}¶
to prevent double initialisation
-
using GFieldCollection_t = muGrid::GlobalFieldCollection¶
-
class FFTWEngine : public muFFT::FFTEngineBase¶
- #include <fftw_engine.hh>
implements the
muFFT::FftEngine_Base
interface using the FFTW libraryPublic Types
-
using Parent = FFTEngineBase¶
base class
Public Functions
-
FFTWEngine() = delete¶
Default constructor.
-
FFTWEngine(const DynCcoord_t &nb_grid_pts, Dim_t nb_dof_per_pixel, Communicator comm = Communicator())¶
Constructor with the domain’s number of grid points in each direciton, the number of components to transform, and the communicator
-
FFTWEngine(const FFTWEngine &other) = delete¶
Copy constructor.
-
FFTWEngine(FFTWEngine &&other) = delete¶
Move constructor.
-
virtual ~FFTWEngine() noexcept¶
Destructor.
-
FFTWEngine &operator=(const FFTWEngine &other) = delete¶
Copy assignment operator.
-
FFTWEngine &operator=(FFTWEngine &&other) = delete¶
Move assignment operator.
-
virtual void initialise(FFT_PlanFlags plan_flags) override¶
compute the plan, etc
-
virtual Workspace_t &fft(Field_t &field) override¶
forward transform
-
using Parent = FFTEngineBase¶
-
class FFTWMPIEngine : public muFFT::FFTEngineBase¶
- #include <fftwmpi_engine.hh>
implements the
muFFT::FFTEngineBase
interface using the FFTW libraryPublic Types
-
using Parent = FFTEngineBase¶
base class
Public Functions
-
FFTWMPIEngine() = delete¶
Default constructor.
-
FFTWMPIEngine(DynCcoord_t nb_grid_pts, Dim_t nb_dof_per_pixel, Communicator comm = Communicator())¶
Constructor with the domain’s number of grid points in each direciton, the number of components to transform, and the communicator
-
FFTWMPIEngine(const FFTWMPIEngine &other) = delete¶
Copy constructor.
-
FFTWMPIEngine(FFTWMPIEngine &&other) = delete¶
Move constructor.
-
virtual ~FFTWMPIEngine() noexcept¶
Destructor.
-
FFTWMPIEngine &operator=(const FFTWMPIEngine &other) = delete¶
Copy assignment operator.
-
FFTWMPIEngine &operator=(FFTWMPIEngine &&other) = delete¶
Move assignment operator.
-
virtual void initialise(FFT_PlanFlags plan_flags) override¶
compute the plan, etc
-
virtual Workspace_t &fft(Field_t &field) override¶
forward transform
-
inline virtual bool is_active() const override¶
return whether this engine is active
Protected Attributes
-
fftw_plan plan_fft = {}¶
holds the plan for forward fourier transform
-
fftw_plan plan_ifft = {}¶
holds the plan for inverse fourier transform
-
ptrdiff_t workspace_size = {}¶
size of workspace buffer returned by planner
-
Real *real_workspace = {}¶
temporary real workspace that is correctly padded
-
bool active = {true}¶
FFTWMPI sometimes assigns zero grid points.
Protected Static Attributes
-
static int nb_engines = {0}¶
number of times this engine has been instatiated
-
using Parent = FFTEngineBase¶
-
class Field¶
- #include <field.hh>
Abstract base class for all fields. A field provides storage discretising a mathematical (scalar, vectorial, tensorial) (real-valued, integer-valued, complex-valued) field on a fixed number of quadrature points per pixel/voxel of a regular grid. Fields defined on the same domains are grouped within
muGrid::FieldCollection
s.Subclassed by muGrid::TypedFieldBase< T >
Public Functions
-
Field() = delete¶
Default constructor.
-
virtual ~Field() = default¶
Destructor.
-
const std::string &get_name() const¶
return the field’s unique name
-
FieldCollection &get_collection() const¶
return a const reference to the field’s collection
-
std::vector<Dim_t> get_shape(Iteration iter_type) const¶
evaluate and return the overall shape of the field (for passing the field to generic multidimensional array objects such as numpy.ndarray)
-
std::vector<Dim_t> get_pixels_shape() const¶
evaluate and return the overall shape of the pixels portion of the field (for passing the field to generic multidimensional array objects such as numpy.ndarray)
-
virtual std::vector<Dim_t> get_components_shape(Iteration iter_type) const¶
evaluate and return the shape of the data contained in a single pixel or quadrature point (for passing the field to generic multidimensional array objects such as numpy.ndarray)
-
Dim_t get_stride(Iteration iter_type) const¶
evaluate and return the number of components in an iterate when iterating over this field
-
virtual const std::type_info &get_stored_typeid() const = 0¶
return the type information of the stored scalar (for compatibility checking)
-
size_t size() const¶
number of entries in the field (= nb_pixel × nb_quad)
-
virtual size_t buffer_size() const = 0¶
size of the internal buffer including the pad region (in scalars)
-
virtual void set_pad_size(size_t pad_size_) = 0¶
add a pad region to the end of the field buffer; required for using this as e.g. an FFT workspace
-
const size_t &get_pad_size() const¶
pad region size
-
virtual void set_zero() = 0¶
initialise field to zero (do more complicated initialisations through fully typed maps)
-
bool is_global() const¶
checks whether this field is registered in a global FieldCollection
Protected Functions
-
Field(const std::string &unique_name, FieldCollection &collection, Dim_t nb_components)¶
Field
s are supposed to only exist in the form ofstd::unique_ptr
s held by a FieldCollection. TheField
constructor is protected to ensure this.- Parameters
unique_name – unique field name (unique within a collection)
nb_components – number of components to store per quadrature point
collection – reference to the holding field collection.
-
virtual void resize(size_t size) = 0¶
resizes the field to the given size
Protected Attributes
- friend FieldCollection
gives field collections the ability to resize() fields
-
size_t current_size = {}¶
maintains a tally of the current size, as it cannot be reliably determined from either
values
oralt_values
alone.
-
const std::string name¶
the field’s unique name
-
FieldCollection &collection¶
reference to the collection this field belongs to
-
const Dim_t nb_components¶
number of components stored per quadrature point (e.g., 3 for a three-dimensional vector, or 9 for a three-dimensional second-rank tensor)
-
size_t pad_size = {}¶
size of padding region at end of buffer
-
Field() = delete¶
-
class FieldCollection¶
- #include <field_collection.hh>
Base class for both
muGrid::GlobalFieldCollection
andmuGrid::LocalFieldCollection
. Manages the a group of fields with the same domain of validitiy (i.e., global fields, or local fields defined on the same pixels).Subclassed by muGrid::GlobalFieldCollection, muGrid::LocalFieldCollection
Public Types
-
enum ValidityDomain¶
domain of validity of the managed fields
Values:
-
enumerator Global¶
-
enumerator Local¶
-
enumerator Global¶
-
using Field_ptr = std::unique_ptr<Field, FieldDestructor<Field>>¶
unique_ptr for holding fields
-
using StateField_ptr = std::unique_ptr<StateField, FieldDestructor<StateField>>¶
unique_ptr for holding state fields
-
using QuadPtIndexIterable = IndexIterable¶
convenience alias
Public Functions
-
FieldCollection() = delete¶
Default constructor.
-
FieldCollection(const FieldCollection &other) = delete¶
Copy constructor.
-
FieldCollection(FieldCollection &&other) = default¶
Move constructor.
-
virtual ~FieldCollection() = default¶
Destructor.
-
FieldCollection &operator=(const FieldCollection &other) = delete¶
Copy assignment operator.
-
FieldCollection &operator=(FieldCollection &&other) = default¶
Move assignment operator.
-
template<typename T>
inline TypedField<T> ®ister_field(const std::string &unique_name, const Dim_t &nb_components)¶ place a new field in the responsibility of this collection (Note, because fields have protected constructors, users can’t create them
Technically, these explicit instantiations are not necessary, as they are implicitly instantiated when the register_<T>field(…) member functions are compiled.
- Parameters
unique_name – unique identifier for this field
nb_components – number of components to be stored per quadrature point (e.g., 4 for a two-dimensional second-rank tensor, or 1 for a scalar field)
-
TypedField<Real> ®ister_real_field(const std::string &unique_name, const Dim_t &nb_components)¶
place a new real-valued field in the responsibility of this collection (Note, because fields have protected constructors, users can’t create them
- Parameters
unique_name – unique identifier for this field
nb_components – number of components to be stored per quadrature point (e.g., 4 for a two-dimensional second-rank tensor, or 1 for a scalar field)
-
TypedField<Complex> ®ister_complex_field(const std::string &unique_name, const Dim_t &nb_components)¶
place a new complex-valued field in the responsibility of this collection (Note, because fields have protected constructors, users can’t create them
- Parameters
unique_name – unique identifier for this field
nb_components – number of components to be stored per quadrature point (e.g., 4 for a two-dimensional second-rank tensor, or 1 for a scalar field)
-
TypedField<Int> ®ister_int_field(const std::string &unique_name, const Dim_t &nb_components)¶
place a new integer-valued field in the responsibility of this collection (Note, because fields have protected constructors, users can’t create them
- Parameters
unique_name – unique identifier for this field
nb_components – number of components to be stored per quadrature point (e.g., 4 for a two-dimensional second-rank tensor, or 1 for a scalar field)
-
TypedField<Uint> ®ister_uint_field(const std::string &unique_name, const Dim_t &nb_components)¶
place a new unsigned integer-valued field in the responsibility of this collection (Note, because fields have protected constructors, users can’t create them
- Parameters
unique_name – unique identifier for this field
nb_components – number of components to be stored per quadrature point (e.g., 4 for a two-dimensional second-rank tensor, or 1 for a scalar field)
-
template<typename T>
inline TypedStateField<T> ®ister_state_field(const std::string &unique_prefix, const Dim_t &nb_memory, const Dim_t &nb_components)¶ place a new state field in the responsibility of this collection (Note, because state fields have protected constructors, users can’t create them
-
TypedStateField<Real> ®ister_real_state_field(const std::string &unique_prefix, const Dim_t &nb_memory, const Dim_t &nb_components)¶
place a new real-valued state field in the responsibility of this collection (Note, because state fields have protected constructors, users can’t create them
- Parameters
unique_prefix – unique idendifier for this state field
nb_memory – number of previous values of this field to store
nb_components – number of scalar components to store per quadrature point
-
TypedStateField<Complex> ®ister_complex_state_field(const std::string &unique_prefix, const Dim_t &nb_memory, const Dim_t &nb_components)¶
place a new complex-valued state field in the responsibility of this collection (Note, because state fields have protected constructors, users can’t create them
- Parameters
unique_prefix – unique idendifier for this state field
nb_memory – number of previous values of this field to store
nb_components – number of scalar components to store per quadrature point
-
TypedStateField<Int> ®ister_int_state_field(const std::string &unique_prefix, const Dim_t &nb_memory, const Dim_t &nb_components)¶
place a new integer-valued state field in the responsibility of this collection (Note, because state fields have protected constructors, users can’t create them
- Parameters
unique_prefix – unique idendifier for this state field
nb_memory – number of previous values of this field to store
nb_components – number of scalar components to store per quadrature point
-
TypedStateField<Uint> ®ister_uint_state_field(const std::string &unique_prefix, const Dim_t &nb_memory, const Dim_t &nb_components)¶
place a new unsigned integer-valued state field in the responsibility of this collection (Note, because state fields have protected constructors, users can’t create them
- Parameters
unique_prefix – unique idendifier for this state field
nb_memory – number of previous values of this field to store
nb_components – number of scalar components to store per quadrature point
-
bool field_exists(const std::string &unique_name) const¶
check whether a field of name ‘unique_name’ has already been registered
-
bool state_field_exists(const std::string &unique_prefix) const¶
check whether a field of name ‘unique_name’ has already been registered
-
const Dim_t &get_nb_entries() const¶
returns the number of entries held by any given field in this collection. This corresponds to nb_pixels × nb_quad_pts, (I.e., a scalar field field and a vector field sharing the the same collection have the same number of entries, even though the vector field has more scalar values.)
-
size_t get_nb_pixels() const¶
returns the number of pixels present in the collection
-
bool has_nb_quad() const¶
check whether the number of quadrature points per pixel/voxel has ben set
-
void set_nb_quad(Dim_t nb_quad_pts_per_pixel)¶
set the number of quadrature points per pixel/voxel. Can only be done once.
-
const Dim_t &get_spatial_dim() const¶
return the spatial dimension of the underlying discretisation grid
-
const ValidityDomain &get_domain() const¶
return the domain of validity (i.e., wher the fields are defined globally (
muGrid::FieldCollection::ValidityDomain::Global
) or locally (muGrid::FieldCollection::ValidityDomain::Local
)
-
bool is_initialised() const¶
whether the collection has been properly initialised (i.e., it knows the number of quadrature points and all its pixels/voxels
-
PixelIndexIterable get_pixel_indices_fast() const¶
return an iterable proxy to the collection which allows to efficiently iterate over the indices fo the collection’s pixels
-
IndexIterable get_pixel_indices() const¶
return an iterable proxy to the collection which allows to iterate over the indices fo the collection’s pixels
-
IndexIterable get_quad_pt_indices() const¶
return an iterable proxy to the collection which allows to iterate over the indices fo the collection’s quadrature points
-
inline std::vector<size_t> get_pixel_ids()¶
-
Field &get_field(const std::string &unique_name)¶
returns a (base-type) reference to the field identified by
unique_name
. Throws amuGrid::FieldCollectionError
if the field does not exist.
-
StateField &get_state_field(const std::string &unique_prefix)¶
returns a (base-type) reference to the state field identified by
unique_prefix
. Throws amuGrid::FieldCollectionError
if the state field does not exist.
-
std::vector<std::string> list_fields() const¶
returns a vector of all field names
preregister a map for latent initialisation
Protected Functions
-
FieldCollection(ValidityDomain domain, const Dim_t &spatial_dimension, const Dim_t &nb_quad_pts)¶
Constructor (not called by user, who constructs either a LocalFieldCollection or a GlobalFieldCollection
- Parameters
domain – Domain of validity, can be global or local
spatial_dimension – spatial dimension of the field (can be muGrid::Unknown, e.g., in the case of the local fields for storing internal material variables)
nb_quad_pts – number of quadrature points per pixel/voxel
-
template<typename T>
TypedField<T> ®ister_field_helper(const std::string &unique_name, const Dim_t &nb_components)¶ internal worker function called by register_<T>_field
-
template<typename T>
TypedStateField<T> ®ister_state_field_helper(const std::string &unique_prefix, const Dim_t &nb_memory, const Dim_t &nb_components)¶ internal worker function called by register_<T>_state_field
-
void allocate_fields()¶
loop through all fields and allocate their memory. Is exclusively called by the daughter classes’
initialise
member function.
-
void initialise_maps()¶
initialise all preregistered maps
Protected Attributes
-
std::map<std::string, StateField_ptr> state_fields = {}¶
storage container for state fields
-
std::vector<std::weak_ptr<std::function<void()>>> init_callbacks = {}¶
Maps registered before initialisation which will need their data_ptr set.
-
ValidityDomain domain¶
domain of validity
-
bool initialised = {false}¶
keeps track of whether the collection has already been initialised
-
std::vector<size_t> pixel_indices = {}¶
Storage for indices of the stored quadrature points in the global field collection. Note that these are not truly global indices, but rather absolute indices within the domain of the local processor. I.e., they are universally valid to address any quadrature point on the local processor, and not for any quadrature point located on anothe processor.
-
enum ValidityDomain¶
-
class FieldCollectionError : public runtime_error¶
- #include <field_collection.hh>
base class for field collection-related exceptions
-
template<class DefaultDestroyable>
struct FieldDestructor¶ - #include <field_collection.hh>
forward declacation of the field’s destructor-functor
Public Functions
-
void operator()(DefaultDestroyable *field)¶
deletes the held field
-
void operator()(DefaultDestroyable *field)¶
-
class FieldError : public runtime_error¶
- #include <field.hh>
base class for field-related exceptions
-
template<typename T, Mapping Mutability>
class FieldMap¶ - #include <field_map.hh>
forward declaration
Dynamically sized field map. Field maps allow iterating over the pixels or quadrature points of a field and to select the shape (in a matrix sense) of the iterate. For example, it allows to iterate in 2×2 matrices over the quadrature points of a strain field for a two-dimensional problem.
Subclassed by muGrid::StaticFieldMap< T, Mutability, MapType, IterationType >
Public Types
-
using Field_t = std::conditional_t<Mutability == Mapping::Const, const TypedFieldBase<T>, TypedFieldBase<T>>¶
const-correct field depending on mapping mutability
-
using Return_t = std::conditional_t<MutVal == Mapping::Const, Eigen::Map<const PlainType>, Eigen::Map<PlainType>>¶
return type for iterators over this- map
-
using EigenRef = Eigen::Ref<const PlainType>¶
Input type for matrix-like values (used for setting uniform values)
-
using PixelEnumeration_t = akantu::containers::ZipContainer<FieldCollection::PixelIndexIterable, FieldMap&>¶
zip-container for iterating over pixel index and stored value simultaneously
-
using Enumeration_t = akantu::containers::ZipContainer<FieldCollection::IndexIterable, FieldMap&>¶
zip-container for iterating over pixel or quadrature point index and stored value simultaneously
Public Functions
-
FieldMap() = delete¶
Default constructor.
-
explicit FieldMap(Field_t &field, Iteration iter_type = Iteration::QuadPt)¶
Constructor from a field. The default case is a map iterating over quadrature points with a matrix of shape (nb_components × 1) per field entry
-
FieldMap(Field_t &field, Dim_t nb_rows, Iteration iter_type = Iteration::QuadPt)¶
Constructor from a field with explicitly chosen shape of iterate. (the number of columns is inferred).
-
virtual ~FieldMap() = default¶
Destructor.
-
FieldMap &operator=(const FieldMap &other) = delete¶
Copy assignment operator (delete because of reference member)
-
FieldMap &operator=(FieldMap &&other) = delete¶
Move assignment operator (delete because of reference member)
-
template<bool IsMutableField = Mutability == Mapping::Mut>
inline std::enable_if_t<IsMutableField, FieldMap> &operator=(const EigenRef &val)¶ Assign a matrix-like value to every entry.
-
template<bool IsMutableField = Mutability == Mapping::Mut>
inline std::enable_if_t<IsMutableField, FieldMap> &operator=(const Scalar &val)¶ Assign a scalar value to every entry.
-
const_iterator cbegin()¶
stl
-
const_iterator cend()¶
stl
-
const_iterator begin() const¶
stl
-
const_iterator end() const¶
stl
-
size_t size() const¶
returns the number of iterates produced by this map (corresponds to the number of field entries if Iteration::Quadpt, or the number of pixels/voxels if Iteration::Pixel);
-
inline Return_t<Mutability> operator[](size_t index)¶
random acces operator
-
void set_data_ptr()¶
query the size from the field’s collection and set data_ptr
-
PixelEnumeration_t enumerate_pixel_indices_fast()¶
return an iterable proxy over pixel indices and stored values simultaneously. Throws a
muGrid::FieldMapError
if the iteration type is over quadrature points
-
Enumeration_t enumerate_indices()¶
return an iterable proxy over pixel/quadrature indices and stored values simultaneously
Public Static Functions
-
static inline constexpr Mapping FieldMutability()¶
determine whether a field is mutably mapped at compile time
-
static inline constexpr bool IsStatic()¶
determine whether a field map is statically sized at compile time
Protected Attributes
-
T *data_ptr = {nullptr}¶
Pointer to mapped data; is also unknown at construction and set in the map’s begin function
-
bool is_initialised = {false}¶
keeps track of whether the map has been initialised.
-
std::shared_ptr<std::function<void()>> callback = {nullptr}¶
shared_ptr used for latent initialisation
-
using Field_t = std::conditional_t<Mutability == Mapping::Const, const TypedFieldBase<T>, TypedFieldBase<T>>¶
-
class FieldMapError : public runtime_error¶
- #include <field_map.hh>
base class for field map-related exceptions
-
template<size_t N>
struct Foreach¶ - #include <iterators.hh>
static for loop
-
template<>
struct Foreach<0>¶ - #include <iterators.hh>
static comparison
-
class FourierDerivative : public muFFT::DerivativeBase¶
- #include <derivative.hh>
Representation of a derivative computed by Fourier interpolation
Public Types
-
using Parent = DerivativeBase¶
base class
Public Functions
-
FourierDerivative() = delete¶
Default constructor.
-
explicit FourierDerivative(Dim_t spatial_dimension, Dim_t direction)¶
Constructor with raw FourierDerivative information.
-
FourierDerivative(const FourierDerivative &other) = default¶
Copy constructor.
-
FourierDerivative(FourierDerivative &&other) = default¶
Move constructor.
-
virtual ~FourierDerivative() = default¶
Destructor.
-
FourierDerivative &operator=(const FourierDerivative &other) = delete¶
Copy assignment operator.
-
FourierDerivative &operator=(FourierDerivative &&other) = delete¶
Move assignment operator.
Protected Attributes
-
Dim_t direction¶
spatial direction in which to perform differentiation
-
using Parent = DerivativeBase¶
-
template<typename Rhs, class CellAdaptor>
struct generic_product_impl<CellAdaptor, Rhs, SparseShape, DenseShape, GemvProduct> : public generic_product_impl_base<CellAdaptor, Rhs, generic_product_impl<CellAdaptor, Rhs>>¶ - #include <cell_adaptor.hh>
Implementation of
muSpectre::CellAdaptor
*Eigen::DenseVector
through a specialization ofEigen::internal::generic_product_impl
:Public Types
-
typedef Product<CellAdaptor, Rhs>::Scalar Scalar¶
undocumented
Public Static Functions
-
template<typename Dest>
static inline void scaleAndAddTo(Dest &dst, const CellAdaptor &lhs, const Rhs &rhs, const Scalar &alpha)¶ undocumented
-
typedef Product<CellAdaptor, Rhs>::Scalar Scalar¶
-
class GlobalFieldCollection : public muGrid::FieldCollection¶
- #include <field_collection_global.hh>
muGrid::GlobalFieldCollection
derives frommuGrid::FieldCollection
and stores global fields that live throughout the whole computational domain, i.e. are defined for every pixel/voxel.Public Types
-
using Parent = FieldCollection¶
alias of base class
-
using DynamicPixels = CcoordOps::DynamicPixels¶
pixel iterator
Public Functions
-
GlobalFieldCollection() = delete¶
Default constructor.
-
GlobalFieldCollection(Dim_t spatial_dimension, Dim_t nb_quad_pts)¶
Constructor
- Parameters
spatial_dimension – number of spatial dimensions, must be 1, 2, 3, or Unknown
nb_quad_pts – number of quadrature points per pixel/voxel
-
GlobalFieldCollection(Dim_t spatial_dimension, Dim_t nb_quad_pts, const DynCcoord_t &nb_grid_pts, const DynCcoord_t &locations = {})¶
Constructor with initialization
- Parameters
spatial_dimension – number of spatial dimensions, must be 1, 2, 3, or Unknown
nb_quad_pts – number of quadrature points per pixel/voxel
-
GlobalFieldCollection(const GlobalFieldCollection &other) = delete¶
Copy constructor.
-
GlobalFieldCollection(GlobalFieldCollection &&other) = default¶
Move constructor.
-
virtual ~GlobalFieldCollection() = default¶
Destructor.
-
GlobalFieldCollection &operator=(const GlobalFieldCollection &other) = delete¶
Copy assignment operator.
-
GlobalFieldCollection &operator=(GlobalFieldCollection &&other) = delete¶
Move assignment operator.
-
const DynamicPixels &get_pixels() const¶
Return the pixels class that allows to iterator over pixels.
-
template<size_t Dim>
inline Dim_t get_index(const Ccoord_t<Dim> &ccoord) const¶ Return index for a ccoord.
-
inline DynCcoord_t get_ccoord(const Dim_t &index) const¶
return coordinates of the i-th pixel
-
void initialise(const DynCcoord_t &nb_grid_pts, const DynCcoord_t &locations = {})¶
freeze the problem size and allocate memory for all fields of the collection. Fields added later on will have their memory allocated upon construction.
-
template<size_t Dim>
inline void initialise(const Ccoord_t<Dim> &nb_grid_pts, const Ccoord_t<Dim> &locations = {})¶ freeze the problem size and allocate memory for all fields of the collection. Fields added later on will have their memory allocated upon construction.
-
void initialise(const DynCcoord_t &nb_grid_pts, const DynCcoord_t &locations, const DynCcoord_t &strides)¶
freeze the problem size and allocate memory for all fields of the collection. Fields added later on will have their memory allocated upon construction.
-
template<size_t Dim>
inline void initialise(const Ccoord_t<Dim> &nb_grid_pts, const Ccoord_t<Dim> &locations, const Ccoord_t<Dim> &strides)¶ freeze the problem size and allocate memory for all fields of the collection. Fields added later on will have their memory allocated upon construction.
-
GlobalFieldCollection get_empty_clone() const¶
obtain a new field collection with the same domain and pixels
Protected Attributes
-
DynamicPixels pixels = {}¶
helper to iterate over the grid
-
using Parent = FieldCollection¶
-
template<Dim_t Dim, class Strain_t, class Tangent_t>
struct Hooke¶ - #include <materials_toolbox.hh>
static inline implementation of Hooke’s law
Public Static Functions
-
static inline constexpr Real compute_lambda(const Real &young, const Real &poisson)¶
compute Lamé’s first constant
- Parameters
young – Young’s modulus
poisson – Poisson’s ratio
-
static inline constexpr Real compute_mu(const Real &young, const Real &poisson)¶
compute Lamé’s second constant (i.e., shear modulus)
- Parameters
young – Young’s modulus
poisson – Poisson’s ratio
-
static inline constexpr Real compute_K(const Real &young, const Real &poisson)¶
compute the bulk modulus
- Parameters
young – Young’s modulus
poisson – Poisson’s ratio
-
static inline Eigen::TensorFixedSize<Real, Eigen::Sizes<Dim, Dim, Dim, Dim>> compute_C(const Real &lambda, const Real &mu)¶
compute the stiffness tensor
- Parameters
lambda – Lamé’s first constant
mu – Lamé’s second constant (i.e., shear modulus)
-
static inline T4Mat<Real, Dim> compute_C_T4(const Real &lambda, const Real &mu)¶
compute the stiffness tensor
- Parameters
lambda – Lamé’s first constant
mu – Lamé’s second constant (i.e., shear modulus)
-
template<class s_t>
static inline auto evaluate_stress(const Real &lambda, const Real &mu, s_t &&E) -> decltype(auto)¶ return stress
- Parameters
lambda – First Lamé’s constant
mu – Second Lamé’s constant (i.e. shear modulus)
E – Green-Lagrange or small strain tensor
-
template<class T_t, class s_t>
static inline auto evaluate_stress(const T_t C, s_t &&E) -> decltype(auto)¶ return stress
- Parameters
C – stiffness tensor (Piola-Kirchhoff 2 (or σ) w.r.t to
E
)E – Green-Lagrange or small strain tensor
-
template<class s_t>
static inline auto evaluate_stress(const Real &lambda, const Real &mu, Tangent_t &&C, s_t &&E) -> decltype(auto)¶ return stress and tangent stiffness
- Parameters
lambda – First Lamé’s constant
mu – Second Lamé’s constant (i.e. shear modulus)
E – Green-Lagrange or small strain tensor
C – stiffness tensor (Piola-Kirchhoff 2 (or σ) w.r.t to
E
)
-
static inline constexpr Real compute_lambda(const Real &young, const Real &poisson)¶
-
class IncompletePixels¶
Public Functions
-
IncompletePixels(const IncompletePixels &other) = default¶
copy constructor
-
IncompletePixels(IncompletePixels &other) = default¶
move constructor
-
virtual ~IncompletePixels() = default¶
-
inline size_t size() const¶
stl conformance
-
IncompletePixels(const IncompletePixels &other) = default¶
-
class IndexIterable¶
- #include <field_collection.hh>
Iterate class for iterating over quadrature point indices of a field collection (i.e. the iterate you get when iterating over the result of
muGrid::FieldCollection::get_quad_pt_indices
).Public Functions
-
IndexIterable() = delete¶
Default constructor.
-
IndexIterable(const IndexIterable &other) = delete¶
Copy constructor.
-
IndexIterable(IndexIterable &&other) = default¶
Move constructor.
-
virtual ~IndexIterable() = default¶
Destructor.
-
IndexIterable &operator=(const IndexIterable &other) = delete¶
Copy assignment operator.
-
IndexIterable &operator=(IndexIterable &&other) = delete¶
Move assignment operator.
-
size_t size() const¶
stl
Protected Functions
-
inline Dim_t get_stride() const¶
evaluate and return the stride with with the fast index of the iterators over the indices of this collection rotate
-
IndexIterable(const FieldCollection &collection, const Iteration &iteration_type)¶
Constructor is protected, because no one ever need to construct this except the fieldcollection
Protected Attributes
- friend FieldCollection
allow the field collection to create
muGrid::FieldCollection::IndexIterable
s
-
const FieldCollection &collection¶
reference back to the proxied collection
-
IndexIterable() = delete¶
-
template<class Derived>
struct is_fixed¶ - #include <eigen_tools.hh>
Helper class to check whether an
Eigen::Array
orEigen::Matrix
is statically sized
-
template<class TestClass>
struct is_matrix¶ - #include <eigen_tools.hh>
Structure to determine whether an expression can be evaluated into a
Eigen::Matrix
,Eigen::Array
, etc. and which helps determine compile-time size
-
template<class T>
struct is_reference_wrapper : public false_type¶
-
template<class Derived>
struct is_square¶ - #include <eigen_tools.hh>
Helper class to check whether an
Eigen::Array
orEigen::Matrix
is a static-size and square.
-
template<class T, Dim_t order>
struct is_tensor¶ - #include <tensor_algebra.hh>
Check whether a given expression represents a Tensor specified order.
-
template<class Strains_t, class Stresses_t, SplitCell is_cell_split = SplitCell::no>
class iterable_proxy¶ - #include <iterable_proxy.hh>
this iterator class is a default for simple laws that just take a strain
Public Types
-
using Strain_t = typename internal::StrainsTComputer<Strains_t>::type¶
expected type for strain values
-
using Stress_t = typename internal::StressesTComputer<Stresses_t>::type¶
expected type for stress values
Public Functions
-
iterable_proxy() = delete¶
Default constructor.
-
template<bool DoNeedTgt = std::tuple_size<Stresses_t>::value == 2, bool DoNeedRate = std::tuple_size<Strain_t>::value == 2>
inline iterable_proxy(MaterialBase &mat, const muGrid::RealField &F, std::enable_if_t<DoNeedRate, const muGrid::RealField> &F_rate, muGrid::RealField
-
using Strain_t = typename internal::StrainsTComputer<Strains_t>::type¶