OGS
ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim > Class Template Reference

Detailed Description

template<typename ShapeFunction, typename IntegrationMethod, int GlobalDim>
class ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >

Definition at line 32 of file HTFEM.h.

#include <HTFEM.h>

Inheritance diagram for ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >:
[legend]
Collaboration diagram for ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >:
[legend]

Public Member Functions

 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
 
- Public Member Functions inherited from ProcessLib::HT::HTLocalAssemblerInterface
 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
 
- Public Member Functions inherited from ProcessLib::LocalAssemblerInterface
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...
 
- Public Member Functions inherited from NumLib::ExtrapolatableElement
virtual ~ExtrapolatableElement ()=default
 

Protected Member Functions

double getHeatEnergyCoefficient (MaterialPropertyLib::VariableArray const &vars, const double porosity, const double fluid_density, const double specific_heat_capacity_fluid, ParameterLib::SpatialPosition const &pos, double const t, double const dt)
 
GlobalDimMatrixType getThermalConductivityDispersivity (MaterialPropertyLib::VariableArray const &vars, const double fluid_density, const double specific_heat_capacity_fluid, const GlobalDimVectorType &velocity, const GlobalDimMatrixType &I, ParameterLib::SpatialPosition const &pos, double const t, double const dt)
 
std::vector< double > const & getIntPtDarcyVelocityLocal (const double t, std::vector< double > const &local_x, std::vector< double > &cache) const
 

Protected Attributes

MeshLib::Element const & _element
 
HTProcessData const & _process_data
 
IntegrationMethod const _integration_method
 
std::vector< IntegrationPointData< NodalRowVectorType, GlobalDimNodalMatrixType >, Eigen::aligned_allocator< IntegrationPointData< NodalRowVectorType, GlobalDimNodalMatrixType > > > _ip_data
 
- Protected Attributes inherited from ProcessLib::HT::HTLocalAssemblerInterface
CoupledSolutionsForStaggeredScheme_coupled_solutions {nullptr}
 

Static Protected Attributes

static const int pressure_index = ShapeFunction::NPOINTS
 
static const int pressure_size = ShapeFunction::NPOINTS
 
static const int temperature_index = 0
 
static const int temperature_size = ShapeFunction::NPOINTS
 

Private Types

using ShapeMatricesType = ShapeMatrixPolicyType< ShapeFunction, GlobalDim >
 
using ShapeMatrices = typename ShapeMatricesType::ShapeMatrices
 
using NodalVectorType = typename ShapeMatricesType::NodalVectorType
 
using NodalRowVectorType = typename ShapeMatricesType::NodalRowVectorType
 
using GlobalDimVectorType = typename ShapeMatricesType::GlobalDimVectorType
 
using GlobalDimNodalMatrixType = typename ShapeMatricesType::GlobalDimNodalMatrixType
 
using GlobalDimMatrixType = typename ShapeMatricesType::GlobalDimMatrixType
 

Member Typedef Documentation

◆ GlobalDimMatrixType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
using ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::GlobalDimMatrixType = typename ShapeMatricesType::GlobalDimMatrixType
private

Definition at line 43 of file HTFEM.h.

◆ GlobalDimNodalMatrixType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
using ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::GlobalDimNodalMatrixType = typename ShapeMatricesType::GlobalDimNodalMatrixType
private

Definition at line 41 of file HTFEM.h.

◆ GlobalDimVectorType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
using ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::GlobalDimVectorType = typename ShapeMatricesType::GlobalDimVectorType
private

Definition at line 40 of file HTFEM.h.

◆ NodalRowVectorType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
using ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::NodalRowVectorType = typename ShapeMatricesType::NodalRowVectorType
private

Definition at line 38 of file HTFEM.h.

◆ NodalVectorType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
using ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::NodalVectorType = typename ShapeMatricesType::NodalVectorType
private

Definition at line 37 of file HTFEM.h.

◆ ShapeMatrices

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
using ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::ShapeMatrices = typename ShapeMatricesType::ShapeMatrices
private

Definition at line 35 of file HTFEM.h.

◆ ShapeMatricesType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
using ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::ShapeMatricesType = ShapeMatrixPolicyType<ShapeFunction, GlobalDim>
private

Definition at line 34 of file HTFEM.h.

Constructor & Destructor Documentation

◆ HTFEM()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::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 
)
inline

Definition at line 46 of file HTFEM.h.

53  _element(element),
54  _process_data(process_data),
55  _integration_method(integration_order)
56  {
57  // This assertion is valid only if all nodal d.o.f. use the same shape
58  // matrices.
59  assert(local_matrix_size == ShapeFunction::NPOINTS * dof_per_node);
60  (void)local_matrix_size;
61  (void)dof_per_node;
62 
63  unsigned const n_integration_points =
64  _integration_method.getNumberOfPoints();
65  _ip_data.reserve(n_integration_points);
66 
67  auto const shape_matrices =
69  GlobalDim>(element, is_axially_symmetric,
71 
72  for (unsigned ip = 0; ip < n_integration_points; ip++)
73  {
74  _ip_data.emplace_back(
75  shape_matrices[ip].N, shape_matrices[ip].dNdx,
76  _integration_method.getWeightedPoint(ip).getWeight() *
77  shape_matrices[ip].integralMeasure *
78  shape_matrices[ip].detJ);
79  }
80  }
IntegrationMethod const _integration_method
Definition: HTFEM.h:163
MeshLib::Element const & _element
Definition: HTFEM.h:160
HTProcessData const & _process_data
Definition: HTFEM.h:161
ShapeMatrixPolicyType< ShapeFunction, GlobalDim > ShapeMatricesType
Definition: HTFEM.h:34
std::vector< IntegrationPointData< NodalRowVectorType, GlobalDimNodalMatrixType >, Eigen::aligned_allocator< IntegrationPointData< NodalRowVectorType, GlobalDimNodalMatrixType > > > _ip_data
Definition: HTFEM.h:168
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::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::_integration_method, ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::_ip_data, and NumLib::initShapeMatrices().

Member Function Documentation

◆ getFlux()

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.

96  {
97  // Eval shape matrices at given point
98  // Note: Axial symmetry is set to false here, because we only need dNdx
99  // here, which is not affected by axial symmetry.
100  auto const shape_matrices =
102  GlobalDim>(
103  _element, false /*is_axially_symmetric*/,
104  std::array{pnt_local_coords})[0];
105 
107  pos.setElementID(this->_element.getID());
108 
110 
111  // local_x contains the local temperature and pressure values
112  double T_int_pt = 0.0;
113  double p_int_pt = 0.0;
114  NumLib::shapeFunctionInterpolate(local_x, shape_matrices.N, T_int_pt,
115  p_int_pt);
116 
117  vars[static_cast<int>(MaterialPropertyLib::Variable::temperature)] =
118  T_int_pt;
119  vars[static_cast<int>(MaterialPropertyLib::Variable::phase_pressure)] =
120  p_int_pt;
121 
122  auto const& medium =
123  *_process_data.media_map->getMedium(_element.getID());
124  auto const& liquid_phase = medium.phase("AqueousLiquid");
125 
126  // TODO (naumov) Temporary value not used by current material models.
127  // Need extension of secondary variables interface.
128  double const dt = std::numeric_limits<double>::quiet_NaN();
129  // fetch permeability, viscosity, density
130  auto const K = MaterialPropertyLib::formEigenTensor<GlobalDim>(
132  .value(vars, pos, t, dt));
133 
134  auto const mu =
135  liquid_phase.property(MaterialPropertyLib::PropertyType::viscosity)
136  .template value<double>(vars, pos, t, dt);
137  GlobalDimMatrixType const K_over_mu = K / mu;
138 
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;
143 
144  if (this->_process_data.has_gravity)
145  {
146  auto const rho_w =
147  liquid_phase
149  .template value<double>(vars, pos, t, dt);
150  auto const b = this->_process_data.specific_body_force;
151  q += K_over_mu * rho_w * b;
152  }
153 
154  Eigen::Vector3d flux;
155  flux.head<GlobalDim>() = q;
156  return flux;
157  }
virtual std::size_t getID() const final
Returns the ID of the element.
Definition: Element.h:82
void setElementID(std::size_t element_id)
typename ShapeMatricesType::GlobalDimVectorType GlobalDimVectorType
Definition: HTFEM.h:40
typename ShapeMatricesType::GlobalDimMatrixType GlobalDimMatrixType
Definition: HTFEM.h:43
std::array< VariableType, static_cast< int >(Variable::number_of_variables)> VariableArray
Definition: VariableType.h:108
static const double q
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)
Definition: Interpolation.h:79
Eigen::VectorXd const specific_body_force
Definition: HTProcessData.h:32
std::unique_ptr< MaterialPropertyLib::MaterialSpatialDistributionMap > media_map
Definition: HTProcessData.h:28

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.

◆ getHeatEnergyCoefficient()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
double ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::getHeatEnergyCoefficient ( MaterialPropertyLib::VariableArray const &  vars,
const double  porosity,
const double  fluid_density,
const double  specific_heat_capacity_fluid,
ParameterLib::SpatialPosition const &  pos,
double const  t,
double const  dt 
)
inlineprotected

Definition at line 170 of file HTFEM.h.

175  {
176  auto const& medium =
177  *_process_data.media_map->getMedium(this->_element.getID());
178  auto const& solid_phase = medium.phase("Solid");
179 
180  auto const specific_heat_capacity_solid =
181  solid_phase
182  .property(
184  .template value<double>(vars, pos, t, dt);
185 
186  auto const solid_density =
187  solid_phase.property(MaterialPropertyLib::PropertyType::density)
188  .template value<double>(vars, pos, t, dt);
189 
190  return solid_density * specific_heat_capacity_solid * (1 - porosity) +
191  fluid_density * specific_heat_capacity_fluid * porosity;
192  }
double fluid_density(const double p, const double T, const double x)

References ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::_process_data, MaterialPropertyLib::density, ProcessLib::TES::fluid_density(), MeshLib::Element::getID(), ProcessLib::HT::HTProcessData::media_map, and MaterialPropertyLib::specific_heat_capacity.

Referenced by ProcessLib::HT::MonolithicHTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::assemble().

◆ getIntPtDarcyVelocityLocal()

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.

242  {
243  std::vector<double> local_p{
244  local_x.data() + pressure_index,
245  local_x.data() + pressure_index + pressure_size};
246  std::vector<double> local_T{
247  local_x.data() + temperature_index,
248  local_x.data() + temperature_index + temperature_size};
249 
250  auto const n_integration_points =
251  _integration_method.getNumberOfPoints();
252 
253  cache.clear();
254  auto cache_mat = MathLib::createZeroedMatrix<
255  Eigen::Matrix<double, GlobalDim, Eigen::Dynamic, Eigen::RowMajor>>(
256  cache, GlobalDim, n_integration_points);
257 
259  pos.setElementID(_element.getID());
260 
262 
263  auto const p_nodal_values = Eigen::Map<const NodalVectorType>(
264  &local_p[0], ShapeFunction::NPOINTS);
265 
266  auto const& medium =
267  *_process_data.media_map->getMedium(_element.getID());
268  auto const& liquid_phase = medium.phase("AqueousLiquid");
269 
270  for (unsigned ip = 0; ip < n_integration_points; ++ip)
271  {
272  auto const& ip_data = _ip_data[ip];
273  auto const& N = ip_data.N;
274  auto const& dNdx = ip_data.dNdx;
275 
276  pos.setIntegrationPoint(ip);
277 
278  double T_int_pt = 0.0;
279  double p_int_pt = 0.0;
280  NumLib::shapeFunctionInterpolate(local_p, N, p_int_pt);
281  NumLib::shapeFunctionInterpolate(local_T, N, T_int_pt);
282 
283  vars[static_cast<int>(MaterialPropertyLib::Variable::temperature)] =
284  T_int_pt;
285  vars[static_cast<int>(
287 
288  // TODO (naumov) Temporary value not used by current material
289  // models. Need extension of secondary variables interface.
290  double const dt = std::numeric_limits<double>::quiet_NaN();
291  auto const K = MaterialPropertyLib::formEigenTensor<GlobalDim>(
293  .value(vars, pos, t, dt));
294 
295  auto const mu =
296  liquid_phase
298  .template value<double>(vars, pos, t, dt);
299  GlobalDimMatrixType const K_over_mu = K / mu;
300 
301  cache_mat.col(ip).noalias() = -K_over_mu * dNdx * p_nodal_values;
302 
304  {
305  auto const rho_w =
306  liquid_phase
308  .template value<double>(vars, pos, t, dt);
309  auto const b = _process_data.specific_body_force;
310  // here it is assumed that the vector b is directed 'downwards'
311  cache_mat.col(ip).noalias() += K_over_mu * rho_w * b;
312  }
313  }
314 
315  return cache;
316  }
static const int pressure_index
Definition: HTFEM.h:319
static const int pressure_size
Definition: HTFEM.h:320
static const int temperature_size
Definition: HTFEM.h:322
static const int temperature_index
Definition: HTFEM.h:321
Eigen::Map< Matrix > createZeroedMatrix(std::vector< double > &data, Eigen::MatrixXd::Index rows, Eigen::MatrixXd::Index cols)
Definition: EigenMapTools.h:32

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().

◆ getShapeMatrix()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
Eigen::Map<const Eigen::RowVectorXd> ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::getShapeMatrix ( const unsigned  integration_point) const
inlineoverridevirtual

Provides the shape matrix at the given integration point.

Implements NumLib::ExtrapolatableElement.

Definition at line 82 of file HTFEM.h.

84  {
85  auto const& N = _ip_data[integration_point].N;
86 
87  // assumes N is stored contiguously in memory
88  return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
89  }

References ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::_ip_data.

◆ getThermalConductivityDispersivity()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
GlobalDimMatrixType ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::getThermalConductivityDispersivity ( MaterialPropertyLib::VariableArray const &  vars,
const double  fluid_density,
const double  specific_heat_capacity_fluid,
const GlobalDimVectorType velocity,
const GlobalDimMatrixType I,
ParameterLib::SpatialPosition const &  pos,
double const  t,
double const  dt 
)
inlineprotected

Definition at line 194 of file HTFEM.h.

200  {
201  auto const& medium =
202  *_process_data.media_map->getMedium(_element.getID());
203 
204  auto const thermal_conductivity =
205  MaterialPropertyLib::formEigenTensor<GlobalDim>(
206  medium
207  .property(
209  .value(vars, pos, t, dt));
210 
211  auto const thermal_dispersivity_longitudinal =
212  medium
215  .template value<double>();
216  auto const thermal_dispersivity_transversal =
217  medium
220  .template value<double>();
221 
222  double const velocity_magnitude = velocity.norm();
223 
224  if (velocity_magnitude < std::numeric_limits<double>::epsilon())
225  {
226  return thermal_conductivity;
227  }
228 
229  GlobalDimMatrixType const thermal_dispersivity =
230  fluid_density * specific_heat_capacity_fluid *
231  (thermal_dispersivity_transversal * velocity_magnitude * I +
232  (thermal_dispersivity_longitudinal -
233  thermal_dispersivity_transversal) /
234  velocity_magnitude * velocity * velocity.transpose());
235 
236  return thermal_conductivity + thermal_dispersivity;
237  }

References ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::_element, ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::_process_data, ProcessLib::TES::fluid_density(), MeshLib::Element::getID(), ProcessLib::HT::HTProcessData::media_map, MaterialPropertyLib::thermal_conductivity, MaterialPropertyLib::thermal_longitudinal_dispersivity, and MaterialPropertyLib::thermal_transversal_dispersivity.

Referenced by ProcessLib::HT::MonolithicHTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::assemble().

Member Data Documentation

◆ _element

◆ _integration_method

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
IntegrationMethod const ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::_integration_method
protected

◆ _ip_data

◆ _process_data

◆ pressure_index

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
const int ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::pressure_index = ShapeFunction::NPOINTS
staticprotected

◆ pressure_size

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
const int ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::pressure_size = ShapeFunction::NPOINTS
staticprotected

◆ temperature_index

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
const int ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::temperature_index = 0
staticprotected

◆ temperature_size

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
const int ProcessLib::HT::HTFEM< ShapeFunction, IntegrationMethod, GlobalDim >::temperature_size = ShapeFunction::NPOINTS
staticprotected

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