31namespace ThermoRichardsFlow
35template <
typename ShapeFunction,
int GlobalDim>
57 bool const is_axially_symmetric,
62 std::string_view
const name,
64 int const integration_order)
override;
68 int const process_id)
override;
71 std::vector<double>
const& local_x,
72 std::vector<double>
const& local_x_prev,
73 std::vector<double>& ,
74 std::vector<double>& ,
75 std::vector<double>& local_rhs_data,
76 std::vector<double>& local_Jac_data)
override;
78 void assemble(
double const t,
double const dt,
79 std::vector<double>
const& local_x,
80 std::vector<double>
const& local_x_prev,
81 std::vector<double>& local_M_data,
82 std::vector<double>& local_K_data,
83 std::vector<double>& local_rhs_data)
override;
87 unsigned const n_integration_points =
90 for (
unsigned ip = 0; ip < n_integration_points; ip++)
93 ip_data.pushBackState();
98 Eigen::VectorXd
const& ,
99 double const ,
double const ,
102 unsigned const n_integration_points =
105 for (
unsigned ip = 0; ip < n_integration_points; ip++)
112 double const t,
double const dt, Eigen::VectorXd
const& local_x,
113 Eigen::VectorXd
const& local_x_prev)
override;
116 const unsigned integration_point)
const override
118 auto const& N =
_ip_data[integration_point].N;
121 return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
126 std::vector<GlobalVector*>
const& x,
127 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
128 std::vector<double>& cache)
const override;
133 std::vector<GlobalVector*>
const& x,
134 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
135 std::vector<double>& cache)
const override;
140 std::vector<GlobalVector*>
const& x,
141 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
142 std::vector<double>& cache)
const override;
146 std::vector<GlobalVector*>
const& x,
147 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
148 std::vector<double>& cache)
const override;
156 std::vector<IpData, Eigen::aligned_allocator<IpData>>
_ip_data;
unsigned getNumberOfPoints() const
void postTimestepConcrete(Eigen::VectorXd const &, Eigen::VectorXd const &, double const, double const, int const) override
static const int temperature_size
void assembleWithJacobian(double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_x_prev, std::vector< double > &, std::vector< double > &, std::vector< double > &local_rhs_data, std::vector< double > &local_Jac_data) override
NumLib::GenericIntegrationMethod const & _integration_method
static const int temperature_index
std::vector< double > const & getIntPtSaturation(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
unsigned getNumberOfIntegrationPoints() const override
std::size_t setIPDataInitialConditions(std::string_view const name, double const *values, int const integration_order) override
MeshLib::Element const & _element
ThermoRichardsFlowLocalAssembler(ThermoRichardsFlowLocalAssembler &&)=delete
ThermoRichardsFlowProcessData & _process_data
std::vector< double > const & getIntPtDryDensitySolid(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
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::vector< double > getSaturation() const override
typename ShapeMatricesType::GlobalDimVectorType GlobalDimVectorType
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_rhs_data) override
ThermoRichardsFlowLocalAssembler(ThermoRichardsFlowLocalAssembler const &)=delete
std::vector< double > const & getIntPtPorosity(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
void setInitialConditionsConcrete(Eigen::VectorXd const local_x, double const t, int const process_id) override
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix(const unsigned integration_point) const override
Provides the shape matrix at the given integration point.
std::vector< double > getPorosity() const override
typename ShapeMatricesType::GlobalDimMatrixType GlobalDimMatrixType
bool const _is_axially_symmetric
static const int pressure_index
void computeSecondaryVariableConcrete(double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev) override
void initializeConcrete() override
std::vector< IpData, Eigen::aligned_allocator< IpData > > _ip_data
static const int pressure_size
MatrixType< GlobalDim, GlobalDim > GlobalDimMatrixType
VectorType< GlobalDim > GlobalDimVectorType