template<typename ShapeFunctionLiquidVelocity, typename ShapeFunctionPressure, typename IntegrationMethod, int GlobalDim>
class ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, IntegrationMethod, GlobalDim >
Definition at line 33 of file StokesFlowFEM.h.
|
| LocalAssemblerData (MeshLib::Element const &element, std::size_t const, bool const is_axially_symmetric, unsigned const integration_order, 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. More...
|
|
void | computeSecondaryVariableConcrete (double const, double const, Eigen::VectorXd const &local_x, Eigen::VectorXd const &) override |
|
| StokesFlowLocalAssemblerInterface ()=default |
|
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 | assembleWithJacobian (double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_xdot, const double dxdot_dx, const double dx_dx, 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_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...
|
|
template<typename ShapeFunctionLiquidVelocity , typename ShapeFunctionPressure , typename IntegrationMethod , int GlobalDim>
void ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, IntegrationMethod, 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 112 of file StokesFlowFEM.h.
120 assert(local_x.size() == local_matrix_size);
122 auto local_p = Eigen::Map<const NodalVectorType>(
126 typename ShapeMatrixTypeLiquidVelocity::template MatrixType<
127 local_matrix_size, local_matrix_size>>(
128 local_K_data, local_matrix_size, local_matrix_size);
130 typename ShapeMatrixTypeLiquidVelocity::template VectorType<
131 local_matrix_size>>(local_b_data, local_matrix_size);
136 local_K.template block<liquid_velocity_size, liquid_velocity_size>(
139 auto Kvp = local_K.template block<liquid_velocity_size, pressure_size>(
142 auto Kpv = local_K.template block<pressure_size, liquid_velocity_size>(
145 auto bv = local_b.template segment<liquid_velocity_size>(
148 unsigned const n_integration_points =
159 assert(GlobalDim == 2);
160 auto const identity_mat =
161 Eigen::Matrix<double, GlobalDim, GlobalDim>::Identity().eval();
162 auto const unit_vector = Eigen::Map<const Eigen::RowVectorXd>(
163 identity_mat.data(), identity_mat.size())
169 auto const& phase = medium.phase(
"AqueousLiquid");
171 for (
unsigned ip(0); ip < n_integration_points; ++ip)
177 auto const& N_p = ip_data.N_p;
178 auto const& N_v_op = ip_data.N_v_op;
180 auto const& dNdx_v = ip_data.dNdx_v;
181 auto const& dNdx_v_op = ip_data.dNdx_v_op;
182 auto const& dNdx_p = ip_data.dNdx_p;
184 auto const& w = ip_data.integration_weight;
186 double p_int_pt = 0.0;
190 vars[
static_cast<int>(
195 .template value<double>(vars, pos, t, dt);
201 auto const K = MaterialPropertyLib::formEigenTensor<GlobalDim>(
203 .value(vars, pos, t, dt));
206 w * N_v_op.transpose() * mu * K.inverse() * N_v_op;
209 for (
int i = 0; i < GlobalDim; ++i)
215 .noalias() += w * dNdx_v.transpose() * mu * dNdx_v;
219 Kvp.noalias() += w * N_v_op.transpose() * dNdx_p;
222 Kpv.noalias() += w * N_p.transpose() * unit_vector * dNdx_v_op;
225 bv.noalias() += w * N_v_op.transpose() * b;
virtual std::size_t getID() const final
Returns the ID of the element.
void setElementID(std::size_t element_id)
void setIntegrationPoint(unsigned integration_point)
static const int liquid_velocity_index
static const int pressure_size
static const int pressure_index
std::array< VariableType, static_cast< int >(Variable::number_of_variables)> VariableArray
Eigen::Map< Matrix > createZeroedMatrix(std::vector< double > &data, Eigen::MatrixXd::Index rows, Eigen::MatrixXd::Index cols)
Eigen::Map< Vector > createZeroedVector(std::vector< double > &data, Eigen::VectorXd::Index size)
void shapeFunctionInterpolate(const NodalValues &nodal_values, const ShapeMatrix &shape_matrix_N, double &interpolated_value, ScalarTypes &... interpolated_values)
Eigen::VectorXd const specific_body_force
an external force that applies in the bulk of the fluid, like gravity.
std::unique_ptr< MaterialPropertyLib::MaterialSpatialDistributionMap > media_map
bool const use_stokes_brinkman_form
References ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::_element, ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::_integration_method, ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::_ip_data, ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::_process_data, MathLib::createZeroedMatrix(), MathLib::createZeroedVector(), MeshLib::Element::getID(), ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::liquid_velocity_index, ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::liquid_velocity_size, ProcessLib::StokesFlow::StokesFlowProcessData::media_map, MaterialPropertyLib::permeability, MaterialPropertyLib::phase_pressure, ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::pressure_index, ProcessLib::StokesFlow::LocalAssemblerData< ShapeFunctionLiquidVelocity, ShapeFunctionPressure, IntegrationMethod, GlobalDim >::pressure_size, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), NumLib::shapeFunctionInterpolate(), ProcessLib::StokesFlow::StokesFlowProcessData::specific_body_force, ProcessLib::StokesFlow::StokesFlowProcessData::use_stokes_brinkman_form, and MaterialPropertyLib::viscosity.