34 template <
typename BMatricesType,
typename ShapeMatrixType,
int DisplacementDim>
46 typename ShapeMatrixType::NodalRowVectorType
N;
47 typename ShapeMatrixType::GlobalDimNodalMatrixType
dNdx;
57 DisplacementDim>::MaterialStateVariables>
72 template <
typename DisplacementVectorType>
76 DisplacementVectorType
const& ,
77 double const degradation)
79 auto linear_elastic_mp =
82 .getMaterialProperties();
84 auto const bulk_modulus = linear_elastic_mp.bulk_modulus(t, x);
85 auto const mu = linear_elastic_mp.mu(t, x);
103 template <
typename ShapeMatrixType>
106 std::vector<ShapeMatrixType, Eigen::aligned_allocator<ShapeMatrixType>>
N;
109 template <
typename ShapeFunction,
typename IntegrationMethod,
116 ShapeFunction::NPOINTS * DisplacementDim;
132 typename ShapeMatricesType::template VectorType<displacement_size>;
134 typename ShapeMatricesType::template VectorType<phasefield_size>;
149 bool const is_axially_symmetric,
150 unsigned const integration_order,
157 unsigned const n_integration_points =
160 _ip_data.reserve(n_integration_points);
163 auto& solid_material =
169 DisplacementDim
> const*>(&solid_material))
172 "For phasefield process only linear elastic solid material "
173 "support is implemented.");
176 auto const shape_matrices =
178 DisplacementDim>(e, is_axially_symmetric,
184 for (
unsigned ip = 0; ip < n_integration_points; ip++)
186 _ip_data.emplace_back(solid_material);
188 ip_data.integration_weight =
190 shape_matrices[ip].integralMeasure * shape_matrices[ip].detJ;
192 static const int kelvin_vector_size =
195 ip_data.eps.setZero(kelvin_vector_size);
196 ip_data.eps_prev.resize(kelvin_vector_size);
197 ip_data.C_tensile.setZero(kelvin_vector_size, kelvin_vector_size);
198 ip_data.C_compressive.setZero(kelvin_vector_size,
200 ip_data.sigma_tensile.setZero(kelvin_vector_size);
201 ip_data.sigma_compressive.setZero(kelvin_vector_size);
202 ip_data.sigma.setZero(kelvin_vector_size);
203 ip_data.strain_energy_tensile = 0.0;
204 ip_data.elastic_energy = 0.0;
206 ip_data.N = shape_matrices[ip].N;
207 ip_data.dNdx = shape_matrices[ip].dNdx;
214 std::vector<double>
const& ,
215 std::vector<double>
const& ,
216 std::vector<double>& ,
217 std::vector<double>& ,
218 std::vector<double>& )
override
221 "PhaseFieldLocalAssembler: assembly without jacobian is not "
226 double const t,
double const dt, Eigen::VectorXd
const& local_x,
227 Eigen::VectorXd
const& local_xdot,
const double dxdot_dx,
228 const double dx_dx,
int const process_id,
229 std::vector<double>& local_M_data, std::vector<double>& local_K_data,
230 std::vector<double>& local_b_data,
231 std::vector<double>& local_Jac_data)
override;
235 unsigned const n_integration_points =
238 for (
unsigned ip = 0; ip < n_integration_points; ip++)
245 double const ,
double const )
override
247 unsigned const n_integration_points =
250 for (
unsigned ip = 0; ip < n_integration_points; ip++)
257 std::size_t mesh_item_id,
259 std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
const&
261 GlobalVector const& x,
double const t,
double& crack_volume,
265 std::size_t mesh_item_id,
267 std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
const&
269 GlobalVector const& x,
double const t,
double& elastic_energy,
270 double& surface_energy,
double& pressure_work,
274 const unsigned integration_point)
const override
279 return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
285 std::vector<GlobalVector*>
const& ,
286 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& ,
287 std::vector<double>& cache)
const override
289 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
295 std::vector<GlobalVector*>
const& ,
296 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& ,
297 std::vector<double>& cache)
const override
299 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
304 double const t,
double const dt, Eigen::VectorXd
const& local_x,
305 std::vector<double>& local_b_data, std::vector<double>& local_Jac_data);
308 double const t,
double const dt, Eigen::VectorXd
const& local_x,
309 std::vector<double>& local_b_data, std::vector<double>& local_Jac_data);
313 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.
void setElementID(std::size_t element_id)
VectorType< _number_of_dof > NodalForceVectorType
Rhs residual.
static constexpr int displacement_index
typename BMatricesType::NodalForceVectorType NodalForceVectorType
std::vector< double > const & getIntPtSigma(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
MeshLib::Element const & _element
typename ShapeMatricesType::template VectorType< displacement_size > DeformationVector
std::vector< double > const & getIntPtEpsilon(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
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
PhaseFieldLocalAssembler(MeshLib::Element const &e, std::size_t const, bool const is_axially_symmetric, unsigned const integration_order, PhaseFieldProcessData< DisplacementDim > &process_data)
typename ShapeMatricesType::NodalMatrixType NodalMatrixType
static const int phase_process_id
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data
typename ShapeMatricesType::template VectorType< phasefield_size > PhaseFieldVector
typename ShapeMatricesType::template MatrixType< phasefield_size, phasefield_size > PhaseFieldMatrix
IntegrationMethod _integration_method
PhaseFieldLocalAssembler(PhaseFieldLocalAssembler const &)=delete
void postTimestepConcrete(Eigen::VectorXd const &, double const, double const) override
ShapeMatrixPolicyType< ShapeFunction, DisplacementDim > ShapeMatricesType
static constexpr int displacement_size
typename ShapeMatricesType::NodalVectorType NodalVectorType
typename ShapeMatricesType::ShapeMatrices ShapeMatrices
static constexpr int phasefield_index
static const int mechanics_process_id
std::vector< IpData, Eigen::aligned_allocator< IpData > > _ip_data
PhaseFieldLocalAssembler(PhaseFieldLocalAssembler &&)=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
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix(const unsigned integration_point) const override
Provides the shape matrix at the given integration point.
void initializeConcrete() override
bool const _is_axially_symmetric
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
void assemble(double const, double const, std::vector< double > const &, std::vector< double > const &, std::vector< double > &, std::vector< double > &, std::vector< double > &) override
PhaseFieldProcessData< DisplacementDim > & _process_data
static constexpr int phasefield_size
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)
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)
std::tuple< MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinVectorType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, MathLib::KelvinVector::KelvinMatrixType< DisplacementDim >, double, double > calculateDegradedStress(double const degradation, double const bulk_modulus, double const mu, MathLib::KelvinVector::KelvinVectorType< DisplacementDim > const &eps)
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
EIGEN_MAKE_ALIGNED_OPERATOR_NEW
ShapeMatrixType::NodalRowVectorType N
BMatricesType::KelvinVectorType eps_prev
double strain_energy_tensile
void updateConstitutiveRelation(double const t, ParameterLib::SpatialPosition const &x, double const, DisplacementVectorType const &, double const degradation)
static constexpr int kelvin_vector_size
BMatricesType::KelvinMatrixType C_compressive
BMatricesType::KelvinVectorType sigma_compressive
ShapeMatrixType::GlobalDimNodalMatrixType dNdx
BMatricesType::KelvinVectorType eps
BMatricesType::KelvinVectorType sigma
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