OGS
ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim > Class Template Reference

Detailed Description

template<typename ShapeFunctionLiquidVelocity, typename ShapeFunctionPressure, int GlobalDim>
class ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >

Definition at line 33 of file StokesFlowFEM.h.

#include <StokesFlowFEM.h>

Inheritance diagram for ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >:
[legend]
Collaboration diagram for ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >:
[legend]

Public Member Functions

 LocalAssemblerData (MeshLib::Element const &element, std::size_t const, NumLib::GenericIntegrationMethod const &integration_method, bool const is_axially_symmetric, StokesFlowProcessData const &process_data)
 
void assemble (double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &, std::vector< double > &, std::vector< double > &local_K_data, std::vector< double > &local_b_data) override
 
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix (const unsigned integration_point) const override
 Provides the shape matrix at the given integration point.
 
void computeSecondaryVariableConcrete (double const, double const, Eigen::VectorXd const &local_x, Eigen::VectorXd const &) override
 
- Public Member Functions inherited from ProcessLib::StokesFlow::StokesFlowLocalAssemblerInterface
 StokesFlowLocalAssemblerInterface ()=default
 
- Public Member Functions inherited from ProcessLib::LocalAssemblerInterface
virtual ~LocalAssemblerInterface ()=default
 
virtual void setInitialConditions (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, double const t, 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_x_prev, int const process_id, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data)
 
virtual void assembleWithJacobian (double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_x_prev, 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 assembleWithJacobianForStaggeredScheme (double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, 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_prev, 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, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, int const process_id)
 
void postNonLinearSolver (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, 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.
 
- Public Member Functions inherited from NumLib::ExtrapolatableElement
virtual ~ExtrapolatableElement ()=default
 

Private Types

using ShapeMatrixTypeLiquidVelocity
 
using ShapeMatrixTypePressure
 
using NodalVectorType = typename ShapeMatrixTypePressure::NodalVectorType
 

Private Attributes

MeshLib::Element const & _element
 
bool const _is_axially_symmetric
 
NumLib::GenericIntegrationMethod const & _integration_method
 
StokesFlowProcessData const & _process_data
 
std::vector< IntegrationPointData< ShapeMatrixTypeLiquidVelocity, ShapeMatrixTypePressure, GlobalDim, ShapeFunctionLiquidVelocity::NPOINTS >, Eigen::aligned_allocator< IntegrationPointData< ShapeMatrixTypeLiquidVelocity, ShapeMatrixTypePressure, GlobalDim, ShapeFunctionLiquidVelocity::NPOINTS > > > _ip_data
 

Static Private Attributes

static const int liquid_velocity_index = 0
 
static const int pressure_index
 
static const int liquid_velocity_size
 
static const int pressure_size = ShapeFunctionPressure::NPOINTS
 

Member Typedef Documentation

◆ NodalVectorType

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
using ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::NodalVectorType = typename ShapeMatrixTypePressure::NodalVectorType
private

Definition at line 48 of file StokesFlowFEM.h.

◆ ShapeMatrixTypeLiquidVelocity

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
using ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::ShapeMatrixTypeLiquidVelocity
private

◆ ShapeMatrixTypePressure

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
using ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::ShapeMatrixTypePressure
private

Constructor & Destructor Documentation

◆ LocalAssemblerData()

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::LocalAssemblerData ( MeshLib::Element const & element,
std::size_t const ,
NumLib::GenericIntegrationMethod const & integration_method,
bool const is_axially_symmetric,
StokesFlowProcessData const & process_data )
inline

Definition at line 51 of file StokesFlowFEM.h.

57 : _element(element),
58 _is_axially_symmetric(is_axially_symmetric),
59 _integration_method(integration_method),
60 _process_data(process_data)
61 {
62 unsigned const n_integration_points =
64 _ip_data.reserve(n_integration_points);
65
66 auto const shape_matrices_v =
67 NumLib::initShapeMatrices<ShapeFunctionLiquidVelocity,
69 element, is_axially_symmetric, _integration_method);
70
71 auto const shape_matrices_p =
72 NumLib::initShapeMatrices<ShapeFunctionPressure,
73 ShapeMatrixTypePressure, GlobalDim>(
74 element, is_axially_symmetric, _integration_method);
75
76 for (unsigned ip = 0; ip < n_integration_points; ip++)
77 {
78 _ip_data.emplace_back(
79 shape_matrices_v[ip].N, shape_matrices_v[ip].dNdx,
80 shape_matrices_p[ip].N, shape_matrices_p[ip].dNdx,
82 shape_matrices_v[ip].integralMeasure *
83 shape_matrices_v[ip].detJ);
84
85 auto& ip_data = _ip_data[ip];
86
87 ip_data.N_v_op = ShapeMatrixTypeLiquidVelocity::template MatrixType<
88 GlobalDim, liquid_velocity_size>::Zero(GlobalDim,
90
91 ip_data.dNdx_v_op =
92 ShapeMatrixTypeLiquidVelocity::template MatrixType<
93 GlobalDim * GlobalDim,
94 liquid_velocity_size>::Zero(GlobalDim * GlobalDim,
96
97 for (int i = 0; i < GlobalDim; ++i)
98 {
99 ip_data.N_v_op
100 .template block<1, liquid_velocity_size / GlobalDim>(
101 i, i * liquid_velocity_size / GlobalDim)
102 .noalias() = shape_matrices_v[ip].N;
103
104 ip_data.dNdx_v_op
105 .template block<GlobalDim,
106 liquid_velocity_size / GlobalDim>(
107 i * GlobalDim, i * liquid_velocity_size / GlobalDim)
108 .noalias() = shape_matrices_v[ip].dNdx;
109 }
110 }
111 }
double getWeight() const
MathLib::WeightedPoint const & getWeightedPoint(unsigned const igp) const
ShapeMatrixPolicyType< ShapeFunctionLiquidVelocity, GlobalDim > ShapeMatrixTypeLiquidVelocity
StokesFlowProcessData const & _process_data
NumLib::GenericIntegrationMethod const & _integration_method
std::vector< IntegrationPointData< ShapeMatrixTypeLiquidVelocity, ShapeMatrixTypePressure, GlobalDim, ShapeFunctionLiquidVelocity::NPOINTS >, Eigen::aligned_allocator< IntegrationPointData< ShapeMatrixTypeLiquidVelocity, ShapeMatrixTypePressure, GlobalDim, ShapeFunctionLiquidVelocity::NPOINTS > > > _ip_data
ShapeMatrixPolicyType< ShapeFunctionPressure, GlobalDim > ShapeMatrixTypePressure
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::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::_integration_method, ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::_ip_data, NumLib::GenericIntegrationMethod::getNumberOfPoints(), MathLib::WeightedPoint::getWeight(), NumLib::GenericIntegrationMethod::getWeightedPoint(), NumLib::initShapeMatrices(), and ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::liquid_velocity_size.

Member Function Documentation

◆ assemble()

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
void ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::assemble ( double const t,
double const dt,
std::vector< double > const & local_x,
std::vector< double > const & ,
std::vector< double > & ,
std::vector< double > & local_K_data,
std::vector< double > & local_b_data )
inlineoverridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 113 of file StokesFlowFEM.h.

119 {
120 auto const local_matrix_size = liquid_velocity_size + pressure_size;
121 assert(local_x.size() == local_matrix_size);
122
123 auto local_p = Eigen::Map<const NodalVectorType>(
124 &local_x[pressure_index], pressure_size);
125
126 auto local_K = MathLib::createZeroedMatrix<
127 typename ShapeMatrixTypeLiquidVelocity::template MatrixType<
128 local_matrix_size, local_matrix_size>>(
129 local_K_data, local_matrix_size, local_matrix_size);
130 auto local_b = MathLib::createZeroedVector<
131 typename ShapeMatrixTypeLiquidVelocity::template VectorType<
132 local_matrix_size>>(local_b_data, local_matrix_size);
133
134 // Get block matrices
135 // Kvv
136 auto Kvv =
137 local_K.template block<liquid_velocity_size, liquid_velocity_size>(
139 // Kvp
140 auto Kvp = local_K.template block<liquid_velocity_size, pressure_size>(
142 // kpv
143 auto Kpv = local_K.template block<pressure_size, liquid_velocity_size>(
145 // bv
146 auto bv = local_b.template segment<liquid_velocity_size>(
148
149 unsigned const n_integration_points =
151
154
155 auto const& b = _process_data.specific_body_force;
156
158
159 // create unit vector
160 assert(GlobalDim == 2);
161 auto const identity_mat =
162 Eigen::Matrix<double, GlobalDim, GlobalDim>::Identity().eval();
163 auto const unit_vector = Eigen::Map<const Eigen::RowVectorXd>(
164 identity_mat.data(), identity_mat.size())
165 .eval();
166
167 // Get material properties
168 auto const& medium =
170 auto const& phase = medium.phase("AqueousLiquid");
171
172 for (unsigned ip(0); ip < n_integration_points; ++ip)
173 {
174 pos.setIntegrationPoint(ip);
175
176 auto const& ip_data = _ip_data[ip];
177
178 auto const& N_p = ip_data.N_p;
179 auto const& N_v_op = ip_data.N_v_op;
180
181 auto const& dNdx_v = ip_data.dNdx_v;
182 auto const& dNdx_v_op = ip_data.dNdx_v_op;
183 auto const& dNdx_p = ip_data.dNdx_p;
184
185 auto const& w = ip_data.integration_weight;
186
187 double p_int_pt = 0.0;
188
189 NumLib::shapeFunctionInterpolate(local_p, N_p, p_int_pt);
190
191 vars.liquid_phase_pressure = p_int_pt;
192
193 // Use the viscosity model to compute the viscosity
195 .template value<double>(vars, pos, t, dt);
196
197 // Kvv
199 {
200 // permeability
201 auto const K = MaterialPropertyLib::formEigenTensor<GlobalDim>(
203 .value(vars, pos, t, dt));
204
205 Kvv.noalias() +=
206 w * N_v_op.transpose() * mu * K.inverse() * N_v_op;
207 }
208
209 for (int i = 0; i < GlobalDim; ++i)
210 {
211 Kvv.template block<liquid_velocity_size / GlobalDim,
212 liquid_velocity_size / GlobalDim>(
213 i * liquid_velocity_size / GlobalDim,
214 i * liquid_velocity_size / GlobalDim)
215 .noalias() += w * dNdx_v.transpose() * mu * dNdx_v;
216 }
217
218 // Kvp
219 Kvp.noalias() += w * N_v_op.transpose() * dNdx_p;
220
221 // Kpv
222 Kpv.noalias() += w * N_p.transpose() * unit_vector * dNdx_v_op;
223
224 // bv
225 bv.noalias() += w * N_v_op.transpose() * b;
226 }
227 }
Phase const & phase(std::size_t index) const
Definition Medium.cpp:33
std::size_t getID() const
Returns the ID of the element.
Definition Element.h:89
void setElementID(std::size_t element_id)
void setIntegrationPoint(unsigned integration_point)
Eigen::Map< Vector > createZeroedVector(std::vector< double > &data, Eigen::VectorXd::Index size)
Eigen::Map< Matrix > createZeroedMatrix(std::vector< double > &data, Eigen::MatrixXd::Index rows, Eigen::MatrixXd::Index cols)
void shapeFunctionInterpolate(const NodalValues &, const ShapeMatrix &)
Eigen::VectorXd const specific_body_force
an external force that applies in the bulk of the fluid, like gravity.
MaterialPropertyLib::MaterialSpatialDistributionMap media_map

References ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::_element, ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::_integration_method, ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::_ip_data, ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::_process_data, MathLib::createZeroedMatrix(), MathLib::createZeroedVector(), MeshLib::Element::getID(), MaterialPropertyLib::MaterialSpatialDistributionMap::getMedium(), NumLib::GenericIntegrationMethod::getNumberOfPoints(), MaterialPropertyLib::VariableArray::liquid_phase_pressure, ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::liquid_velocity_index, ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::liquid_velocity_size, ProcessLib::StokesFlow::StokesFlowProcessData::media_map, MaterialPropertyLib::permeability, MaterialPropertyLib::Medium::phase(), ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::pressure_index, ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::pressure_size, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), NumLib::detail::shapeFunctionInterpolate(), ProcessLib::StokesFlow::StokesFlowProcessData::specific_body_force, ProcessLib::StokesFlow::StokesFlowProcessData::use_stokes_brinkman_form, and MaterialPropertyLib::viscosity.

◆ computeSecondaryVariableConcrete()

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
void ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::computeSecondaryVariableConcrete ( double const ,
double const ,
Eigen::VectorXd const & local_x,
Eigen::VectorXd const &  )
inlineoverridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 238 of file StokesFlowFEM.h.

243 {
244 auto const local_p =
245 local_x.template segment<pressure_size>(pressure_index);
246
248 ShapeFunctionPressure,
249 typename ShapeFunctionLiquidVelocity::MeshElement, GlobalDim>(
252 }
void interpolateToHigherOrderNodes(MeshLib::Element const &element, bool const is_axially_symmetric, Eigen::MatrixBase< EigenMatrixType > const &node_values, MeshLib::PropertyVector< double > &interpolated_values_global_vector)
MeshLib::PropertyVector< double > * pressure_interpolated

References ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::_element, ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::_is_axially_symmetric, ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::_process_data, NumLib::interpolateToHigherOrderNodes(), ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::pressure_index, and ProcessLib::StokesFlow::StokesFlowProcessData::pressure_interpolated.

◆ getShapeMatrix()

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
Eigen::Map< const Eigen::RowVectorXd > ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::getShapeMatrix ( const unsigned integration_point) const
inlineoverridevirtual

Provides the shape matrix at the given integration point.

Implements NumLib::ExtrapolatableElement.

Definition at line 229 of file StokesFlowFEM.h.

231 {
232 auto const& N_p = _ip_data[integration_point].N_p;
233
234 // assumes N_p is stored contiguously in memory
235 return Eigen::Map<const Eigen::RowVectorXd>(N_p.data(), N_p.size());
236 }

References ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::_ip_data.

Member Data Documentation

◆ _element

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
MeshLib::Element const& ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::_element
private

◆ _integration_method

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
NumLib::GenericIntegrationMethod const& ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::_integration_method
private

◆ _ip_data

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
std::vector<IntegrationPointData<ShapeMatrixTypeLiquidVelocity, ShapeMatrixTypePressure, GlobalDim, ShapeFunctionLiquidVelocity::NPOINTS>, Eigen::aligned_allocator<IntegrationPointData< ShapeMatrixTypeLiquidVelocity, ShapeMatrixTypePressure, GlobalDim, ShapeFunctionLiquidVelocity::NPOINTS> > > ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::_ip_data
private

◆ _is_axially_symmetric

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
bool const ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::_is_axially_symmetric
private

◆ _process_data

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
StokesFlowProcessData const& ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::_process_data
private

◆ liquid_velocity_index

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
const int ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::liquid_velocity_index = 0
staticprivate

◆ liquid_velocity_size

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
const int ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::liquid_velocity_size
staticprivate

◆ pressure_index

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
const int ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::pressure_index
staticprivate

◆ pressure_size

template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , int GlobalDim>
const int ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, GlobalDim >::pressure_size = ShapeFunctionPressure::NPOINTS
staticprivate

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