OGS
ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim > Class Template Reference

Detailed Description

template<typename ShapeFunction, int DisplacementDim>
class ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >

Definition at line 81 of file LargeDeformationFEM.h.

#include <LargeDeformationFEM.h>

Inheritance diagram for ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >:
[legend]
Collaboration diagram for ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >:
[legend]

Public Types

using ShapeMatricesType
 
using NodalMatrixType = typename ShapeMatricesType::NodalMatrixType
 
using GlobalDimNodalMatrixType
 
using NodalVectorType = typename ShapeMatricesType::NodalVectorType
 
using ShapeMatrices = typename ShapeMatricesType::ShapeMatrices
 
using BMatricesType = BMatrixPolicyType<ShapeFunction, DisplacementDim>
 
using BMatrixType = typename BMatricesType::BMatrixType
 
using StiffnessMatrixType = typename BMatricesType::StiffnessMatrixType
 
using NodalForceVectorType = typename BMatricesType::NodalForceVectorType
 
using NodalDisplacementVectorType
 
using GMatricesType = GMatrixPolicyType<ShapeFunction, DisplacementDim>
 
using GradientVectorType = typename GMatricesType::GradientVectorType
 
using GradientMatrixType = typename GMatricesType::GradientMatrixType
 
using VectorTypeForFbar = typename BMatricesType::NodalForceVectorType
 
using IpData
 

Public Member Functions

 LargeDeformationLocalAssembler (LargeDeformationLocalAssembler const &)=delete
 
 LargeDeformationLocalAssembler (LargeDeformationLocalAssembler &&)=delete
 
 LargeDeformationLocalAssembler (MeshLib::Element const &e, std::size_t const, NumLib::GenericIntegrationMethod const &integration_method, bool const is_axially_symmetric, LargeDeformationProcessData< DisplacementDim > &process_data)
 
void initializeConcrete () override
 
double computeOutputStrainData (double const detF0, BMatrixType const &B, GradientVectorType const &grad_u, Eigen::Ref< Eigen::VectorXd const > const &u, typename ConstitutiveRelations::OutputData< DisplacementDim > &output_data) const
 
ConstitutiveRelations::ConstitutiveData< DisplacementDim > updateConstitutiveRelations (double const alpha, GradientMatrixType const &G, Eigen::Ref< Eigen::VectorXd const > const &u_prev, ParameterLib::SpatialPosition const &x_position, double const t, double const dt, typename ConstitutiveRelations::ConstitutiveSetting< DisplacementDim > &CS, MaterialPropertyLib::Medium const &medium, typename ConstitutiveRelations::StatefulData< DisplacementDim > &current_state, typename ConstitutiveRelations::StatefulDataPrev< DisplacementDim > const &prev_state, MaterialStateData< DisplacementDim > &material_state, typename ConstitutiveRelations::OutputData< DisplacementDim > &output_data) const
 
void assemble (double const, double const, std::vector< double > const &, std::vector< double > const &, std::vector< double > &, std::vector< double > &, std::vector< double > &) override
 
std::optional< std::tuple< double, VectorTypeForFbar > > computeFBarVariables (bool const compute_detF0_only, Eigen::VectorXd const &local_x) const
 
void assembleWithJacobian (double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_x_prev, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data) override
 
void postTimestepConcrete (Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, double const t, double const dt, int const) override
 
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix (const unsigned integration_point) const override
 Provides the shape matrix at the given integration point.
 
- Public Member Functions inherited from ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >
 LargeDeformationLocalAssemblerInterface (MeshLib::Element const &e, NumLib::GenericIntegrationMethod const &integration_method, bool const is_axially_symmetric, LargeDeformationProcessData< DisplacementDim > &process_data)
 
std::size_t setIPDataInitialConditions (std::string_view name, double const *values, int const integration_order)
 Returns number of read integration points.
 
unsigned getNumberOfIntegrationPoints () const
 
int getMaterialID () const
 
std::vector< double > getMaterialStateVariableInternalState (std::function< std::span< double >(typename MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables &)> const &get_values_span, int const &n_components) const
 
MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables const & getMaterialStateVariablesAt (unsigned integration_point) const
 
- 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.
 
- Public Member Functions inherited from NumLib::ExtrapolatableElement
virtual ~ExtrapolatableElement ()=default
 

Static Public Attributes

static constexpr auto & N_u_op
 

Static Private Member Functions

static constexpr auto localDOF (std::vector< double > const &x)
 

Private Attributes

std::vector< IpData, Eigen::aligned_allocator< IpData > > _ip_data
 
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data
 

Static Private Attributes

static const int displacement_size
 

Additional Inherited Members

- Static Public Member Functions inherited from ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >
static auto getReflectionDataForOutput ()
 
- Protected Attributes inherited from ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >
LargeDeformationProcessData< DisplacementDim > & process_data_
 
std::vector< MaterialStateData< DisplacementDim > > material_states_
 
std::vector< typename ConstitutiveRelations::StatefulData< DisplacementDim > > current_states_
 
std::vector< typename ConstitutiveRelations::StatefulDataPrev< DisplacementDim > > prev_states_
 
std::vector< typename ConstitutiveRelations::OutputData< DisplacementDim > > output_data_
 
NumLib::GenericIntegrationMethod const & integration_method_
 
MeshLib::Element const & element_
 
bool const is_axially_symmetric_
 
ConstitutiveRelations::SolidConstitutiveRelation< DisplacementDim > const & solid_material_
 

Member Typedef Documentation

◆ BMatricesType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::BMatricesType = BMatrixPolicyType<ShapeFunction, DisplacementDim>

Definition at line 92 of file LargeDeformationFEM.h.

◆ BMatrixType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::BMatrixType = typename BMatricesType::BMatrixType

Definition at line 94 of file LargeDeformationFEM.h.

◆ GlobalDimNodalMatrixType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::GlobalDimNodalMatrixType
Initial value:
MatrixType< GlobalDim, ShapeFunction::NPOINTS > GlobalDimNodalMatrixType

Definition at line 88 of file LargeDeformationFEM.h.

◆ GMatricesType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::GMatricesType = GMatrixPolicyType<ShapeFunction, DisplacementDim>

Definition at line 100 of file LargeDeformationFEM.h.

◆ GradientMatrixType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::GradientMatrixType = typename GMatricesType::GradientMatrixType

Definition at line 102 of file LargeDeformationFEM.h.

◆ GradientVectorType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::GradientVectorType = typename GMatricesType::GradientVectorType

Definition at line 101 of file LargeDeformationFEM.h.

◆ IpData

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::IpData
Initial value:
IntegrationPointData<BMatricesType, ShapeMatricesType, DisplacementDim>

Definition at line 106 of file LargeDeformationFEM.h.

◆ NodalDisplacementVectorType

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

Definition at line 97 of file LargeDeformationFEM.h.

◆ NodalForceVectorType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::NodalForceVectorType = typename BMatricesType::NodalForceVectorType

Definition at line 96 of file LargeDeformationFEM.h.

◆ NodalMatrixType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::NodalMatrixType = typename ShapeMatricesType::NodalMatrixType

Definition at line 87 of file LargeDeformationFEM.h.

◆ NodalVectorType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::NodalVectorType = typename ShapeMatricesType::NodalVectorType

Definition at line 90 of file LargeDeformationFEM.h.

◆ ShapeMatrices

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::ShapeMatrices = typename ShapeMatricesType::ShapeMatrices

Definition at line 91 of file LargeDeformationFEM.h.

◆ ShapeMatricesType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::ShapeMatricesType

◆ StiffnessMatrixType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::StiffnessMatrixType = typename BMatricesType::StiffnessMatrixType

Definition at line 95 of file LargeDeformationFEM.h.

◆ VectorTypeForFbar

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::VectorTypeForFbar = typename BMatricesType::NodalForceVectorType

Definition at line 104 of file LargeDeformationFEM.h.

Constructor & Destructor Documentation

◆ LargeDeformationLocalAssembler() [1/3]

template<typename ShapeFunction , int DisplacementDim>
ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::LargeDeformationLocalAssembler ( LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim > const & )
delete

◆ LargeDeformationLocalAssembler() [2/3]

template<typename ShapeFunction , int DisplacementDim>
ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::LargeDeformationLocalAssembler ( LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim > && )
delete

◆ LargeDeformationLocalAssembler() [3/3]

template<typename ShapeFunction , int DisplacementDim>
ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::LargeDeformationLocalAssembler ( MeshLib::Element const & e,
std::size_t const ,
NumLib::GenericIntegrationMethod const & integration_method,
bool const is_axially_symmetric,
LargeDeformationProcessData< DisplacementDim > & process_data )
inline

Definition at line 116 of file LargeDeformationFEM.h.

123 e, integration_method, is_axially_symmetric, process_data)
124 {
125 unsigned const n_integration_points =
127
128 _ip_data.resize(n_integration_points);
129 _secondary_data.N.resize(n_integration_points);
130
131 auto const shape_matrices =
133 DisplacementDim>(
134 e, is_axially_symmetric, this->integration_method_);
135
136 for (unsigned ip = 0; ip < n_integration_points; ip++)
137 {
138 auto& ip_data = _ip_data[ip];
139 auto const& sm = shape_matrices[ip];
140 _ip_data[ip].integration_weight =
142 sm.integralMeasure * sm.detJ;
143
144 ip_data.N_u = sm.N;
145 ip_data.dNdx_u = sm.dNdx;
146
147 _secondary_data.N[ip] = shape_matrices[ip].N;
148 }
149 }
constexpr double getWeight() const
MathLib::WeightedPoint const & getWeightedPoint(unsigned const igp) const
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data
std::vector< IpData, Eigen::aligned_allocator< IpData > > _ip_data
ShapeMatrixPolicyType< ShapeFunction, DisplacementDim > ShapeMatricesType
std::vector< typename ShapeMatricesType::ShapeMatrices, Eigen::aligned_allocator< typename ShapeMatricesType::ShapeMatrices > > initShapeMatrices(MeshLib::Element const &e, bool const is_axially_symmetric, IntegrationMethod const &integration_method)
LargeDeformationLocalAssemblerInterface(MeshLib::Element const &e, NumLib::GenericIntegrationMethod const &integration_method, bool const is_axially_symmetric, LargeDeformationProcessData< DisplacementDim > &process_data)
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N

References ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::_ip_data, ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::_secondary_data, NumLib::GenericIntegrationMethod::getNumberOfPoints(), MathLib::WeightedPoint::getWeight(), NumLib::GenericIntegrationMethod::getWeightedPoint(), NumLib::initShapeMatrices(), ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::integration_method_, and ProcessLib::LargeDeformation::SecondaryData< ShapeMatrixType >::N.

Member Function Documentation

◆ assemble()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< 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 278 of file LargeDeformationFEM.h.

284 {
285 OGS_FATAL(
286 "LargeDeformationLocalAssembler: assembly without jacobian is not "
287 "implemented.");
288 }
#define OGS_FATAL(...)
Definition Error.h:26

References OGS_FATAL.

◆ assembleWithJacobian()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::assembleWithJacobian ( double const t,
double const dt,
std::vector< double > const & local_x,
std::vector< double > const & local_x_prev,
std::vector< double > & local_b_data,
std::vector< double > & local_Jac_data )
inlineoverridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 317 of file LargeDeformationFEM.h.

322 {
323 auto const local_matrix_size = local_x.size();
324
326 local_Jac_data, local_matrix_size, local_matrix_size);
327
329 local_b_data, local_matrix_size);
330
331 auto [u] = localDOF(local_x);
332 auto [u_prev] = localDOF(local_x_prev);
333
334 unsigned const n_integration_points =
336
338 x_position.setElementID(this->element_.getID());
339
340 typename ConstitutiveRelations::ConstitutiveSetting<DisplacementDim>
341 constitutive_setting;
342 auto const& medium =
343 *this->process_data_.media_map.getMedium(this->element_.getID());
344
345 auto const f_bar_variables =
346 computeFBarVariables(/*compute_detF0_only?*/ false, u);
347
348 for (unsigned ip = 0; ip < n_integration_points; ip++)
349 {
350 auto const& w = _ip_data[ip].integration_weight;
351 auto const& N = _ip_data[ip].N_u;
352 auto const& dNdx = _ip_data[ip].dNdx_u;
353
354 auto const x_coord =
355 NumLib::interpolateXCoordinate<ShapeFunction,
357 this->element_, N);
358
359 // For the 2D case the 33-component is needed (and the four entries
360 // of the non-symmetric matrix); In 3d there are nine entries.
361 GradientMatrixType G(DisplacementDim * DisplacementDim +
362 (DisplacementDim == 2 ? 1 : 0),
363 ShapeFunction::NPOINTS * DisplacementDim);
364 Deformation::computeGMatrix<DisplacementDim,
365 ShapeFunction::NPOINTS>(
366 dNdx, G, this->is_axially_symmetric_, N, x_coord);
367
368 GradientVectorType const grad_u = G * u;
369
370 auto const B_0 = LinearBMatrix::computeBMatrix<
371 DisplacementDim, ShapeFunction::NPOINTS,
373 dNdx, N, x_coord, this->is_axially_symmetric_);
374
375 auto const B_N = NonLinearBMatrix::computeBMatrix<
376 DisplacementDim, ShapeFunction::NPOINTS,
378 dNdx, N, grad_u, x_coord, this->is_axially_symmetric_);
379
380 double const detF0 =
381 f_bar_variables ? std::get<double>(*f_bar_variables) : 0.0;
382 double const alpha = computeOutputStrainData(
383 detF0, B_0 + 0.5 * B_N, grad_u, u, this->output_data_[ip]);
384
385 auto const CD = updateConstitutiveRelations(
386 alpha, G, u_prev, x_position, t, dt, constitutive_setting,
387 medium, this->current_states_[ip], this->prev_states_[ip],
388 this->material_states_[ip], this->output_data_[ip]);
389
390 BMatrixType B = B_0 + B_N;
391
392 auto const& sigma = this->current_states_[ip].stress_data.sigma;
393 auto const& b = *CD.volumetric_body_force;
394 auto const& C = CD.s_mech_data.stiffness_tensor;
395
396 auto const sigma_geom =
399
400 if (!f_bar_variables)
401 {
402 local_b.noalias() -=
403 (B.transpose() * sigma - N_u_op(N).transpose() * b) * w;
404
405 local_Jac.noalias() += G.transpose() * sigma_geom * G * w;
406
407 local_Jac.noalias() += B.transpose() * C * B * w;
408 }
409 else // with F bar
410 {
411 double const alpha_p2 = alpha * alpha;
412 local_Jac.noalias() +=
413 alpha_p2 * G.transpose() * sigma_geom * G * w;
414
415 auto const q0 = std::get<VectorTypeForFbar>(*f_bar_variables);
416
417 auto const& F =
418 this->output_data_[ip]
419 .deformation_gradient_data.deformation_gradient;
421 DisplacementDim, ShapeFunction::NPOINTS, VectorTypeForFbar,
424 dNdx, F, this->is_axially_symmetric_);
425
426 auto const q0_q = q0 - q;
427
428 auto const S_q = sigma * (q0_q).transpose();
429
430 local_Jac.noalias() +=
431 2.0 * alpha_p2 * S_q.transpose() * B * w / DisplacementDim;
432
433 auto const twoEplsI =
435 alpha_p2, this->output_data_[ip].eps_data.eps,
436 this->is_axially_symmetric_);
437
438 // B bar:
439 BMatrixType const Bbar =
440 alpha_p2 *
441 (B + twoEplsI * (q0_q).transpose() / DisplacementDim);
442
443 local_Jac.noalias() +=
444 2.0 * Bbar.transpose() * S_q * w / DisplacementDim;
445
446 local_Jac.noalias() += Bbar.transpose() * C * Bbar * w;
447
448 local_Jac.noalias() -=
449 alpha_p2 * sigma.dot(twoEplsI) *
450 (q0 * q0.transpose() - q * q.transpose()) * w /
451 DisplacementDim;
452
453 local_b.noalias() -=
454 (Bbar.transpose() * sigma - N_u_op(N).transpose() * b) * w;
455 }
456 }
457 }
std::size_t getID() const
Returns the ID of the element.
Definition Element.h:89
void setElementID(std::size_t element_id)
MatrixType< _kelvin_vector_size, _number_of_dof > BMatrixType
std::optional< std::tuple< double, VectorTypeForFbar > > computeFBarVariables(bool const compute_detF0_only, Eigen::VectorXd const &local_x) const
typename GMatricesType::GradientVectorType GradientVectorType
double computeOutputStrainData(double const detF0, BMatrixType const &B, GradientVectorType const &grad_u, Eigen::Ref< Eigen::VectorXd const > const &u, typename ConstitutiveRelations::OutputData< DisplacementDim > &output_data) const
static constexpr auto localDOF(std::vector< double > const &x)
typename BMatricesType::NodalForceVectorType VectorTypeForFbar
ConstitutiveRelations::ConstitutiveData< DisplacementDim > updateConstitutiveRelations(double const alpha, GradientMatrixType const &G, Eigen::Ref< Eigen::VectorXd const > const &u_prev, ParameterLib::SpatialPosition const &x_position, double const t, double const dt, typename ConstitutiveRelations::ConstitutiveSetting< DisplacementDim > &CS, MaterialPropertyLib::Medium const &medium, typename ConstitutiveRelations::StatefulData< DisplacementDim > &current_state, typename ConstitutiveRelations::StatefulDataPrev< DisplacementDim > const &prev_state, MaterialStateData< DisplacementDim > &material_state, typename ConstitutiveRelations::OutputData< DisplacementDim > &output_data) const
typename GMatricesType::GradientMatrixType GradientMatrixType
Eigen::Matrix< double, 3, 3 > kelvinVectorToTensor(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)
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)
void computeGMatrix(DNDX_Type const &dNdx, GMatrixType &g_matrix, const bool is_axially_symmetric, N_Type const &N, const double radius)
Fills a G-matrix based on given shape function dN/dx values.
Definition GMatrix.h:25
Eigen::Matrix< double, MPL::tensorSize(DisplacementDim), MPL::tensorSize(DisplacementDim)> computeSigmaGeom(Eigen::Matrix3d const &sigma_tensor)
BMatrixType computeBMatrix(DNDX_Type const &dNdx, N_Type const &N, const double radius, const bool is_axially_symmetric)
Fills a B-matrix based on given shape function dN/dx values.
BMatrixType computeBMatrix(DNDX_Type const &dNdx, N_Type const &N, GradientVectorType const &grad_u, const double radius, const bool is_axially_symmetric)
VectorTypeForFbar computeQVector(DNDX_Type const &dNdx, GradientVectorType const &F, bool const)
MathLib::KelvinVector::KelvinVectorType< DisplacementDim > compute2EPlusI(double const alpha_p2, MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &eps_bar, bool const is_axially_symmetric)
std::vector< MaterialStateData< DisplacementDim > > material_states_
std::vector< typename ConstitutiveRelations::StatefulDataPrev< DisplacementDim > > prev_states_
std::vector< typename ConstitutiveRelations::StatefulData< DisplacementDim > > current_states_
std::vector< typename ConstitutiveRelations::OutputData< DisplacementDim > > output_data_

References ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::_ip_data, ProcessLib::NonLinearFbar::compute2EPlusI(), ProcessLib::LinearBMatrix::computeBMatrix(), ProcessLib::NonLinearBMatrix::computeBMatrix(), ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::computeFBarVariables(), ProcessLib::Deformation::computeGMatrix(), ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::computeOutputStrainData(), ProcessLib::NonLinearFbar::computeQVector(), ProcessLib::LargeDeformation::computeSigmaGeom(), MathLib::createZeroedMatrix(), MathLib::createZeroedVector(), ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::current_states_, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::element_, MeshLib::Element::getID(), NumLib::GenericIntegrationMethod::getNumberOfPoints(), ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::integration_method_, NumLib::interpolateXCoordinate(), ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::is_axially_symmetric_, MathLib::KelvinVector::kelvinVectorToTensor(), ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::localDOF(), ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::material_states_, ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::N_u_op, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::output_data_, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::prev_states_, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::process_data_, ParameterLib::SpatialPosition::setElementID(), and ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::updateConstitutiveRelations().

◆ computeFBarVariables()

template<typename ShapeFunction , int DisplacementDim>
std::optional< std::tuple< double, VectorTypeForFbar > > ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::computeFBarVariables ( bool const compute_detF0_only,
Eigen::VectorXd const & local_x ) const
inline

Definition at line 290 of file LargeDeformationFEM.h.

292 {
293 if (this->process_data_.bar_det_f_type ==
295 {
296 return std::nullopt;
297 }
298
299 if (this->process_data_.bar_det_f_type ==
301 {
303 DisplacementDim, GradientVectorType, VectorTypeForFbar,
304 NodalVectorType, ShapeFunction, ShapeMatricesType, IpData>(
305 _ip_data, compute_detF0_only, local_x,
306 this->integration_method_, this->element_,
308 }
309
311 DisplacementDim, GradientVectorType, GradientMatrixType,
312 VectorTypeForFbar, ShapeFunction, ShapeMatricesType>(
313 compute_detF0_only, local_x, this->element_,
315 }
IntegrationPointData< BMatricesType, ShapeMatricesType, DisplacementDim > IpData
typename ShapeMatricesType::NodalVectorType NodalVectorType
std::tuple< double, VectorTypeForFbar > computeFBarInitialVariablesElementCenter(bool const compute_detF0_only, Eigen::VectorXd const &u, MeshLib::Element const &element, bool const is_axially_symmetric)
std::tuple< double, VectorTypeForFbar > computeFBarInitialVariablesAverage(std::vector< IpData, Eigen::aligned_allocator< IpData > > const &ip_data, bool const compute_detF0_only, Eigen::VectorXd const &u, NumLib::GenericIntegrationMethod const &integration_method, MeshLib::Element const &element, bool const is_axially_symmetric)

References ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::_ip_data, ProcessLib::NonLinearFbar::computeFBarInitialVariablesAverage(), ProcessLib::NonLinearFbar::computeFBarInitialVariablesElementCenter(), ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::element_, ProcessLib::NonLinearFbar::ELEMENT_AVERAGE, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::integration_method_, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::is_axially_symmetric_, ProcessLib::NonLinearFbar::NONE, and ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::process_data_.

Referenced by ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::assembleWithJacobian(), and ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::postTimestepConcrete().

◆ computeOutputStrainData()

template<typename ShapeFunction , int DisplacementDim>
double ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::computeOutputStrainData ( double const detF0,
BMatrixType const & B,
GradientVectorType const & grad_u,
Eigen::Ref< Eigen::VectorXd const > const & u,
typename ConstitutiveRelations::OutputData< DisplacementDim > & output_data ) const
inline

Definition at line 187 of file LargeDeformationFEM.h.

195 {
196 // Note: Here B=B_{linear}+B_{nonlinear}/2 For Green-Lagrange strain.
197 output_data.eps_data.eps = B * u;
198 output_data.deformation_gradient_data.deformation_gradient =
200 output_data.deformation_gradient_data.volume_ratio =
202 output_data.deformation_gradient_data.deformation_gradient);
203
204 if (this->process_data_.bar_det_f_type ==
206 {
207 return 1.0;
208 }
209
210 double const detF = output_data.deformation_gradient_data.volume_ratio;
211 double const J_ratio = detF0 / detF;
212
213 if (J_ratio < .0)
214 {
215 OGS_FATAL(
216 "det(F0)/det(F) is negative with det(F0) = {:g}, and det(F) = "
217 "{:g} ",
218 detF0, detF);
219 }
220
221 double const alpha =
222 DisplacementDim == 3 ? std::cbrt(J_ratio) : std::sqrt(J_ratio);
223
224 output_data.deformation_gradient_data.deformation_gradient
225 .template segment<DisplacementDim * DisplacementDim>(0) *= alpha;
226 output_data.deformation_gradient_data.volume_ratio *=
227 std::pow(alpha, DisplacementDim);
228
229 double const alpha_p2 = alpha * alpha;
230 output_data.eps_data.eps =
231 alpha_p2 * output_data.eps_data.eps +
232 0.5 * (alpha_p2 - 1) *
235 return alpha;
236 }
double determinant(Eigen::MatrixBase< Derived > const &tensor)
Computes determinant of a vectorized tensor.
MathLib::KelvinVector::KelvinVectorType< DisplacementDim > identityForF(bool const is_axially_symmetric)

References ProcessLib::LargeDeformation::ConstitutiveRelations::OutputData< DisplacementDim >::deformation_gradient_data, MathLib::VectorizedTensor::determinant(), ProcessLib::LargeDeformation::ConstitutiveRelations::OutputData< DisplacementDim >::eps_data, MathLib::VectorizedTensor::identity(), ProcessLib::NonLinearFbar::identityForF(), ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::is_axially_symmetric_, ProcessLib::NonLinearFbar::NONE, OGS_FATAL, and ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::process_data_.

Referenced by ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::assembleWithJacobian(), and ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::postTimestepConcrete().

◆ getShapeMatrix()

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

Provides the shape matrix at the given integration point.

Implements NumLib::ExtrapolatableElement.

Definition at line 530 of file LargeDeformationFEM.h.

532 {
533 auto const& N = _secondary_data.N[integration_point];
534
535 // assumes N is stored contiguously in memory
536 return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
537 }

References ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::_secondary_data, and ProcessLib::LargeDeformation::SecondaryData< ShapeMatrixType >::N.

◆ initializeConcrete()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::initializeConcrete ( )
inlineoverridevirtual

Set initial stress from parameter.

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 151 of file LargeDeformationFEM.h.

152 {
153 unsigned const n_integration_points =
155 for (unsigned ip = 0; ip < n_integration_points; ip++)
156 {
157 auto const& ip_data = _ip_data[ip];
158
159 ParameterLib::SpatialPosition const x_position{
160 std::nullopt, this->element_.getID(),
162 NumLib::interpolateCoordinates<ShapeFunction,
164 this->element_, ip_data.N_u))};
165
167 if (this->process_data_.initial_stress != nullptr)
168 {
169 this->current_states_[ip].stress_data.sigma.noalias() =
171 DisplacementDim>((*this->process_data_.initial_stress)(
172 std::numeric_limits<
173 double>::quiet_NaN() /* time independent */,
174 x_position));
175 }
176
177 double const t = 0; // TODO (naumov) pass t from top
178 auto& material_state = this->material_states_[ip];
179 this->solid_material_.initializeInternalStateVariables(
180 t, x_position, *material_state.material_state_variables);
181
182 material_state.pushBackState();
183 this->prev_states_[ip] = this->current_states_[ip];
184 }
185 }
Eigen::Matrix< double, Eigen::MatrixBase< Derived >::RowsAtCompileTime, 1 > symmetricTensorToKelvinVector(Eigen::MatrixBase< Derived > const &v)
std::array< double, 3 > interpolateCoordinates(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)
ConstitutiveRelations::SolidConstitutiveRelation< DisplacementDim > const & solid_material_

References ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::_ip_data, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::current_states_, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::element_, MeshLib::Element::getID(), NumLib::GenericIntegrationMethod::getNumberOfPoints(), ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::integration_method_, NumLib::interpolateCoordinates(), ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::material_states_, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::prev_states_, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::process_data_, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::solid_material_, and MathLib::KelvinVector::symmetricTensorToKelvinVector().

◆ localDOF()

template<typename ShapeFunction , int DisplacementDim>
static constexpr auto ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::localDOF ( std::vector< double > const & x)
inlinestaticconstexprprivate

◆ postTimestepConcrete()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::postTimestepConcrete ( Eigen::VectorXd const & local_x,
Eigen::VectorXd const & local_x_prev,
double const t,
double const dt,
int const  )
inlineoverridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 459 of file LargeDeformationFEM.h.

463 {
464 unsigned const n_integration_points =
466
468 x_position.setElementID(this->element_.getID());
469
470 typename ConstitutiveRelations::ConstitutiveSetting<DisplacementDim>
471 constitutive_setting;
472
473 auto& medium =
474 *this->process_data_.media_map.getMedium(this->element_.getID());
475
476 auto const f_bar_variables =
477 computeFBarVariables(/*compute_detF0_only?*/ true, local_x);
478
479 for (unsigned ip = 0; ip < n_integration_points; ip++)
480 {
481 auto const& N = _ip_data[ip].N_u;
482 auto const& dNdx = _ip_data[ip].dNdx_u;
483 auto const x_coord =
484 NumLib::interpolateXCoordinate<ShapeFunction,
486 this->element_, N);
487
488 // For the 2D case the 33-component is needed (and the four entries
489 // of the non-symmetric matrix); In 3d there are nine entries.
490 GradientMatrixType G(DisplacementDim * DisplacementDim +
491 (DisplacementDim == 2 ? 1 : 0),
492 ShapeFunction::NPOINTS * DisplacementDim);
493 Deformation::computeGMatrix<DisplacementDim,
494 ShapeFunction::NPOINTS>(
495 dNdx, G, this->is_axially_symmetric_, N, x_coord);
496
497 GradientVectorType const grad_u = G * local_x;
498
499 // B for Green-Lagrange strain.
501 DisplacementDim, ShapeFunction::NPOINTS,
503 dNdx, N, x_coord, this->is_axially_symmetric_);
504
505 B.noalias() += 0.5 * NonLinearBMatrix::computeBMatrix<
506 DisplacementDim, ShapeFunction::NPOINTS,
508 dNdx, N, grad_u, x_coord,
510
511 double const detF0 =
512 f_bar_variables ? std::get<double>(*f_bar_variables) : 0.0;
513 double const alpha = computeOutputStrainData(
514 detF0, B, grad_u, local_x, this->output_data_[ip]);
515
517 alpha, G, local_x_prev, x_position, t, dt, constitutive_setting,
518 medium, this->current_states_[ip], this->prev_states_[ip],
519 this->material_states_[ip], this->output_data_[ip]);
520
521 this->material_states_[ip].pushBackState();
522 }
523
524 for (unsigned ip = 0; ip < n_integration_points; ip++)
525 {
526 this->prev_states_[ip] = this->current_states_[ip];
527 }
528 }

References ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::_ip_data, ProcessLib::LinearBMatrix::computeBMatrix(), ProcessLib::NonLinearBMatrix::computeBMatrix(), ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::computeFBarVariables(), ProcessLib::Deformation::computeGMatrix(), ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::computeOutputStrainData(), ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::current_states_, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::element_, MeshLib::Element::getID(), NumLib::GenericIntegrationMethod::getNumberOfPoints(), ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::integration_method_, NumLib::interpolateXCoordinate(), ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::is_axially_symmetric_, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::material_states_, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::output_data_, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::prev_states_, ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::process_data_, ParameterLib::SpatialPosition::setElementID(), and ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::updateConstitutiveRelations().

◆ updateConstitutiveRelations()

template<typename ShapeFunction , int DisplacementDim>
ConstitutiveRelations::ConstitutiveData< DisplacementDim > ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::updateConstitutiveRelations ( double const alpha,
GradientMatrixType const & G,
Eigen::Ref< Eigen::VectorXd const > const & u_prev,
ParameterLib::SpatialPosition const & x_position,
double const t,
double const dt,
typename ConstitutiveRelations::ConstitutiveSetting< DisplacementDim > & CS,
MaterialPropertyLib::Medium const & medium,
typename ConstitutiveRelations::StatefulData< DisplacementDim > & current_state,
typename ConstitutiveRelations::StatefulDataPrev< DisplacementDim > const & prev_state,
MaterialStateData< DisplacementDim > & material_state,
typename ConstitutiveRelations::OutputData< DisplacementDim > & output_data ) const
inline

Definition at line 239 of file LargeDeformationFEM.h.

254 {
255 double const T_ref =
256 this->process_data_.reference_temperature
257 ? (*this->process_data_.reference_temperature)(t, x_position)[0]
258 : std::numeric_limits<double>::quiet_NaN();
259
260 typename ConstitutiveRelations::ConstitutiveModels<DisplacementDim>
261 models(this->process_data_, this->solid_material_);
262 typename ConstitutiveRelations::ConstitutiveTempData<DisplacementDim>
263 tmp;
264 typename ConstitutiveRelations::ConstitutiveData<DisplacementDim> CD;
265
266 CS.eval(
267 models, t, dt, x_position, //
268 medium, //
269 T_ref, //
270 output_data.deformation_gradient_data, //
271 alpha * (G * u_prev +
273 current_state, prev_state, material_state, tmp, CD);
274
275 return CD;
276 }

References ProcessLib::LargeDeformation::ConstitutiveRelations::OutputData< DisplacementDim >::deformation_gradient_data, ProcessLib::LargeDeformation::ConstitutiveRelations::ConstitutiveSetting< DisplacementDim >::eval(), MathLib::VectorizedTensor::identity(), ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::process_data_, and ProcessLib::LargeDeformation::LargeDeformationLocalAssemblerInterface< DisplacementDim >::solid_material_.

Referenced by ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::assembleWithJacobian(), and ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::postTimestepConcrete().

Member Data Documentation

◆ _ip_data

◆ _secondary_data

◆ displacement_size

template<typename ShapeFunction , int DisplacementDim>
const int ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::displacement_size
staticprivate
Initial value:
=
ShapeFunction::NPOINTS * DisplacementDim

Definition at line 551 of file LargeDeformationFEM.h.

◆ N_u_op

template<typename ShapeFunction , int DisplacementDim>
auto& ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::N_u_op
staticconstexpr
Initial value:
DisplacementDim, typename ShapeMatricesType::NodalRowVectorType>
constexpr Eigen::CwiseNullaryOp< EigenBlockMatrixViewFunctor< D, M >, typename EigenBlockMatrixViewFunctor< D, M >::Matrix > eigenBlockMatrixView(const Eigen::MatrixBase< M > &matrix)
RowVectorType< ShapeFunction::NPOINTS > NodalRowVectorType

Definition at line 109 of file LargeDeformationFEM.h.

Referenced by ProcessLib::LargeDeformation::LargeDeformationLocalAssembler< ShapeFunction, DisplacementDim >::assembleWithJacobian().


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