35 template <
typename BMatricesType,
typename ShapeMatrixType,
int DisplacementDim>
47 typename ShapeMatrixType::NodalRowVectorType
N;
48 typename ShapeMatrixType::GlobalDimNodalMatrixType
dNdx;
58 DisplacementDim>::MaterialStateVariables>
73 template <
typename DisplacementVectorType>
77 DisplacementVectorType
const& ,
78 double const degradation,
81 auto linear_elastic_mp =
84 .getMaterialProperties();
86 auto const bulk_modulus = linear_elastic_mp.bulk_modulus(
t, x);
87 auto const mu = linear_elastic_mp.mu(
t, x);
89 switch (energy_split_model)
95 MaterialLib::Solids::Phasefield::
96 calculateIsotropicDegradedStress<DisplacementDim>(
104 MaterialLib::Solids::Phasefield::
105 calculateVolDevDegradedStress<DisplacementDim>(
114 calculateIsotropicDegradedStressWithRankineEnergy<
127 template <
typename ShapeMatrixType>
130 std::vector<ShapeMatrixType, Eigen::aligned_allocator<ShapeMatrixType>>
N;
133 template <
typename ShapeFunction,
int DisplacementDim>
139 ShapeFunction::NPOINTS * DisplacementDim;
155 typename ShapeMatricesType::template VectorType<displacement_size>;
157 typename ShapeMatricesType::template VectorType<phasefield_size>;
173 bool const is_axially_symmetric,
180 unsigned const n_integration_points =
183 _ip_data.reserve(n_integration_points);
186 auto& solid_material =
192 DisplacementDim
> const*>(&solid_material))
195 "For phasefield process only linear elastic solid material "
196 "support is implemented.");
199 auto const shape_matrices =
201 DisplacementDim>(e, is_axially_symmetric,
207 for (
unsigned ip = 0; ip < n_integration_points; ip++)
209 _ip_data.emplace_back(solid_material);
211 ip_data.integration_weight =
213 shape_matrices[ip].integralMeasure * shape_matrices[ip].detJ;
215 static const int kelvin_vector_size =
218 ip_data.eps.setZero(kelvin_vector_size);
219 ip_data.eps_prev.resize(kelvin_vector_size);
220 ip_data.C_tensile.setZero(kelvin_vector_size, kelvin_vector_size);
221 ip_data.C_compressive.setZero(kelvin_vector_size,
223 ip_data.sigma_tensile.setZero(kelvin_vector_size);
224 ip_data.sigma_compressive.setZero(kelvin_vector_size);
225 ip_data.sigma.setZero(kelvin_vector_size);
226 ip_data.strain_energy_tensile = 0.0;
227 ip_data.elastic_energy = 0.0;
229 ip_data.N = shape_matrices[ip].N;
230 ip_data.dNdx = shape_matrices[ip].dNdx;
237 std::vector<double>
const& ,
238 std::vector<double>
const& ,
239 std::vector<double>& ,
240 std::vector<double>& ,
241 std::vector<double>& )
override
244 "PhaseFieldLocalAssembler: assembly without jacobian is not "
249 double const t,
double const dt, Eigen::VectorXd
const& local_x,
250 Eigen::VectorXd
const& local_xdot,
int const process_id,
251 std::vector<double>& local_M_data, std::vector<double>& local_K_data,
252 std::vector<double>& local_b_data,
253 std::vector<double>& local_Jac_data)
override;
257 unsigned const n_integration_points =
260 for (
unsigned ip = 0; ip < n_integration_points; ip++)
268 double const )
override
270 unsigned const n_integration_points =
273 for (
unsigned ip = 0; ip < n_integration_points; ip++)
280 std::size_t mesh_item_id,
282 std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
const&
288 std::size_t mesh_item_id,
290 std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
const&
292 GlobalVector const& x,
double const t,
double& elastic_energy,
293 double& surface_energy,
double& pressure_work,
297 const unsigned integration_point)
const override
302 return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
308 std::vector<GlobalVector*>
const& ,
309 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& ,
310 std::vector<double>& cache)
const override
312 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
318 std::vector<GlobalVector*>
const& ,
319 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& ,
320 std::vector<double>& cache)
const override
322 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
327 double const t,
double const dt, Eigen::VectorXd
const& local_x,
328 std::vector<double>& local_b_data, std::vector<double>& local_Jac_data);
331 double const t,
double const dt, Eigen::VectorXd
const& local_x,
332 std::vector<double>& local_b_data, std::vector<double>& local_Jac_data);
336 std::vector<IpData, Eigen::aligned_allocator<IpData>>
_ip_data;
Global vector based on Eigen vector.
virtual std::size_t getID() const final
Returns the ID of the element.
MathLib::WeightedPoint const & getWeightedPoint(unsigned const igp) const
unsigned getNumberOfPoints() const
void setElementID(std::size_t element_id)
VectorType< _number_of_dof > NodalForceVectorType
Rhs residual.
typename ShapeMatricesType::ShapeMatrices ShapeMatrices
void postTimestepConcrete(Eigen::VectorXd const &, double const, double const) override
static constexpr int phasefield_size
static const int phase_process_id
typename ShapeMatricesType::template VectorType< phasefield_size > PhaseFieldVector
PhaseFieldProcessData< DisplacementDim > & _process_data
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 constexpr int displacement_index
void assembleWithJacobianForDeformationEquations(double const t, double const dt, Eigen::VectorXd const &local_x, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
void assembleWithJacobianPhaseFieldEquations(double const t, double const dt, Eigen::VectorXd const &local_x, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
MeshLib::Element const & _element
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data
typename BMatricesType::NodalForceVectorType NodalForceVectorType
NumLib::GenericIntegrationMethod const & _integration_method
typename ShapeMatricesType::template VectorType< displacement_size > DeformationVector
static const int mechanics_process_id
static constexpr int displacement_size
ShapeMatrixPolicyType< ShapeFunction, DisplacementDim > ShapeMatricesType
PhaseFieldLocalAssembler(PhaseFieldLocalAssembler const &)=delete
void assembleWithJacobianForStaggeredScheme(double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_xdot, 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
void initializeConcrete() override
std::vector< IpData, Eigen::aligned_allocator< IpData > > _ip_data
static constexpr int phasefield_index
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix(const unsigned integration_point) const override
Provides the shape matrix at the given integration point.
void computeEnergy(std::size_t mesh_item_id, std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap >> const &dof_tables, GlobalVector const &x, double const t, double &elastic_energy, double &surface_energy, double &pressure_work, CoupledSolutionsForStaggeredScheme const *const cpl_xs) override
typename ShapeMatricesType::template MatrixType< phasefield_size, phasefield_size > PhaseFieldMatrix
bool const _is_axially_symmetric
PhaseFieldLocalAssembler(MeshLib::Element const &e, std::size_t const, NumLib::GenericIntegrationMethod const &integration_method, bool const is_axially_symmetric, PhaseFieldProcessData< DisplacementDim > &process_data)
PhaseFieldLocalAssembler(PhaseFieldLocalAssembler &&)=delete
void computeCrackIntegral(std::size_t mesh_item_id, std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap >> const &dof_tables, GlobalVector const &x, double const t, double &crack_volume, CoupledSolutionsForStaggeredScheme const *const cpl_xs) override
std::vector< double > const & getIntPtEpsilon(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
typename ShapeMatricesType::NodalMatrixType NodalMatrixType
typename ShapeMatricesType::NodalVectorType NodalVectorType
std::vector< double > const & getIntPtSigma(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
MechanicsBase< DisplacementDim > & selectSolidConstitutiveRelation(std::map< int, std::unique_ptr< MechanicsBase< DisplacementDim >>> const &constitutive_relations, MeshLib::PropertyVector< int > const *const material_ids, std::size_t const element_id)
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), kelvin_vector_dimensions(DisplacementDim), Eigen::RowMajor > KelvinMatrixType
std::vector< typename ShapeMatricesType::ShapeMatrices, Eigen::aligned_allocator< typename ShapeMatricesType::ShapeMatrices > > initShapeMatrices(MeshLib::Element const &e, bool const is_axially_symmetric, IntegrationMethod const &integration_method)
MatrixType< ShapeFunction::NPOINTS, ShapeFunction::NPOINTS > NodalMatrixType
NumLib::ShapeMatrices< NodalRowVectorType, DimNodalMatrixType, DimMatrixType, GlobalDimNodalMatrixType > ShapeMatrices
VectorType< ShapeFunction::NPOINTS > NodalVectorType
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N
BMatricesType::KelvinVectorType sigma_tensile
MaterialLib::Solids::MechanicsBase< DisplacementDim > const & solid_material
double integration_weight
ShapeMatrixType::NodalRowVectorType N
BMatricesType::KelvinVectorType eps_prev
double strain_energy_tensile
BMatricesType::KelvinMatrixType C_compressive
BMatricesType::KelvinVectorType sigma_compressive
ShapeMatrixType::GlobalDimNodalMatrixType dNdx
BMatricesType::KelvinVectorType eps
BMatricesType::KelvinVectorType sigma
void updateConstitutiveRelation(double const t, ParameterLib::SpatialPosition const &x, double const, DisplacementVectorType const &, double const degradation, EnergySplitModel const energy_split_model)
IntegrationPointData(MaterialLib::Solids::MechanicsBase< DisplacementDim > const &solid_material)
BMatricesType::KelvinMatrixType C_tensile
double history_variable_prev
std::unique_ptr< typename MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables > material_state_variables
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N