OGS
HydroMechanicsLocalAssemblerMatrix.h
Go to the documentation of this file.
1 
11 #pragma once
12 
13 #include <vector>
14 
21 
22 namespace ProcessLib
23 {
24 namespace LIE
25 {
26 namespace HydroMechanics
27 {
28 namespace MPL = MaterialPropertyLib;
29 
30 template <typename ShapeFunctionDisplacement, typename ShapeFunctionPressure,
31  int GlobalDim>
34 {
35 public:
37  HydroMechanicsLocalAssemblerMatrix const&) = delete;
39  delete;
40 
42  MeshLib::Element const& e,
43  std::size_t const n_variables,
44  std::size_t const local_matrix_size,
45  std::vector<unsigned> const& dofIndex_to_localIndex,
46  NumLib::GenericIntegrationMethod const& integration_method,
47  bool const is_axially_symmetric,
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 
60 protected:
61  void assembleWithJacobianConcrete(double const t, double const dt,
62  Eigen::VectorXd const& local_x,
63  Eigen::VectorXd const& local_x_dot,
64  Eigen::VectorXd& local_rhs,
65  Eigen::MatrixXd& local_Jac) override;
66 
68  double const t, double const dt,
69  Eigen::Ref<const Eigen::VectorXd> const& p,
70  Eigen::Ref<const Eigen::VectorXd> const& p_dot,
71  Eigen::Ref<const Eigen::VectorXd> const& u,
72  Eigen::Ref<const Eigen::VectorXd> const& u_dot,
73  Eigen::Ref<Eigen::VectorXd> rhs_p, Eigen::Ref<Eigen::VectorXd> rhs_u,
74  Eigen::Ref<Eigen::MatrixXd> J_pp, Eigen::Ref<Eigen::MatrixXd> J_pu,
75  Eigen::Ref<Eigen::MatrixXd> J_uu, Eigen::Ref<Eigen::MatrixXd> J_up);
76 
78  double const t, double const dt,
79  Eigen::VectorXd const& local_x) override;
80 
82  double const t, double const dt,
83  Eigen::Ref<const Eigen::VectorXd> const& p,
84  Eigen::Ref<const Eigen::VectorXd> const& u);
85 
87  double const t, Eigen::Ref<Eigen::VectorXd> p);
88  void setPressureDotOfInactiveNodes(Eigen::Ref<Eigen::VectorXd> p_dot);
89 
90  // Types for displacement.
93  using BMatricesType =
95 
96  // Types for pressure.
99 
104  GlobalDim,
105  ShapeFunctionDisplacement::NPOINTS>;
106 
107  using GlobalDimVector = Eigen::Matrix<double, GlobalDim, 1>;
108 
110 
111  std::vector<IntegrationPointDataType,
112  Eigen::aligned_allocator<IntegrationPointDataType>>
114 
115  static const int pressure_index = 0;
116  static const int pressure_size = ShapeFunctionPressure::NPOINTS;
117  static const int displacement_index = ShapeFunctionPressure::NPOINTS;
118  static const int displacement_size =
119  ShapeFunctionDisplacement::NPOINTS * GlobalDim;
120  static const int kelvin_vector_size =
122 };
123 
124 } // namespace HydroMechanics
125 } // namespace LIE
126 } // namespace ProcessLib
127 
ShapeMatrixPolicyType< ShapeFunctionDisplacement, GlobalDim > ShapeMatricesTypeDisplacement
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 postTimestepConcreteWithBlockVectors(double const t, double const dt, Eigen::Ref< const Eigen::VectorXd > const &p, Eigen::Ref< const Eigen::VectorXd > const &u)
std::vector< IntegrationPointDataType, Eigen::aligned_allocator< IntegrationPointDataType > > _ip_data
IntegrationPointDataMatrix< BMatricesType, ShapeMatricesTypeDisplacement, ShapeMatricesTypePressure, GlobalDim, ShapeFunctionDisplacement::NPOINTS > IntegrationPointDataType
BMatrixPolicyType< ShapeFunctionDisplacement, GlobalDim > BMatricesType
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)
HydroMechanicsLocalAssemblerMatrix(HydroMechanicsLocalAssemblerMatrix const &)=delete
void postTimestepConcreteWithVector(double const t, double const dt, Eigen::VectorXd const &local_x) override
HydroMechanicsLocalAssemblerMatrix(HydroMechanicsLocalAssemblerMatrix &&)=delete
void preTimestepConcrete(std::vector< double > const &, double const, double const) override
ShapeMatrixPolicyType< ShapeFunctionPressure, GlobalDim > ShapeMatricesTypePressure
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
Definition: KelvinVector.h:24
static const double u
static const double t