template<typename ShapeFunction, typename IntegrationMethod, int GlobalDim>
class ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >
Definition at line 32 of file HTFEM.h.
|
| HTFEM (MeshLib::Element const &element, std::size_t const local_matrix_size, bool const is_axially_symmetric, unsigned const integration_order, HTProcessData const &process_data, const unsigned dof_per_node) |
|
Eigen::Map< const Eigen::RowVectorXd > | getShapeMatrix (const unsigned integration_point) const override |
| Provides the shape matrix at the given integration point. More...
|
|
Eigen::Vector3d | getFlux (MathLib::Point3d const &pnt_local_coords, double const t, std::vector< double > const &local_x) const override |
|
| HTLocalAssemblerInterface ()=default |
|
void | setStaggeredCoupledSolutions (std::size_t const, CoupledSolutionsForStaggeredScheme *const coupling_term) |
|
virtual std::vector< double > const & | getIntPtDarcyVelocity (const double, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &) const =0 |
|
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 | assemble (double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_xdot, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data) |
|
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< std::vector< double >> const &) const |
| Fits to staggered scheme. More...
|
|
template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
Eigen::Vector3d ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::getFlux |
( |
MathLib::Point3d const & |
pnt_local_coords, |
|
|
double const |
t, |
|
|
std::vector< double > const & |
local_x |
|
) |
| const |
|
inlineoverridevirtual |
Computes the flux in the point pnt_local_coords
that is given in local coordinates using the values from local_x
.
Implements ProcessLib::HT::HTLocalAssemblerInterface.
Definition at line 93 of file HTFEM.h.
100 auto const shape_matrices =
104 std::array{pnt_local_coords})[0];
112 double T_int_pt = 0.0;
113 double p_int_pt = 0.0;
124 auto const& liquid_phase = medium.phase(
"AqueousLiquid");
128 double const dt = std::numeric_limits<double>::quiet_NaN();
130 auto const K = MaterialPropertyLib::formEigenTensor<GlobalDim>(
132 .value(vars, pos, t, dt));
136 .template value<double>(vars, pos, t, dt);
139 auto const p_nodal_values = Eigen::Map<const NodalVectorType>(
140 &local_x[local_x.size() / 2], ShapeFunction::NPOINTS);
142 -K_over_mu * shape_matrices.dNdx * p_nodal_values;
149 .template value<double>(vars, pos, t, dt);
151 q += K_over_mu * rho_w * b;
154 Eigen::Vector3d flux;
155 flux.head<GlobalDim>() =
q;
virtual std::size_t getID() const final
Returns the ID of the element.
void setElementID(std::size_t element_id)
typename ShapeMatricesType::GlobalDimVectorType GlobalDimVectorType
typename ShapeMatricesType::GlobalDimMatrixType GlobalDimMatrixType
std::array< VariableType, static_cast< int >(Variable::number_of_variables)> VariableArray
std::vector< typename ShapeMatricesType::ShapeMatrices, Eigen::aligned_allocator< typename ShapeMatricesType::ShapeMatrices > > computeShapeMatrices(MeshLib::Element const &e, bool const is_axially_symmetric, PointContainer const &points)
void shapeFunctionInterpolate(const NodalValues &nodal_values, const ShapeMatrix &shape_matrix_N, double &interpolated_value, ScalarTypes &... interpolated_values)
Eigen::VectorXd const specific_body_force
std::unique_ptr< MaterialPropertyLib::MaterialSpatialDistributionMap > media_map
References ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::_element, ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::_process_data, NumLib::computeShapeMatrices(), MaterialPropertyLib::density, MeshLib::Element::getID(), ProcessLib::HT::HTProcessData::has_gravity, ProcessLib::HT::HTProcessData::media_map, MaterialPropertyLib::permeability, MaterialPropertyLib::phase_pressure, MathLib::q, ParameterLib::SpatialPosition::setElementID(), NumLib::shapeFunctionInterpolate(), ProcessLib::HT::HTProcessData::specific_body_force, MaterialPropertyLib::temperature, and MaterialPropertyLib::viscosity.
template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
std::vector<double> const& ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::getIntPtDarcyVelocityLocal |
( |
const double |
t, |
|
|
std::vector< double > const & |
local_x, |
|
|
std::vector< double > & |
cache |
|
) |
| const |
|
inlineprotected |
Definition at line 239 of file HTFEM.h.
243 std::vector<double> local_p{
246 std::vector<double> local_T{
250 auto const n_integration_points =
255 Eigen::Matrix<double, GlobalDim, Eigen::Dynamic, Eigen::RowMajor>>(
256 cache, GlobalDim, n_integration_points);
263 auto const p_nodal_values = Eigen::Map<const NodalVectorType>(
264 &local_p[0], ShapeFunction::NPOINTS);
268 auto const& liquid_phase = medium.phase(
"AqueousLiquid");
270 for (
unsigned ip = 0; ip < n_integration_points; ++ip)
273 auto const&
N = ip_data.N;
274 auto const& dNdx = ip_data.dNdx;
276 pos.setIntegrationPoint(ip);
278 double T_int_pt = 0.0;
279 double p_int_pt = 0.0;
285 vars[
static_cast<int>(
290 double const dt = std::numeric_limits<double>::quiet_NaN();
291 auto const K = MaterialPropertyLib::formEigenTensor<GlobalDim>(
293 .value(vars, pos, t, dt));
298 .template value<double>(vars, pos, t, dt);
301 cache_mat.col(ip).noalias() = -K_over_mu * dNdx * p_nodal_values;
308 .template value<double>(vars, pos, t, dt);
311 cache_mat.col(ip).noalias() += K_over_mu * rho_w * b;
static const int pressure_index
static const int pressure_size
static const int temperature_size
static const int temperature_index
Eigen::Map< Matrix > createZeroedMatrix(std::vector< double > &data, Eigen::MatrixXd::Index rows, Eigen::MatrixXd::Index cols)
References ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::_element, ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::_integration_method, ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::_ip_data, ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::_process_data, MathLib::createZeroedMatrix(), MaterialPropertyLib::density, MeshLib::Element::getID(), ProcessLib::HT::HTProcessData::has_gravity, ProcessLib::HT::HTProcessData::media_map, MaterialPropertyLib::permeability, MaterialPropertyLib::phase_pressure, ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::pressure_index, ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::pressure_size, ParameterLib::SpatialPosition::setElementID(), NumLib::shapeFunctionInterpolate(), ProcessLib::HT::HTProcessData::specific_body_force, MaterialPropertyLib::temperature, ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::temperature_index, ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::temperature_size, and MaterialPropertyLib::viscosity.
Referenced by ProcessLib::HT::MonolithicHTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::getIntPtDarcyVelocity().