33 namespace HydroMechanics
37 template <
typename BMatricesType,
typename ShapeMatricesTypeDisplacement,
38 typename ShapeMatricesTypePressure,
int DisplacementDim,
int NPoints>
50 typename ShapeMatricesTypeDisplacement::template MatrixType<
51 DisplacementDim, NPoints * DisplacementDim>
57 typename ShapeMatricesTypeDisplacement::NodalRowVectorType
N_u;
58 typename ShapeMatricesTypeDisplacement::GlobalDimNodalMatrixType
dNdx_u;
60 typename ShapeMatricesTypePressure::NodalRowVectorType
N_p;
61 typename ShapeMatricesTypePressure::GlobalDimNodalMatrixType
dNdx_p;
65 DisplacementDim>::MaterialStateVariables>
85 double const temperature)
92 auto const null_state =
solid_material.createMaterialStateVariables();
96 variable_array[
static_cast<int>(MPL::Variable::stress)].emplace<KV>(
98 variable_array[
static_cast<int>(MPL::Variable::mechanical_strain)]
99 .emplace<KV>(KV::Zero());
101 .emplace<double>(temperature);
103 variable_array_prev[
static_cast<int>(MPL::Variable::stress)]
104 .emplace<KV>(KV::Zero());
105 variable_array_prev[
static_cast<int>(MPL::Variable::mechanical_strain)]
106 .emplace<KV>(KV::Zero());
108 .emplace<double>(temperature);
111 solid_material.integrateStress(variable_array_prev, variable_array,
112 t, x_position, dt, *null_state);
116 OGS_FATAL(
"Computation of elastic tangent stiffness failed.");
120 std::move(std::get<2>(*solution));
125 template <
typename DisplacementVectorType>
131 DisplacementVectorType
const& ,
135 variable_array_prev[
static_cast<int>(
143 variable_array_prev[
static_cast<int>(
148 variable_array_prev, variable_array, t, x_position, dt,
153 OGS_FATAL(
"Computation of local constitutive relation failed.");
167 template <
typename ShapeMatrixType>
170 std::vector<ShapeMatrixType, Eigen::aligned_allocator<ShapeMatrixType>>
N_u;
173 template <
typename ShapeFunctionDisplacement,
typename ShapeFunctionPressure,
174 typename IntegrationMethod,
int DisplacementDim>
198 bool const is_axially_symmetric,
199 unsigned const integration_order,
204 std::string
const&
name,
205 double const* values,
206 int const integration_order)
override;
209 std::vector<double>
const& ,
210 std::vector<double>
const& ,
211 std::vector<double>& ,
212 std::vector<double>& ,
213 std::vector<double>& )
override
216 "HydroMechanicsLocalAssembler: assembly without jacobian is not "
221 std::vector<double>
const& local_x,
222 std::vector<double>
const& local_xdot,
223 const double ,
const double ,
224 std::vector<double>& ,
225 std::vector<double>& ,
226 std::vector<double>& local_rhs_data,
227 std::vector<double>& local_Jac_data)
override;
230 const double t,
double const dt, Eigen::VectorXd
const& local_x,
231 Eigen::VectorXd
const& local_xdot,
const double dxdot_dx,
232 const double dx_dx,
int const process_id,
233 std::vector<double>& local_M_data, std::vector<double>& local_K_data,
234 std::vector<double>& local_b_data,
235 std::vector<double>& local_Jac_data)
override;
239 unsigned const n_integration_points =
242 for (
unsigned ip = 0; ip < n_integration_points; ip++)
252 ShapeFunctionDisplacement,
260 double>::quiet_NaN() ,
264 ip_data.pushBackState();
270 double const )
override
272 unsigned const n_integration_points =
275 for (
unsigned ip = 0; ip < n_integration_points; ip++)
282 double const t,
double const dt, Eigen::VectorXd
const& local_xs,
283 Eigen::VectorXd
const& local_x_dot)
override;
286 std::vector<double>
const& local_xdot,
287 double const t,
double const dt,
288 bool const use_monolithic_scheme,
289 int const process_id)
override;
293 bool const use_monolithic_scheme,
294 int const process_id)
override;
297 const unsigned integration_point)
const override
302 return Eigen::Map<const Eigen::RowVectorXd>(N_u.data(), N_u.size());
308 std::vector<double>
getSigma()
const override;
310 std::vector<double>
getEpsilon()
const override;
314 std::vector<GlobalVector*>
const& x,
315 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
316 std::vector<double>& cache)
const override;
320 std::vector<GlobalVector*>
const& ,
321 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& ,
322 std::vector<double>& cache)
const override
324 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
330 std::vector<GlobalVector*>
const& ,
331 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& ,
332 std::vector<double>& cache)
const override
334 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
358 const double t,
double const dt, Eigen::VectorXd
const& local_x,
359 std::vector<double>& local_b_data, std::vector<double>& local_Jac_data);
385 const double t,
double const dt, Eigen::VectorXd
const& local_x,
386 Eigen::VectorXd
const& local_xdot, std::vector<double>& local_b_data,
387 std::vector<double>& local_Jac_data);
392 DisplacementDim>::MaterialStateVariables
const&
403 ShapeFunctionDisplacement::NPOINTS>;
404 std::vector<IpData, Eigen::aligned_allocator<IpData>>
_ip_data;
417 ShapeFunctionDisplacement::NPOINTS * DisplacementDim;
virtual std::size_t getID() const final
Returns the ID of the element.
unsigned getNumberOfIntegrationPoints() const override
void postNonLinearSolverConcrete(std::vector< double > const &local_x, std::vector< double > const &local_xdot, double const t, double const dt, bool const use_monolithic_scheme, int const process_id) override
void computeSecondaryVariableConcrete(double const t, double const dt, Eigen::VectorXd const &local_xs, Eigen::VectorXd const &local_x_dot) override
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix(const unsigned integration_point) const override
Provides the shape matrix at the given integration point.
void initializeConcrete() override
void postTimestepConcrete(Eigen::VectorXd const &, double const, double const) override
HydroMechanicsProcessData< DisplacementDim > & _process_data
static int const KelvinVectorSize
HydroMechanicsLocalAssembler(HydroMechanicsLocalAssembler &&)=delete
void assembleWithJacobianForDeformationEquations(const double t, double const dt, Eigen::VectorXd const &local_x, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
std::vector< double > getEpsilon() const override
static const int displacement_index
std::vector< IpData, Eigen::aligned_allocator< IpData > > _ip_data
ShapeMatrixPolicyType< ShapeFunctionPressure, DisplacementDim > ShapeMatricesTypePressure
Eigen::Matrix< double, KelvinVectorSize, 1 > SymmetricTensor
std::size_t setIPDataInitialConditions(std::string const &name, double const *values, int const integration_order) override
Returns number of read integration points.
bool const _is_axially_symmetric
std::vector< double > getSigma() const override
std::vector< double > const & getIntPtEpsilon(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
void assembleWithJacobianForStaggeredScheme(const double t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_xdot, const double dxdot_dx, const double dx_dx, int const process_id, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data) override
BMatrixPolicyType< ShapeFunctionDisplacement, DisplacementDim > BMatricesType
void assembleWithJacobianForPressureEquations(const double t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_xdot, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
static const int pressure_size
SecondaryData< typename ShapeMatricesTypeDisplacement::ShapeMatrices::ShapeType > _secondary_data
HydroMechanicsLocalAssembler(HydroMechanicsLocalAssembler const &)=delete
void assembleWithJacobian(double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_xdot, const double, const double, std::vector< double > &, std::vector< double > &, std::vector< double > &local_rhs_data, std::vector< double > &local_Jac_data) override
ShapeMatrixPolicyType< ShapeFunctionDisplacement, DisplacementDim > ShapeMatricesTypeDisplacement
static const int displacement_size
std::vector< double > const & getIntPtSigma(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
MeshLib::Element const & _element
void setInitialConditionsConcrete(std::vector< double > const &local_x, double const t, bool const use_monolithic_scheme, int const process_id) override
static const int pressure_index
IntegrationMethod _integration_method
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::vector< double > const & getIntPtDarcyVelocity(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables const & getMaterialStateVariablesAt(unsigned integration_point) const override
std::array< VariableType, static_cast< int >(Variable::number_of_variables)> VariableArray
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), kelvin_vector_dimensions(DisplacementDim), Eigen::RowMajor > KelvinMatrixType
Eigen::Matrix< double, Eigen::MatrixBase< Derived >::RowsAtCompileTime, 1 > symmetricTensorToKelvinVector(Eigen::MatrixBase< Derived > const &v)
MathLib::TemplatePoint< double, 3 > Point3d
std::array< double, 3 > interpolateCoordinates(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)
double integration_weight
EIGEN_MAKE_ALIGNED_OPERATOR_NEW
BMatricesType::KelvinMatrixType updateConstitutiveRelation(MPL::VariableArray const &variable_array, double const t, ParameterLib::SpatialPosition const &x_position, double const dt, DisplacementVectorType const &, double const T)
ShapeMatricesTypeDisplacement::GlobalDimNodalMatrixType dNdx_u
std::unique_ptr< typename MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables > material_state_variables
ShapeMatricesTypeDisplacement::NodalRowVectorType N_u
BMatricesType::KelvinVectorType eps
ShapeMatricesTypeDisplacement::template MatrixType< DisplacementDim, NPoints *DisplacementDim > N_u_op
ShapeMatricesTypePressure::GlobalDimNodalMatrixType dNdx_p
MathLib::KelvinVector::KelvinMatrixType< DisplacementDim > computeElasticTangentStiffness(double const t, ParameterLib::SpatialPosition const &x_position, double const dt, double const temperature)
ShapeMatricesTypePressure::NodalRowVectorType N_p
BMatricesType::KelvinVectorType sigma_eff_prev
BMatricesType::KelvinVectorType eps_prev
MaterialLib::Solids::MechanicsBase< DisplacementDim > const & solid_material
IntegrationPointData(MaterialLib::Solids::MechanicsBase< DisplacementDim > const &solid_material)
BMatricesType::KelvinVectorType sigma_eff
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N_u