OGS
HydroMechanicsLocalAssemblerFracture.h
Go to the documentation of this file.
1 
11 #pragma once
12 
13 #include <vector>
14 
16 
19 
22 
23 namespace ProcessLib
24 {
25 namespace LIE
26 {
27 namespace HydroMechanics
28 {
29 template <typename ShapeFunctionDisplacement,
30  typename ShapeFunctionPressure,
31  typename IntegrationMethod,
32  int GlobalDim>
35 {
36 public:
38  HydroMechanicsLocalAssemblerFracture const&) = delete;
41 
43  MeshLib::Element const& e,
44  std::size_t const local_matrix_size,
45  std::vector<unsigned> const& dofIndex_to_localIndex,
46  bool const is_axially_symmetric,
47  unsigned const integration_order,
49 
50  void preTimestepConcrete(std::vector<double> const& /*local_x*/,
51  double const /*t*/,
52  double const /*delta_t*/) override
53  {
54  for (auto& data : _ip_data)
55  {
56  data.pushBackState();
57  }
58  }
59 
61  const double t, double const dt,
62  Eigen::VectorXd const& local_x) override;
63 
64 private:
65  void assembleWithJacobianConcrete(double const t, double const dt,
66  Eigen::VectorXd const& local_x,
67  Eigen::VectorXd const& local_xdot,
68  Eigen::VectorXd& local_b,
69  Eigen::MatrixXd& local_J) override;
70 
72  double const t, double const dt,
73  Eigen::Ref<const Eigen::VectorXd> const& p,
74  Eigen::Ref<const Eigen::VectorXd> const& p_dot,
75  Eigen::Ref<const Eigen::VectorXd> const& g,
76  Eigen::Ref<const Eigen::VectorXd> const& g_dot,
77  Eigen::Ref<Eigen::VectorXd> rhs_p, Eigen::Ref<Eigen::VectorXd> rhs_g,
78  Eigen::Ref<Eigen::MatrixXd> J_pp, Eigen::Ref<Eigen::MatrixXd> J_pg,
79  Eigen::Ref<Eigen::MatrixXd> J_gg, Eigen::Ref<Eigen::MatrixXd> J_gp);
80 
81  // Types for displacement.
83  ShapeMatrixPolicyType<ShapeFunctionDisplacement,
84  ShapeFunctionDisplacement::DIM>;
85  using HMatricesType =
88 
89  // Types for pressure.
91  ShapeMatrixPolicyType<ShapeFunctionPressure,
92  ShapeFunctionPressure::DIM>;
93 
94  // Types for the integration point data
99  GlobalDim>;
100 
101  using GlobalDimVector = Eigen::Matrix<double, GlobalDim, 1>;
103 
105 
106  std::vector<IntegrationPointDataType,
107  Eigen::aligned_allocator<IntegrationPointDataType>>
109 
110  static const int pressure_index = 0;
111  static const int pressure_size = ShapeFunctionPressure::NPOINTS;
112  static const int displacement_index = ShapeFunctionPressure::NPOINTS;
113  static const int displacement_size =
114  ShapeFunctionDisplacement::NPOINTS * GlobalDim;
115 };
116 
117 } // namespace HydroMechanics
118 } // namespace LIE
119 } // namespace ProcessLib
120 
MatrixType< DisplacementDim, _number_of_dof > HMatrixType
Definition: HMatrixUtils.h:42
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
HydroMechanicsLocalAssemblerFracture(HydroMechanicsLocalAssemblerFracture &&)=delete
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
ShapeMatrixPolicyType< ShapeFunctionPressure, ShapeFunctionPressure::DIM > ShapeMatricesTypePressure
std::vector< IntegrationPointDataType, Eigen::aligned_allocator< IntegrationPointDataType > > _ip_data
HydroMechanicsLocalAssemblerFracture(HydroMechanicsLocalAssemblerFracture const &)=delete
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)
ShapeMatrixPolicyType< ShapeFunctionDisplacement, ShapeFunctionDisplacement::DIM > ShapeMatricesTypeDisplacement
IntegrationPointDataFracture< HMatricesType, ShapeMatricesTypeDisplacement, ShapeMatricesTypePressure, GlobalDim > IntegrationPointDataType
VectorType< ShapeFunction::DIM > DimVectorType