35 namespace ThermoHydroMechanics
41 template <
typename ShapeMatrixType>
44 std::vector<ShapeMatrixType, Eigen::aligned_allocator<ShapeMatrixType>>
N_u;
47 template <
typename ShapeFunctionDisplacement,
typename ShapeFunctionPressure,
48 typename IntegrationMethod,
int DisplacementDim>
79 bool const is_axially_symmetric,
80 unsigned const integration_order,
85 std::string
const&
name,
87 int const integration_order)
override;
90 std::vector<double>
const& ,
91 std::vector<double>
const& ,
92 std::vector<double>& ,
93 std::vector<double>& ,
94 std::vector<double>& )
override
97 "ThermoHydroMechanicsLocalAssembler: assembly without Jacobian is "
102 std::vector<double>
const& local_x,
103 std::vector<double>
const& local_xdot,
104 const double ,
const double ,
105 std::vector<double>& ,
106 std::vector<double>& ,
107 std::vector<double>& local_rhs_data,
108 std::vector<double>& local_Jac_data)
override;
112 unsigned const n_integration_points =
115 for (
unsigned ip = 0; ip < n_integration_points; ip++)
125 ShapeFunctionDisplacement,
132 std::numeric_limits<double>::quiet_NaN()
139 ip_data.pushBackState();
145 double const )
override
147 unsigned const n_integration_points =
150 for (
unsigned ip = 0; ip < n_integration_points; ip++)
157 double const t,
double const dt, Eigen::VectorXd
const& local_x,
158 Eigen::VectorXd
const& local_x_dot)
override;
161 std::vector<double>
const& local_xdot,
162 double const t,
double const dt,
163 bool const use_monolithic_scheme,
164 int const process_id)
override;
167 const unsigned integration_point)
const override
172 return Eigen::Map<const Eigen::RowVectorXd>(N_u.data(), N_u.size());
177 std::vector<GlobalVector*>
const& x,
178 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
179 std::vector<double>& cache)
const override;
186 constexpr
int kelvin_vector_size =
189 return transposeInPlace<kelvin_vector_size>(
190 [
this](std::vector<double>& values)
197 return ProcessLib::setIntegrationPointKelvinVectorData<DisplacementDim>(
203 std::vector<GlobalVector*>
const& ,
204 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& ,
205 std::vector<double>& cache)
const override
207 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
213 constexpr
int kelvin_vector_size =
216 return transposeInPlace<kelvin_vector_size>(
217 [
this](std::vector<double>& values)
223 std::vector<GlobalVector*>
const& ,
224 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& ,
225 std::vector<double>& cache)
const override
227 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
239 ShapeFunctionDisplacement::NPOINTS>;
240 std::vector<IpData, Eigen::aligned_allocator<IpData>>
_ip_data;
257 ShapeFunctionDisplacement::NPOINTS * DisplacementDim;
virtual std::size_t getID() const final
Returns the ID of the element.
std::size_t setSigma(double const *values)
void assemble(double const, double const, std::vector< double > const &, std::vector< double > const &, std::vector< double > &, std::vector< double > &, std::vector< double > &) override
static const int temperature_size
void initializeConcrete() override
std::vector< double > getEpsilon() const override
ThermoHydroMechanicsLocalAssembler(ThermoHydroMechanicsLocalAssembler &&)=delete
static const int pressure_index
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix(const unsigned integration_point) const override
Provides the shape matrix at the given integration point.
void postNonLinearSolverConcrete(std::vector< double > const &local_x, std::vector< double > const &local_xdot, double const t, double const dt, bool const use_monolithic_scheme, int const process_id) override
void postTimestepConcrete(Eigen::VectorXd const &, double const, double const) override
ThermoHydroMechanicsLocalAssembler(ThermoHydroMechanicsLocalAssembler const &)=delete
static const int pressure_size
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 const int temperature_index
typename ShapeMatricesTypePressure::GlobalDimMatrixType GlobalDimMatrixType
typename ShapeMatricesTypePressure::GlobalDimVectorType GlobalDimVectorType
std::vector< double > getSigma() const override
ShapeMatrixPolicyType< ShapeFunctionDisplacement, DisplacementDim > ShapeMatricesTypeDisplacement
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
virtual std::vector< double > const & getIntPtEpsilon(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
std::vector< IpData, Eigen::aligned_allocator< IpData > > _ip_data
ThermoHydroMechanicsProcessData< DisplacementDim > & _process_data
static const int displacement_size
void computeSecondaryVariableConcrete(double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_dot) override
bool const _is_axially_symmetric
MeshLib::Element const & _element
ShapeMatrixPolicyType< ShapeFunctionPressure, DisplacementDim > ShapeMatricesTypePressure
std::vector< double > const & getIntPtSigma(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
SecondaryData< typename ShapeMatricesTypeDisplacement::ShapeMatrices::ShapeType > _secondary_data
IntegrationMethod _integration_method
static int const KelvinVectorSize
BMatrixPolicyType< ShapeFunctionDisplacement, DisplacementDim > BMatricesType
Eigen::Matrix< double, KelvinVectorSize, 1 > SymmetricTensor
std::size_t setIPDataInitialConditions(std::string const &name, double const *values, int const integration_order) override
Returns number of read integration points.
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
BMatricesType::KelvinVectorType sigma_eff
BMatricesType::KelvinVectorType eps
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N_u