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, typename ShapeFunctionPressure,
30 
31  int GlobalDim>
34 {
35 public:
37  HydroMechanicsLocalAssemblerFracture const&) = delete;
40 
42  MeshLib::Element const& e,
43  std::size_t const local_matrix_size,
44  std::vector<unsigned> const& dofIndex_to_localIndex,
45  NumLib::GenericIntegrationMethod const& integration_method,
46  bool const is_axially_symmetric,
48 
49  void preTimestepConcrete(std::vector<double> const& /*local_x*/,
50  double const /*t*/,
51  double const /*delta_t*/) override
52  {
53  for (auto& data : _ip_data)
54  {
55  data.pushBackState();
56  }
57  }
58 
60  const double t, double const dt,
61  Eigen::VectorXd const& local_x) override;
62 
63 private:
64  void assembleWithJacobianConcrete(double const t, double const dt,
65  Eigen::VectorXd const& local_x,
66  Eigen::VectorXd const& local_xdot,
67  Eigen::VectorXd& local_b,
68  Eigen::MatrixXd& local_J) override;
69 
71  double const t, double const dt,
72  Eigen::Ref<const Eigen::VectorXd> const& p,
73  Eigen::Ref<const Eigen::VectorXd> const& p_dot,
74  Eigen::Ref<const Eigen::VectorXd> const& g,
75  Eigen::Ref<const Eigen::VectorXd> const& g_dot,
76  Eigen::Ref<Eigen::VectorXd> rhs_p, Eigen::Ref<Eigen::VectorXd> rhs_g,
77  Eigen::Ref<Eigen::MatrixXd> J_pp, Eigen::Ref<Eigen::MatrixXd> J_pg,
78  Eigen::Ref<Eigen::MatrixXd> J_gg, Eigen::Ref<Eigen::MatrixXd> J_gp);
79 
80  // Types for displacement.
83  using HMatricesType =
86 
87  // Types for pressure.
90 
91  // Types for the integration point data
96  GlobalDim>;
97 
98  using GlobalDimVectorType = Eigen::Matrix<double, GlobalDim, 1>;
99 
101 
102  std::vector<IntegrationPointDataType,
103  Eigen::aligned_allocator<IntegrationPointDataType>>
105 
106  static const int pressure_index = 0;
107  static const int pressure_size = ShapeFunctionPressure::NPOINTS;
108  static const int displacement_index = ShapeFunctionPressure::NPOINTS;
109  static const int displacement_size =
110  ShapeFunctionDisplacement::NPOINTS * GlobalDim;
111 };
112 
113 } // namespace HydroMechanics
114 } // namespace LIE
115 } // namespace ProcessLib
116 
MatrixType< DisplacementDim, _number_of_dof > HMatrixType
Definition: HMatrixUtils.h:42
ShapeMatrixPolicyType< ShapeFunctionPressure, GlobalDim > 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, GlobalDim > ShapeMatricesTypeDisplacement
void preTimestepConcrete(std::vector< double > const &, double const, double const) override
IntegrationPointDataFracture< HMatricesType, ShapeMatricesTypeDisplacement, ShapeMatricesTypePressure, GlobalDim > IntegrationPointDataType
void postTimestepConcreteWithVector(const double t, double const dt, Eigen::VectorXd const &local_x) override
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
HydroMechanicsLocalAssemblerFracture(HydroMechanicsLocalAssemblerFracture &&)=delete
static const double t