template<typename ShapeFunction, int GlobalDim>
class ProcessLib::HT::HTFEM< ShapeFunction, GlobalDim >
Definition at line 34 of file HTFEM.h.
|
| HTFEM (MeshLib::Element const &element, std::size_t const local_matrix_size, NumLib::GenericIntegrationMethod const &integration_method, bool const is_axially_symmetric, 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.
|
|
Eigen::Vector3d | getFlux (MathLib::Point3d const &pnt_local_coords, double const t, std::vector< double > const &local_x) const override |
|
| HTLocalAssemblerInterface ()=default |
|
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 |
|
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 | assemble (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) |
|
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_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_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< std::vector< double > > const &) const |
| Fits to staggered scheme.
|
|
template<typename ShapeFunction , int GlobalDim>
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 101 of file HTFEM.h.
104 {
105
106
107
108 auto const shape_matrices =
110 GlobalDim>(
112 std::array{pnt_local_coords})[0];
113
116
118
119
120 double T_int_pt = 0.0;
121 double p_int_pt = 0.0;
123 p_int_pt);
124
127
128 auto const& medium =
130 auto const& liquid_phase = medium.
phase(
"AqueousLiquid");
131
132
133
134 double const dt = std::numeric_limits<double>::quiet_NaN();
135
138 .value(vars, pos, t, dt));
139
140 auto const mu =
142 .template value<double>(vars, pos, t, dt);
144
145 auto const p_nodal_values = Eigen::Map<const NodalVectorType>(
146 &local_x[local_x.size() / 2], ShapeFunction::NPOINTS);
148 -K_over_mu * shape_matrices.dNdx * p_nodal_values;
149
151 {
152 auto const rho_w =
153 liquid_phase
155 .template value<double>(vars, pos, t, dt);
156 auto const b =
158 [this->_element.getID()];
159 q += K_over_mu * rho_w * b;
160 }
161
162 Eigen::Vector3d flux;
163 flux.head<GlobalDim>() = q;
164 return flux;
165 }
Medium * getMedium(std::size_t element_id)
Phase const & phase(std::size_t index) const
double liquid_phase_pressure
typename ShapeMatricesType::GlobalDimMatrixType GlobalDimMatrixType
typename ShapeMatricesType::GlobalDimVectorType GlobalDimVectorType
Eigen::Matrix< double, GlobalDim, GlobalDim > formEigenTensor(MaterialPropertyLib::PropertyDataType const &values)
void shapeFunctionInterpolate(const NodalValues &, const ShapeMatrix &)
std::vector< typename ShapeMatricesType::ShapeMatrices, Eigen::aligned_allocator< typename ShapeMatricesType::ShapeMatrices > > computeShapeMatrices(MeshLib::Element const &e, bool const is_axially_symmetric, PointContainer const &points)
MaterialPropertyLib::MaterialSpatialDistributionMap media_map
std::vector< Eigen::VectorXd > const projected_specific_body_force_vectors
Projected specific body force vector: R * R^T * b.
References ProcessLib::HT::HTFEM< ShapeFunction, GlobalDim >::_element, ProcessLib::HT::HTFEM< ShapeFunction, GlobalDim >::_process_data, NumLib::computeShapeMatrices(), MaterialPropertyLib::density, MaterialPropertyLib::formEigenTensor(), MeshLib::Element::getID(), MaterialPropertyLib::MaterialSpatialDistributionMap::getMedium(), ProcessLib::HT::HTProcessData::has_gravity, MaterialPropertyLib::VariableArray::liquid_phase_pressure, ProcessLib::HT::HTProcessData::media_map, MaterialPropertyLib::permeability, MaterialPropertyLib::Medium::phase(), ProcessLib::HT::HTProcessData::projected_specific_body_force_vectors, ParameterLib::SpatialPosition::setElementID(), NumLib::detail::shapeFunctionInterpolate(), MaterialPropertyLib::VariableArray::temperature, and MaterialPropertyLib::viscosity.
template<typename ShapeFunction , int GlobalDim>
std::vector< double > const & ProcessLib::HT::HTFEM< ShapeFunction, 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.
242 {
243 std::vector<double> local_p{
246 std::vector<double> local_T{
249
250 auto const n_integration_points =
252
253 cache.clear();
255 Eigen::Matrix<double, GlobalDim, Eigen::Dynamic, Eigen::RowMajor>>(
256 cache, GlobalDim, n_integration_points);
257
260
262
263 auto const p_nodal_values = Eigen::Map<const NodalVectorType>(
264 &local_p[0], ShapeFunction::NPOINTS);
265
266 auto const& medium =
268 auto const& liquid_phase = medium.
phase(
"AqueousLiquid");
269
270 auto const& Ns =
273
274 for (unsigned ip = 0; ip < n_integration_points; ++ip)
275 {
277 auto const& dNdx = ip_data.dNdx;
278 auto const&
N = Ns[ip];
279
280 pos.setIntegrationPoint(ip);
281
282 double T_int_pt = 0.0;
283 double p_int_pt = 0.0;
286
289
290
291
292 double const dt = std::numeric_limits<double>::quiet_NaN();
295 .value(vars, pos, t, dt));
296
297 auto const mu =
298 liquid_phase
300 .template value<double>(vars, pos, t, dt);
302
303 cache_mat.col(ip).noalias() = -K_over_mu * dNdx * p_nodal_values;
304
306 {
307 auto const rho_w =
308 liquid_phase
310 .template value<double>(vars, pos, t, dt);
311 auto const b =
314
315 cache_mat.col(ip).noalias() += K_over_mu * rho_w * b;
316 }
317 }
318
319 return cache;
320 }
auto const & NsHigherOrder() const
static const int temperature_index
static const int temperature_size
static const int pressure_size
static const int pressure_index
Eigen::Map< Matrix > createZeroedMatrix(std::vector< double > &data, Eigen::MatrixXd::Index rows, Eigen::MatrixXd::Index cols)
NumLib::ShapeMatrixCache shape_matrix_cache
caches for each mesh element type the shape matrix
References ProcessLib::HT::HTFEM< ShapeFunction, GlobalDim >::_element, ProcessLib::HT::HTFEM< ShapeFunction, GlobalDim >::_integration_method, ProcessLib::HT::HTFEM< ShapeFunction, GlobalDim >::_ip_data, ProcessLib::HT::HTFEM< ShapeFunction, GlobalDim >::_process_data, MathLib::createZeroedMatrix(), MaterialPropertyLib::density, MaterialPropertyLib::formEigenTensor(), MeshLib::Element::getID(), MaterialPropertyLib::MaterialSpatialDistributionMap::getMedium(), NumLib::GenericIntegrationMethod::getNumberOfPoints(), ProcessLib::HT::HTProcessData::has_gravity, MaterialPropertyLib::VariableArray::liquid_phase_pressure, ProcessLib::HT::HTProcessData::media_map, NumLib::ShapeMatrixCache::NsHigherOrder(), MaterialPropertyLib::permeability, MaterialPropertyLib::Medium::phase(), ProcessLib::HT::HTFEM< ShapeFunction, GlobalDim >::pressure_index, ProcessLib::HT::HTFEM< ShapeFunction, GlobalDim >::pressure_size, ProcessLib::HT::HTProcessData::projected_specific_body_force_vectors, ParameterLib::SpatialPosition::setElementID(), ProcessLib::HT::HTProcessData::shape_matrix_cache, NumLib::detail::shapeFunctionInterpolate(), MaterialPropertyLib::VariableArray::temperature, ProcessLib::HT::HTFEM< ShapeFunction, GlobalDim >::temperature_index, ProcessLib::HT::HTFEM< ShapeFunction, GlobalDim >::temperature_size, and MaterialPropertyLib::viscosity.
Referenced by ProcessLib::HT::MonolithicHTFEM< ShapeFunction, GlobalDim >::getIntPtDarcyVelocity().