32 namespace ThermoMechanicalPhaseField
34 template <
typename BMatricesType,
typename ShapeMatrixType,
int DisplacementDim>
46 typename ShapeMatrixType::NodalRowVectorType
N;
47 typename ShapeMatrixType::GlobalDimNodalMatrixType
dNdx;
55 typename ShapeMatrixType::GlobalDimVectorType
heatflux;
59 DisplacementDim>::MaterialStateVariables>
74 template <
typename DisplacementVectorType>
78 DisplacementVectorType
const& ,
81 double const degradation)
85 auto linear_elastic_mp =
88 .getMaterialProperties();
90 auto const bulk_modulus = linear_elastic_mp.bulk_modulus(t, x);
91 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,
120 ShapeFunction::NPOINTS * DisplacementDim;
138 typename ShapeMatricesType::template VectorType<temperature_size>;
140 typename ShapeMatricesType::template VectorType<displacement_size>;
142 typename ShapeMatricesType::template VectorType<phasefield_size>;
154 bool const is_axially_symmetric,
155 unsigned const integration_order,
157 int const mechanics_related_process_id,
158 int const phase_field_process_id,
159 int const heat_conduction_process_id)
168 unsigned const n_integration_points =
171 _ip_data.reserve(n_integration_points);
174 auto& solid_material =
180 DisplacementDim
> const*>(&solid_material))
183 "For phasefield process only linear elastic solid material "
184 "support is implemented.");
187 auto const shape_matrices =
189 DisplacementDim>(e, is_axially_symmetric,
195 for (
unsigned ip = 0; ip < n_integration_points; ip++)
197 _ip_data.emplace_back(solid_material);
199 ip_data.integration_weight =
201 shape_matrices[ip].integralMeasure * shape_matrices[ip].detJ;
203 static const int kelvin_vector_size =
206 ip_data.eps.setZero(kelvin_vector_size);
207 ip_data.eps_prev.resize(kelvin_vector_size);
208 ip_data.eps_m.setZero(kelvin_vector_size);
209 ip_data.C_tensile.setZero(kelvin_vector_size, kelvin_vector_size);
210 ip_data.C_compressive.setZero(kelvin_vector_size,
212 ip_data.sigma_tensile.setZero(kelvin_vector_size);
213 ip_data.sigma_compressive.setZero(kelvin_vector_size);
214 ip_data.heatflux.setZero(DisplacementDim);
215 ip_data.sigma.setZero(kelvin_vector_size);
216 ip_data.strain_energy_tensile = 0.0;
217 ip_data.elastic_energy = 0.0;
219 ip_data.N = shape_matrices[ip].N;
220 ip_data.dNdx = shape_matrices[ip].dNdx;
227 std::vector<double>
const& ,
228 std::vector<double>
const& ,
229 std::vector<double>& ,
230 std::vector<double>& ,
231 std::vector<double>& )
override
234 "ThermoMechanicalPhaseFieldLocalAssembler: assembly without "
235 "Jacobian is not implemented.");
239 double const t,
double const dt, Eigen::VectorXd
const& local_x,
240 Eigen::VectorXd
const& local_xdot,
const double dxdot_dx,
241 const double dx_dx,
int const process_id,
242 std::vector<double>& local_M_data, std::vector<double>& local_K_data,
243 std::vector<double>& local_b_data,
244 std::vector<double>& local_Jac_data)
override;
248 unsigned const n_integration_points =
251 for (
unsigned ip = 0; ip < n_integration_points; ip++)
258 double const ,
double const )
override
260 unsigned const n_integration_points =
263 for (
unsigned ip = 0; ip < n_integration_points; ip++)
270 const unsigned integration_point)
const override
275 return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
281 std::vector<GlobalVector*>
const& ,
282 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& ,
283 std::vector<double>& cache)
const override
285 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
291 std::vector<GlobalVector*>
const& ,
292 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& ,
293 std::vector<double>& cache)
const override
295 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
301 std::vector<GlobalVector*>
const& ,
302 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& ,
303 std::vector<double>& cache)
const override
307 auto const n_integration_points =
_ip_data.size();
311 double, DisplacementDim, Eigen::Dynamic, Eigen::RowMajor>>(
312 cache, DisplacementDim, n_integration_points);
314 for (
unsigned ip = 0; ip < n_integration_points; ++ip)
316 auto const& heatflux =
_ip_data[ip].heatflux;
318 for (
typename KelvinVectorType::Index component = 0;
319 component < DisplacementDim;
322 cache_mat(component, ip) = heatflux[component];
330 double const t,
double const dt, Eigen::VectorXd
const& local_x,
331 std::vector<double>& local_b_data, std::vector<double>& local_Jac_data);
334 double const t,
double const dt, Eigen::VectorXd
const& local_x,
335 Eigen::VectorXd
const& local_xdot, std::vector<double>& local_b_data,
336 std::vector<double>& local_Jac_data);
339 double const t, Eigen::VectorXd
const& local_x,
340 std::vector<double>& local_b_data, std::vector<double>& local_Jac_data);
344 std::vector<IpData, Eigen::aligned_allocator<IpData>>
_ip_data;
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.
VectorType< _kelvin_vector_size > KelvinVectorType
static constexpr int phasefield_size
void assemble(double const, double const, std::vector< double > const &, std::vector< double > const &, std::vector< double > &, std::vector< double > &, std::vector< double > &) override
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
static constexpr int phasefield_index
std::vector< double > const & getIntPtEpsilon(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
int const _mechanics_related_process_id
ID of the processes that contains mechanical process.
ThermoMechanicalPhaseFieldLocalAssembler(MeshLib::Element const &e, std::size_t const, bool const is_axially_symmetric, unsigned const integration_order, ThermoMechanicalPhaseFieldProcessData< DisplacementDim > &process_data, int const mechanics_related_process_id, int const phase_field_process_id, int const heat_conduction_process_id)
void assembleWithJacobianForPhaseFieldEquations(double const t, Eigen::VectorXd const &local_x, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
bool const _is_axially_symmetric
IntegrationMethod _integration_method
MeshLib::Element const & _element
typename ShapeMatricesType::template VectorType< temperature_size > TemperatureVector
typename ShapeMatricesType::template VectorType< displacement_size > DeformationVector
typename ShapeMatricesType::template VectorType< phasefield_size > PhaseFieldVector
std::vector< double > const & getIntPtSigma(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
void initializeConcrete() override
void assembleWithJacobianForHeatConductionEquations(double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_xdot, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
ThermoMechanicalPhaseFieldLocalAssembler(ThermoMechanicalPhaseFieldLocalAssembler const &)=delete
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data
ThermoMechanicalPhaseFieldProcessData< DisplacementDim > & _process_data
typename BMatricesType::NodalForceVectorType NodalForceVectorType
typename ShapeMatricesType::GlobalDimVectorType GlobalDimVectorType
static constexpr int displacement_index
int const _phase_field_process_id
ID of phase field process.
std::vector< IpData, Eigen::aligned_allocator< IpData > > _ip_data
typename ShapeMatricesType::ShapeMatrices ShapeMatrices
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix(const unsigned integration_point) const override
Provides the shape matrix at the given integration point.
std::vector< double > const & getIntPtHeatFlux(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
static constexpr int temperature_size
void postTimestepConcrete(Eigen::VectorXd const &, double const, double const) override
ShapeMatrixPolicyType< ShapeFunction, DisplacementDim > ShapeMatricesType
ThermoMechanicalPhaseFieldLocalAssembler(ThermoMechanicalPhaseFieldLocalAssembler &&)=delete
static constexpr int displacement_size
int const _heat_conduction_process_id
ID of heat conduction process.
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)
static constexpr int temperature_index
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
Eigen::Map< Matrix > createZeroedMatrix(std::vector< double > &data, Eigen::MatrixXd::Index rows, Eigen::MatrixXd::Index cols)
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)
NumLib::ShapeMatrices< NodalRowVectorType, DimNodalMatrixType, DimMatrixType, GlobalDimNodalMatrixType > ShapeMatrices
VectorType< GlobalDim > GlobalDimVectorType
static Eigen::Matrix< double, KelvinVectorSize, 1 > const identity2
Kelvin mapping of 2nd order identity tensor.
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N
BMatricesType::KelvinVectorType eps_prev
void updateConstitutiveRelation(double const t, ParameterLib::SpatialPosition const &x, double const, DisplacementVectorType const &, double const alpha, double const delta_T, double const degradation)
ShapeMatrixType::GlobalDimVectorType heatflux
std::unique_ptr< typename MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables > material_state_variables
BMatricesType::KelvinVectorType eps
BMatricesType::KelvinVectorType eps_m
BMatricesType::KelvinVectorType sigma_tensile
BMatricesType::KelvinMatrixType C_tensile
ShapeMatrixType::NodalRowVectorType N
double strain_energy_tensile
EIGEN_MAKE_ALIGNED_OPERATOR_NEW
BMatricesType::KelvinVectorType sigma
IntegrationPointData(MaterialLib::Solids::MechanicsBase< DisplacementDim > const &solid_material)
BMatricesType::KelvinMatrixType C_compressive
BMatricesType::KelvinVectorType sigma_compressive
ShapeMatrixType::GlobalDimNodalMatrixType dNdx
MaterialLib::Solids::MechanicsBase< DisplacementDim > const & solid_material
double integration_weight
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N