32 namespace ThermoRichardsMechanics
38 template <
typename ShapeMatrixType>
41 std::vector<ShapeMatrixType, Eigen::aligned_allocator<ShapeMatrixType>>
N_u;
44 template <
typename ShapeFunctionDisplacement,
typename ShapeFunction,
45 typename IntegrationMethod,
int DisplacementDim>
67 ShapeFunctionDisplacement::NPOINTS>;
83 bool const is_axially_symmetric,
84 unsigned const integration_order,
89 std::string
const&
name,
91 int const integration_order)
override;
95 bool const use_monolithic_scheme,
96 int const process_id)
override;
99 std::vector<double>
const& local_x,
100 std::vector<double>
const& local_xdot,
101 const double ,
const double ,
102 std::vector<double>& ,
103 std::vector<double>& ,
104 std::vector<double>& local_rhs_data,
105 std::vector<double>& local_Jac_data)
override;
109 unsigned const n_integration_points =
112 for (
unsigned ip = 0; ip < n_integration_points; ip++)
122 ShapeFunctionDisplacement,
130 double>::quiet_NaN() ,
134 ip_data.pushBackState();
140 double const )
override
142 unsigned const n_integration_points =
145 for (
unsigned ip = 0; ip < n_integration_points; ip++)
152 double const t,
double const dt, Eigen::VectorXd
const& local_x,
153 Eigen::VectorXd
const& local_x_dot)
override;
156 const unsigned integration_point)
const override
161 return Eigen::Map<const Eigen::RowVectorXd>(N_u.data(), N_u.size());
164 std::vector<double>
getSigma()
const override;
168 std::vector<GlobalVector*>
const& x,
169 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
170 std::vector<double>& cache)
const override;
175 std::vector<GlobalVector*>
const& x,
176 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
177 std::vector<double>& cache)
const override;
182 std::vector<GlobalVector*>
const& x,
183 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
184 std::vector<double>& cache)
const override;
189 std::vector<GlobalVector*>
const& x,
190 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
191 std::vector<double>& cache)
const override;
195 std::vector<GlobalVector*>
const& x,
196 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
197 std::vector<double>& cache)
const override;
202 std::vector<GlobalVector*>
const& x,
203 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
204 std::vector<double>& cache)
const override;
206 std::vector<double>
getEpsilon()
const override;
209 std::vector<GlobalVector*>
const& x,
210 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
211 std::vector<double>& cache)
const override;
215 std::vector<GlobalVector*>
const& x,
216 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
217 std::vector<double>& cache)
const override;
223 DisplacementDim>::MaterialStateVariables
const&
229 std::vector<IpData, Eigen::aligned_allocator<IpData>>
ip_data_;
244 ShapeFunctionDisplacement::NPOINTS * DisplacementDim;
virtual std::size_t getID() const final
Returns the ID of the element.
VectorType< _kelvin_vector_size > KelvinVectorType
BMatrixPolicyType< ShapeFunctionDisplacement, DisplacementDim > BMatricesType
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
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix(const unsigned integration_point) const override
Provides the shape matrix at the given integration point.
ShapeMatrixPolicyType< ShapeFunctionDisplacement, DisplacementDim > ShapeMatricesTypeDisplacement
std::vector< double > const & getIntPtSigma(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
MeshLib::Element const & element_
static const int pressure_index
ThermoRichardsMechanicsProcessData< DisplacementDim > & process_data_
void initializeConcrete() override
void computeSecondaryVariableConcrete(double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_dot) override
typename ShapeMatricesType::GlobalDimVectorType GlobalDimVectorType
std::vector< double > getPorosity() const override
typename BMatricesType::KelvinVectorType KelvinVectorType
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
std::vector< double > getEpsilon() 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
static const int displacement_index
SecondaryData< typename ShapeMatricesTypeDisplacement::ShapeMatrices::ShapeType > secondary_data_
static const int pressure_size
unsigned getNumberOfIntegrationPoints() const override
Eigen::Matrix< double, KelvinVectorSize, 1 > SymmetricTensor
static const int temperature_size
static const int displacement_size
std::vector< double > getSwellingStress() const override
void setInitialConditionsConcrete(std::vector< double > const &local_x, double const t, bool const use_monolithic_scheme, int const process_id) override
typename ShapeMatricesType::GlobalDimMatrixType GlobalDimMatrixType
std::vector< double > const & getIntPtEpsilon(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
std::vector< IpData, Eigen::aligned_allocator< IpData > > ip_data_
static const int temperature_index
std::vector< double > getTransportPorosity() const override
ThermoRichardsMechanicsLocalAssembler(ThermoRichardsMechanicsLocalAssembler const &)=delete
std::size_t setIPDataInitialConditions(std::string const &name, double const *values, int const integration_order) override
void assembleWithJacobian(double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_xdot, const double, const double, std::vector< double > &, std::vector< double > &, std::vector< double > &local_rhs_data, std::vector< double > &local_Jac_data) override
static int const KelvinVectorSize
ThermoRichardsMechanicsLocalAssembler(ThermoRichardsMechanicsLocalAssembler &&)=delete
std::vector< double > const & getIntPtTransportPorosity(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 & getIntPtSwellingStress(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
void postTimestepConcrete(Eigen::VectorXd const &, double const, double const) override
IntegrationMethod integration_method_
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
bool const is_axially_symmetric_
MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables const & getMaterialStateVariablesAt(unsigned integration_point) const override
std::vector< double > getSaturation() const override
std::vector< double > getSigma() const override
ShapeMatrixPolicyType< ShapeFunction, DisplacementDim > ShapeMatricesType
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
Eigen::Matrix< double, Eigen::MatrixBase< Derived >::RowsAtCompileTime, 1 > symmetricTensorToKelvinVector(Eigen::MatrixBase< Derived > const &v)
MathLib::TemplatePoint< double, 3 > Point3d
std::array< double, 3 > interpolateCoordinates(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)
MatrixType< GlobalDim, GlobalDim > GlobalDimMatrixType
VectorType< GlobalDim > GlobalDimVectorType
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N_u
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N_u