OGS
ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim > Class Template Reference

Detailed Description

template<typename ShapeFunctionDisplacement, typename ShapeFunctionPressure, int GlobalDim>
class ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >

Definition at line 32 of file HydroMechanicsLocalAssemblerMatrix.h.

#include <HydroMechanicsLocalAssemblerMatrix.h>

Inheritance diagram for ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >:
[legend]
Collaboration diagram for ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >:
[legend]

Public Member Functions

 HydroMechanicsLocalAssemblerMatrix (HydroMechanicsLocalAssemblerMatrix const &)=delete
 
 HydroMechanicsLocalAssemblerMatrix (HydroMechanicsLocalAssemblerMatrix &&)=delete
 
 HydroMechanicsLocalAssemblerMatrix (MeshLib::Element const &e, std::size_t const n_variables, std::size_t const local_matrix_size, std::vector< unsigned > const &dofIndex_to_localIndex, NumLib::GenericIntegrationMethod const &integration_method, bool const is_axially_symmetric, HydroMechanicsProcessData< GlobalDim > &process_data)
 
void preTimestepConcrete (std::vector< double > const &, double const, double const) override
 
- Public Member Functions inherited from ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerInterface
 HydroMechanicsLocalAssemblerInterface (MeshLib::Element const &element, bool const is_axially_symmetric, std::size_t n_local_size, std::vector< unsigned > dofIndex_to_localIndex)
 
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 &local_xdot_, std::vector< double > &, std::vector< double > &, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data) override
 
void postTimestepConcrete (Eigen::VectorXd const &local_x_, const double t, double const dt) override
 
- Public Member Functions inherited from ProcessLib::LocalAssemblerInterface
virtual ~LocalAssemblerInterface ()=default
 
void setInitialConditions (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, double const t, bool const use_monolithic_scheme, 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_xdot, 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_xdot, 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)
 
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_dot, 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, double const t, double const dt)
 
void postNonLinearSolver (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, GlobalVector const &xdot, double const t, double const dt, bool const use_monolithic_scheme, 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. More...
 

Protected Types

using ShapeMatricesTypeDisplacement = ShapeMatrixPolicyType< ShapeFunctionDisplacement, GlobalDim >
 
using BMatricesType = BMatrixPolicyType< ShapeFunctionDisplacement, GlobalDim >
 
using ShapeMatricesTypePressure = ShapeMatrixPolicyType< ShapeFunctionPressure, GlobalDim >
 
using IntegrationPointDataType = IntegrationPointDataMatrix< BMatricesType, ShapeMatricesTypeDisplacement, ShapeMatricesTypePressure, GlobalDim, ShapeFunctionDisplacement::NPOINTS >
 
using GlobalDimVector = Eigen::Matrix< double, GlobalDim, 1 >
 

Protected Member Functions

void assembleWithJacobianConcrete (double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_dot, Eigen::VectorXd &local_rhs, Eigen::MatrixXd &local_Jac) override
 
void assembleBlockMatricesWithJacobian (double const t, double const dt, Eigen::Ref< const Eigen::VectorXd > const &p, Eigen::Ref< const Eigen::VectorXd > const &p_dot, Eigen::Ref< const Eigen::VectorXd > const &u, Eigen::Ref< const Eigen::VectorXd > const &u_dot, Eigen::Ref< Eigen::VectorXd > rhs_p, Eigen::Ref< Eigen::VectorXd > rhs_u, Eigen::Ref< Eigen::MatrixXd > J_pp, Eigen::Ref< Eigen::MatrixXd > J_pu, Eigen::Ref< Eigen::MatrixXd > J_uu, Eigen::Ref< Eigen::MatrixXd > J_up)
 
void postTimestepConcreteWithVector (double const t, double const dt, Eigen::VectorXd const &local_x) override
 
void postTimestepConcreteWithBlockVectors (double const t, double const dt, Eigen::Ref< const Eigen::VectorXd > const &p, Eigen::Ref< const Eigen::VectorXd > const &u)
 
void setPressureOfInactiveNodes (double const t, Eigen::Ref< Eigen::VectorXd > p)
 
void setPressureDotOfInactiveNodes (Eigen::Ref< Eigen::VectorXd > p_dot)
 

Protected Attributes

HydroMechanicsProcessData< GlobalDim > & _process_data
 
std::vector< IntegrationPointDataType, Eigen::aligned_allocator< IntegrationPointDataType > > _ip_data
 
- Protected Attributes inherited from ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerInterface
MeshLib::Element const & _element
 
bool const _is_axially_symmetric
 

Static Protected Attributes

static const int pressure_index = 0
 
static const int pressure_size = ShapeFunctionPressure::NPOINTS
 
static const int displacement_index = ShapeFunctionPressure::NPOINTS
 
static const int displacement_size
 
static const int kelvin_vector_size
 

Member Typedef Documentation

◆ BMatricesType

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
using ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::BMatricesType = BMatrixPolicyType<ShapeFunctionDisplacement, GlobalDim>
protected

Definition at line 93 of file HydroMechanicsLocalAssemblerMatrix.h.

◆ GlobalDimVector

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
using ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::GlobalDimVector = Eigen::Matrix<double, GlobalDim, 1>
protected

Definition at line 107 of file HydroMechanicsLocalAssemblerMatrix.h.

◆ IntegrationPointDataType

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
using ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::IntegrationPointDataType = IntegrationPointDataMatrix<BMatricesType, ShapeMatricesTypeDisplacement, ShapeMatricesTypePressure, GlobalDim, ShapeFunctionDisplacement::NPOINTS>
protected

Definition at line 100 of file HydroMechanicsLocalAssemblerMatrix.h.

◆ ShapeMatricesTypeDisplacement

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
using ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::ShapeMatricesTypeDisplacement = ShapeMatrixPolicyType<ShapeFunctionDisplacement, GlobalDim>
protected

Definition at line 91 of file HydroMechanicsLocalAssemblerMatrix.h.

◆ ShapeMatricesTypePressure

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
using ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::ShapeMatricesTypePressure = ShapeMatrixPolicyType<ShapeFunctionPressure, GlobalDim>
protected

Definition at line 97 of file HydroMechanicsLocalAssemblerMatrix.h.

Constructor & Destructor Documentation

◆ HydroMechanicsLocalAssemblerMatrix() [1/3]

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::HydroMechanicsLocalAssemblerMatrix ( HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim > const &  )
delete

◆ HydroMechanicsLocalAssemblerMatrix() [2/3]

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::HydroMechanicsLocalAssemblerMatrix ( HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim > &&  )
delete

◆ HydroMechanicsLocalAssemblerMatrix() [3/3]

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::HydroMechanicsLocalAssemblerMatrix ( MeshLib::Element const &  e,
std::size_t const  n_variables,
std::size_t const  local_matrix_size,
std::vector< unsigned > const &  dofIndex_to_localIndex,
NumLib::GenericIntegrationMethod const &  integration_method,
bool const  is_axially_symmetric,
HydroMechanicsProcessData< GlobalDim > &  process_data 
)

Definition at line 31 of file HydroMechanicsLocalAssemblerMatrix-impl.h.

41  e, is_axially_symmetric,
42  (n_variables - 1) * ShapeFunctionDisplacement::NPOINTS * GlobalDim +
43  ShapeFunctionPressure::NPOINTS,
44  dofIndex_to_localIndex),
45  _process_data(process_data)
46 {
47  unsigned const n_integration_points =
48  integration_method.getNumberOfPoints();
49 
50  _ip_data.reserve(n_integration_points);
51 
52  auto const shape_matrices_u =
53  NumLib::initShapeMatrices<ShapeFunctionDisplacement,
55  e, is_axially_symmetric, integration_method);
56 
57  auto const shape_matrices_p =
58  NumLib::initShapeMatrices<ShapeFunctionPressure,
59  ShapeMatricesTypePressure, GlobalDim>(
60  e, is_axially_symmetric, integration_method);
61 
63  _process_data.solid_materials, _process_data.material_ids, e.getID());
64 
66  x_position.setElementID(e.getID());
67  for (unsigned ip = 0; ip < n_integration_points; ip++)
68  {
69  x_position.setIntegrationPoint(ip);
70 
71  _ip_data.emplace_back(solid_material);
72  auto& ip_data = _ip_data[ip];
73  auto const& sm_u = shape_matrices_u[ip];
74  auto const& sm_p = shape_matrices_p[ip];
75  ip_data.integration_weight =
76  sm_u.detJ * sm_u.integralMeasure *
77  integration_method.getWeightedPoint(ip).getWeight();
78  ip_data.darcy_velocity.setZero();
79 
80  ip_data.N_u = sm_u.N;
81  ip_data.dNdx_u = sm_u.dNdx;
82  ip_data.H_u.setZero(GlobalDim, displacement_size);
83  for (int i = 0; i < GlobalDim; ++i)
84  {
85  ip_data.H_u
86  .template block<1, displacement_size / GlobalDim>(
87  i, i * displacement_size / GlobalDim)
88  .noalias() = ip_data.N_u;
89  }
90 
91  ip_data.N_p = sm_p.N;
92  ip_data.dNdx_p = sm_p.dNdx;
93 
94  ip_data.sigma_eff.setZero(kelvin_vector_size);
95  ip_data.eps.setZero(kelvin_vector_size);
96 
97  ip_data.sigma_eff_prev.resize(kelvin_vector_size);
98  ip_data.eps_prev.resize(kelvin_vector_size);
99  ip_data.C.resize(kelvin_vector_size, kelvin_vector_size);
100 
101  auto const initial_effective_stress =
102  _process_data.initial_effective_stress(0, x_position);
103  for (unsigned i = 0; i < kelvin_vector_size; i++)
104  {
105  ip_data.sigma_eff[i] = initial_effective_stress[i];
106  ip_data.sigma_eff_prev[i] = initial_effective_stress[i];
107  }
108  }
109 }
void setElementID(std::size_t element_id)
void setIntegrationPoint(unsigned integration_point)
HydroMechanicsLocalAssemblerInterface(MeshLib::Element const &element, bool const is_axially_symmetric, std::size_t n_local_size, std::vector< unsigned > dofIndex_to_localIndex)
ShapeMatrixPolicyType< ShapeFunctionDisplacement, GlobalDim > ShapeMatricesTypeDisplacement
std::vector< IntegrationPointDataType, Eigen::aligned_allocator< IntegrationPointDataType > > _ip_data
ShapeMatrixPolicyType< ShapeFunctionPressure, GlobalDim > ShapeMatricesTypePressure
MechanicsBase< DisplacementDim > & selectSolidConstitutiveRelation(std::map< int, std::unique_ptr< MechanicsBase< DisplacementDim >>> const &constitutive_relations, MeshLib::PropertyVector< int > const *const material_ids, std::size_t const element_id)
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)

References ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::_ip_data, ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::_process_data, ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::displacement_size, MeshLib::Element::getID(), NumLib::GenericIntegrationMethod::getNumberOfPoints(), MathLib::WeightedPoint::getWeight(), NumLib::GenericIntegrationMethod::getWeightedPoint(), NumLib::initShapeMatrices(), ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::kelvin_vector_size, MaterialLib::Solids::selectSolidConstitutiveRelation(), ParameterLib::SpatialPosition::setElementID(), and ParameterLib::SpatialPosition::setIntegrationPoint().

Member Function Documentation

◆ assembleBlockMatricesWithJacobian()

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
void ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::assembleBlockMatricesWithJacobian ( double const  t,
double const  dt,
Eigen::Ref< const Eigen::VectorXd > const &  p,
Eigen::Ref< const Eigen::VectorXd > const &  p_dot,
Eigen::Ref< const Eigen::VectorXd > const &  u,
Eigen::Ref< const Eigen::VectorXd > const &  u_dot,
Eigen::Ref< Eigen::VectorXd >  rhs_p,
Eigen::Ref< Eigen::VectorXd >  rhs_u,
Eigen::Ref< Eigen::MatrixXd >  J_pp,
Eigen::Ref< Eigen::MatrixXd >  J_pu,
Eigen::Ref< Eigen::MatrixXd >  J_uu,
Eigen::Ref< Eigen::MatrixXd >  J_up 
)
protected

Definition at line 155 of file HydroMechanicsLocalAssemblerMatrix-impl.h.

165 {
166  assert(this->_element.getDimension() == GlobalDim);
167 
168  typename ShapeMatricesTypePressure::NodalMatrixType laplace_p =
169  ShapeMatricesTypePressure::NodalMatrixType::Zero(pressure_size,
170  pressure_size);
171 
172  typename ShapeMatricesTypePressure::NodalMatrixType storage_p =
173  ShapeMatricesTypePressure::NodalMatrixType::Zero(pressure_size,
174  pressure_size);
175 
176  typename ShapeMatricesTypeDisplacement::template MatrixType<
178  Kup = ShapeMatricesTypeDisplacement::template MatrixType<
180  pressure_size);
181 
182  auto const& gravity_vec = _process_data.specific_body_force;
183 
184  MPL::VariableArray variables;
185  MPL::VariableArray variables_prev;
187  x_position.setElementID(_element.getID());
188 
189  unsigned const n_integration_points = _ip_data.size();
190  for (unsigned ip = 0; ip < n_integration_points; ip++)
191  {
192  x_position.setIntegrationPoint(ip);
193 
194  auto& ip_data = _ip_data[ip];
195  auto const& ip_w = ip_data.integration_weight;
196  auto const& N_u = ip_data.N_u;
197  auto const& dNdx_u = ip_data.dNdx_u;
198  auto const& N_p = ip_data.N_p;
199  auto const& dNdx_p = ip_data.dNdx_p;
200  auto const& H_u = ip_data.H_u;
201 
202  auto const x_coord =
203  NumLib::interpolateXCoordinate<ShapeFunctionDisplacement,
205  _element, N_u);
206  auto const B =
208  ShapeFunctionDisplacement::NPOINTS,
209  typename BMatricesType::BMatrixType>(
210  dNdx_u, N_u, x_coord, _is_axially_symmetric);
211 
212  auto const& eps_prev = ip_data.eps_prev;
213  auto const& sigma_eff_prev = ip_data.sigma_eff_prev;
214  auto& sigma_eff = ip_data.sigma_eff;
215 
216  auto& eps = ip_data.eps;
217  auto& state = ip_data.material_state_variables;
218 
219  auto const alpha = _process_data.biot_coefficient(t, x_position)[0];
220  auto const rho_sr = _process_data.solid_density(t, x_position)[0];
221  auto const rho_fr = _process_data.fluid_density(t, x_position)[0];
222  auto const porosity = _process_data.porosity(t, x_position)[0];
223 
224  double const rho = rho_sr * (1 - porosity) + porosity * rho_fr;
225  auto const& identity2 =
227 
228  eps.noalias() = B * u;
229 
230  variables.mechanical_strain
232 
233  variables_prev.stress
235  sigma_eff_prev);
236  variables_prev.mechanical_strain
238  eps_prev);
239  variables_prev.temperature = _process_data.reference_temperature;
240 
241  auto&& solution = _ip_data[ip].solid_material.integrateStress(
242  variables_prev, variables, t, x_position, dt, *state);
243 
244  if (!solution)
245  {
246  OGS_FATAL("Computation of local constitutive relation failed.");
247  }
248 
250  std::tie(sigma_eff, state, C) = std::move(*solution);
251 
252  J_uu.noalias() += B.transpose() * C * B * ip_w;
253 
254  rhs_u.noalias() -= B.transpose() * sigma_eff * ip_w;
255  rhs_u.noalias() -= -H_u.transpose() * rho * gravity_vec * ip_w;
256 
257  //
258  // pressure equation, pressure part and displacement equation, pressure
259  // part
260  //
261  if (!_process_data.deactivate_matrix_in_flow) // Only for hydraulically
262  // active matrix
263  {
264  Kup.noalias() += B.transpose() * alpha * identity2 * N_p * ip_w;
265 
266  double const k_over_mu =
267  _process_data.intrinsic_permeability(t, x_position)[0] /
268  _process_data.fluid_viscosity(t, x_position)[0];
269  double const S = _process_data.specific_storage(t, x_position)[0];
270 
271  auto q = ip_data.darcy_velocity.head(GlobalDim);
272  q.noalias() = -k_over_mu * (dNdx_p * p + rho_fr * gravity_vec);
273 
274  laplace_p.noalias() +=
275  dNdx_p.transpose() * k_over_mu * dNdx_p * ip_w;
276  storage_p.noalias() += N_p.transpose() * S * N_p * ip_w;
277 
278  rhs_p.noalias() +=
279  dNdx_p.transpose() * rho_fr * k_over_mu * gravity_vec * ip_w;
280  }
281  }
282 
283  // displacement equation, pressure part
284  J_up.noalias() -= Kup;
285 
286  // pressure equation, pressure part.
287  J_pp.noalias() += laplace_p + storage_p / dt;
288 
289  // pressure equation, displacement part.
290  J_pu.noalias() += Kup.transpose() / dt;
291 
292  // pressure equation
293  rhs_p.noalias() -=
294  laplace_p * p + storage_p * p_dot + Kup.transpose() * u_dot;
295 
296  // displacement equation
297  rhs_u.noalias() -= -Kup * p;
298 }
#define OGS_FATAL(...)
Definition: Error.h:26
std::variant< std::monostate, Eigen::Matrix< double, 4, 1 >, Eigen::Matrix< double, 6, 1 > > mechanical_strain
Definition: VariableType.h:165
std::variant< std::monostate, Eigen::Matrix< double, 4, 1 >, Eigen::Matrix< double, 6, 1 > > stress
Definition: VariableType.h:175
virtual constexpr unsigned getDimension() const =0
Get dimension of the mesh element.
virtual std::size_t getID() const final
Returns the ID of the element.
Definition: Element.h:89
MatrixType< _kelvin_vector_size, _number_of_dof > BMatrixType
Definition: BMatrixPolicy.h:55
@ rho
density. For some models, rho substitutes p
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
Definition: KelvinVector.h:49
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), kelvin_vector_dimensions(DisplacementDim), Eigen::RowMajor > KelvinMatrixType
Definition: KelvinVector.h:57
static const double q
static const double u
static const double t
double interpolateXCoordinate(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)
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.
Definition: LinearBMatrix.h:42
MatrixType< ShapeFunction::NPOINTS, ShapeFunction::NPOINTS > NodalMatrixType

References MaterialPropertyLib::alpha, ProcessLib::LinearBMatrix::computeBMatrix(), NumLib::interpolateXCoordinate(), MaterialPropertyLib::VariableArray::mechanical_strain, OGS_FATAL, MathLib::q, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), MaterialPropertyLib::VariableArray::stress, MathLib::t, MaterialPropertyLib::VariableArray::temperature, and MathLib::u.

◆ assembleWithJacobianConcrete()

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
void ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::assembleWithJacobianConcrete ( double const  t,
double const  dt,
Eigen::VectorXd const &  local_x,
Eigen::VectorXd const &  local_x_dot,
Eigen::VectorXd &  local_rhs,
Eigen::MatrixXd &  local_Jac 
)
overrideprotectedvirtual

Implements ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerInterface.

Reimplemented in ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrixNearFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >.

Definition at line 115 of file HydroMechanicsLocalAssemblerMatrix-impl.h.

120 {
121  auto p = const_cast<Eigen::VectorXd&>(local_x).segment(pressure_index,
122  pressure_size);
123  auto p_dot = const_cast<Eigen::VectorXd&>(local_x_dot)
124  .segment(pressure_index, pressure_size);
125 
126  if (_process_data.deactivate_matrix_in_flow)
127  {
130  }
131 
132  auto u = local_x.segment(displacement_index, displacement_size);
133  auto u_dot = local_x_dot.segment(displacement_index, displacement_size);
134 
135  auto rhs_p = local_rhs.template segment<pressure_size>(pressure_index);
136  auto rhs_u =
137  local_rhs.template segment<displacement_size>(displacement_index);
138 
139  auto J_pp = local_Jac.template block<pressure_size, pressure_size>(
141  auto J_pu = local_Jac.template block<pressure_size, displacement_size>(
143  auto J_uu = local_Jac.template block<displacement_size, displacement_size>(
145  auto J_up = local_Jac.template block<displacement_size, pressure_size>(
147 
148  assembleBlockMatricesWithJacobian(t, dt, p, p_dot, u, u_dot, rhs_p, rhs_u,
149  J_pp, J_pu, J_uu, J_up);
150 }
void assembleBlockMatricesWithJacobian(double const t, double const dt, Eigen::Ref< const Eigen::VectorXd > const &p, Eigen::Ref< const Eigen::VectorXd > const &p_dot, Eigen::Ref< const Eigen::VectorXd > const &u, Eigen::Ref< const Eigen::VectorXd > const &u_dot, Eigen::Ref< Eigen::VectorXd > rhs_p, Eigen::Ref< Eigen::VectorXd > rhs_u, Eigen::Ref< Eigen::MatrixXd > J_pp, Eigen::Ref< Eigen::MatrixXd > J_pu, Eigen::Ref< Eigen::MatrixXd > J_uu, Eigen::Ref< Eigen::MatrixXd > J_up)

References MathLib::t, and MathLib::u.

◆ postTimestepConcreteWithBlockVectors()

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
void ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::postTimestepConcreteWithBlockVectors ( double const  t,
double const  dt,
Eigen::Ref< const Eigen::VectorXd > const &  p,
Eigen::Ref< const Eigen::VectorXd > const &  u 
)
protected

Definition at line 321 of file HydroMechanicsLocalAssemblerMatrix-impl.h.

326 {
327  MPL::VariableArray variables;
328  MPL::VariableArray variables_prev;
330  x_position.setElementID(_element.getID());
331 
332  unsigned const n_integration_points = _ip_data.size();
333  for (unsigned ip = 0; ip < n_integration_points; ip++)
334  {
335  x_position.setIntegrationPoint(ip);
336 
337  auto& ip_data = _ip_data[ip];
338 
339  auto const& eps_prev = ip_data.eps_prev;
340  auto const& sigma_eff_prev = ip_data.sigma_eff_prev;
341 
342  auto& eps = ip_data.eps;
343  auto& sigma_eff = ip_data.sigma_eff;
344  auto& state = ip_data.material_state_variables;
345 
346  auto const& N_u = ip_data.N_u;
347  auto const& dNdx_u = ip_data.dNdx_u;
348 
349  auto const x_coord =
350  NumLib::interpolateXCoordinate<ShapeFunctionDisplacement,
352  _element, N_u);
353  auto const B =
355  ShapeFunctionDisplacement::NPOINTS,
356  typename BMatricesType::BMatrixType>(
357  dNdx_u, N_u, x_coord, _is_axially_symmetric);
358 
359  eps.noalias() = B * u;
360 
361  variables.mechanical_strain
363 
364  variables_prev.stress
366  sigma_eff_prev);
367  variables_prev.mechanical_strain
369  eps_prev);
370  variables_prev.temperature = _process_data.reference_temperature;
371 
372  auto&& solution = _ip_data[ip].solid_material.integrateStress(
373  variables_prev, variables, t, x_position, dt, *state);
374 
375  if (!solution)
376  {
377  OGS_FATAL("Computation of local constitutive relation failed.");
378  }
379 
381  std::tie(sigma_eff, state, C) = std::move(*solution);
382 
383  if (!_process_data.deactivate_matrix_in_flow) // Only for hydraulically
384  // active matrix
385  {
386  double const k_over_mu =
387  _process_data.intrinsic_permeability(t, x_position)[0] /
388  _process_data.fluid_viscosity(t, x_position)[0];
389  auto const rho_fr = _process_data.fluid_density(t, x_position)[0];
390  auto const& gravity_vec = _process_data.specific_body_force;
391  auto const& dNdx_p = ip_data.dNdx_p;
392 
393  ip_data.darcy_velocity.head(GlobalDim).noalias() =
394  -k_over_mu * (dNdx_p * p + rho_fr * gravity_vec);
395  }
396  }
397 
398  int n = GlobalDim == 2 ? 4 : 6;
399  Eigen::VectorXd ele_stress = Eigen::VectorXd::Zero(n);
400  Eigen::VectorXd ele_strain = Eigen::VectorXd::Zero(n);
401  GlobalDimVector ele_velocity = GlobalDimVector::Zero();
402 
403  for (auto const& ip_data : _ip_data)
404  {
405  ele_stress += ip_data.sigma_eff;
406  ele_strain += ip_data.eps;
407  ele_velocity += ip_data.darcy_velocity;
408  }
409 
410  ele_stress /= static_cast<double>(n_integration_points);
411  ele_strain /= static_cast<double>(n_integration_points);
412  ele_velocity /= static_cast<double>(n_integration_points);
413 
414  auto const element_id = _element.getID();
415  (*_process_data.mesh_prop_stress_xx)[_element.getID()] = ele_stress[0];
416  (*_process_data.mesh_prop_stress_yy)[_element.getID()] = ele_stress[1];
417  (*_process_data.mesh_prop_stress_zz)[_element.getID()] = ele_stress[2];
418  (*_process_data.mesh_prop_stress_xy)[_element.getID()] = ele_stress[3];
419  if (GlobalDim == 3)
420  {
421  (*_process_data.mesh_prop_stress_yz)[_element.getID()] = ele_stress[4];
422  (*_process_data.mesh_prop_stress_xz)[_element.getID()] = ele_stress[5];
423  }
424 
425  (*_process_data.mesh_prop_strain_xx)[_element.getID()] = ele_strain[0];
426  (*_process_data.mesh_prop_strain_yy)[_element.getID()] = ele_strain[1];
427  (*_process_data.mesh_prop_strain_zz)[_element.getID()] = ele_strain[2];
428  (*_process_data.mesh_prop_strain_xy)[_element.getID()] = ele_strain[3];
429  if (GlobalDim == 3)
430  {
431  (*_process_data.mesh_prop_strain_yz)[_element.getID()] = ele_strain[4];
432  (*_process_data.mesh_prop_strain_xz)[_element.getID()] = ele_strain[5];
433  }
434 
435  for (unsigned i = 0; i < GlobalDim; i++)
436  {
437  (*_process_data.mesh_prop_velocity)[element_id * GlobalDim + i] =
438  ele_velocity[i];
439  }
440 
442  ShapeFunctionPressure, typename ShapeFunctionDisplacement::MeshElement,
443  GlobalDim>(_element, _is_axially_symmetric, p,
444  *_process_data.mesh_prop_nodal_p);
445 }
void interpolateToHigherOrderNodes(MeshLib::Element const &element, bool const is_axially_symmetric, Eigen::MatrixBase< EigenMatrixType > const &node_values, MeshLib::PropertyVector< double > &interpolated_values_global_vector)

References ProcessLib::LinearBMatrix::computeBMatrix(), NumLib::interpolateToHigherOrderNodes(), NumLib::interpolateXCoordinate(), MaterialPropertyLib::VariableArray::mechanical_strain, OGS_FATAL, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), MaterialPropertyLib::VariableArray::stress, MathLib::t, MaterialPropertyLib::VariableArray::temperature, and MathLib::u.

◆ postTimestepConcreteWithVector()

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
void ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::postTimestepConcreteWithVector ( double const  t,
double const  dt,
Eigen::VectorXd const &  local_x 
)
overrideprotectedvirtual

Implements ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerInterface.

Reimplemented in ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrixNearFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >.

Definition at line 304 of file HydroMechanicsLocalAssemblerMatrix-impl.h.

306 {
307  auto p = const_cast<Eigen::VectorXd&>(local_x).segment(pressure_index,
308  pressure_size);
309  if (_process_data.deactivate_matrix_in_flow)
310  {
312  }
313  auto u = local_x.segment(displacement_index, displacement_size);
314 
316 }
void postTimestepConcreteWithBlockVectors(double const t, double const dt, Eigen::Ref< const Eigen::VectorXd > const &p, Eigen::Ref< const Eigen::VectorXd > const &u)

References MathLib::t, and MathLib::u.

◆ preTimestepConcrete()

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
void ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::preTimestepConcrete ( std::vector< double > const &  ,
double const  ,
double const   
)
inlineoverridevirtual

◆ setPressureDotOfInactiveNodes()

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
void ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::setPressureDotOfInactiveNodes ( Eigen::Ref< Eigen::VectorXd >  p_dot)
protected

Definition at line 472 of file HydroMechanicsLocalAssemblerMatrix-impl.h.

473 {
474  for (unsigned i = 0; i < pressure_size; i++)
475  {
476  // only inactive nodes
477  if (_process_data.p_element_status->isActiveNode(_element.getNode(i)))
478  {
479  continue;
480  }
481  p_dot[i] = 0;
482  }
483 }
virtual const Node * getNode(unsigned idx) const =0

◆ setPressureOfInactiveNodes()

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
void ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::setPressureOfInactiveNodes ( double const  t,
Eigen::Ref< Eigen::VectorXd >  p 
)
protected

Definition at line 450 of file HydroMechanicsLocalAssemblerMatrix-impl.h.

452 {
454  x_position.setElementID(_element.getID());
455  for (unsigned i = 0; i < pressure_size; i++)
456  {
457  // only inactive nodes
458  if (_process_data.p_element_status->isActiveNode(_element.getNode(i)))
459  {
460  continue;
461  }
462  x_position.setNodeID(getNodeIndex(_element, i));
463  auto const p0 = (*_process_data.p0)(t, x_position)[0];
464  p[i] = p0;
465  }
466 }
void setNodeID(std::size_t node_id)
std::size_t getNodeIndex(Element const &element, unsigned const idx)
Definition: Element.cpp:219

References MeshLib::getNodeIndex(), ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setNodeID(), and MathLib::t.

Member Data Documentation

◆ _ip_data

◆ _process_data

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
HydroMechanicsProcessData<GlobalDim>& ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::_process_data
protected

◆ displacement_index

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
const int ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::displacement_index = ShapeFunctionPressure::NPOINTS
staticprotected

Definition at line 117 of file HydroMechanicsLocalAssemblerMatrix.h.

◆ displacement_size

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
const int ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::displacement_size
staticprotected

◆ kelvin_vector_size

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
const int ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::kelvin_vector_size
staticprotected
Initial value:
=
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
Definition: KelvinVector.h:24

Definition at line 120 of file HydroMechanicsLocalAssemblerMatrix.h.

Referenced by ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::HydroMechanicsLocalAssemblerMatrix().

◆ pressure_index

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
const int ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::pressure_index = 0
staticprotected

Definition at line 115 of file HydroMechanicsLocalAssemblerMatrix.h.

◆ pressure_size

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , int GlobalDim>
const int ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerMatrix< ShapeFunctionDisplacement, ShapeFunctionPressure, GlobalDim >::pressure_size = ShapeFunctionPressure::NPOINTS
staticprotected

Definition at line 116 of file HydroMechanicsLocalAssemblerMatrix.h.


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