14 #include <Eigen/Dense>
34 template <
typename ShapeFunction,
typename IntegrationMethod,
int GlobalDim>
36 :
public HTFEM<ShapeFunction, IntegrationMethod, GlobalDim>
45 typename ShapeMatricesType::template VectorType<
NUM_NODAL_DOF *
46 ShapeFunction::NPOINTS>;
56 std::size_t
const local_matrix_size,
57 bool is_axially_symmetric,
58 unsigned const integration_order,
60 :
HTFEM<ShapeFunction, IntegrationMethod, GlobalDim>(
61 element, local_matrix_size, is_axially_symmetric,
67 std::vector<double>
const& local_x,
68 std::vector<double>
const& ,
69 std::vector<double>& local_M_data,
70 std::vector<double>& local_K_data,
71 std::vector<double>& local_b_data)
override
73 auto const local_matrix_size = local_x.size();
76 assert(local_matrix_size == ShapeFunction::NPOINTS *
NUM_NODAL_DOF);
78 auto local_M = MathLib::createZeroedMatrix<LocalMatrixType>(
79 local_M_data, local_matrix_size, local_matrix_size);
80 auto local_K = MathLib::createZeroedMatrix<LocalMatrixType>(
81 local_K_data, local_matrix_size, local_matrix_size);
82 auto local_b = MathLib::createZeroedVector<LocalVectorType>(
83 local_b_data, local_matrix_size);
85 auto KTT = local_K.template block<temperature_size, temperature_size>(
87 auto MTT = local_M.template block<temperature_size, temperature_size>(
89 auto Kpp = local_K.template block<pressure_size, pressure_size>(
91 auto Mpp = local_M.template block<pressure_size, pressure_size>(
93 auto Bp = local_b.template block<pressure_size, 1>(
pressure_index, 0);
98 auto p_nodal_values = Eigen::Map<const NodalVectorType>(
104 auto const& liquid_phase = medium.phase(
"AqueousLiquid");
105 auto const& solid_phase = medium.phase(
"Solid");
107 auto const& b = process_data.specific_body_force;
110 GlobalDimMatrixType::Identity(GlobalDim, GlobalDim));
114 unsigned const n_integration_points =
117 for (
unsigned ip(0); ip < n_integration_points; ip++)
121 auto const& ip_data = this->
_ip_data[ip];
122 auto const& N = ip_data.N;
123 auto const& dNdx = ip_data.dNdx;
124 auto const& w = ip_data.integration_weight;
126 double T_int_pt = 0.0;
127 double p_int_pt = 0.0;
133 vars[
static_cast<int>(
136 vars[
static_cast<int>(
140 auto const specific_storage =
142 .template value<double>(vars, pos, t, dt);
144 auto const porosity =
146 .template value<double>(vars, pos, t, dt);
150 auto const intrinsic_permeability =
151 MaterialPropertyLib::formEigenTensor<GlobalDim>(
155 .value(vars, pos, t, dt));
157 auto const specific_heat_capacity_fluid =
160 .template value<double>(vars, pos, t, dt);
166 .template value<double>(vars, pos, t, dt);
172 .template value<double>(vars, pos, t, dt);
176 process_data.has_gravity
187 (dNdx.transpose() * thermal_conductivity_dispersivity * dNdx +
188 N.transpose() * velocity.transpose() * dNdx *
fluid_density *
189 specific_heat_capacity_fluid) *
191 Kpp.noalias() += w * dNdx.transpose() * K_over_mu * dNdx;
195 specific_heat_capacity_fluid, pos, t, dt) *
197 Mpp.noalias() += w * N.transpose() * specific_storage * N;
198 if (process_data.has_gravity)
209 std::vector<GlobalVector*>
const& x,
210 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
211 std::vector<double>& cache)
const override
213 int const process_id = 0;
216 assert(!indices.empty());
217 auto const& local_x = x[process_id]->get(indices);
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)
IntegrationMethod const _integration_method
typename ShapeMatricesType::NodalVectorType NodalVectorType
static const int pressure_index
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)
MeshLib::Element const & _element
static const int pressure_size
typename ShapeMatricesType::GlobalDimVectorType GlobalDimVectorType
HTProcessData const & _process_data
typename ShapeMatricesType::ShapeMatrices ShapeMatrices
static const int temperature_index
std::vector< IntegrationPointData< NodalRowVectorType, GlobalDimNodalMatrixType >, Eigen::aligned_allocator< IntegrationPointData< NodalRowVectorType, GlobalDimNodalMatrixType > > > _ip_data
typename ShapeMatricesType::GlobalDimMatrixType GlobalDimMatrixType
typename ShapeMatricesType::NodalRowVectorType NodalRowVectorType
std::vector< double > const & getIntPtDarcyVelocityLocal(const double t, std::vector< double > const &local_x, std::vector< double > &cache) const
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)
MonolithicHTFEM(MeshLib::Element const &element, std::size_t const local_matrix_size, bool is_axially_symmetric, unsigned const integration_order, HTProcessData const &process_data)
typename ShapeMatricesType::template MatrixType< NUM_NODAL_DOF *ShapeFunction::NPOINTS, NUM_NODAL_DOF *ShapeFunction::NPOINTS > LocalMatrixType
typename ShapeMatricesType::template VectorType< NUM_NODAL_DOF *ShapeFunction::NPOINTS > LocalVectorType
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 ShapeMatricesType::GlobalDimVectorType GlobalDimVectorType
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 override
std::array< VariableType, static_cast< int >(Variable::number_of_variables)> VariableArray
void shapeFunctionInterpolate(const NodalValues &nodal_values, const ShapeMatrix &shape_matrix_N, double &interpolated_value, ScalarTypes &... interpolated_values)
std::vector< GlobalIndexType > getIndices(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)
const unsigned NUM_NODAL_DOF
double fluid_density(const double p, const double T, const double x)
NumLib::ShapeMatrices< NodalRowVectorType, DimNodalMatrixType, DimMatrixType, GlobalDimNodalMatrixType > ShapeMatrices
MatrixType< GlobalDim, GlobalDim > GlobalDimMatrixType
VectorType< GlobalDim > GlobalDimVectorType
VectorType< ShapeFunction::NPOINTS > NodalVectorType
RowVectorType< ShapeFunction::NPOINTS > NodalRowVectorType
std::unique_ptr< MaterialPropertyLib::MaterialSpatialDistributionMap > media_map