33 namespace RichardsMechanics
39 template <
typename ShapeMatrixType>
42 std::vector<ShapeMatrixType, Eigen::aligned_allocator<ShapeMatrixType>>
N_u;
45 template <
typename ShapeFunctionDisplacement,
typename ShapeFunctionPressure,
46 typename IntegrationMethod,
int DisplacementDim>
66 ShapeFunctionDisplacement::NPOINTS>;
81 bool const is_axially_symmetric,
82 unsigned const integration_order,
87 std::string
const&
name,
89 int const integration_order)
override;
93 bool const use_monolithic_scheme,
94 int const process_id)
override;
96 void assemble(
double const t,
double const dt,
97 std::vector<double>
const& local_x,
98 std::vector<double>
const& local_xdot,
99 std::vector<double>& local_M_data,
100 std::vector<double>& local_K_data,
101 std::vector<double>& local_rhs_data)
override;
104 std::vector<double>
const& local_x,
105 std::vector<double>
const& local_xdot,
106 const double ,
const double ,
107 std::vector<double>& ,
108 std::vector<double>& ,
109 std::vector<double>& local_rhs_data,
110 std::vector<double>& local_Jac_data)
override;
113 double const t,
double const dt, Eigen::VectorXd
const& local_x,
114 Eigen::VectorXd
const& local_xdot,
const double dxdot_dx,
115 const double dx_dx,
int const process_id,
116 std::vector<double>& local_M_data, std::vector<double>& local_K_data,
117 std::vector<double>& local_b_data,
118 std::vector<double>& local_Jac_data)
override;
122 unsigned const n_integration_points =
125 for (
unsigned ip = 0; ip < n_integration_points; ip++)
135 ShapeFunctionDisplacement,
143 double>::quiet_NaN() ,
147 ip_data.pushBackState();
153 double const )
override
155 unsigned const n_integration_points =
158 for (
unsigned ip = 0; ip < n_integration_points; ip++)
165 double const t,
double const dt, Eigen::VectorXd
const& local_x,
166 Eigen::VectorXd
const& local_x_dot)
override;
169 const unsigned integration_point)
const override
174 return Eigen::Map<const Eigen::RowVectorXd>(N_u.data(), N_u.size());
177 std::vector<double>
getSigma()
const override;
181 std::vector<GlobalVector*>
const& x,
182 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
183 std::vector<double>& cache)
const override;
188 std::vector<GlobalVector*>
const& x,
189 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
190 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;
209 std::vector<GlobalVector*>
const& x,
210 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
211 std::vector<double>& cache)
const override;
216 std::vector<GlobalVector*>
const& x,
217 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
218 std::vector<double>& cache)
const override;
222 std::vector<GlobalVector*>
const& x,
223 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
224 std::vector<double>& cache)
const override;
229 std::vector<GlobalVector*>
const& x,
230 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
231 std::vector<double>& cache)
const override;
233 std::vector<double>
getEpsilon()
const override;
236 std::vector<GlobalVector*>
const& x,
237 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
238 std::vector<double>& cache)
const override;
243 MaterialStateVariables&)>
const& get_values_span,
244 int const& n_components)
const override;
248 std::vector<GlobalVector*>
const& x,
249 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
250 std::vector<double>& cache)
const override;
279 double const t,
double const dt, Eigen::VectorXd
const& local_x,
280 Eigen::VectorXd
const& local_xdot,
const double dxdot_dx,
281 const double dx_dx, std::vector<double>& local_M_data,
282 std::vector<double>& local_K_data, std::vector<double>& local_b_data,
283 std::vector<double>& local_Jac_data);
314 double const t,
double const dt, Eigen::VectorXd
const& local_x,
315 Eigen::VectorXd
const& local_xdot,
const double dxdot_dx,
316 const double dx_dx, std::vector<double>& local_M_data,
317 std::vector<double>& local_K_data, std::vector<double>& local_b_data,
318 std::vector<double>& local_Jac_data);
323 DisplacementDim>::MaterialStateVariables
const&
329 std::vector<IpData, Eigen::aligned_allocator<IpData>>
_ip_data;
342 ShapeFunctionDisplacement::NPOINTS * DisplacementDim;
virtual std::size_t getID() const final
Returns the ID of the element.
VectorType< _kelvin_vector_size > KelvinVectorType
void postTimestepConcrete(Eigen::VectorXd const &, double const, double const) override
RichardsMechanicsProcessData< DisplacementDim > & _process_data
MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables const & getMaterialStateVariablesAt(unsigned integration_point) const override
bool const _is_axially_symmetric
ShapeMatrixPolicyType< ShapeFunctionDisplacement, DisplacementDim > ShapeMatricesTypeDisplacement
ShapeMatrixPolicyType< ShapeFunctionPressure, DisplacementDim > ShapeMatricesTypePressure
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
std::vector< double > getMicroPressure() const override
static const int displacement_size
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
std::vector< double > getSwellingStress() const override
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
BMatrixPolicyType< ShapeFunctionDisplacement, DisplacementDim > BMatricesType
typename BMatricesType::KelvinVectorType KelvinVectorType
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
std::vector< IpData, Eigen::aligned_allocator< IpData > > _ip_data
void assembleWithJacobianForDeformationEquations(double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd 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)
void computeSecondaryVariableConcrete(double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_dot) override
std::vector< double > getMaterialStateVariableInternalState(std::function< BaseLib::DynamicSpan< double >(typename MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables &)> const &get_values_span, int const &n_components) const override
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
unsigned getNumberOfIntegrationPoints() const override
IntegrationMethod _integration_method
static int const KelvinVectorSize
std::vector< double > getTransportPorosity() const override
std::vector< double > getMicroSaturation() const override
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
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
void setInitialConditionsConcrete(std::vector< double > const &local_x, double const t, bool const use_monolithic_scheme, int const process_id) override
std::vector< double > const & getIntPtMicroPressure(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 > getSigma() const override
std::vector< double > getEpsilon() const override
std::vector< double > getPorosity() const override
static const int pressure_size
void initializeConcrete() override
RichardsMechanicsLocalAssembler(RichardsMechanicsLocalAssembler &&)=delete
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) override
std::vector< double > getSaturation() const override
SecondaryData< typename ShapeMatricesTypeDisplacement::ShapeMatrices::ShapeType > _secondary_data
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix(const unsigned integration_point) const override
Provides the shape matrix at the given integration point.
MeshLib::Element const & _element
std::vector< double > const & getIntPtMicroSaturation(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
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_rhs_data) override
void assembleWithJacobianForPressureEquations(double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd 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)
static const int pressure_index
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
std::size_t setIPDataInitialConditions(std::string const &name, double const *values, int const integration_order) override
typename ShapeMatricesTypePressure::GlobalDimMatrixType GlobalDimMatrixType
RichardsMechanicsLocalAssembler(RichardsMechanicsLocalAssembler const &)=delete
Eigen::Matrix< double, KelvinVectorSize, 1 > SymmetricTensor
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
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
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N_u
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N_u