15 #include <Eigen/Dense>
35 template <
typename NodalRowVectorType,
typename GlobalDimNodalMatrixType>
39 GlobalDimNodalMatrixType
const& dNdx_,
40 double const& integration_weight_)
45 NodalRowVectorType
const N;
46 GlobalDimNodalMatrixType
const dNdx;
61 std::vector<GlobalVector*>
const& x,
62 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
63 std::vector<double>& cache)
const = 0;
66 template <
typename ShapeFunction,
typename IntegrationMethod,
int GlobalDim>
88 bool const is_axially_symmetric,
89 unsigned const integration_order,
95 unsigned const n_integration_points =
97 _ip_data.reserve(n_integration_points);
99 auto const& shape_matrices =
101 GlobalDim>(element, is_axially_symmetric,
104 for (
unsigned ip = 0; ip < n_integration_points; ip++)
107 shape_matrices[ip].N, shape_matrices[ip].dNdx,
109 shape_matrices[ip].integralMeasure *
110 shape_matrices[ip].detJ);
114 void assemble(
double const t,
double const dt,
115 std::vector<double>
const& local_x,
116 std::vector<double>
const& ,
117 std::vector<double>& local_M_data,
118 std::vector<double>& local_K_data,
119 std::vector<double>& local_b_data)
override;
125 std::vector<double>
const& local_x)
const override;
128 const unsigned integration_point)
const override
130 auto const& N =
_ip_data[integration_point].N;
133 return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
138 std::vector<GlobalVector*>
const& x,
139 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
140 std::vector<double>& velocity_cache)
const override;
148 Eigen::aligned_allocator<
159 Eigen::Map<NodalMatrixType>& local_K,
160 Eigen::Map<NodalVectorType>& local_b,
164 double const rho_L,
DimVectorType const& specific_body_force,
165 bool const has_gravity);
168 Eigen::Map<const NodalVectorType>
const& local_p,
172 double const rho_L,
DimVectorType const& specific_body_force,
173 bool const has_gravity);
183 Eigen::Map<NodalMatrixType>& local_K,
184 Eigen::Map<NodalVectorType>& local_b,
188 double const rho_L,
DimVectorType const& specific_body_force,
189 bool const has_gravity);
192 Eigen::Map<const NodalVectorType>
const& local_p,
196 double const rho_L,
DimVectorType const& specific_body_force,
197 bool const has_gravity);
200 template <
typename LaplacianGravityVelocityCalculator>
202 std::vector<double>
const& local_x,
203 std::vector<double>& local_M_data,
204 std::vector<double>& local_K_data,
205 std::vector<double>& local_b_data);
207 template <
typename LaplacianGravityVelocityCalculator,
208 typename VelocityCacheType>
210 const double t,
const double dt, std::vector<double>
const& local_x,
212 VelocityCacheType& darcy_velocity_at_ips)
const;
214 template <
typename VelocityCacheType>
217 std::vector<double>
const& local_x,
219 VelocityCacheType& darcy_velocity_at_ips)
const;
virtual std::vector< double > const & getIntPtDarcyVelocity(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
std::vector< IntegrationPointData< NodalRowVectorType, GlobalDimNodalMatrixType >, Eigen::aligned_allocator< IntegrationPointData< NodalRowVectorType, GlobalDimNodalMatrixType > > > _ip_data
void assembleMatrixAndVector(double const t, double const dt, std::vector< double > const &local_x, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data)
void assemble(double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data) override
typename LocalAssemblerTraits::LocalMatrix NodalMatrixType
const LiquidFlowData & _process_data
typename ShapeMatricesType::DimVectorType DimVectorType
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix(const unsigned integration_point) const override
Provides the shape matrix at the given integration point.
std::vector< double > const & getIntPtDarcyVelocity(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &velocity_cache) const override
typename ShapeMatricesType::NodalRowVectorType NodalRowVectorType
typename ShapeMatricesType::ShapeMatrices ShapeMatrices
ShapeMatrixPolicyType< ShapeFunction > ShapeMatricesType
typename ShapeMatricesType::GlobalDimNodalMatrixType GlobalDimNodalMatrixType
Eigen::Vector3d getFlux(MathLib::Point3d const &p_local_coords, double const t, std::vector< double > const &local_x) const override
LiquidFlowLocalAssembler(MeshLib::Element const &element, std::size_t const, bool const is_axially_symmetric, unsigned const integration_order, LiquidFlowData const &process_data)
void computeDarcyVelocity(bool const is_scalar_permeability, const double t, const double dt, std::vector< double > const &local_x, ParameterLib::SpatialPosition const &pos, VelocityCacheType &darcy_velocity_at_ips) const
typename ShapeMatricesType::DimMatrixType DimMatrixType
IntegrationMethod const _integration_method
typename LocalAssemblerTraits::LocalVector NodalVectorType
void computeProjectedDarcyVelocity(const double t, const double dt, std::vector< double > const &local_x, ParameterLib::SpatialPosition const &pos, VelocityCacheType &darcy_velocity_at_ips) const
MeshLib::Element const & _element
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)
const unsigned NUM_NODAL_DOF
VectorType< ShapeFunction::DIM > DimVectorType
NumLib::ShapeMatrices< NodalRowVectorType, DimNodalMatrixType, DimMatrixType, GlobalDimNodalMatrixType > ShapeMatrices
MatrixType< ShapeFunction::DIM, ShapeFunction::DIM > DimMatrixType
MatrixType< GlobalDim, ShapeFunction::NPOINTS > GlobalDimNodalMatrixType
RowVectorType< ShapeFunction::NPOINTS > NodalRowVectorType
double const integration_weight
IntegrationPointData(NodalRowVectorType const &N_, GlobalDimNodalMatrixType const &dNdx_, double const &integration_weight_)
GlobalDimNodalMatrixType const dNdx
EIGEN_MAKE_ALIGNED_OPERATOR_NEW
NodalRowVectorType const N
static void calculateLaplacianAndGravityTerm(Eigen::Map< NodalMatrixType > &local_K, Eigen::Map< NodalVectorType > &local_b, IntegrationPointData< NodalRowVectorType, GlobalDimNodalMatrixType > const &ip_data, DimMatrixType const &permeability, double const mu, double const rho_L, DimVectorType const &specific_body_force, bool const has_gravity)
static Eigen::Matrix< double, ShapeFunction::DIM, 1 > calculateVelocity(Eigen::Map< const NodalVectorType > const &local_p, IntegrationPointData< NodalRowVectorType, GlobalDimNodalMatrixType > const &ip_data, DimMatrixType const &permeability, double const mu, double const rho_L, DimVectorType const &specific_body_force, bool const has_gravity)
static void calculateLaplacianAndGravityTerm(Eigen::Map< NodalMatrixType > &local_K, Eigen::Map< NodalVectorType > &local_b, IntegrationPointData< NodalRowVectorType, GlobalDimNodalMatrixType > const &ip_data, DimMatrixType const &permeability, double const mu, double const rho_L, DimVectorType const &specific_body_force, bool const has_gravity)
static Eigen::Matrix< double, ShapeFunction::DIM, 1 > calculateVelocity(Eigen::Map< const NodalVectorType > const &local_p, IntegrationPointData< NodalRowVectorType, GlobalDimNodalMatrixType > const &ip_data, DimMatrixType const &permeability, double const mu, double const rho_L, DimVectorType const &specific_body_force, bool const has_gravity)
Vector< NNodes *NodalDOF > LocalVector
Matrix< NNodes *NodalDOF, NNodes *NodalDOF > LocalMatrix