Reference¶
Common¶
-
void
akantu
::
initialize
(const std::string &input_file, int &argc, char **&argv)¶ initialize the static part of akantu and read the global input_file
-
void
akantu
::
initialize
(int &argc, char **&argv)¶ initialize the static part of akantu
-
using
akantu
::
UInt
= unsigned int¶
-
using
akantu
::
Int
= int¶
-
using
akantu
::
Real
= double¶
-
enum
akantu
::
ElementType
¶ type of elements
Values:
-
enumerator
_not_defined
¶
-
enumerator
_cohesive_1d_2
¶
-
enumerator
_cohesive_2d_4
¶
-
enumerator
_cohesive_2d_6
¶
-
enumerator
_cohesive_3d_12
¶
-
enumerator
_cohesive_3d_16
¶
-
enumerator
_cohesive_3d_6
¶
-
enumerator
_cohesive_3d_8
¶
-
enumerator
_point_1
¶
-
enumerator
_segment_2
¶
-
enumerator
_segment_3
¶
-
enumerator
_triangle_3
¶
-
enumerator
_triangle_6
¶
-
enumerator
_quadrangle_4
¶
-
enumerator
_quadrangle_8
¶
-
enumerator
_tetrahedron_4
¶
-
enumerator
_tetrahedron_10
¶
-
enumerator
_pentahedron_6
¶
-
enumerator
_pentahedron_15
¶
-
enumerator
_hexahedron_8
¶
-
enumerator
_hexahedron_20
¶
-
enumerator
_bernoulli_beam_2
¶
-
enumerator
_bernoulli_beam_3
¶
-
enumerator
_discrete_kirchhoff_triangle_18
¶
-
enumerator
_max_element_type
¶
-
enumerator
-
enum
akantu
::
ModelType
¶ Values:
-
enumerator
model
¶
-
enumerator
solid_mechanics_model
¶
-
enumerator
solid_mechanics_model_cohesive
¶
-
enumerator
heat_transfer_model
¶
-
enumerator
structural_mechanics_model
¶
-
enumerator
embedded_model
¶
-
enumerator
-
enum
akantu
::
AnalysisMethod
¶ enum AnalysisMethod type of solving method used to solve the equation of motion
Values:
-
enumerator
_static
= 0¶
-
enumerator
_implicit_dynamic
= 1¶
-
enumerator
_explicit_lumped_mass
= 2¶
-
enumerator
_explicit_lumped_capacity
= 2¶
-
enumerator
_explicit_consistent_mass
= 3¶
-
enumerator
-
enum
akantu
::
SolveConvergenceCriteria
¶ enum SolveConvergenceCriteria different convergence criteria
Values:
-
enumerator
_residual
¶ Use residual to test the convergence.
-
enumerator
_solution
¶ Use solution to test the convergence.
-
enumerator
_residual_mass_wgh
¶ Use residual weighted by inv. nodal mass to testb
-
enumerator
-
class
akantu
::
ArrayBase
¶ class that afford to store vectors in static memory
Subclassed by akantu::ArrayDataLayer< Element >, akantu::ArrayDataLayer< Entity >, akantu::ArrayDataLayer< Int >, akantu::ArrayDataLayer< char >, akantu::ArrayDataLayer< bool >, akantu::ArrayDataLayer< NodeFlag >, akantu::ArrayDataLayer< UInt >, akantu::ArrayDataLayer< Real >, akantu::ArrayDataLayer< akantu::Element >, akantu::ArrayDataLayer< T, allocation_trait >
Public Functions
-
ArrayBase
(const ID &id = "")¶
-
~ArrayBase
() = default¶
-
bool
empty
() const __attribute__((warn_unused_result))¶
-
void
printself
(std::ostream &stream, int indent = 0) const = 0¶ function to print the containt of the class
-
const ID &
getID
() const¶ Get the name of th array.
-
void
setID
(const ID &id)¶ Set the name of th array.
-
-
template<typename
T
, ArrayAllocationTypeallocation_trait
= ArrayAllocationTrait<T>::value>
classakantu
::
ArrayDataLayer
: public akantu::ArrayBase¶ Subclassed by akantu::Array< T >, akantu::Array< T, is_scal >
Public Types
-
using
reference
= value_type&¶
-
using
pointer_type
= value_type*¶
-
using
const_reference
= const value_type&¶
Public Functions
-
~ArrayDataLayer
() override = default¶
-
ArrayDataLayer
(UInt size = 0, UInt nb_component = 1, const ID &id = "")¶ Allocation of a new vector.
-
ArrayDataLayer
(UInt size, UInt nb_component, const_reference value, const ID &id = "")¶ Allocation of a new vector with a default value.
-
ArrayDataLayer
(const ArrayDataLayer &vect, const ID &id = "")¶ Copy constructor (deep copy)
-
ArrayDataLayer
(const std::vector<value_type> &vect)¶ Copy constructor (deep copy)
-
ArrayDataLayer &
operator=
(const ArrayDataLayer &other)¶
-
ArrayDataLayer
(ArrayDataLayer &&other) noexcept = default¶
-
ArrayDataLayer &
operator=
(ArrayDataLayer &&other) noexcept = default¶
-
void
push_back
(const_reference value)¶ append a tuple of size nb_component containing value
-
template<template<typename> class
C
, typename = std::enable_if_t<aka::is_tensor<C<T>>::value or aka::is_tensor_proxy<C<T>>::value>>
voidpush_back
(const C<T> &new_elem)¶ append a vector
-
using
-
template<typename
T
, boolis_scal
>
classakantu
::
Array
: public akantu::ArrayDataLayer<T>¶ Public Types
-
using
value_type
= typename parent::value_type¶
-
using
reference
= typename parent::reference¶
-
using
pointer_type
= typename parent::pointer_type¶
-
using
const_reference
= typename parent::const_reference¶
-
using
const_scalar_iterator
= const_iterator<T>¶ const_iterator for Array of nb_component = 1
-
using
vector_iterator
= iterator<Vector<T>>¶ iterator returning Vectors of size n on entries of Array with nb_component = n
-
using
const_vector_iterator
= const_iterator<Vector<T>>¶ const_iterator returning Vectors of n size on entries of Array with nb_component = n
-
using
matrix_iterator
= iterator<Matrix<T>>¶ iterator returning Matrices of size (m, n) on entries of Array with nb_component = m*n
-
using
const_matrix_iterator
= const_iterator<Matrix<T>>¶ const iterator returning Matrices of size (m, n) on entries of Array with nb_component = m*n
-
using
tensor3_iterator
= iterator<Tensor3<T>>¶ iterator returning Tensor3 of size (m, n, k) on entries of Array with nb_component = m*n*k
-
using
const_tensor3_iterator
= const_iterator<Tensor3<T>>¶ const iterator returning Tensor3 of size (m, n, k) on entries of Array with nb_component = m*n*k
Public Functions
-
~Array
() override¶
-
Array
()¶
-
Array
(UInt size, UInt nb_component, const_reference value, const ID &id = "")¶ Allocation of a new vector with a default value.
-
UInt
find
(const_reference elem) const¶ search elem in the vector, return the position of the first occurrence or -1 if not found
-
void
push_back
(const_reference value)¶ append a value to the end of the Array
-
template<template<typename> class
C
, typename = std::enable_if_t<aka::is_tensor<C<T>>::value or aka::is_tensor_proxy<C<T>>::value>>
voidpush_back
(const C<T> &new_elem)¶
-
template<typename
Ret
>
voidpush_back
(const iterator<Ret> &it)¶ append the content of the iterator at the end of the Array
-
template<template<typename> class
C
, typename = std::enable_if_t<aka::is_tensor<C<T>>::value or aka::is_tensor_proxy<C<T>>::value>>
UIntfind
(const C<T> &elem)¶
-
void
set
(T t)¶ set all entries of the array to the value t
- Parameters
t
: value to fill the array with
-
void
zero
()¶ set the array to T{}
-
void
clear
()¶ resize the array to 0
-
template<template<typename> class
C
, typename = std::enable_if_t<aka::is_tensor<C<T>>::value or aka::is_tensor_proxy<C<T>>::value>>
voidset
(const C<T> &vm)¶ set all tuples of the array to a given vector or matrix
-
void
copy
(const Array<T, is_scal> &other, bool no_sanity_check = false)¶ copy another Array in the current Array, the no_sanity_check allows you to force the copy in cases where you know what you do with two non matching Arrays in terms of n
-
void
printself
(std::ostream &stream, int indent = 0) const override¶ function to print the containt of the class
-
const_reference
operator()
(UInt i, UInt j = 0) const¶ return a const reference to the j-th entry of the i-th tuple
-
const_reference
operator[]
(UInt i) const¶ return a const reference to the ith component of the 1D array
-
Array<ElementType> &
operator*=
(const ElementType&)¶
-
Array<ElementType> &
operator-=
(const Array<ElementType>&)¶
-
Array<ElementType> &
operator+=
(const Array<ElementType>&)¶
-
using
-
template<typename
T
, typenameSupportType
>
classakantu
::
ElementTypeMapArray
: public akantu::ElementTypeMap<std::unique_ptr<Array<T>>, SupportType>¶ Public Types
-
using
type_iterator
= typename parent::type_iterator¶
Public Functions
-
void
operator=
(const ElementTypeMapArray&) = delete¶ standard assigment (copy) operator
-
ElementTypeMapArray
(const ElementTypeMapArray&)¶
-
void
copy
(const ElementTypeMapArray &other)¶ explicit copy
-
ElementTypeMapArray
(const ID &id = "by_element_type_array", const ID &parent_id = "no_parent")¶ Constructor
- Parameters
id
: optional: identifier (string)parent_id
: optional: parent identifier. for organizational purposes only
-
Array<T> &
alloc
(UInt size, UInt nb_component, const SupportType &type, GhostType ghost_type, const T &default_value = T())¶ allocate memory for a new array
- Return
a reference to the allocated array
- Parameters
size
: number of tuples of the new arraynb_component
: tuple sizetype
: the type under which the array is indexed in the mapghost_type
: whether to add the field to the data map or the ghost_data mapdefault_value
: the default value to use to fill the array
-
void
alloc
(UInt size, UInt nb_component, const SupportType &type, const T &default_value = T())¶ allocate memory for a new array in both the data and the ghost_data map
- Parameters
size
: number of tuples of the new arraynb_component
: tuple sizetype
: the type under which the array is indexed in the mapdefault_value
: the default value to use to fill the array
-
const Array<T> &
operator()
(const SupportType &type, GhostType ghost_type = _not_ghost) const¶
-
const T &
operator()
(const Element &element, UInt component = 0) const¶ access the data of an element, this combine the map and array accessor
-
T &
operator()
(const Element &element, UInt component = 0)¶ access the data of an element, this combine the map and array accessor
-
decltype(auto)
get
(const Element &element)¶ access the data of an element, this combine the map and array accessor
-
decltype(auto)
get
(const Element &element) const¶
-
Array<T> &
operator()
(const SupportType &type, GhostType ghost_type = _not_ghost)¶
-
void
setArray
(const SupportType &type, GhostType ghost_type, const Array<T> &vect)¶ insert data of a new type (not yet present) into the map.
- Return
stored data corresponding to type.
- Parameters
type
: type of data (if this type is already present in the map, an exception is thrown).ghost_type
: optional: by default, the data map for non-ghost elements is searchedvect
: the vector to include into the map
-
void
free
()¶ frees all memory related to the data
-
void
clear
()¶
-
bool
empty
() const __attribute__((warn_unused_result))¶
-
void
zero
()¶ set all values in the ElementTypeMap to zero
-
void
onElementsRemoved
(const ElementTypeMapArray<UInt> &new_numbering)¶ deletes and reorders entries in the stored arrays
- Parameters
new_numbering
: a ElementTypeMapArray of new indices. UInt(-1) indicates deleted entries.
-
void
printself
(std::ostream &stream, int indent = 0) const override¶ text output helper
-
void
setID
(const ID &id)¶ set the id
- Parameters
id
: the new name
-
auto
getID
() const -> ID¶ return the id
-
ElementTypeMap<UInt>
getNbComponents
(UInt dim = _all_dimensions, GhostType requested_ghost_type = _not_ghost, ElementKind kind = _ek_not_defined) const¶
-
template<class
Func
>
voidinitialize
(const Func &f, const T &default_value, bool do_not_default)¶ initialize the arrays in accordance to a functor
-
template<typename ...
pack
>
voidinitialize
(const Mesh &mesh, pack&&... _pack)¶ initialize with sizes and number of components in accordance of a mesh content
-
template<typename ...
pack
>
voidinitialize
(const FEEngine &fe_engine, pack&&... _pack)¶ initialize with sizes and number of components in accordance of a fe engine content (aka integration points)
-
ID
getName
() const¶ get the name of the internal field
-
template<typename ...
pack
>
UIntsize
(pack&&... _pack) const¶ get the size of the ElementTypeMapArray<T>
- Parameters
[in] _pack
:optional arguments can be any of:
_spatial_dimension
the dimension to consider (default: _all_dimensions)_ghost_type
(default: _not_ghost)_element_kind
(default: _ek_not_defined)_all_ghost_types
(default: false)
-
bool
isNodal
() const¶
-
void
isNodal
(bool is_nodal)¶
-
using
-
template<typename
T
>
classakantu
::
Vector
: public akantu::TensorStorage<T, 1, Vector<T>>¶ Public Types
-
using
value_type
= typename parent::value_type¶
Public Functions
-
Vector
()¶
-
const_iterator
begin
() const¶
-
const_iterator
end
() const¶
-
~Vector
() override = default¶
-
template<typename
Func
, typenameAcc
>
decltype(auto)accumulate
(const Vector<T> &v, Acc &&accumulator, Func &&func) const¶
-
void
printself
(std::ostream &stream, int indent = 0) const¶ function to print the containt of the class
-
using
Mesh¶
-
class
akantu
::
Mesh
: public akantu::EventHandlerManager<MeshEventHandler>, public akantu::GroupManager, public akantu::MeshData, public akantu::Dumpable¶ This class contaisn the coordinates of the nodes in the Mesh.nodes akant::Array, and the connectivity. The connectivity are stored in by element types.
In order to loop on all element you have to loop on all types like this :
for(auto & type : mesh.elementTypes()) { UInt nb_element = mesh.getNbElement(type); const Array<UInt> & conn = mesh.getConnectivity(type); for(UInt e = 0; e < nb_element; ++e) { ... } } or for_each_element(mesh, [](Element & element) { std::cout << element << std::endl });
Public Types
-
using
type_iterator
= ElementTypeMapArray<UInt, ElementType>::type_iterator¶
-
using
ElementTypesIteratorHelper
= ElementTypeMapArray<UInt, ElementType>::ElementTypesIteratorHelper¶
Public Functions
-
Mesh
(UInt spatial_dimension, const ID &id = "mesh")¶ constructor that create nodes coordinates array
-
Mesh
(UInt spatial_dimension, Communicator &communicator, const ID &id = "mesh")¶ mesh not distributed and not using the default communicator
constructor that use an existing nodes coordinates array, by getting the vector of coordinates
-
~Mesh
() override¶
-
void
read
(const std::string &filename, const MeshIOType &mesh_io_type = _miot_auto)¶ read the mesh from a file
-
void
write
(const std::string &filename, const MeshIOType &mesh_io_type = _miot_auto)¶ write the mesh to a file
-
template<typename ...
pack
>
std::enable_if_t<are_named_argument<pack...>::value>distribute
(pack&&... _pack)¶ with the arguments to pass to the partitionner
-
bool
isDistributed
() const¶ defines is the mesh is distributed or not
-
void
makePeriodic
(const SpatialDirection &direction)¶ set the periodicity in a given direction
-
void
makePeriodic
(const SpatialDirection &direction, const ID &list_1, const ID &list_2)¶
-
bool
isPeriodic
() const¶ defines if the mesh is periodic or not
-
bool
isPeriodic
(const SpatialDirection&) const¶
-
UInt
getPeriodicMaster
(UInt slave) const¶ get the master node for a given slave nodes, except if node not a slave
-
decltype(auto)
getPeriodicSlaves
(UInt master) const¶ get an iterable list of slaves for a given master node
-
void
printself
(std::ostream &stream, int indent = 0) const override¶ function to print the containt of the class
-
template<typename
T
>
voidextractNodalValuesFromElement
(const Array<T> &nodal_values, T *local_coord, const UInt *connectivity, UInt n_nodes, UInt nb_degree_of_freedom) const¶ extract coordinates of nodes from an element
-
void
addConnectivityType
(ElementType type, GhostType ghost_type = _not_ghost)¶ add a Array of connectivity for the given ElementType and GhostType .
-
void
initNormals
()¶ initialize normals
-
void
getGlobalConnectivity
(ElementTypeMapArray<UInt> &global_connectivity)¶ get global connectivity array
-
const ID &
getID
() const¶ get the id of the mesh
-
UInt
getSpatialDimension
() const¶ get the spatial dimension of the mesh = number of component of the coordinates
-
Array<Real> &
getNormals
(const ElementType &el_type, GhostType ghost_type = _not_ghost)¶ get the normals for the elements
-
const Array<UInt> &
getGlobalNodesIds
() const¶ get the Array of global ids of the nodes (only used in parallel)
-
const BBox &
getBBox
() const¶
-
const BBox &
getLocalBBox
() const¶
-
const Array<UInt> &
getConnectivity
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶ get the connectivity Array for a given type
-
Array<UInt> &
getConnectivity
(const ElementType &el_type, GhostType ghost_type = _not_ghost)¶
-
const ElementTypeMapArray<UInt> &
getConnectivities
() const¶
-
UInt
getNbElement
(ElementType type, GhostType ghost_type = _not_ghost) const¶ get the number of element of a type in the mesh
-
UInt
getNbElement
(UInt spatial_dimension = _all_dimensions, GhostType ghost_type = _not_ghost, ElementKind kind = _ek_not_defined) const¶ get the number of element for a given ghost_type and a given dimension
-
void
getBarycenter
(const Element &element, Vector<Real> &barycenter) const¶ compute the barycenter of a given element
-
void
getBarycenters
(Array<Real> &barycenter, ElementType type, GhostType ghost_type) const¶
-
const auto &
getElementToSubelement
() const¶ get the element connected to a subelement (element of lower dimension)
-
const auto &
getElementToSubelement
(ElementType el_type, GhostType ghost_type = _not_ghost) const¶ get the element connected to a subelement
-
const auto &
getElementToSubelement
(const Element &element) const¶ get the elements connected to a subelement
-
const auto &
getSubelementToElement
() const¶ get the subelement (element of lower dimension) connected to a element
-
const auto &
getSubelementToElement
(ElementType el_type, GhostType ghost_type = _not_ghost) const¶ get the subelement connected to an element
-
VectorProxy<Element>
getSubelementToElement
(const Element &element) const¶ get the subelement (element of lower dimension) connected to a element
-
VectorProxy<UInt>
getConnectivity
(const Element &element) const¶ get connectivity of a given element
-
template<typename
T
>
const Array<T> &getData
(const ID &data_name, ElementType el_type, GhostType ghost_type = _not_ghost) const¶ get a name field associated to the mesh
-
template<typename
T
>
Array<T> &getData
(const ID &data_name, ElementType el_type, GhostType ghost_type = _not_ghost)¶ get a name field associated to the mesh
-
template<typename
T
>
const ElementTypeMapArray<T> &getData
(const ID &data_name) const¶ get a name field associated to the mesh
-
template<typename
T
>
ElementTypeMapArray<T> &getData
(const ID &data_name)¶ get a name field associated to the mesh
-
template<typename
T
>
ElementTypeMap<UInt>getNbDataPerElem
(ElementTypeMapArray<T> &array)¶
-
template<typename
T
>
Array<T> &getDataPointer
(const std::string &data_name, ElementType el_type, GhostType ghost_type = _not_ghost, UInt nb_component = 1, bool size_to_nb_element = true, bool resize_with_parent = false)¶ templated getter returning the pointer to data in MeshData (modifiable)
-
template<typename
T
>
Array<T> &getDataPointer
(const ID &data_name, ElementType el_type, GhostType ghost_type, UInt nb_component, bool size_to_nb_element, bool resize_with_parent, const T &defaul_)¶
-
auto
hasMeshFacets
() const¶
-
bool
isMeshFacets
() const¶
-
DumperIOHelper &
getGroupDumper
(const std::string &dumper_name, const std::string &group_name)¶ return the dumper from a group and and a dumper name
-
auto
getFacetConnectivity
(const Element &element, UInt t = 0) const¶ get connectivity of facets for a given element
-
template<typename ...
pack
>
ElementTypesIteratorHelperelementTypes
(pack&&... _pack) const¶
-
decltype(auto)
firstType
(UInt dim = _all_dimensions, GhostType ghost_type = _not_ghost, ElementKind kind = _ek_regular) const¶
-
decltype(auto)
lastType
(UInt dim = _all_dimensions, GhostType ghost_type = _not_ghost, ElementKind kind = _ek_regular) const¶
-
const ElementSynchronizer &
getElementSynchronizer
() const¶
-
ElementSynchronizer &
getElementSynchronizer
()¶
-
const NodeSynchronizer &
getNodeSynchronizer
() const¶
-
NodeSynchronizer &
getNodeSynchronizer
()¶
-
const PeriodicNodeSynchronizer &
getPeriodicNodeSynchronizer
() const¶
-
PeriodicNodeSynchronizer &
getPeriodicNodeSynchronizer
()¶
-
const auto &
getCommunicator
() const¶
-
auto &
getCommunicator
()¶
-
const auto &
getPeriodicMasterSlaves
() const¶
-
template<>
voidsendEvent
(NewElementsEvent &event)¶
-
template<>
voidsendEvent
(NewNodesEvent &event)¶
-
template<>
voidsendEvent
(RemovedElementsEvent &event)¶
-
template<>
voidsendEvent
(RemovedNodesEvent &event)¶
Public Static Functions
-
UInt
getNbNodesPerElement
(ElementType type)¶ get the number of nodes per element for a given element type
-
ElementType
getP1ElementType
(ElementType type)¶ get the number of nodes per element for a given element type considered as a first order element
-
ElementKind
getKind
(ElementType type)¶ get the kind of the element type
-
UInt
getSpatialDimension
(ElementType type)¶ get spatial dimension of a type of element
-
UInt
getNbFacetsPerElement
(ElementType type)¶ get number of facets of a given element type
-
UInt
getNbFacetsPerElement
(ElementType type, UInt t)¶ get number of facets of a given element type
-
auto
getFacetLocalConnectivity
(ElementType type, UInt t = 0)¶ get local connectivity of a facet for a given facet type
-
UInt
getNbFacetTypes
(ElementType type, UInt t = 0)¶ get the number of type of the surface element associated to a given element type
-
constexpr auto
getFacetType
(ElementType type, UInt t = 0)¶ get the type of the surface element associated to a given element
-
constexpr auto
getAllFacetTypes
(ElementType type)¶ get all the type of the surface element associated to a given element
-
class
PeriodicSlaves
¶ Public Functions
-
PeriodicSlaves
(const PeriodicSlaves &other) = default¶
-
PeriodicSlaves
(PeriodicSlaves &&other) = default¶
-
PeriodicSlaves &
operator=
(const PeriodicSlaves &other) = default¶
-
auto
begin
() const¶
-
auto
end
() const¶
-
class
const_iterator
¶ Public Functions
-
const_iterator
(internal_iterator it)¶
-
const_iterator
operator++
()¶
-
bool
operator!=
(const const_iterator &other)¶
-
auto
operator*
()¶
-
-
-
using
-
class
akantu
::
FEEngine
: public akantu::MeshEventHandler¶ The generic FEEngine class derived in a FEEngineTemplate class containing the shape functions and the integration method
Public Types
-
using
ElementTypesIteratorHelper
= ElementTypeMapArray<Real, ElementType>::ElementTypesIteratorHelper¶
Public Functions
-
~FEEngine
() override¶
-
void
initShapeFunctions
(GhostType ghost_type = _not_ghost) = 0¶ pre-compute all the shape functions, their derivatives and the jacobians
-
void
integrate
(const Array<Real> &f, Array<Real> &intf, UInt nb_degree_of_freedom, ElementType type, GhostType ghost_type = _not_ghost, const Array<UInt> &filter_elements = empty_filter) const = 0¶ integrate f for all elements of type “type”
-
Real
integrate
(const Array<Real> &f, ElementType type, GhostType ghost_type = _not_ghost, const Array<UInt> &filter_elements = empty_filter) const = 0¶ integrate a scalar value f on all elements of type “type”
-
void
integrateOnIntegrationPoints
(const Array<Real> &f, Array<Real> &intf, UInt nb_degree_of_freedom, ElementType type, GhostType ghost_type = _not_ghost, const Array<UInt> &filter_elements = empty_filter) const = 0¶ integrate f for all integration points of type “type” but don’t sum over all integration points
-
Real
integrate
(const Vector<Real> &f, ElementType type, UInt index, GhostType ghost_type = _not_ghost) const = 0¶ integrate one element scalar value on all elements of type “type”
-
UInt
getNbIntegrationPoints
(ElementType type, GhostType ghost_type = _not_ghost) const = 0¶ get the number of integration points
-
const Array<Real> &
getShapes
(ElementType type, GhostType ghost_type = _not_ghost, UInt id = 0) const = 0¶ get the precomputed shapes
-
const Array<Real> &
getShapesDerivatives
(ElementType type, GhostType ghost_type = _not_ghost, UInt id = 0) const = 0¶ get the derivatives of shapes
-
const Matrix<Real> &
getIntegrationPoints
(ElementType type, GhostType ghost_type = _not_ghost) const = 0¶ get integration points
-
void
gradientOnIntegrationPoints
(const Array<Real> &u, Array<Real> &nablauq, UInt nb_degree_of_freedom, ElementType type, GhostType ghost_type = _not_ghost, const Array<UInt> &filter_elements = empty_filter) const = 0¶ Compute the gradient nablauq on the integration points of an element type from nodal values u
-
void
interpolateOnIntegrationPoints
(const Array<Real> &u, Array<Real> &uq, UInt nb_degree_of_freedom, ElementType type, GhostType ghost_type = _not_ghost, const Array<UInt> &filter_elements = empty_filter) const = 0¶ Interpolate a nodal field u at the integration points of an element type -> uq
-
void
interpolateOnIntegrationPoints
(const Array<Real> &u, ElementTypeMapArray<Real> &uq, const ElementTypeMapArray<UInt> *filter_elements = nullptr) const = 0¶ Interpolate a nodal field u at the integration points of many element types -> uq
-
void
computeBtD
(const Array<Real> &Ds, Array<Real> &BtDs, ElementType type, GhostType ghost_type = _not_ghost, const Array<UInt> &filter_elements = empty_filter) const = 0¶ pre multiplies a tensor by the shapes derivaties
-
void
computeBtDB
(const Array<Real> &Ds, Array<Real> &BtDBs, UInt order_d, ElementType type, GhostType ghost_type = _not_ghost, const Array<UInt> &filter_elements = empty_filter) const = 0¶ left and right multiplies a tensor by the shapes derivaties
-
void
computeNtb
(const Array<Real> &bs, Array<Real> &Ntbs, ElementType type, GhostType ghost_type = _not_ghost, const Array<UInt> &filter_elements = empty_filter) const = 0¶ left multiples a vector by the shape functions
-
void
computeIntegrationPointsCoordinates
(ElementTypeMapArray<Real> &integration_points_coordinates, const ElementTypeMapArray<UInt> *filter_elements = nullptr) const = 0¶ Compute the interpolation point position in the global coordinates for many element types
-
void
computeIntegrationPointsCoordinates
(Array<Real> &integration_points_coordinates, ElementType type, GhostType ghost_type = _not_ghost, const Array<UInt> &filter_elements = empty_filter) const = 0¶ Compute the interpolation point position in the global coordinates for an element type
-
void
initElementalFieldInterpolationFromIntegrationPoints
(const ElementTypeMapArray<Real> &interpolation_points_coordinates, ElementTypeMapArray<Real> &interpolation_points_coordinates_matrices, ElementTypeMapArray<Real> &integration_points_coordinates_inv_matrices, const ElementTypeMapArray<UInt> *element_filter) const = 0¶ Build pre-computed matrices for interpolation of field form integration points at other given positions (interpolation_points)
-
void
interpolateElementalFieldFromIntegrationPoints
(const ElementTypeMapArray<Real> &field, const ElementTypeMapArray<Real> &interpolation_points_coordinates, ElementTypeMapArray<Real> &result, GhostType ghost_type, const ElementTypeMapArray<UInt> *element_filter) const = 0¶ interpolate field at given position (interpolation_points) from given values of this field at integration points (field)
-
void
interpolateElementalFieldFromIntegrationPoints
(const ElementTypeMapArray<Real> &field, const ElementTypeMapArray<Real> &interpolation_points_coordinates_matrices, const ElementTypeMapArray<Real> &integration_points_coordinates_inv_matrices, ElementTypeMapArray<Real> &result, GhostType ghost_type, const ElementTypeMapArray<UInt> *element_filter) const = 0¶ Interpolate field at given position from given values of this field at integration points (field) using matrices precomputed with initElementalFieldInterplationFromIntegrationPoints
-
void
interpolate
(const Vector<Real> &real_coords, const Matrix<Real> &nodal_values, Vector<Real> &interpolated, const Element &element) const = 0¶ interpolate on a phyiscal point inside an element
-
void
computeShapes
(const Vector<Real> &real_coords, UInt elem, ElementType type, Vector<Real> &shapes, GhostType ghost_type = _not_ghost) const = 0¶ compute the shape on a provided point
-
void
computeShapeDerivatives
(const Vector<Real> &real_coords, UInt element, ElementType type, Matrix<Real> &shape_derivatives, GhostType ghost_type = _not_ghost) const = 0¶ compute the shape derivatives on a provided point
-
void
assembleFieldLumped
(const std::function<void(Matrix<Real>&, const Element&)> &field_funct, const ID &matrix_id, const ID &dof_id, DOFManager &dof_manager, ElementType type, GhostType ghost_type = _not_ghost, ) const = 0¶ assembles the lumped version of
\[ \int N^t rho N \]
-
void
assembleFieldMatrix
(const std::function<void(Matrix<Real>&, const Element&)> &field_funct, const ID &matrix_id, const ID &dof_id, DOFManager &dof_manager, ElementType type, GhostType ghost_type = _not_ghost, ) const = 0¶ assembles the matrix
\[ \int N^t rho N \]
-
void
computeNormalsOnIntegrationPoints
(GhostType ghost_type = _not_ghost) = 0¶ pre-compute normals on integration points
-
void
computeNormalsOnIntegrationPoints
(const Array<Real>&, GhostType = _not_ghost)¶ pre-compute normals on integration points
-
void
computeNormalsOnIntegrationPoints
(const Array<Real>&, Array<Real>&, ElementType, GhostType = _not_ghost) const¶ pre-compute normals on integration points
-
void
printself
(std::ostream &stream, int indent = 0) const¶ function to print the containt of the class
-
ElementTypesIteratorHelper
elementTypes
(UInt dim = _all_dimensions, GhostType ghost_type = _not_ghost, ElementKind kind = _ek_regular) const¶
-
UInt
getElementDimension
() const¶ get the dimension of the element handeled by this fe_engine object
-
const Array<Real> &
getNormalsOnIntegrationPoints
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶ get the normals on integration points
-
const ShapeFunctions &
getShapeFunctionsInterface
() const = 0¶ get the shape function class (probably useless: see getShapeFunction in fe_engine_template.hh)
-
const Integrator &
getIntegratorInterface
() const = 0¶ get the integrator class (probably useless: see getIntegrator in fe_engine_template.hh)
-
ID
getID
() const¶
Public Static Functions
-
template<typename
T
>
voidextractNodalToElementField
(const Mesh &mesh, const Array<T> &nodal_f, Array<T> &elemental_f, ElementType type, GhostType ghost_type = _not_ghost, const Array<UInt> &filter_elements = empty_filter)¶ extract the nodal values and store them per element
-
template<typename
T
>
voidfilterElementalData
(const Mesh &mesh, const Array<T> &elem_f, Array<T> &filtered_f, ElementType type, GhostType ghost_type = _not_ghost, const Array<UInt> &filter_elements = empty_filter)¶ filter a field
-
Real
getElementInradius
(const Matrix<Real> &coord, ElementType type)¶ get the in-radius of an element
-
ElementType
getCohesiveElementType
(ElementType type_facet)¶ get cohesive element type for a given facet type
-
Vector<ElementType>
getIGFEMElementTypes
(ElementType type)¶ get igfem element type for a given regular type
-
InterpolationType
getInterpolationType
(ElementType el_type)¶ get the interpolation element associated to an element type
-
using
Models¶
Common¶
-
class
akantu::BC::Dirichlet
::
FixedValue
: public akantu::BC::Dirichlet::DirichletFunctor¶
-
class
akantu::BC::Dirichlet
::
FlagOnly
: public akantu::BC::Dirichlet::DirichletFunctor¶
-
class
akantu::BC::Dirichlet
::
IncrementValue
: public akantu::BC::Dirichlet::DirichletFunctor¶
Warning
doxygenclass: Cannot find class “akantu::BC::Neumann::FromStress” in doxygen xml output for project “Akantu” from directory: /builds/akantu/akantu/build/doc/dev-doc/xml
Warning
doxygenclass: Cannot find class “akantu::BC::Neumann::FromTraction” in doxygen xml output for project “Akantu” from directory: /builds/akantu/akantu/build/doc/dev-doc/xml
-
template<class
ModelType
>
classakantu
::
BoundaryCondition
¶ Public Functions
-
BoundaryCondition
()¶
-
void
initBC
(ModelType &model, Array<Real> &primal, Array<Real> &dual)¶ Initialize the boundary conditions.
-
void
initBC
(ModelType &model, Array<Real> &primal, Array<Real> &primal_increment, Array<Real> &dual)¶
-
template<typename
FunctorType
>
voidapplyBC
(const FunctorType &func)¶ Apply the boundary conditions.
-
template<class
FunctorType
>
voidapplyBC
(const FunctorType &func, const std::string &group_name)¶
-
template<class
FunctorType
>
voidapplyBC
(const FunctorType &func, const ElementGroup &element_group)¶
-
template<class
FunctorType
, BC::Functor::Typetype
= FunctorType::type>
structTemplateFunctionWrapper
¶
-
template<typename FunctorType> _dirichlet >
Public Static Functions
-
void
applyBC
(const FunctorType &func, const ElementGroup &group, BoundaryCondition<ModelType> &bc_instance)¶
-
void
-
template<typename FunctorType> _neumann >
Public Static Functions
-
void
applyBC
(const FunctorType &func, const ElementGroup &group, BoundaryCondition<ModelType> &bc_instance)
-
void
applyBC
(const FunctorType &func, const ElementGroup &group, BoundaryCondition<ModelType> &bc_instance, GhostType ghost_type)¶
-
void
-
Warning
doxygenclass: Cannot find class “akantu::BoundaryConditionFunctor” in doxygen xml output for project “Akantu” from directory: /builds/akantu/akantu/build/doc/dev-doc/xml
-
template<class
EventHandler
>
classakantu
::
EventHandlerManager
¶ Public Functions
-
~EventHandlerManager
() = default¶
-
void
registerEventHandler
(EventHandler &event_handler, EventHandlerPriority priority = _ehp_highest)¶ register a new EventHandler to the Manager. The register object will then be informed about the events the manager observes.
-
void
unregisterEventHandler
(EventHandler &event_handler)¶ unregister a EventHandler object. This object will not be notified anymore about the events this manager observes.
-
-
class
akantu
::
Model
: public akantu::ModelSolver, public akantu::MeshEventHandler¶ Subclassed by akantu::HeatTransferModel, akantu::SolidMechanicsModel, akantu::StructuralMechanicsModel
Public Functions
-
~Model
() override¶
-
template<typename ...
pack
>
std::enable_if_t<are_named_argument<pack...>::value>initFull
(pack&&... _pack)¶
-
template<typename ...
pack
>
std::enable_if_t<not are_named_argument<pack...>::value>initFull
(pack&&... _pack)¶
-
void
initNewSolver
(const AnalysisMethod &method)¶ initialize a new solver if needed
-
void
dumpGroup
(const std::string &group_name)¶ Dump the data for a given group.
-
void
dumpGroup
(const std::string &group_name, const std::string &dumper_name)¶
-
void
dumpGroup
()¶ Dump the data for all boundaries.
-
void
setGroupDirectory
(const std::string &directory, const std::string &group_name)¶ Set the directory for a given group.
-
void
setGroupDirectory
(const std::string &directory)¶ Set the directory for all boundaries.
-
void
setGroupBaseName
(const std::string &basename, const std::string &group_name)¶ Set the base name for a given group.
-
DumperIOHelper &
getGroupDumper
(const std::string &group_name)¶ Get the internal dumper of a given group.
-
void
updateDataForNonLocalCriterion
(__attribute__((unused)) ElementTypeMapReal &criterion)¶
-
const ID &
getID
() const¶ get id of model
-
void
synchronizeBoundaries
()¶ synchronize the boundary in case of parallel run
-
FEEngine &
getFEEngine
(const ID &name = "") const¶ return the fem object associated with a provided name
-
FEEngine &
getFEEngineBoundary
(const ID &name = "")¶ return the fem boundary object associated with a provided name
-
template<typename
FEEngineClass
>
voidregisterFEEngineObject
(const std::string &name, Mesh &mesh, UInt spatial_dimension)¶ register a fem object associated with name
-
void
unRegisterFEEngineObject
(const std::string &name)¶ unregister a fem object associated with name
-
SynchronizerRegistry &
getSynchronizerRegistry
()¶ return the synchronizer registry
-
template<typename
FEEngineClass
>
FEEngineClass &getFEEngineClass
(std::string name = "") const¶ return the fem object associated with a provided name
-
template<typename
FEEngineClass
>
FEEngineClass &getFEEngineClassBoundary
(std::string name = "")¶ return the fem boundary object associated with a provided name
-
AnalysisMethod
getAnalysisMethod
() const¶ Get the type of analysis method used.
-
void
setTextModeToDumper
()¶
-
void
addDumpField
(const std::string &field_id)¶
-
void
addDumpFieldVector
(const std::string &field_id)¶
-
void
addDumpFieldToDumper
(const std::string &dumper_name, const std::string &field_id)¶
-
void
addDumpFieldVectorToDumper
(const std::string &dumper_name, const std::string &field_id)¶
-
void
addDumpFieldTensorToDumper
(const std::string &dumper_name, const std::string &field_id)¶
-
void
addDumpFieldTensor
(const std::string &field_id)¶
-
void
setBaseName
(const std::string &field_id)¶
-
void
setBaseNameToDumper
(const std::string &dumper_name, const std::string &basename)¶
-
void
addDumpGroupField
(const std::string &field_id, const std::string &group_name)¶
-
void
addDumpGroupFieldToDumper
(const std::string &dumper_name, const std::string &field_id, const std::string &group_name, ElementKind element_kind, bool padding_flag)¶
-
void
addDumpGroupFieldToDumper
(const std::string &dumper_name, const std::string &field_id, const std::string &group_name, UInt spatial_dimension, ElementKind element_kind, bool padding_flag)¶
-
void
removeDumpGroupField
(const std::string &field_id, const std::string &group_name)¶
-
void
removeDumpGroupFieldFromDumper
(const std::string &dumper_name, const std::string &field_id, const std::string &group_name)¶
-
void
addDumpGroupFieldVector
(const std::string &field_id, const std::string &group_name)¶
-
void
addDumpGroupFieldVectorToDumper
(const std::string &dumper_name, const std::string &field_id, const std::string &group_name)¶
-
std::shared_ptr<dumpers::Field>
createNodalFieldReal
(const std::string&, const std::string&, bool)¶
-
std::shared_ptr<dumpers::Field>
createNodalFieldUInt
(const std::string&, const std::string&, bool)¶
-
std::shared_ptr<dumpers::Field>
createNodalFieldBool
(const std::string&, const std::string&, bool)¶
-
std::shared_ptr<dumpers::Field>
createElementalField
(const std::string&, const std::string&, bool, UInt, ElementKind)¶
-
void
setDirectory
(const std::string &directory)¶
-
void
setDirectoryToDumper
(const std::string &dumper_name, const std::string &directory)¶
-
void
dump
(const std::string &dumper_name)¶
-
void
dump
()¶
-
-
class
akantu
::
NonLocalManagerCallback
¶ Subclassed by akantu::SolidMechanicsModel
Public Functions
-
void
initializeNonLocal
()¶
-
void
insertIntegrationPointsInNeighborhoods
(GhostType ghost_type) = 0¶
-
void
computeNonLocalStresses
(GhostType ghost_type) = 0¶
-
void
updateLocalInternal
(ElementTypeMapReal &internal_flat, GhostType ghost_type, ElementKind kind) = 0¶ update the values of the non local internal
-
void
updateNonLocalInternal
(ElementTypeMapReal &internal_flat, GhostType ghost_type, ElementKind kind) = 0¶ copy the results of the averaging in the materials
-
void
Solvers¶
-
class
akantu
::
ModelSolver
: public akantu::Parsable, public akantu::SolverCallback, public akantu::SynchronizerRegistry¶ Subclassed by akantu::Model
Public Functions
-
~ModelSolver
() override¶
-
void
initDOFManager
()¶ initialize the dof manager based on solver type passed in the input file
-
void
initDOFManager
(const ID &solver_type)¶ initialize the dof manager based on the used chosen solver type
-
void
solveStep
(const ID &solver_id = "")¶ solve a step using a given pre instantiated time step solver and nondynamic linear solver
-
void
solveStep
(SolverCallback &callback, const ID &solver_id = "")¶ solve a step using a given pre instantiated time step solver and non linear solver with a user defined callback instead of the model itself /!\ This can mess up everything
-
void
getNewSolver
(const ID &solver_id, TimeStepSolverType time_step_solver_type, NonLinearSolverType non_linear_solver_type = NonLinearSolverType::_auto)¶ Initialize a time solver that can be used afterwards with its id.
-
void
setIntegrationScheme
(const ID &solver_id, const ID &dof_id, const IntegrationSchemeType &integration_scheme_type, IntegrationScheme::SolutionType solution_type = IntegrationScheme::_not_defined)¶ set an integration scheme for a given dof and a given solver
-
void
setIntegrationScheme
(const ID &solver_id, const ID &dof_id, IntegrationScheme &integration_scheme, IntegrationScheme::SolutionType solution_type = IntegrationScheme::_not_defined)¶ set an externally instantiated integration scheme
-
void
predictor
() override¶ Predictor interface for the callback.
-
void
corrector
() override¶ Corrector interface for the callback.
-
TimeStepSolverType
getDefaultSolverType
() const¶ Default time step solver to instantiate for this model.
-
ModelSolverOptions
getDefaultSolverOptions
(const TimeStepSolverType &type) const¶ Default configurations for a given time step solver.
-
DOFManager &
getDOFManager
()¶ get access to the internal dof manager
-
bool
hasSolver
(const ID &solver_id) const¶ set the parameter ‘param’ of the solver ‘solver_id’
get the parameter ‘param’ of the solver ‘solver_id’ answer to the question “does the solver exists ?”
-
void
setDefaultSolver
(const ID &solver_id)¶ changes the current default solver
-
bool
hasDefaultSolver
() const¶ is a default solver defined
-
bool
hasIntegrationScheme
(const ID &solver_id, const ID &dof_id) const¶ is an integration scheme set for a given solver and a given dof
-
TimeStepSolver &
getTimeStepSolver
(const ID &solver_id = "")¶
-
NonLinearSolver &
getNonLinearSolver
(const ID &solver_id = "")¶
-
const TimeStepSolver &
getTimeStepSolver
(const ID &solver_id = "") const¶
-
const NonLinearSolver &
getNonLinearSolver
(const ID &solver_id = "") const¶
-
-
class
akantu
::
DOFManager
: protected akantu::MeshEventHandler¶ Subclassed by akantu::DOFManagerDefault
Public Functions
-
DOFManager
(const ID &id = "dof_manager")¶
-
~DOFManager
() override¶
-
void
registerDOFs
(const ID &dof_id, Array<Real> &dofs_array, const DOFSupportType &support_type)¶ register an array of degree of freedom
-
void
registerDOFs
(const ID &dof_id, Array<Real> &dofs_array, const ID &support_group)¶ the dof as an implied type of _dst_nodal and is defined only on a subset of nodes
-
void
registerDOFsPrevious
(const ID &dof_id, Array<Real> &dofs_array)¶ register an array of previous values of the degree of freedom
-
void
registerDOFsIncrement
(const ID &dof_id, Array<Real> &dofs_array)¶ register an array of increment of degree of freedom
-
void
registerDOFsDerivative
(const ID &dof_id, UInt order, Array<Real> &dofs_derivative)¶ register an array of derivatives for a particular dof array
-
void
registerBlockedDOFs
(const ID &dof_id, Array<bool> &blocked_dofs)¶ register array representing the blocked degree of freedoms
-
void
assembleToResidual
(const ID &dof_id, Array<Real> &array_to_assemble, Real scale_factor = 1.)¶ Assemble an array to the global residual array.
-
void
assembleToLumpedMatrix
(const ID &dof_id, Array<Real> &array_to_assemble, const ID &lumped_mtx, Real scale_factor = 1.)¶ Assemble an array to the global lumped matrix array.
-
void
assembleElementalArrayLocalArray
(const Array<Real> &elementary_vect, Array<Real> &array_assembeled, ElementType type, GhostType ghost_type, Real scale_factor = 1., const Array<UInt> &filter_elements = empty_filter)¶ Assemble elementary values to a local array of the size nb_nodes * nb_dof_per_node. The dof number is implicitly considered as conn(el, n) * nb_nodes_per_element + d. With 0 < n < nb_nodes_per_element and 0 < d < nb_dof_per_node
-
void
assembleElementalArrayToResidual
(const ID &dof_id, const Array<Real> &elementary_vect, ElementType type, GhostType ghost_type, Real scale_factor = 1., const Array<UInt> &filter_elements = empty_filter)¶ Assemble elementary values to the global residual array. The dof number is implicitly considered as conn(el, n) * nb_nodes_per_element + d. With 0 < n < nb_nodes_per_element and 0 < d < nb_dof_per_node
-
void
assembleElementalArrayToLumpedMatrix
(const ID &dof_id, const Array<Real> &elementary_vect, const ID &lumped_mtx, ElementType type, GhostType ghost_type, Real scale_factor = 1., const Array<UInt> &filter_elements = empty_filter)¶ Assemble elementary values to a global array corresponding to a lumped matrix
-
void
assembleElementalMatricesToMatrix
(const ID &matrix_id, const ID &dof_id, const Array<Real> &elementary_mat, ElementType type, GhostType ghost_type = _not_ghost, const MatrixType &elemental_matrix_type = _symmetric, const Array<UInt> &filter_elements = empty_filter) = 0¶ Assemble elementary values to the global residual array. The dof number is implicitly considered as conn(el, n) * nb_nodes_per_element + d. With 0 < n < nb_nodes_per_element and 0 < d < nb_dof_per_node
-
void
assembleMatMulVectToArray
(const ID &dof_id, const ID &A_id, const Array<Real> &x, Array<Real> &array, Real scale_factor = 1) = 0¶ multiply a vector by a matrix and assemble the result to the residual
-
void
assembleLumpedMatMulVectToResidual
(const ID &dof_id, const ID &A_id, const Array<Real> &x, Real scale_factor = 1) = 0¶ multiply a vector by a lumped matrix and assemble the result to the residual
-
void
assemblePreassembledMatrix
(const ID &dof_id_m, const ID &dof_id_n, const ID &matrix_id, const TermsToAssemble &terms) = 0¶ assemble coupling terms between to dofs
-
void
assembleMatMulVectToResidual
(const ID &dof_id, const ID &A_id, const Array<Real> &x, Real scale_factor = 1)¶ multiply a vector by a matrix and assemble the result to the residual
-
void
assembleMatMulDOFsToResidual
(const ID &A_id, Real scale_factor = 1)¶ multiply the dofs by a matrix and assemble the result to the residual
-
void
updateGlobalBlockedDofs
()¶ updates the global blocked_dofs array
-
void
zeroResidual
()¶ sets the residual to 0
-
void
zeroMatrix
(const ID &mtx)¶ sets the matrix to 0
-
void
zeroLumpedMatrix
(const ID &mtx)¶ sets the lumped matrix to 0
-
void
applyBoundary
(const ID &matrix_id = "J")¶
-
void
getLumpedMatrixPerDOFs
(const ID &dof_id, const ID &lumped_mtx, Array<Real> &lumped)¶ extract a lumped matrix part corresponding to a given dof
-
void
splitSolutionPerDOFs
()¶ splits the solution storage from a global view to the per dof storages
-
Int
globalToLocalEquationNumber
(Int global) const¶ return the local index of the global equation number
-
Int
localToGlobalEquationNumber
(Int local) const¶ converts local equation numbers to global equation numbers;
-
NodeFlag
getDOFFlag
(Int local_id) const¶ get the array of dof types (use only if you know what you do…)
-
std::vector<ID>
getDOFIDs
() const¶ Retrieve all the registered DOFs.
-
Array<Real> &
getDOFs
(const ID &dofs_id)¶ Get a reference to the registered dof array for a given id.
-
DOFSupportType
getSupportType
(const ID &dofs_id) const¶ Get the support type of a given dof.
-
bool
hasDOFs
(const ID &dof_id) const¶ are the dofs registered
-
Array<Real> &
getDOFsDerivatives
(const ID &dofs_id, UInt order)¶ Get a reference to the registered dof derivatives array for a given id.
-
const Array<bool> &
getBlockedDOFs
(const ID &dofs_id) const¶ Get a reference to the blocked dofs array registered for the given id.
-
bool
hasBlockedDOFs
(const ID &dofs_id) const¶ Does the dof has a blocked array.
-
Array<Real> &
getDOFsIncrement
(const ID &dofs_id)¶ Get a reference to the registered dof increment array for a given id.
-
bool
hasDOFsIncrement
(const ID &dofs_id) const¶ Does the dof has a increment array.
-
bool
hasPreviousDOFs
(const ID &dofs_id) const¶ Get a reference to the registered dof array for previous step values a given id
-
void
savePreviousDOFs
(const ID &dofs_id)¶ saves the values from dofs to previous dofs
-
const Array<Real> &
getSolution
(const ID &dofs_id) const¶ Get a reference to the solution array registered for the given id.
-
Array<Real> &
getSolution
(const ID &dofs_id)¶ Get a reference to the solution array registered for the given id.
-
SparseMatrix &
getNewMatrix
(const ID &matrix_id, const MatrixType &matrix_type) = 0¶ Get an instance of a new SparseMatrix.
-
SparseMatrix &
getNewMatrix
(const ID &matrix_id, const ID &matrix_to_copy_id) = 0¶ Get an instance of a new SparseMatrix as a copy of the SparseMatrix matrix_to_copy_id
-
const Array<Int> &
getLocalEquationsNumbers
(const ID &dof_id) const¶ Get the equation numbers corresponding to a dof_id. This might be used to access the matrix.
-
void
getArrayPerDOFs
(const ID &dof_id, const SolverVector &global, Array<Real> &local) = 0¶ extract degrees of freedom (identified by ID) from a global solver array
-
SparseMatrix &
getMatrix
(const ID &matrix_id)¶ Get the reference of an existing matrix.
-
bool
hasMatrix
(const ID &matrix_id) const¶ check if the given matrix exists
-
SolverVector &
getNewLumpedMatrix
(const ID &matrix_id) = 0¶ Get an instance of a new lumped matrix.
-
const SolverVector &
getLumpedMatrix
(const ID &matrix_id) const¶ Get the lumped version of a given matrix.
-
SolverVector &
getLumpedMatrix
(const ID &matrix_id)¶ Get the lumped version of a given matrix.
-
bool
hasLumpedMatrix
(const ID &matrix_id) const¶ check if the given matrix exists
-
NonLinearSolver &
getNewNonLinearSolver
(const ID &nls_solver_id, const NonLinearSolverType &_non_linear_solver_type) = 0¶ Get instance of a non linear solver.
-
NonLinearSolver &
getNonLinearSolver
(const ID &nls_solver_id)¶ get instance of a non linear solver
-
bool
hasNonLinearSolver
(const ID &solver_id) const¶ check if the given solver exists
-
TimeStepSolver &
getNewTimeStepSolver
(const ID &time_step_solver_id, const TimeStepSolverType &type, NonLinearSolver &non_linear_solver, SolverCallback &solver_callback) = 0¶ Get instance of a time step solver.
-
TimeStepSolver &
getTimeStepSolver
(const ID &time_step_solver_id)¶ get instance of a time step solver
-
bool
hasTimeStepSolver
(const ID &solver_id) const¶ check if the given solver exists
-
const auto &
getCommunicator
() const¶
-
auto &
getCommunicator
()¶
-
const auto &
getSolution
() const¶
-
auto &
getSolution
()¶
-
const auto &
getResidual
() const¶
-
auto &
getResidual
()¶
-
void
onNodesAdded
(const Array<UInt> &nodes_list, const NewNodesEvent &event) override¶ function to implement to react on akantu::NewNodesEvent
-
void
onNodesRemoved
(const Array<UInt> &nodes_list, const Array<UInt> &new_numbering, const RemovedNodesEvent &event) override¶ function to implement to react on akantu::RemovedNodesEvent
-
void
onElementsAdded
(const Array<Element> &elements_list, const NewElementsEvent &event) override¶ function to implement to react on akantu::NewElementsEvent
-
void
onElementsRemoved
(const Array<Element> &elements_list, const ElementTypeMapArray<UInt> &new_numbering, const RemovedElementsEvent &event) override¶ function to implement to react on akantu::RemovedElementsEvent
-
void
onElementsChanged
(const Array<Element> &old_elements_list, const Array<Element> &new_elements_list, const ElementTypeMapArray<UInt> &new_numbering, const ChangedElementsEvent &event) override¶ function to implement to react on akantu::ChangedElementsEvent
-
-
class
akantu
::
NonLinearSolver
: public akantu::Parsable¶ Subclassed by akantu::NonLinearSolverLinear, akantu::NonLinearSolverLumped, akantu::NonLinearSolverNewtonRaphson
Public Functions
-
NonLinearSolver
(DOFManager &dof_manager, const NonLinearSolverType &non_linear_solver_type, const ID &id = "non_linear_solver")¶
-
~NonLinearSolver
() override¶
-
void
solve
(SolverCallback &callback) = 0¶ solve the system described by the jacobian matrix, and rhs contained in the dof manager
-
-
class
akantu
::
NonLinearSolverNewtonRaphson
: public akantu::NonLinearSolver¶ Public Functions
-
NonLinearSolverNewtonRaphson
(DOFManagerDefault &dof_manager, const NonLinearSolverType &non_linear_solver_type, const ID &id = "non_linear_solver_newton_raphson")¶
-
~NonLinearSolverNewtonRaphson
() override¶
-
void
solve
(SolverCallback &solver_callback) override¶ Function that solve the non linear system described by the dof manager and the solver callback functions
-
SparseSolverMumps &
getSolver
()¶
-
const SparseSolverMumps &
getSolver
() const¶
-
Solid Mechanics Model¶
-
class
akantu
::
SolidMechanicsModel
: public akantu::Model, public akantu::DataAccessor<Element>, public akantu::DataAccessor<UInt>, public akantu::BoundaryCondition<SolidMechanicsModel>, public akantu::NonLocalManagerCallback, public akantu::EventHandlerManager<SolidMechanicsModelEventHandler>¶ Subclassed by akantu::SolidMechanicsModelCohesive
Public Types
-
using
MyFEEngineType
= FEEngineTemplate<IntegratorGauss, ShapeLagrange>¶
Public Functions
-
SolidMechanicsModel
(Mesh &mesh, UInt dim = _all_dimensions, const ID &id = "solid_mechanics_model", ModelType model_type = ModelType::_solid_mechanics_model)¶ A solid mechanics model need a mesh and a dimension to be created. the model by it self can not do a lot, the good init functions should be called in order to configure the model depending on what we want to do.
- Parameters
mesh
: mesh representing the model we want to simulatedim
: spatial dimension of the problem, if dim = 0 (default value) the dimension of the problem is assumed to be the on of the meshid
: an id to identify the modelmodel_type
: this is an internal parameter for inheritance purposes
-
~SolidMechanicsModel
() override¶
-
void
initMaterials
()¶ initialize all internal arrays for materials
-
void
assembleStiffnessMatrix
()¶ assembles the stiffness matrix,
-
void
assembleInternalForces
()¶ assembles the internal forces in the array internal_forces
This function computes the internal forces as \(F_{int} = \int_{\Omega} N \sigma d\Omega@\)
-
bool
isDefaultSolverExplicit
()¶
-
Material &
registerNewMaterial
(const ID &mat_name, const ID &mat_type, const ID &opt_param)¶ register an empty material of a given type
-
void
reassignMaterial
()¶ reassigns materials depending on the material selector
-
void
applyEigenGradU
(const Matrix<Real> &prescribed_eigen_grad_u, const ID &material_name, GhostType ghost_type = _not_ghost)¶ apply a constant eigen_grad_u on all quadrature points of a given material
-
void
assembleMassLumped
()¶ assemble the lumped mass matrix
-
void
assembleMass
()¶ assemble the mass matrix for consistent mass resolutions
-
void
packData
(CommunicationBuffer &buffer, const Array<Element> &elements, const SynchronizationTag &tag) const override¶
-
void
unpackData
(CommunicationBuffer &buffer, const Array<Element> &elements, const SynchronizationTag &tag) override¶
-
void
packData
(CommunicationBuffer &buffer, const Array<UInt> &dofs, const SynchronizationTag &tag) const override¶
-
void
unpackData
(CommunicationBuffer &buffer, const Array<UInt> &dofs, const SynchronizationTag &tag) override¶
-
void
onDump
()¶
-
bool
isInternal
(const std::string &field_name, ElementKind element_kind)¶ decide wether a field is a material internal or not
-
ElementTypeMap<UInt>
getInternalDataPerElem
(const std::string &field_name, ElementKind kind)¶ give the amount of data per element
-
ElementTypeMapArray<Real> &
flattenInternal
(const std::string &field_name, ElementKind kind, GhostType ghost_type = _not_ghost)¶ flatten a given material internal field
-
void
flattenAllRegisteredInternals
(ElementKind kind)¶ flatten all the registered material internals
-
std::shared_ptr<dumpers::Field>
createNodalFieldReal
(const std::string &field_name, const std::string &group_name, bool padding_flag) override¶
-
std::shared_ptr<dumpers::Field>
createNodalFieldBool
(const std::string &field_name, const std::string &group_name, bool padding_flag) override¶
-
std::shared_ptr<dumpers::Field>
createElementalField
(const std::string &field_name, const std::string &group_name, bool padding_flag, UInt spatial_dimension, ElementKind kind) override¶
-
void
dump
(const std::string &dumper_name) override¶
-
void
dump
() override¶
-
auto &
getDisplacement
()¶ get the SolidMechanicsModel::displacement array
-
const auto &
getDisplacement
() const¶ get the SolidMechanicsModel::displacement array
-
const auto &
getPreviousDisplacement
() const¶ get the SolidMechanicsModel::previous_displacement array
-
const auto &
getIncrement
() const¶ get the SolidMechanicsModel::displacement_increment array
-
auto &
getIncrement
()¶ get the SolidMechanicsModel::displacement_increment array
-
const auto &
getMass
() const¶ get the lumped SolidMechanicsModel::mass array
-
auto &
getVelocity
()¶ get the SolidMechanicsModel::velocity array
-
const auto &
getVelocity
() const¶ get the SolidMechanicsModel::velocity array
-
auto &
getAcceleration
()¶ get the SolidMechanicsModel::acceleration array
-
const auto &
getAcceleration
() const¶ get the SolidMechanicsModel::acceleration array
-
auto &
getExternalForce
()¶ get the SolidMechanicsModel::external_force array
-
const auto &
getExternalForce
() const¶ get the SolidMechanicsModel::external_force array
-
auto &
getInternalForce
()¶ get the SolidMechanicsModel::internal_force array (internal forces)
-
const auto &
getInternalForce
() const¶ get the SolidMechanicsModel::internal_force array (internal forces)
-
auto &
getBlockedDOFs
()¶ get the SolidMechanicsModel::blocked_dofs array
-
const auto &
getBlockedDOFs
() const¶ get the SolidMechanicsModel::blocked_dofs array
-
decltype(auto)
getMaterials
()¶ get an iterable on the materials
-
decltype(auto)
getMaterials
() const¶ get an iterable on the materials
-
const Material &
getMaterial
(UInt mat_index) const¶ get a particular material (by numerical material index)
-
const Material &
getMaterial
(const std::string &name) const¶ get a particular material (by material name)
-
Real
getEnergy
(const std::string &energy_id, ElementType type, UInt index)¶ compute the energy for an element
-
Real
getEnergy
(const std::string &energy_id, const Element &element)¶ compute the energy for an element
-
const ElementTypeMapArray<UInt> &
getMaterialByElement
() const¶
-
const ElementTypeMapArray<UInt> &
getMaterialLocalNumbering
() const¶
-
const Array<UInt> &
getMaterialByElement
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶ vectors containing local material element index for each global element index
-
const Array<UInt> &
getMaterialLocalNumbering
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶
-
MaterialSelector &
getMaterialSelector
()¶
-
NonLocalManager &
getNonLocalManager
() const¶ Access the non_local_manager interface.
-
class
NewMaterialElementsEvent
: public akantu::NewElementsEvent¶
-
using
Warning
doxygenclass: Cannot find class “akantu::SolidMechanicsModelOptions” in doxygen xml output for project “Akantu” from directory: /builds/akantu/akantu/build/doc/dev-doc/xml
-
class
akantu
::
MaterialSelector
: public std::enable_shared_from_this<MaterialSelector>¶ main class to assign same or different materials for different elements
Subclassed by akantu::DefaultMaterialCohesiveSelector, akantu::DefaultMaterialSelector, akantu::ElementDataMaterialSelector< T >, akantu::MaterialCohesiveRulesSelector, akantu::MeshDataMaterialCohesiveSelector
-
template<typename
T
>
classakantu
::
MeshDataMaterialSelector
: public akantu::ElementDataMaterialSelector<T>¶ class to use mesh data information to assign different materials where name is the tag value: tag_0, tag_1
Public Functions
-
MeshDataMaterialSelector
(const std::string &name, const SolidMechanicsModel &model, UInt first_index = 1)¶
-
-
class
akantu
::
Material
: public akantu::DataAccessor<Element>, public akantu::Parsable, public akantu::MeshEventHandler, protected akantu::SolidMechanicsModelEventHandler¶ Interface of all materials Prerequisites for a new material
inherit from this class
implement the following methods:
virtual Real getStableTimeStep(Real h, const Element & element = ElementNull); virtual void computeStress(ElementType el_type, GhostType ghost_type = _not_ghost); virtual void computeTangentStiffness(ElementType el_type, Array<Real> & tangent_matrix, GhostType ghost_type = _not_ghost);
Subclassed by akantu::PlaneStressToolbox< spatial_dimension >, akantu::MaterialCohesive, akantu::MaterialElasticLinearAnisotropic< Dim >, akantu::MaterialThermal< spatial_dimension >, akantu::PlaneStressToolbox< dim, ParentMaterial >
Public Functions
-
Material
(SolidMechanicsModel &model, const ID &id = "")¶ Initialize material with defaults.
-
Material
(SolidMechanicsModel &model, UInt dim, const Mesh &mesh, FEEngine &fe_engine, const ID &id = "")¶ Initialize material with custom mesh & fe_engine.
-
~Material
() override¶ Destructor.
-
void
extrapolateInternal
(const ID &id, const Element &element, const Matrix<Real> &points, Matrix<Real> &extrapolated)¶ extrapolate internal values
-
Real
getCelerity
(const Element &element) const¶ get a material celerity to compute the stable time step (default: is the push wave speed)
-
template<typename
T
>
voidregisterInternal
(InternalField<T>&)¶
-
template<typename
T
>
voidunregisterInternal
(InternalField<T>&)¶
-
void
initMaterial
()¶ initialize the material computed parameter
-
void
assembleInternalForces
(GhostType ghost_type)¶ compute the residual for this material
assemble the residual for this material
Compute the internal forces by assembling \(\int_{e} \sigma_e \frac{\partial \varphi}{\partial X} dX \)
- Parameters
[in] ghost_type
: compute the internal forces for _ghost or _not_ghost element
-
void
savePreviousState
()¶ save the stress in the previous_stress if needed
-
void
restorePreviousState
()¶ restore the stress from previous_stress if needed
-
void
computeAllStresses
(GhostType ghost_type = _not_ghost)¶ compute the stresses for this material
Compute the stress from the gradu
- Parameters
[in] ghost_type
: compute the residual for _ghost or _not_ghost element
-
void
computeAllCauchyStresses
(GhostType ghost_type = _not_ghost)¶
-
void
setToSteadyState
(GhostType ghost_type = _not_ghost)¶ set material to steady state
-
void
assembleStiffnessMatrix
(GhostType ghost_type)¶ compute the stiffness matrix
Compute the stiffness matrix by assembling \(\int_{\omega} B^t \times D \times B d\omega \)
- Parameters
[in] ghost_type
: compute the residual for _ghost or _not_ghost element
-
UInt
addElement
(ElementType type, UInt element, GhostType ghost_type)¶ add an element to the local mesh filter
-
void
printself
(std::ostream &stream, int indent = 0) const override¶ function to print the contain of the class
-
void
interpolateStress
(ElementTypeMapArray<Real> &result, GhostType ghost_type = _not_ghost)¶ interpolate stress on given positions for each element by means of a geometrical interpolation on quadrature points
-
void
interpolateStressOnFacets
(ElementTypeMapArray<Real> &result, ElementTypeMapArray<Real> &by_elem_result, GhostType ghost_type = _not_ghost)¶ interpolate stress on given positions for each element by means of a geometrical interpolation on quadrature points and store the results per facet
-
void
initElementalFieldInterpolation
(const ElementTypeMapArray<Real> &interpolation_points_coordinates)¶ function to initialize the elemental field interpolation function by inverting the quadrature points’ coordinates
-
template<UInt
dim
>
voidStoCauchy
(ElementType el_type, GhostType ghost_type = _not_ghost)¶ Computation of Cauchy stress tensor in the case of finite deformation from the 2nd Piola-Kirchhoff for a given element type
-
template<UInt
dim
>
voidStoCauchy
(const Matrix<Real> &F, const Matrix<Real> &S, Matrix<Real> &sigma, const Real &C33 = 1.0) const¶ Computation the Cauchy stress the 2nd Piola-Kirchhoff and the deformation gradient
-
void
packData
(CommunicationBuffer &buffer, const Array<Element> &elements, const SynchronizationTag &tag) const override¶
-
void
unpackData
(CommunicationBuffer &buffer, const Array<Element> &elements, const SynchronizationTag &tag) override¶
-
template<typename
T
>
voidpackElementDataHelper
(const ElementTypeMapArray<T> &data_to_pack, CommunicationBuffer &buffer, const Array<Element> &elements, const ID &fem_id = ID()) const¶
-
template<typename
T
>
voidunpackElementDataHelper
(ElementTypeMapArray<T> &data_to_unpack, CommunicationBuffer &buffer, const Array<Element> &elements, const ID &fem_id = ID())¶
-
void
onNodesAdded
(const Array<UInt>&, const NewNodesEvent&) override¶ function to implement to react on akantu::NewNodesEvent
-
void
onNodesRemoved
(const Array<UInt>&, const Array<UInt>&, const RemovedNodesEvent&) override¶ function to implement to react on akantu::RemovedNodesEvent
-
void
onElementsAdded
(const Array<Element> &element_list, const NewElementsEvent &event) override¶ function to implement to react on akantu::NewElementsEvent
-
void
onElementsRemoved
(const Array<Element> &element_list, const ElementTypeMapArray<UInt> &new_numbering, const RemovedElementsEvent &event) override¶ function to implement to react on akantu::RemovedElementsEvent
-
void
onElementsChanged
(const Array<Element>&, const Array<Element>&, const ElementTypeMapArray<UInt>&, const ChangedElementsEvent&) override¶ function to implement to react on akantu::ChangedElementsEvent
-
void
beforeSolveStep
()¶
-
void
afterSolveStep
(bool converged = true)¶
-
void
onDamageIteration
() override¶ function to implement to react on akantu::BeginningOfDamageIterationEvent
-
void
onDamageUpdate
() override¶ function to implement to react on akantu::AfterDamageEvent
-
void
onDump
() override¶ function to implement to react on akantu::BeforeDumpEvent
-
const std::string &
getName
() const¶
-
const SolidMechanicsModel &
getModel
() const¶
-
const ID &
getID
() const¶
-
Real
getPotentialEnergy
()¶ return the potential energy for the subset of elements contained by the material
-
Real
getPotentialEnergy
(ElementType &type, UInt index)¶ return the potential energy for the provided element
-
Real
getEnergy
(const std::string &type)¶ return the energy (identified by id) for the subset of elements contained by the material
-
Real
getEnergy
(const std::string &energy_id, ElementType type, UInt index)¶ return the energy (identified by id) for the provided element
-
const Array<UInt> &
getElementFilter
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶
-
const Array<Real> &
getGradU
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶
-
const Array<Real> &
getStress
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶
-
const Array<Real> &
getPotentialEnergy
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶
-
const ElementTypeMapArray<Real> &
getGradU
() const¶
-
const ElementTypeMapArray<Real> &
getStress
() const¶
-
const ElementTypeMapArray<UInt> &
getElementFilter
() const¶
-
bool
isNonLocal
() const¶
-
template<typename
T
>
const Array<T> &getArray
(const ID &id, ElementType type, GhostType ghost_type = _not_ghost) const¶
-
template<typename
T
>
Array<T> &getArray
(const ID &id, ElementType type, GhostType ghost_type = _not_ghost)¶
-
template<typename
T
>
const InternalField<T> &getInternal
(const ID &id) const¶
-
template<typename
T
>
InternalField<T> &getInternal
(const ID &id)¶
-
template<typename
T
>
boolisInternal
(const ID &id, ElementKind element_kind) const¶
-
template<typename
T
>
ElementTypeMap<UInt>getInternalDataPerElem
(const ID &id, ElementKind element_kind) const¶
-
bool
isFiniteDeformation
() const¶
-
bool
isInelasticDeformation
() const¶
-
const Parameter &
getParam
(const ID ¶m) const¶
-
template<typename
T
>
voidflattenInternal
(const std::string &field_id, ElementTypeMapArray<T> &internal_flat, GhostType ghost_type = _not_ghost, ElementKind element_kind = _ek_not_defined) const¶
-
void
applyEigenGradU
(const Matrix<Real> &prescribed_eigen_grad_u, GhostType = _not_ghost)¶ apply a constant eigengrad_u everywhere in the material
-
bool
hasMatrixChanged
(const ID &id)¶
-
MatrixType
getMatrixType
(const ID &id)¶
-
bool
hasStiffnessMatrixChanged
()¶ specify if the matrix need to be recomputed for this material
-
MatrixType
getTangentType
()¶ specify the type of matrix, if not overloaded the material is not valid for static or implicit computations
-
template<typename
T
>
const InternalField<T> &getInternal
([[gnu::unused]] const ID &int_id) const
-
template<typename
T
>
InternalField<T> &getInternal
([[gnu::unused]] const ID &int_id)
-
template<>
const InternalField<Real> &getInternal
(const ID &int_id) const¶
-
template<>
InternalField<Real> &getInternal
(const ID &int_id)¶
Public Static Functions
-
UInt
getCauchyStressMatrixSize
(UInt dim)¶ Size of the Stress matrix for the case of finite deformation see: Bathe et al, IJNME, Vol 9, 353-386, 1975
-
template<UInt
dim
>
voidsetCauchyStressMatrix
(const Matrix<Real> &S_t, Matrix<Real> &sigma)¶ Sets the stress matrix according to Bathe et al, IJNME, Vol 9, 353-386, 1975
-
template<UInt
dim
>
decltype(auto)stressToVoigt
(const Matrix<Real> &stress)¶ write the stress tensor in the Voigt notation.
-
template<UInt
dim
>
decltype(auto)strainToVoigt
(const Matrix<Real> &strain)¶ write the strain tensor in the Voigt notation.
-
template<UInt
dim
>
voidvoigtToStress
(const Vector<Real> &voigt, Matrix<Real> &stress)¶ write a voigt vector to stress
-
MaterialFactory &
getFactory
()¶ static method to reteive the material factory
-
template<typename
T
>
classakantu
::
InternalField
: public akantu::ElementTypeMapArray<T>¶ class for the internal fields of materials to store values for each quadrature
Subclassed by akantu::CohesiveInternalField< bool >, akantu::CohesiveInternalField< UInt >, akantu::CohesiveInternalField< T >, akantu::FacetInternalField< T >, akantu::RandomInternalField< T, BaseField, Generator >
Public Functions
-
~InternalField
() override¶
-
InternalField
(const ID &id, Material &material, FEEngine &fem, const ElementTypeMapArray<UInt> &element_filter)¶ This constructor is only here to let cohesive elements compile.
-
InternalField
(const ID &id, Material &material, UInt dim, FEEngine &fem, const ElementTypeMapArray<UInt> &element_filter)¶ More general constructor.
-
InternalField
(const ID &id, const InternalField<T> &other)¶
-
InternalField
operator=
(const InternalField&) = delete¶
-
void
setElementKind
(ElementKind element_kind)¶ function to reset the element kind for the internal
-
void
initializeHistory
()¶ activate the history of this field
-
void
resize
()¶ resize the arrays and set the new element to 0
-
void
reset
()¶ reset all the fields to the default value
-
void
saveCurrentValues
()¶ save the current values in the history
-
void
restorePreviousValues
()¶ restore the previous values from the history
-
void
removeIntegrationPoints
(const ElementTypeMapArray<UInt> &new_numbering)¶ remove the quadrature points corresponding to suppressed elements
-
void
printself
(std::ostream &stream, int = 0) const override¶ print the content
-
operator T
() const¶ get the default value
-
decltype(auto)
elementTypes
(GhostType ghost_type = _not_ghost) const¶ get filter types for range loop
-
decltype(auto)
filterTypes
(GhostType ghost_type = _not_ghost) const¶ get filter types for range loop
-
const Array<UInt> &
getFilter
(ElementType type, GhostType ghost_type = _not_ghost) const¶ get the array for a given type of the element_filter
-
Array<T> &
operator()
(ElementType type, GhostType ghost_type = _not_ghost)¶ get the Array corresponding to the type en ghost_type specified
-
const Array<T> &
operator()
(ElementType type, GhostType ghost_type = _not_ghost) const¶
-
Array<T> &
previous
(ElementType type, GhostType ghost_type = _not_ghost)¶
-
const Array<T> &
previous
(ElementType type, GhostType ghost_type = _not_ghost) const¶
-
InternalField<T> &
previous
()¶
-
const InternalField<T> &
previous
() const¶
-
bool
hasHistory
() const¶ check if the history is used or not
-
ElementKind
getElementKind
() const¶ get the kind treated by the internal
-
Solid Mechanics Model Cohesive¶
-
class
akantu
::
SolidMechanicsModelCohesive
: public akantu::SolidMechanicsModel, public akantu::SolidMechanicsModelEventHandler¶ Public Types
-
using
MyFEEngineCohesiveType
= FEEngineTemplate<IntegratorGauss, ShapeLagrange, _ek_cohesive>¶
-
using
MyFEEngineFacetType
= FEEngineTemplate<IntegratorGauss, ShapeLagrange, _ek_regular, FacetsCohesiveIntegrationOrderFunctor>¶
Public Functions
-
SolidMechanicsModelCohesive
(Mesh &mesh, UInt dim = _all_dimensions, const ID &id = "solid_mechanics_model_cohesive")¶
-
~SolidMechanicsModelCohesive
() override¶
-
void
assembleInternalForces
() override¶ assemble the residual for the explicit scheme
-
UInt
checkCohesiveStress
()¶ function to perform a stress check on each facet and insert cohesive elements if needed (returns the number of new cohesive elements)
-
void
interpolateStress
()¶ interpolate stress on facets
-
void
updateAutomaticInsertion
()¶ update automatic insertion after a change in the element inserter
-
void
insertIntrinsicElements
()¶ insert intrinsic cohesive elements
-
void
afterSolveStep
(bool converged = true) override¶
-
void
onDump
() override¶
-
void
addDumpGroupFieldToDumper
(const std::string &dumper_name, const std::string &field_id, const std::string &group_name, ElementKind element_kind, bool padding_flag) override¶
-
void
packData
(CommunicationBuffer &buffer, const Array<Element> &elements, const SynchronizationTag &tag) const override¶
-
void
unpackData
(CommunicationBuffer &buffer, const Array<Element> &elements, const SynchronizationTag &tag) override¶
-
const Array<Real> &
getStressOnFacets
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶ get stress on facets vector
-
Array<UInt> &
getFacetMaterial
(const ElementType &el_type, GhostType ghost_type = _not_ghost)¶ get facet material
-
const Array<UInt> &
getFacetMaterial
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶ get facet material
-
const ElementTypeMapArray<UInt> &
getFacetMaterial
() const¶ get facet material
-
const Array<Real> &
getTangents
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶ - Todo:
THIS HAS TO BE CHANGED
-
CohesiveElementInserter &
getElementInserter
()¶ get element inserter
-
bool
getIsExtrinsic
() const¶ get is_extrinsic boolean
-
ElementSynchronizer &
getCohesiveSynchronizer
()¶ get cohesive elements synchronizer
-
class
NewCohesiveNodesEvent
: public akantu::NewNodesEvent¶
-
using
-
class
akantu
::
FragmentManager
: public akantu::GroupManager¶ Public Functions
-
FragmentManager
(SolidMechanicsModelCohesive &model, bool dump_data = true, const ID &id = "fragment_manager")¶
-
void
buildFragments
(Real damage_limit = 1.)¶ build fragment list (cohesive elements are considered broken if damage >= damage_limit)
-
void
computeCenterOfMass
()¶ compute fragments’ center of mass
-
void
computeVelocity
()¶ compute fragments’ velocity
-
void
computeInertiaMoments
()¶ computes principal moments of inertia with respect to the center of mass of each fragment
Given the distance \( \mathbf{r} \) between a quadrature point and its center of mass, the moment of inertia is computed as
\[ I_\mathrm{CM} = \mathrm{tr}(\mathbf{r}\mathbf{r}^\mathrm{T}) \mathbf{I} - \mathbf{r}\mathbf{r}^\mathrm{T} \]for more information check Wikipedia (http://en.wikipedia.org/wiki/Moment_of_inertia#Identities_for_a_skew-symmetric_matrix)
-
void
computeNbElementsPerFragment
()¶ compute number of elements per fragment
-
Heat Transfer Model¶
-
class
akantu
::
HeatTransferModel
: public akantu::Model, public akantu::DataAccessor<Element>, public akantu::DataAccessor<UInt>¶ Public Types
-
using
FEEngineType
= FEEngineTemplate<IntegratorGauss, ShapeLagrange>¶
Public Functions
-
~HeatTransferModel
() override¶
-
void
assembleCapacityLumped
()¶ calculate the lumped capacity vector for heat transfer problem
-
void
assembleConductivityMatrix
()¶ assemble the conductivity matrix
-
void
assembleCapacity
()¶ assemble the conductivity matrix
-
void
computeRho
(Array<Real> &rho, ElementType type, GhostType ghost_type)¶ compute the capacity on quadrature points
-
void
packData
(CommunicationBuffer &buffer, const Array<Element> &elements, const SynchronizationTag &tag) const override¶
-
void
unpackData
(CommunicationBuffer &buffer, const Array<Element> &elements, const SynchronizationTag &tag) override¶
-
void
packData
(CommunicationBuffer &buffer, const Array<UInt> &indexes, const SynchronizationTag &tag) const override¶
-
void
unpackData
(CommunicationBuffer &buffer, const Array<UInt> &indexes, const SynchronizationTag &tag) override¶
-
std::shared_ptr<dumpers::Field>
createNodalFieldReal
(const std::string &field_name, const std::string &group_name, bool padding_flag) override¶
-
std::shared_ptr<dumpers::Field>
createNodalFieldBool
(const std::string &field_name, const std::string &group_name, bool padding_flag) override¶
-
std::shared_ptr<dumpers::Field>
createElementalField
(const std::string &field_name, const std::string &group_name, bool padding_flag, UInt spatial_dimension, ElementKind kind) override¶
-
const Array<Real> &
getTemperatureGradient
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶ get the temperature gradient
-
const Array<Real> &
getConductivityOnQpoints
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶ get the conductivity on q points
-
const Array<Real> &
getTemperatureOnQpoints
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶ get the conductivity on q points
-
const Array<Real> &
getKgradT
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶ internal variables
-
Real
getEnergy
(const std::string &energy_id, ElementType type, UInt index)¶ get the energy denominated by thermal
-
Real
getThermalEnergy
(ElementType type, UInt index)¶ get the thermal energy for a given element
-
using
Structural Mechanics Model¶
Warning
doxygenclass: Cannot find class “akantu::StructuralMaterial” in doxygen xml output for project “Akantu” from directory: /builds/akantu/akantu/build/doc/dev-doc/xml
-
class
akantu
::
StructuralMechanicsModel
: public akantu::Model¶ Public Types
-
using
MyFEEngineType
= FEEngineTemplate<IntegratorGauss, ShapeStructural, _ek_structural>¶
Public Functions
-
StructuralMechanicsModel
(Mesh &mesh, UInt dim = _all_dimensions, const ID &id = "structural_mechanics_model")¶
-
~StructuralMechanicsModel
() override¶
-
void
initFullImpl
(const ModelOptions &options) override¶ Init full model.
-
MatrixType
getMatrixType
(const ID &matrix_id) override¶ get the type of matrix needed
-
void
assembleResidual
() override¶ callback to assemble the residual (rhs)
callback to assemble the residual StructuralMechanicsModel::(rhs)
-
void
assembleResidual
(const ID &residual_part) override¶ callback to assemble the rhs parts, (e.g. internal_forces + external_forces)
-
bool
canSplitResidual
() override¶ tells if the residual can be computed in separated parts
-
std::shared_ptr<dumpers::Field>
createNodalFieldReal
(const std::string &field_name, const std::string &group_name, bool padding_flag) override¶
-
std::shared_ptr<dumpers::Field>
createNodalFieldBool
(const std::string &field_name, const std::string &group_name, bool padding_flag) override¶
-
std::shared_ptr<dumpers::Field>
createElementalField
(const std::string &field_name, const std::string &group_name, bool padding_flag, UInt spatial_dimension, ElementKind kind) override¶
-
Array<Real> &
getAcceleration
() const¶ get the StructuralMechanicsModel::acceleration vector, updated by StructuralMechanicsModel::updateAcceleration
-
Array<Real> &
getInternalForce
() const¶ get the StructuralMechanicsModel::internal_force vector (boundary forces)
-
const Array<Real> &
getRotationMatrix
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶
-
const Array<Real> &
getStress
(const ElementType &el_type, GhostType ghost_type = _not_ghost) const¶
-
Array<UInt> &
getElementMaterial
(const ElementType &el_type, GhostType ghost_type = _not_ghost)¶
-
Array<UInt> &
getSet_ID
(const ElementType &el_type, GhostType ghost_type = _not_ghost)¶
-
UInt
addMaterial
(StructuralMaterial &material, const ID &name = "")¶ This function adds the
StructuralMaterial
material to the list of materials managed by *this.It is important that this function might invalidate all references to structural materials, that were previously optained by
getMaterial()
.- Return
The ID of the material that was added.
- Note
The return type is is new.
- Parameters
material
: The new material.
-
const StructuralMaterial &
getMaterialByElement
(const Element &element) const¶
-
const StructuralMaterial &
getMaterial
(UInt material_index) const¶ Returns the ith material of *this.
- Parameters
i
: The ith material
-
const StructuralMaterial &
getMaterial
(const ID &name) const¶
-
void
computeForcesByGlobalTractionArray
(const Array<Real> &traction_global, ElementType type)¶ Compute Linear load function set in global axis.
-
void
computeForcesByLocalTractionArray
(const Array<Real> &tractions, ElementType type)¶ Compute Linear load function set in local axis.
-
using
Synchronizers¶
-
template<class
T
>
classakantu
::
DataAccessor
: public virtual akantu::DataAccessorBase¶ Public Functions
-
DataAccessor
() = default¶
-
~DataAccessor
() override = default¶
-
UInt
getNbData
(const Array<T> &elements, const SynchronizationTag &tag) const = 0¶ get the number of data to exchange for a given array of T (elements or dofs) and a given akantu::SynchronizationTag
-
Input/Output¶
-
class
akantu
::
Dumpable
¶ Subclassed by akantu::ElementGroup, akantu::Mesh, akantu::NodeGroup
Public Functions
-
Dumpable
()¶
-
~Dumpable
()¶
-
template<class
T
>
voidregisterDumper
(const std::string &dumper_name, const std::string &file_name = "", bool is_default = false)¶ create a new dumper (of templated type T) and register it under dumper_name. file_name is used for construction of T. is default states if this dumper is the default dumper.
register an externally created dumper
-
void
addDumpMesh
(const Mesh &mesh, UInt spatial_dimension = _all_dimensions, GhostType ghost_type = _not_ghost, ElementKind element_kind = _ek_not_defined)¶ register a mesh to the default dumper
-
void
addDumpMeshToDumper
(const std::string &dumper_name, const Mesh &mesh, UInt spatial_dimension = _all_dimensions, GhostType ghost_type = _not_ghost, ElementKind element_kind = _ek_not_defined)¶ register a mesh to the default identified by its name
-
void
addDumpFilteredMesh
(const Mesh &mesh, const ElementTypeMapArray<UInt> &elements_filter, const Array<UInt> &nodes_filter, UInt spatial_dimension = _all_dimensions, GhostType ghost_type = _not_ghost, ElementKind element_kind = _ek_not_defined)¶ register a filtered mesh as the default dumper
-
void
addDumpFilteredMeshToDumper
(const std::string &dumper_name, const Mesh &mesh, const ElementTypeMapArray<UInt> &elements_filter, const Array<UInt> &nodes_filter, UInt spatial_dimension = _all_dimensions, GhostType ghost_type = _not_ghost, ElementKind element_kind = _ek_not_defined)¶ register a filtered mesh and provides a name
-
void
addDumpField
(const std::string &field_id)¶ to implement
-
void
addDumpFieldToDumper
(const std::string &dumper_name, const std::string &field_id)¶ to implement
add a field
-
template<typename
T
>
voidaddDumpFieldExternalToDumper
(const std::string &dumper_name, const std::string &field_id, const Array<T> &field)¶
-
template<typename
T
>
voidaddDumpFieldExternal
(const std::string &field_id, const ElementTypeMapArray<T> &field, UInt spatial_dimension = _all_dimensions, GhostType ghost_type = _not_ghost, ElementKind element_kind = _ek_not_defined)¶
-
template<typename
T
>
voidaddDumpFieldExternalToDumper
(const std::string &dumper_name, const std::string &field_id, const ElementTypeMapArray<T> &field, UInt spatial_dimension = _all_dimensions, GhostType ghost_type = _not_ghost, ElementKind element_kind = _ek_not_defined)¶
-
void
removeDumpField
(const std::string &field_id)¶
-
void
removeDumpFieldFromDumper
(const std::string &dumper_name, const std::string &field_id)¶
-
void
addDumpFieldVector
(const std::string &field_id)¶
-
void
addDumpFieldVectorToDumper
(const std::string &dumper_name, const std::string &field_id)¶
-
void
addDumpFieldTensor
(const std::string &field_id)¶
-
void
addDumpFieldTensorToDumper
(const std::string &dumper_name, const std::string &field_id)¶
-
void
setDirectory
(const std::string &directory)¶
-
void
setDirectoryToDumper
(const std::string &dumper_name, const std::string &directory)¶
-
void
setBaseName
(const std::string &basename)¶
-
void
setBaseNameToDumper
(const std::string &dumper_name, const std::string &basename)¶
-
void
setTextModeToDumper
(const std::string &dumper_name)¶
-
void
setTextModeToDumper
()¶
-
void
dump
()¶
-
void
dump
(const std::string &dumper_name)¶
-
DumperIOHelper &
getDumper
()¶
-
DumperIOHelper &
getDumper
(const std::string &dumper_name)¶
-
std::string
getDefaultDumperName
() const¶
-
-
class
akantu
::
DumperIOHelper
: public std::enable_shared_from_this<DumperIOHelper>¶ Subclassed by akantu::DumperParaview, akantu::DumperText
Public Functions
-
DumperIOHelper
()¶
-
~DumperIOHelper
()¶
-
void
registerMesh
(const Mesh &mesh, UInt spatial_dimension = _all_dimensions, GhostType ghost_type = _not_ghost, ElementKind element_kind = _ek_not_defined)¶ register a given Mesh for the current dumper
-
void
registerFilteredMesh
(const Mesh &mesh, const ElementTypeMapArray<UInt> &elements_filter, const Array<UInt> &nodes_filter, UInt spatial_dimension = _all_dimensions, GhostType ghost_type = _not_ghost, ElementKind element_kind = _ek_not_defined)¶ register a filtered Mesh (provided filter lists) for the current dumper
register a Field object identified by name and provided by pointer
-
void
unRegisterField
(const std::string &field_id)¶ remove the Field identified by name from managed fields
register a VariableBase object identified by name and provided by pointer
-
void
unRegisterVariable
(const std::string &variable_id)¶ remove a VariableBase identified by name from managed fields
-
void
dump
()¶ request dump: this calls IOHelper dump routine
-
void
dump
(UInt step)¶ request dump: this first set the current step and then calls IOHelper dump routine
-
void
dump
(Real current_time, UInt step)¶ request dump: this first set the current step and current time and then calls IOHelper dump routine
-
void
setParallelContext
(bool is_parallel)¶ set the parallel context for IOHeper
-
void
setDirectory
(const std::string &directory)¶ set the directory where to generate the dumped files
-
void
setBaseName
(const std::string &basename)¶ set the base name (needed by most IOHelper dumpers)
-
iohelper::Dumper &
getDumper
() const¶ direct access to the iohelper::Dumper object
-
-
class
akantu
::
DumperParaview
: public akantu::DumperIOHelper¶
-
class
akantu
::
DumperText
: public akantu::DumperIOHelper¶ Public Functions
-
DumperText
(const std::string &basename = "dumper_text", iohelper::TextDumpMode mode = iohelper::_tdm_space, bool parallel = true)¶
-
~DumperText
() override = default¶
-
void
registerMesh
(const Mesh &mesh, UInt spatial_dimension = _all_dimensions, GhostType ghost_type = _not_ghost, ElementKind element_kind = _ek_not_defined) override¶ register a given Mesh for the current dumper
-
void
registerFilteredMesh
(const Mesh &mesh, const ElementTypeMapArray<UInt> &elements_filter, const Array<UInt> &nodes_filter, UInt spatial_dimension = _all_dimensions, GhostType ghost_type = _not_ghost, ElementKind element_kind = _ek_not_defined) override¶ register a filtered Mesh (provided filter lists) for the current dumper
-
void
setBaseName
(const std::string &basename) override¶ set the base name (needed by most IOHelper dumpers)
-
Warning
doxygenclass: Cannot find class “akantu::Field” in doxygen xml output for project “Akantu” from directory: /builds/akantu/akantu/build/doc/dev-doc/xml
-
class
akantu
::
Parser
: public akantu::ParserSection¶ Root of parsing tree, represents the global ParserSection.