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

Detailed Description

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

Definition at line 33 of file HydroMechanicsLocalAssemblerFracture.h.

#include <HydroMechanicsLocalAssemblerFracture.h>

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

Public Member Functions

 HydroMechanicsLocalAssemblerFracture (HydroMechanicsLocalAssemblerFracture const &)=delete
 
 HydroMechanicsLocalAssemblerFracture (HydroMechanicsLocalAssemblerFracture &&)=delete
 
 HydroMechanicsLocalAssemblerFracture (MeshLib::Element const &e, std::size_t const local_matrix_size, std::vector< unsigned > const &dofIndex_to_localIndex, bool const is_axially_symmetric, unsigned const integration_order, HydroMechanicsProcessData< GlobalDim > &process_data)
 
void preTimestepConcrete (std::vector< double > const &, double const, double const) override
 
void postTimestepConcreteWithVector (const double t, double const dt, Eigen::VectorXd const &local_x) 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_, const double, const double, 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, 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)
 
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 std::vector< double > interpolateNodalValuesToIntegrationPoints (std::vector< double > const &)
 
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...
 

Private Types

using ShapeMatricesTypeDisplacement = ShapeMatrixPolicyType< ShapeFunctionDisplacement, ShapeFunctionDisplacement::DIM >
 
using HMatricesType = HMatrixPolicyType< ShapeFunctionDisplacement, GlobalDim >
 
using HMatrixType = typename HMatricesType::HMatrixType
 
using ShapeMatricesTypePressure = ShapeMatrixPolicyType< ShapeFunctionPressure, ShapeFunctionPressure::DIM >
 
using IntegrationPointDataType = IntegrationPointDataFracture< HMatricesType, ShapeMatricesTypeDisplacement, ShapeMatricesTypePressure, GlobalDim >
 
using GlobalDimVector = Eigen::Matrix< double, GlobalDim, 1 >
 
using DimVectorType = typename ShapeMatricesTypePressure::DimVectorType
 

Private Member Functions

void assembleWithJacobianConcrete (double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_xdot, Eigen::VectorXd &local_b, Eigen::MatrixXd &local_J) 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 &g, Eigen::Ref< const Eigen::VectorXd > const &g_dot, Eigen::Ref< Eigen::VectorXd > rhs_p, Eigen::Ref< Eigen::VectorXd > rhs_g, Eigen::Ref< Eigen::MatrixXd > J_pp, Eigen::Ref< Eigen::MatrixXd > J_pg, Eigen::Ref< Eigen::MatrixXd > J_gg, Eigen::Ref< Eigen::MatrixXd > J_gp)
 

Private Attributes

HydroMechanicsProcessData< GlobalDim > & _process_data
 
std::vector< IntegrationPointDataType, Eigen::aligned_allocator< IntegrationPointDataType > > _ip_data
 

Static Private 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
 

Additional Inherited Members

- Protected Attributes inherited from ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerInterface
MeshLib::Element const & _element
 
bool const _is_axially_symmetric
 

Member Typedef Documentation

◆ DimVectorType

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , typename IntegrationMethod , int GlobalDim>
using ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::DimVectorType = typename ShapeMatricesTypePressure::DimVectorType
private

Definition at line 102 of file HydroMechanicsLocalAssemblerFracture.h.

◆ GlobalDimVector

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

Definition at line 101 of file HydroMechanicsLocalAssemblerFracture.h.

◆ HMatricesType

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , typename IntegrationMethod , int GlobalDim>
using ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::HMatricesType = HMatrixPolicyType<ShapeFunctionDisplacement, GlobalDim>
private

Definition at line 85 of file HydroMechanicsLocalAssemblerFracture.h.

◆ HMatrixType

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , typename IntegrationMethod , int GlobalDim>
using ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::HMatrixType = typename HMatricesType::HMatrixType
private

Definition at line 87 of file HydroMechanicsLocalAssemblerFracture.h.

◆ IntegrationPointDataType

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , typename IntegrationMethod , int GlobalDim>
using ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::IntegrationPointDataType = IntegrationPointDataFracture<HMatricesType, ShapeMatricesTypeDisplacement, ShapeMatricesTypePressure, GlobalDim>
private

Definition at line 95 of file HydroMechanicsLocalAssemblerFracture.h.

◆ ShapeMatricesTypeDisplacement

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , typename IntegrationMethod , int GlobalDim>
using ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::ShapeMatricesTypeDisplacement = ShapeMatrixPolicyType<ShapeFunctionDisplacement, ShapeFunctionDisplacement::DIM>
private

Definition at line 82 of file HydroMechanicsLocalAssemblerFracture.h.

◆ ShapeMatricesTypePressure

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , typename IntegrationMethod , int GlobalDim>
using ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::ShapeMatricesTypePressure = ShapeMatrixPolicyType<ShapeFunctionPressure, ShapeFunctionPressure::DIM>
private

Definition at line 90 of file HydroMechanicsLocalAssemblerFracture.h.

Constructor & Destructor Documentation

◆ HydroMechanicsLocalAssemblerFracture() [1/3]

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

◆ HydroMechanicsLocalAssemblerFracture() [2/3]

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

◆ HydroMechanicsLocalAssemblerFracture() [3/3]

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , typename IntegrationMethod , int GlobalDim>
ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::HydroMechanicsLocalAssemblerFracture ( MeshLib::Element const &  e,
std::size_t const  local_matrix_size,
std::vector< unsigned > const &  dofIndex_to_localIndex,
bool const  is_axially_symmetric,
unsigned const  integration_order,
HydroMechanicsProcessData< GlobalDim > &  process_data 
)

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

47  e, is_axially_symmetric,
48  ShapeFunctionDisplacement::NPOINTS * GlobalDim +
49  ShapeFunctionPressure::NPOINTS,
50  dofIndex_to_localIndex),
51  _process_data(process_data)
52 {
53  assert(e.getDimension() == GlobalDim - 1);
54 
55  IntegrationMethod integration_method(integration_order);
56  unsigned const n_integration_points =
57  integration_method.getNumberOfPoints();
58 
59  _ip_data.reserve(n_integration_points);
60 
61  auto const shape_matrices_u =
62  NumLib::initShapeMatrices<ShapeFunctionDisplacement,
64  e, is_axially_symmetric, integration_method);
65 
66  auto const shape_matrices_p =
67  NumLib::initShapeMatrices<ShapeFunctionPressure,
68  ShapeMatricesTypePressure, GlobalDim>(
69  e, is_axially_symmetric, integration_method);
70 
71  auto const& frac_prop = *_process_data.fracture_property;
72 
73  // Get element nodes for aperture0 interpolation from nodes to integration
74  // point. The aperture0 parameter is time-independent.
76  aperture0_node_values = frac_prop.aperture0.getNodalValuesOnElement(
77  e, /*time independent*/ 0);
78 
80  x_position.setElementID(e.getID());
81  for (unsigned ip = 0; ip < n_integration_points; ip++)
82  {
83  x_position.setIntegrationPoint(ip);
84 
85  _ip_data.emplace_back(*_process_data.fracture_model);
86  auto const& sm_u = shape_matrices_u[ip];
87  auto const& sm_p = shape_matrices_p[ip];
88  auto& ip_data = _ip_data[ip];
89  ip_data.integration_weight =
90  sm_u.detJ * sm_u.integralMeasure *
91  integration_method.getWeightedPoint(ip).getWeight();
92 
93  ip_data.H_u.setZero(GlobalDim,
94  ShapeFunctionDisplacement::NPOINTS * GlobalDim);
96  GlobalDim, ShapeFunctionDisplacement::NPOINTS,
98  HMatrixType>(sm_u.N, ip_data.H_u);
99  ip_data.N_p = sm_p.N;
100  ip_data.dNdx_p = sm_p.dNdx;
101 
102  // Initialize current time step values
103  ip_data.w.setZero(GlobalDim);
104  ip_data.sigma_eff.setZero(GlobalDim);
105 
106  // Previous time step values are not initialized and are set later.
107  ip_data.w_prev.resize(GlobalDim);
108  ip_data.sigma_eff_prev.resize(GlobalDim);
109 
110  ip_data.C.resize(GlobalDim, GlobalDim);
111 
112  ip_data.aperture0 = aperture0_node_values.dot(sm_u.N);
113  ip_data.aperture = ip_data.aperture0;
114 
115  ip_data.permeability_state =
116  frac_prop.permeability_model->getNewState();
117 
118  auto const initial_effective_stress =
119  _process_data.initial_fracture_effective_stress(0, x_position);
120  for (int i = 0; i < GlobalDim; i++)
121  {
122  ip_data.sigma_eff[i] = initial_effective_stress[i];
123  ip_data.sigma_eff_prev[i] = initial_effective_stress[i];
124  }
125  }
126 }
void setElementID(std::size_t element_id)
void setIntegrationPoint(unsigned integration_point)
ShapeMatrixPolicyType< ShapeFunctionPressure, ShapeFunctionPressure::DIM > ShapeMatricesTypePressure
std::vector< IntegrationPointDataType, Eigen::aligned_allocator< IntegrationPointDataType > > _ip_data
ShapeMatrixPolicyType< ShapeFunctionDisplacement, ShapeFunctionDisplacement::DIM > ShapeMatricesTypeDisplacement
HydroMechanicsLocalAssemblerInterface(MeshLib::Element const &element, bool const is_axially_symmetric, std::size_t n_local_size, std::vector< unsigned > dofIndex_to_localIndex)
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)
void computeHMatrix(N_Type const &N, HMatrixType &H)
Fills a H-matrix based on given shape function.
Definition: HMatrixUtils.h:53
VectorType< ShapeFunction::NPOINTS > NodalVectorType
RowVectorType< ShapeFunction::NPOINTS > NodalRowVectorType

References ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::_ip_data, ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::_process_data, ProcessLib::computeHMatrix(), MeshLib::Element::getDimension(), MeshLib::Element::getID(), NumLib::initShapeMatrices(), ParameterLib::SpatialPosition::setElementID(), and ParameterLib::SpatialPosition::setIntegrationPoint().

Member Function Documentation

◆ assembleBlockMatricesWithJacobian()

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , typename IntegrationMethod , int GlobalDim>
void ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, 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 &  g,
Eigen::Ref< const Eigen::VectorXd > const &  g_dot,
Eigen::Ref< Eigen::VectorXd >  rhs_p,
Eigen::Ref< Eigen::VectorXd >  rhs_g,
Eigen::Ref< Eigen::MatrixXd >  J_pp,
Eigen::Ref< Eigen::MatrixXd >  J_pg,
Eigen::Ref< Eigen::MatrixXd >  J_gg,
Eigen::Ref< Eigen::MatrixXd >  J_gp 
)
private

Definition at line 163 of file HydroMechanicsLocalAssemblerFracture-impl.h.

173 {
174  auto const& frac_prop = *_process_data.fracture_property;
175  auto const& R = frac_prop.R;
176 
177  // the index of a normal (normal to a fracture plane) component
178  // in a displacement vector
179  auto constexpr index_normal = GlobalDim - 1;
180 
181  typename ShapeMatricesTypePressure::NodalMatrixType laplace_p =
182  ShapeMatricesTypePressure::NodalMatrixType::Zero(pressure_size,
183  pressure_size);
184 
185  typename ShapeMatricesTypePressure::NodalMatrixType storage_p =
186  ShapeMatricesTypePressure::NodalMatrixType::Zero(pressure_size,
187  pressure_size);
188 
189  typename ShapeMatricesTypeDisplacement::template MatrixType<
191  Kgp = ShapeMatricesTypeDisplacement::template MatrixType<
193  pressure_size);
194 
195  // Projection of the body force vector at the element.
196  Eigen::MatrixXd const global2local_rotation =
197  R.template topLeftCorner<ShapeFunctionPressure::DIM, GlobalDim>();
198 
199  DimVectorType const gravity_vec =
200  global2local_rotation * _process_data.specific_body_force;
201 
203  x_position.setElementID(_element.getID());
204 
205  unsigned const n_integration_points = _ip_data.size();
206  for (unsigned ip = 0; ip < n_integration_points; ip++)
207  {
208  x_position.setIntegrationPoint(ip);
209 
210  auto& ip_data = _ip_data[ip];
211  auto const& ip_w = ip_data.integration_weight;
212  auto const& N_p = ip_data.N_p;
213  auto const& dNdx_p = ip_data.dNdx_p;
214  auto const& H_g = ip_data.H_u;
215  auto const& identity2 =
217 
218  auto& mat = ip_data.fracture_material;
219  auto& effective_stress = ip_data.sigma_eff;
220  auto const& effective_stress_prev = ip_data.sigma_eff_prev;
221  auto& w = ip_data.w;
222  auto const& w_prev = ip_data.w_prev;
223  auto& C = ip_data.C;
224  auto& state = *ip_data.material_state_variables;
225  auto& b_m = ip_data.aperture;
226 
227  double const S = frac_prop.specific_storage(t, x_position)[0];
228  double const mu = _process_data.fluid_viscosity(t, x_position)[0];
229  auto const alpha = frac_prop.biot_coefficient(t, x_position)[0];
230  auto const rho_fr = _process_data.fluid_density(t, x_position)[0];
231 
232  // displacement jumps in local coordinates
233  w.noalias() = R * H_g * g;
234 
235  // aperture
236  b_m = ip_data.aperture0 + w[index_normal];
237  if (b_m < 0.0)
238  {
239  DBUG(
240  "Element {:d}, gp {:d}: Fracture aperture is {:g}, but it must "
241  "be "
242  "non-negative. Setting it to zero.",
243  _element.getID(), ip, b_m);
244  b_m = 0;
245  }
246 
247  auto const initial_effective_stress =
248  _process_data.initial_fracture_effective_stress(0, x_position);
249 
250  Eigen::Map<typename HMatricesType::ForceVectorType const> const stress0(
251  initial_effective_stress.data(), initial_effective_stress.size());
252 
253  // local C, local stress
254  mat.computeConstitutiveRelation(
255  t, x_position, ip_data.aperture0, stress0, w_prev, w,
256  effective_stress_prev, effective_stress, C, state);
257 
258  auto& k = ip_data.permeability;
259  k = frac_prop.permeability_model->permeability(
260  ip_data.permeability_state.get(), ip_data.aperture0, b_m);
261 
262  // derivative of permeability respect to aperture
263  double const local_dk_db =
264  frac_prop.permeability_model->dpermeability_daperture(
265  ip_data.permeability_state.get(), ip_data.aperture0, b_m);
266 
267  //
268  // displacement equation, displacement jump part
269  //
270  rhs_g.noalias() -=
271  H_g.transpose() * R.transpose() * effective_stress * ip_w;
272  J_gg.noalias() += H_g.transpose() * R.transpose() * C * R * H_g * ip_w;
273 
274  //
275  // displacement equation, pressure part
276  //
277  Kgp.noalias() +=
278  H_g.transpose() * R.transpose() * alpha * identity2 * N_p * ip_w;
279 
280  //
281  // pressure equation, pressure part.
282  //
283  storage_p.noalias() += N_p.transpose() * b_m * S * N_p * ip_w;
284  laplace_p.noalias() +=
285  dNdx_p.transpose() * b_m * k / mu * dNdx_p * ip_w;
286  rhs_p.noalias() +=
287  dNdx_p.transpose() * b_m * k / mu * rho_fr * gravity_vec * ip_w;
288 
289  //
290  // pressure equation, displacement jump part.
291  //
292  auto const grad_head_over_mu =
293  ((dNdx_p * p + rho_fr * gravity_vec) / mu).eval();
294  Eigen::Matrix<double, 1, displacement_size> const mT_R_Hg =
295  identity2.transpose() * R * H_g;
296  // velocity in global coordinates
297  ip_data.darcy_velocity =
298  -global2local_rotation.transpose() * k * grad_head_over_mu;
299  J_pg.noalias() += N_p.transpose() * S * N_p * p_dot * mT_R_Hg * ip_w;
300  J_pg.noalias() +=
301  dNdx_p.transpose() * k * grad_head_over_mu * mT_R_Hg * ip_w;
302  J_pg.noalias() += dNdx_p.transpose() * b_m * local_dk_db *
303  grad_head_over_mu * mT_R_Hg * ip_w;
304  }
305 
306  // displacement equation, pressure part
307  J_gp.noalias() -= Kgp;
308 
309  // pressure equation, pressure part.
310  J_pp.noalias() += laplace_p + storage_p / dt;
311 
312  // pressure equation, displacement jump part.
313  J_pg.noalias() += Kgp.transpose() / dt;
314 
315  // pressure equation
316  rhs_p.noalias() -=
317  laplace_p * p + storage_p * p_dot + Kgp.transpose() * g_dot;
318 
319  // displacement equation
320  rhs_g.noalias() -= -Kgp * p;
321 }
void DBUG(char const *fmt, Args const &... args)
Definition: Logging.h:27
virtual std::size_t getID() const final
Returns the ID of the element.
Definition: Element.h:82
MatrixType< ShapeFunction::NPOINTS, ShapeFunction::NPOINTS > NodalMatrixType

References MaterialPropertyLib::alpha, DBUG(), MaterialPropertyLib::effective_stress, ParameterLib::SpatialPosition::setElementID(), and ParameterLib::SpatialPosition::setIntegrationPoint().

◆ assembleWithJacobianConcrete()

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , typename IntegrationMethod , int GlobalDim>
void ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::assembleWithJacobianConcrete ( double const  t,
double const  dt,
Eigen::VectorXd const &  local_x,
Eigen::VectorXd const &  local_xdot,
Eigen::VectorXd &  local_b,
Eigen::MatrixXd &  local_J 
)
overrideprivatevirtual

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

Definition at line 132 of file HydroMechanicsLocalAssemblerFracture-impl.h.

137 {
138  auto const p = local_x.segment(pressure_index, pressure_size);
139  auto const p_dot = local_xdot.segment(pressure_index, pressure_size);
140  auto const g = local_x.segment(displacement_index, displacement_size);
141  auto const g_dot =
142  local_xdot.segment(displacement_index, displacement_size);
143 
144  auto rhs_p = local_b.segment(pressure_index, pressure_size);
145  auto rhs_g = local_b.segment(displacement_index, displacement_size);
146  auto J_pp = local_J.block(pressure_index, pressure_index, pressure_size,
147  pressure_size);
148  auto J_pg = local_J.block(pressure_index, displacement_index, pressure_size,
150  auto J_gp = local_J.block(displacement_index, pressure_index,
152  auto J_gg = local_J.block(displacement_index, displacement_index,
154 
155  assembleBlockMatricesWithJacobian(t, dt, p, p_dot, g, g_dot, rhs_p, rhs_g,
156  J_pp, J_pg, J_gg, J_gp);
157 }
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 &g, Eigen::Ref< const Eigen::VectorXd > const &g_dot, Eigen::Ref< Eigen::VectorXd > rhs_p, Eigen::Ref< Eigen::VectorXd > rhs_g, Eigen::Ref< Eigen::MatrixXd > J_pp, Eigen::Ref< Eigen::MatrixXd > J_pg, Eigen::Ref< Eigen::MatrixXd > J_gg, Eigen::Ref< Eigen::MatrixXd > J_gp)

◆ postTimestepConcreteWithVector()

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

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

Definition at line 327 of file HydroMechanicsLocalAssemblerFracture-impl.h.

330 {
331  auto const nodal_g = local_x.segment(displacement_index, displacement_size);
332 
333  auto const& frac_prop = *_process_data.fracture_property;
334  auto const& R = frac_prop.R;
335  // the index of a normal (normal to a fracture plane) component
336  // in a displacement vector
337  auto constexpr index_normal = GlobalDim - 1;
338 
340  x_position.setElementID(_element.getID());
341 
342  unsigned const n_integration_points = _ip_data.size();
343  for (unsigned ip = 0; ip < n_integration_points; ip++)
344  {
345  x_position.setIntegrationPoint(ip);
346 
347  auto& ip_data = _ip_data[ip];
348  auto const& H_g = ip_data.H_u;
349  auto& mat = ip_data.fracture_material;
350  auto& effective_stress = ip_data.sigma_eff;
351  auto const& effective_stress_prev = ip_data.sigma_eff_prev;
352  auto& w = ip_data.w;
353  auto const& w_prev = ip_data.w_prev;
354  auto& C = ip_data.C;
355  auto& state = *ip_data.material_state_variables;
356  auto& b_m = ip_data.aperture;
357 
358  // displacement jumps in local coordinates
359  w.noalias() = R * H_g * nodal_g;
360 
361  // aperture
362  b_m = ip_data.aperture0 + w[index_normal];
363  if (b_m < 0.0)
364  {
365  DBUG(
366  "Element {:d}, gp {:d}: Fracture aperture is {:g}, but it is "
367  "expected to be non-negative.",
368  _element.getID(), ip, b_m);
369  }
370 
371  auto const initial_effective_stress =
372  _process_data.initial_fracture_effective_stress(0, x_position);
373 
374  Eigen::Map<typename HMatricesType::ForceVectorType const> const stress0(
375  initial_effective_stress.data(), initial_effective_stress.size());
376 
377  // local C, local stress
378  mat.computeConstitutiveRelation(
379  t, x_position, ip_data.aperture0, stress0, w_prev, w,
380  effective_stress_prev, effective_stress, C, state);
381  }
382 
383  double ele_b = 0;
384  double ele_k = 0;
385  typename HMatricesType::ForceVectorType ele_sigma_eff =
386  HMatricesType::ForceVectorType::Zero(GlobalDim);
387  typename HMatricesType::ForceVectorType ele_w =
388  HMatricesType::ForceVectorType::Zero(GlobalDim);
389  double ele_Fs = -std::numeric_limits<double>::max();
390  GlobalDimVector ele_velocity = GlobalDimVector::Zero();
391  for (auto const& ip : _ip_data)
392  {
393  ele_b += ip.aperture;
394  ele_k += ip.permeability;
395  ele_w += ip.w;
396  ele_sigma_eff += ip.sigma_eff;
397  ele_Fs = std::max(
398  ele_Fs, ip.material_state_variables->getShearYieldFunctionValue());
399  ele_velocity += ip.darcy_velocity;
400  }
401  ele_b /= static_cast<double>(n_integration_points);
402  ele_k /= static_cast<double>(n_integration_points);
403  ele_w /= static_cast<double>(n_integration_points);
404  ele_sigma_eff /= static_cast<double>(n_integration_points);
405  ele_velocity /= static_cast<double>(n_integration_points);
406  auto const element_id = _element.getID();
407  (*_process_data.mesh_prop_b)[element_id] = ele_b;
408  (*_process_data.mesh_prop_k_f)[element_id] = ele_k;
409  (*_process_data.mesh_prop_w_n)[element_id] = ele_w[index_normal];
410  (*_process_data.mesh_prop_w_s)[element_id] = ele_w[0];
411  (*_process_data.mesh_prop_fracture_stress_normal)[element_id] =
412  ele_sigma_eff[index_normal];
413  (*_process_data.mesh_prop_fracture_stress_shear)[element_id] =
414  ele_sigma_eff[0];
415  (*_process_data.mesh_prop_fracture_shear_failure)[element_id] = ele_Fs;
416 
417  if (GlobalDim == 3)
418  {
419  (*_process_data.mesh_prop_w_s2)[element_id] = ele_w[1];
420  (*_process_data.mesh_prop_fracture_stress_shear2)[element_id] =
421  ele_sigma_eff[1];
422  }
423 
424  for (unsigned i = 0; i < GlobalDim; i++)
425  {
426  (*_process_data.mesh_prop_velocity)[element_id * GlobalDim + i] =
427  ele_velocity[i];
428  }
429 }
VectorType< DisplacementDim > ForceVectorType
Definition: HMatrixUtils.h:45

References DBUG(), MaterialPropertyLib::effective_stress, ParameterLib::SpatialPosition::setElementID(), and ParameterLib::SpatialPosition::setIntegrationPoint().

◆ preTimestepConcrete()

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

Member Data Documentation

◆ _ip_data

◆ _process_data

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , typename IntegrationMethod , int GlobalDim>
HydroMechanicsProcessData<GlobalDim>& ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::_process_data
private

◆ displacement_index

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , typename IntegrationMethod , int GlobalDim>
const int ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::displacement_index = ShapeFunctionPressure::NPOINTS
staticprivate

Definition at line 112 of file HydroMechanicsLocalAssemblerFracture.h.

◆ displacement_size

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , typename IntegrationMethod , int GlobalDim>
const int ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::displacement_size
staticprivate
Initial value:
=
ShapeFunctionDisplacement::NPOINTS * GlobalDim

Definition at line 113 of file HydroMechanicsLocalAssemblerFracture.h.

◆ pressure_index

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , typename IntegrationMethod , int GlobalDim>
const int ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::pressure_index = 0
staticprivate

Definition at line 110 of file HydroMechanicsLocalAssemblerFracture.h.

◆ pressure_size

template<typename ShapeFunctionDisplacement , typename ShapeFunctionPressure , typename IntegrationMethod , int GlobalDim>
const int ProcessLib::LIE::HydroMechanics::HydroMechanicsLocalAssemblerFracture< ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::pressure_size = ShapeFunctionPressure::NPOINTS
staticprivate

Definition at line 111 of file HydroMechanicsLocalAssemblerFracture.h.


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