OGS
ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim > Class Template Reference

Detailed Description

template<typename ShapeFunction, int DisplacementDim>
class ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >

Definition at line 35 of file SmallDeformationLocalAssemblerMatrix.h.

#include <SmallDeformationLocalAssemblerMatrix.h>

Inheritance diagram for ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >:
[legend]
Collaboration diagram for ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >:
[legend]

Public Types

using ShapeMatricesType
using NodalMatrixType = typename ShapeMatricesType::NodalMatrixType
using NodalVectorType = typename ShapeMatricesType::NodalVectorType
using ShapeMatrices = typename ShapeMatricesType::ShapeMatrices
using BMatricesType = BMatrixPolicyType<ShapeFunction, DisplacementDim>
using BBarMatrixType = typename BMatricesType::BBarMatrixType
using StiffnessMatrixType = typename BMatricesType::StiffnessMatrixType
using NodalForceVectorType = typename BMatricesType::NodalForceVectorType
using NodalDisplacementVectorType

Public Member Functions

 SmallDeformationLocalAssemblerMatrix (SmallDeformationLocalAssemblerMatrix const &)=delete
 SmallDeformationLocalAssemblerMatrix (SmallDeformationLocalAssemblerMatrix &&)=delete
 SmallDeformationLocalAssemblerMatrix (MeshLib::Element const &e, std::size_t const local_matrix_size, NumLib::GenericIntegrationMethod const &integration_method, bool const is_axially_symmetric, SmallDeformationProcessData< DisplacementDim > &process_data)
void assemble (double const, double const, std::vector< double > const &, std::vector< double > const &, std::vector< double > &, std::vector< double > &, std::vector< double > &) override
void assembleWithJacobian (double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data) override
void preTimestepConcrete (std::vector< double > const &, double const, double const) override
void computeSecondaryVariableConcreteWithVector (double const t, Eigen::VectorXd const &local_x) override
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix (const unsigned integration_point) const override
 Provides the shape matrix at the given integration point.
std::vector< double > const & getIntPtSigma (const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
std::vector< double > const & getIntPtEpsilon (const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
std::vector< double > const & getIntPtFractureStress (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
std::vector< double > const & getIntPtFractureAperture (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Public Member Functions inherited from ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerInterface
 SmallDeformationLocalAssemblerInterface ()
 SmallDeformationLocalAssemblerInterface (std::size_t n_local_size, std::vector< unsigned > dofIndex_to_localIndex)
virtual void assembleWithJacobian (double const, double const, Eigen::VectorXd const &, Eigen::VectorXd &, Eigen::MatrixXd &)
void computeSecondaryVariableConcrete (double const t, double const, Eigen::VectorXd const &local_x, Eigen::VectorXd const &) override
Public Member Functions inherited from ProcessLib::LocalAssemblerInterface
virtual ~LocalAssemblerInterface ()=default
virtual void setInitialConditions (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, double const t, int const process_id)
virtual void initialize (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)
virtual void preAssemble (double const, double const, std::vector< double > const &)
virtual void assembleForStaggeredScheme (double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, int const process_id, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data)
virtual void assembleWithJacobianForStaggeredScheme (double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, int const process_id, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
virtual void computeSecondaryVariable (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_prev, int const process_id)
virtual void preTimestep (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, double const t, double const delta_t)
virtual void postTimestep (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, int const process_id)
void postNonLinearSolver (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, int const process_id)
virtual Eigen::Vector3d getFlux (MathLib::Point3d const &, double const, std::vector< double > const &) const
virtual Eigen::Vector3d getFlux (MathLib::Point3d const &, double const, std::vector< std::vector< double > > const &) const
 Fits to staggered scheme.
virtual std::optional< VectorSegmentgetVectorDeformationSegment () const
Public Member Functions inherited from NumLib::ExtrapolatableElement
virtual ~ExtrapolatableElement ()=default

Private Types

using IntegrationPointDataType

Private Member Functions

std::optional< BBarMatrixTypegetDilatationalBBarMatrix () const

Private Attributes

SmallDeformationProcessData< DisplacementDim > & _process_data
std::vector< IntegrationPointDataType, Eigen::aligned_allocator< IntegrationPointDataType > > _ip_data
NumLib::GenericIntegrationMethod const & _integration_method
MeshLib::Element const & _element
bool const _is_axially_symmetric
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data

Member Typedef Documentation

◆ BBarMatrixType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::BBarMatrixType = typename BMatricesType::BBarMatrixType

Definition at line 45 of file SmallDeformationLocalAssemblerMatrix.h.

◆ BMatricesType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::BMatricesType = BMatrixPolicyType<ShapeFunction, DisplacementDim>

Definition at line 44 of file SmallDeformationLocalAssemblerMatrix.h.

◆ IntegrationPointDataType

◆ NodalDisplacementVectorType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::NodalDisplacementVectorType
Initial value:
VectorType< _number_of_dof > NodalForceVectorType
Rhs residual.

Definition at line 49 of file SmallDeformationLocalAssemblerMatrix.h.

◆ NodalForceVectorType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::NodalForceVectorType = typename BMatricesType::NodalForceVectorType

Definition at line 48 of file SmallDeformationLocalAssemblerMatrix.h.

◆ NodalMatrixType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::NodalMatrixType = typename ShapeMatricesType::NodalMatrixType

Definition at line 41 of file SmallDeformationLocalAssemblerMatrix.h.

◆ NodalVectorType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::NodalVectorType = typename ShapeMatricesType::NodalVectorType

Definition at line 42 of file SmallDeformationLocalAssemblerMatrix.h.

◆ ShapeMatrices

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::ShapeMatrices = typename ShapeMatricesType::ShapeMatrices

Definition at line 43 of file SmallDeformationLocalAssemblerMatrix.h.

◆ ShapeMatricesType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::ShapeMatricesType
Initial value:
EigenFixedShapeMatrixPolicy< ShapeFunction, GlobalDim > ShapeMatrixPolicyType

Definition at line 39 of file SmallDeformationLocalAssemblerMatrix.h.

◆ StiffnessMatrixType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::StiffnessMatrixType = typename BMatricesType::StiffnessMatrixType

Definition at line 47 of file SmallDeformationLocalAssemblerMatrix.h.

Constructor & Destructor Documentation

◆ SmallDeformationLocalAssemblerMatrix() [1/3]

template<typename ShapeFunction, int DisplacementDim>
ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::SmallDeformationLocalAssemblerMatrix ( SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim > const & )
delete

◆ SmallDeformationLocalAssemblerMatrix() [2/3]

template<typename ShapeFunction, int DisplacementDim>
ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::SmallDeformationLocalAssemblerMatrix ( SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim > && )
delete

◆ SmallDeformationLocalAssemblerMatrix() [3/3]

template<typename ShapeFunction, int DisplacementDim>
ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::SmallDeformationLocalAssemblerMatrix ( MeshLib::Element const & e,
std::size_t const local_matrix_size,
NumLib::GenericIntegrationMethod const & integration_method,
bool const is_axially_symmetric,
SmallDeformationProcessData< DisplacementDim > & process_data )

Definition at line 38 of file SmallDeformationLocalAssemblerMatrix-impl.h.

47 _element(e),
49{
50 unsigned const n_integration_points =
51 _integration_method.getNumberOfPoints();
52
55
56 auto const shape_matrices =
60
62 _process_data.solid_materials, _process_data.material_ids, e.getID());
63
64 for (unsigned ip = 0; ip < n_integration_points; ip++)
65 {
66 _ip_data.emplace_back(solid_material);
67 auto& ip_data = _ip_data[ip];
68 auto const& sm = shape_matrices[ip];
69 ip_data.N_u = sm.N;
70 ip_data.dNdx_u = sm.dNdx;
71 ip_data.integration_weight =
72 _integration_method.getWeightedPoint(ip).getWeight() *
73 sm.integralMeasure * sm.detJ;
74
75 // Initialize current time step values
76 static const int kelvin_vector_size =
78 ip_data._sigma.setZero(kelvin_vector_size);
79 ip_data._eps.setZero(kelvin_vector_size);
80
81 // Previous time step values are not initialized and are set later.
82 ip_data._sigma_prev.resize(kelvin_vector_size);
83 ip_data._eps_prev.resize(kelvin_vector_size);
84
86
87 _secondary_data.N[ip] = sm.N;
88 }
89}
std::vector< IntegrationPointDataType, Eigen::aligned_allocator< IntegrationPointDataType > > _ip_data
auto & selectSolidConstitutiveRelation(SolidMaterialsMap const &constitutive_relations, MeshLib::PropertyVector< int > const *const material_ids, std::size_t const element_id)
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.

References _element, _integration_method, _ip_data, _is_axially_symmetric, _process_data, _secondary_data, MeshLib::Element::getID(), NumLib::initShapeMatrices(), MathLib::KelvinVector::kelvin_vector_dimensions(), and MaterialLib::Solids::selectSolidConstitutiveRelation().

Member Function Documentation

◆ assemble()

template<typename ShapeFunction, int DisplacementDim>
void ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::assemble ( double const ,
double const ,
std::vector< double > const & ,
std::vector< double > const & ,
std::vector< double > & ,
std::vector< double > & ,
std::vector< double > &  )
inlineoverridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 64 of file SmallDeformationLocalAssemblerMatrix.h.

70 {
72 "SmallDeformationLocalAssembler: assembly without jacobian is not "
73 "implemented.");
74 }
#define OGS_FATAL(...)
Definition Error.h:26

References OGS_FATAL.

◆ assembleWithJacobian()

template<typename ShapeFunction, int DisplacementDim>
void ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::assembleWithJacobian ( double const t,
double const dt,
std::vector< double > const & local_x,
std::vector< double > const & ,
std::vector< double > & local_b_data,
std::vector< double > & local_Jac_data )
overridevirtual

Reimplemented from ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerInterface.

Definition at line 92 of file SmallDeformationLocalAssemblerMatrix-impl.h.

98{
99 assert(_element.getDimension() == DisplacementDim);
100
101 auto const local_matrix_size = local_x.size();
102
105
108
109 unsigned const n_integration_points =
110 _integration_method.getNumberOfPoints();
111
115 x_position.setElementID(_element.getID());
116
118
119 for (unsigned ip = 0; ip < n_integration_points; ip++)
120 {
121 auto const& w = _ip_data[ip].integration_weight;
122
123 auto const& N = _ip_data[ip].N_u;
124 auto const& dNdx = _ip_data[ip].dNdx_u;
125 auto const x_coord =
127 _element, N);
128
133
134 auto const& eps_prev = _ip_data[ip]._eps_prev;
135 auto const& sigma_prev = _ip_data[ip]._sigma_prev;
136
137 auto& eps = _ip_data[ip]._eps;
138 auto& sigma = _ip_data[ip]._sigma;
139 auto& state = _ip_data[ip]._material_state_variables;
140
141 eps.noalias() =
144
145 variables.mechanical_strain
147 eps);
148 variables_prev.stress
150 sigma_prev);
151 variables_prev.mechanical_strain
153 eps_prev);
154 variables_prev.temperature = _process_data.reference_temperature;
155
156 auto&& solution = _ip_data[ip]._solid_material.integrateStress(
158
159 if (!solution)
160 {
161 OGS_FATAL("Computation of local constitutive relation failed.");
162 }
163
166
167 local_b.noalias() -= B.transpose() * sigma * w;
168 local_Jac.noalias() += B.transpose() * C * B * w;
169 }
170}
MatrixType< _kelvin_vector_size, _number_of_dof > BMatrixType
Eigen::Map< Vector > createZeroedVector(std::vector< double > &data, Eigen::VectorXd::Index size)
Eigen::Map< Matrix > createZeroedMatrix(std::vector< double > &data, Eigen::MatrixXd::Index rows, Eigen::MatrixXd::Index cols)
double interpolateXCoordinate(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)

References _element, _integration_method, _ip_data, _is_axially_symmetric, _process_data, ProcessLib::LinearBMatrix::computeBMatrixPossiblyWithBbar(), MathLib::createZeroedMatrix(), MathLib::createZeroedVector(), getDilatationalBBarMatrix(), NumLib::interpolateXCoordinate(), MaterialPropertyLib::VariableArray::mechanical_strain, OGS_FATAL, ParameterLib::SpatialPosition::setElementID(), MaterialPropertyLib::VariableArray::stress, and MaterialPropertyLib::VariableArray::temperature.

◆ computeSecondaryVariableConcreteWithVector()

template<typename ShapeFunction, int DisplacementDim>
void ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::computeSecondaryVariableConcreteWithVector ( double const t,
Eigen::VectorXd const & local_x )
overridevirtual

Implements ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerInterface.

Definition at line 173 of file SmallDeformationLocalAssemblerMatrix-impl.h.

176{
177 // Compute average value per element
180 auto const e_id = _element.getID();
181
182 unsigned const n_integration_points =
183 _integration_method.getNumberOfPoints();
184 for (unsigned ip = 0; ip < n_integration_points; ip++)
185 {
186 sigma_avg += _ip_data[ip]._sigma;
187 }
189
191 &(*_process_data.element_stresses)[e_id * KV::RowsAtCompileTime]) =
193}
Eigen::Matrix< double, 4, 1 > kelvinVectorToSymmetricTensor(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)

References _element, _integration_method, _ip_data, _process_data, and MathLib::KelvinVector::kelvinVectorToSymmetricTensor().

◆ getDilatationalBBarMatrix()

template<typename ShapeFunction, int DisplacementDim>
std::optional< BBarMatrixType > ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::getDilatationalBBarMatrix ( ) const
inlineprivate

◆ getIntPtEpsilon()

template<typename ShapeFunction, int DisplacementDim>
std::vector< double > const & ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::getIntPtEpsilon ( const double t,
std::vector< GlobalVector * > const & x,
std::vector< NumLib::LocalToGlobalIndexMap const * > const & dof_table,
std::vector< double > & cache ) const
overridevirtual

◆ getIntPtFractureAperture()

template<typename ShapeFunction, int DisplacementDim>
std::vector< double > const & ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::getIntPtFractureAperture ( const double ,
std::vector< GlobalVector * > const & ,
std::vector< NumLib::LocalToGlobalIndexMap const * > const & ,
std::vector< double > & cache ) const
inlineoverridevirtual

Implements ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerInterface.

Definition at line 129 of file SmallDeformationLocalAssemblerMatrix.h.

134 {
135 cache.resize(0);
136 return cache;
137 }

◆ getIntPtFractureStress()

template<typename ShapeFunction, int DisplacementDim>
std::vector< double > const & ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::getIntPtFractureStress ( const double ,
std::vector< GlobalVector * > const & ,
std::vector< NumLib::LocalToGlobalIndexMap const * > const & ,
std::vector< double > & cache ) const
inlineoverridevirtual

Implements ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerInterface.

Definition at line 119 of file SmallDeformationLocalAssemblerMatrix.h.

124 {
125 cache.resize(0);
126 return cache;
127 }

◆ getIntPtSigma()

◆ getShapeMatrix()

template<typename ShapeFunction, int DisplacementDim>
Eigen::Map< const Eigen::RowVectorXd > ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::getShapeMatrix ( const unsigned integration_point) const
inlineoverridevirtual

Provides the shape matrix at the given integration point.

Implements NumLib::ExtrapolatableElement.

Definition at line 98 of file SmallDeformationLocalAssemblerMatrix.h.

100 {
101 auto const& N = _secondary_data.N[integration_point];
102
103 // assumes N is stored contiguously in memory
104 return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
105 }

References _secondary_data.

◆ preTimestepConcrete()

template<typename ShapeFunction, int DisplacementDim>
void ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::preTimestepConcrete ( std::vector< double > const & ,
double const ,
double const  )
inlineoverridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 82 of file SmallDeformationLocalAssemblerMatrix.h.

85 {
86 unsigned const n_integration_points =
87 _integration_method.getNumberOfPoints();
88
89 for (unsigned ip = 0; ip < n_integration_points; ip++)
90 {
91 _ip_data[ip].pushBackState();
92 }
93 }

References _integration_method, and _ip_data.

Member Data Documentation

◆ _element

◆ _integration_method

◆ _ip_data

◆ _is_axially_symmetric

template<typename ShapeFunction, int DisplacementDim>
bool const ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::_is_axially_symmetric
private

◆ _process_data

◆ _secondary_data

template<typename ShapeFunction, int DisplacementDim>
SecondaryData<typename ShapeMatrices::ShapeType> ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrix< ShapeFunction, DisplacementDim >::_secondary_data
private

The documentation for this class was generated from the following files: