25template <
int DisplacementDim>
29 std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
30 std::vector<std::unique_ptr<ParameterLib::ParameterBase>>
const& parameters,
31 unsigned const integration_order,
32 std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>&&
36 bool const use_monolithic_scheme)
37 :
Process(std::move(name), mesh, std::move(jacobian_assembler), parameters,
38 integration_order, std::move(process_variables),
39 std::move(secondary_variables), use_monolithic_scheme),
40 _process_data(std::move(process_data))
42 if (use_monolithic_scheme)
45 "Monolithic scheme is not implemented for the PhaseField process.");
52template <
int DisplacementDim>
58template <
int DisplacementDim>
61 const int process_id)
const
66 auto const& l = *_local_to_global_index_map;
67 return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
68 &l.getGhostIndices(), &this->_sparsity_pattern};
72 auto const& l = *_local_to_global_index_map_single_component;
73 return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
74 &l.getGhostIndices(), &_sparsity_pattern_with_single_component};
77template <
int DisplacementDim>
84 return *_local_to_global_index_map;
88 return *_local_to_global_index_map_single_component;
91template <
int DisplacementDim>
95 const int mechanics_process_id = 0;
96 constructDofTableOfSpecifiedProcessStaggeredScheme(mechanics_process_id);
100 std::vector<MeshLib::MeshSubset> all_mesh_subsets_single_component{
101 *_mesh_subset_all_nodes};
102 _local_to_global_index_map_single_component =
103 std::make_unique<NumLib::LocalToGlobalIndexMap>(
104 std::move(all_mesh_subsets_single_component),
108 assert(_local_to_global_index_map_single_component);
112 *_local_to_global_index_map_single_component, _mesh);
115template <
int DisplacementDim>
119 unsigned const integration_order)
127 _secondary_variables.addSecondaryVariable(
130 DisplacementDim>::RowsAtCompileTime,
131 getExtrapolator(), _local_assemblers,
132 &LocalAssemblerInterface::getIntPtSigma));
134 _secondary_variables.addSecondaryVariable(
137 DisplacementDim>::RowsAtCompileTime,
138 getExtrapolator(), _local_assemblers,
139 &LocalAssemblerInterface::getIntPtEpsilon));
141 _secondary_variables.addSecondaryVariable(
144 DisplacementDim>::RowsAtCompileTime,
145 getExtrapolator(), _local_assemblers,
146 &LocalAssemblerInterface::getIntPtSigmaTensile));
148 _secondary_variables.addSecondaryVariable(
151 DisplacementDim>::RowsAtCompileTime,
152 getExtrapolator(), _local_assemblers,
153 &LocalAssemblerInterface::getIntPtSigmaCompressive));
155 _secondary_variables.addSecondaryVariable(
158 DisplacementDim>::RowsAtCompileTime,
159 getExtrapolator(), _local_assemblers,
160 &LocalAssemblerInterface::getIntPtEpsilonTensile));
165 *_local_to_global_index_map);
168template <
int DisplacementDim>
170 std::map<
int, std::shared_ptr<MaterialPropertyLib::Medium>>
const& media)
174 const int mechanical_process_id = 0;
175 initializeProcessBoundaryConditionsAndSourceTerms(
176 *_local_to_global_index_map, mechanical_process_id, media);
178 const int phasefield_process_id = 1;
179 initializeProcessBoundaryConditionsAndSourceTerms(
180 *_local_to_global_index_map_single_component, phasefield_process_id,
184template <
int DisplacementDim>
186 const double t,
double const dt, std::vector<GlobalVector*>
const& x,
187 std::vector<GlobalVector*>
const& x_prev,
int const process_id,
190 DBUG(
"Assemble PhaseFieldProcess.");
192 std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
198 "Assemble the equations of phase field in "
199 "PhaseFieldProcess for the staggered scheme.");
204 "Assemble the equations of deformation in "
205 "PhaseFieldProcess for the staggered scheme.");
207 dof_tables.emplace_back(_local_to_global_index_map_single_component.get());
208 dof_tables.emplace_back(_local_to_global_index_map.get());
219template <
int DisplacementDim>
221 const double t,
double const dt, std::vector<GlobalVector*>
const& x,
222 std::vector<GlobalVector*>
const& x_prev,
int const process_id,
225 std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
231 "Assemble the Jacobian equations of phase field in "
232 "PhaseFieldProcess for the staggered scheme.");
237 "Assemble the Jacobian equations of deformation in "
238 "PhaseFieldProcess for the staggered scheme.");
240 dof_tables.emplace_back(_local_to_global_index_map.get());
241 dof_tables.emplace_back(_local_to_global_index_map_single_component.get());
249 x_prev, process_id, M, K, b, Jac);
254 std::transform(_nodal_forces->begin(), _nodal_forces->end(),
255 _nodal_forces->begin(), [](
double val) { return -val; });
259template <
int DisplacementDim>
261 std::vector<GlobalVector*>
const& x,
double const t,
double const dt,
262 const int process_id)
264 DBUG(
"PreTimestep PhaseFieldProcess {:d}.", process_id);
266 _process_data.injected_volume = t;
268 _x_previous_timestep =
279template <
int DisplacementDim>
281 std::vector<GlobalVector*>
const& x,
282 std::vector<GlobalVector*>
const& ,
const double t,
283 const double ,
int const process_id)
285 if (isPhaseFieldProcess(process_id))
287 DBUG(
"PostTimestep PhaseFieldProcess.");
289 _process_data.elastic_energy = 0.0;
290 _process_data.surface_energy = 0.0;
291 _process_data.pressure_work = 0.0;
293 std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
295 dof_tables.emplace_back(_local_to_global_index_map.get());
296 dof_tables.emplace_back(
297 _local_to_global_index_map_single_component.get());
300 getProcessVariables(process_id)[0];
303 &LocalAssemblerInterface::computeEnergy, _local_assemblers,
305 _process_data.elastic_energy, _process_data.surface_energy,
306 _process_data.pressure_work);
309 double const elastic_energy = _process_data.elastic_energy;
310 MPI_Allreduce(&elastic_energy, &_process_data.elastic_energy, 1,
311 MPI_DOUBLE, MPI_SUM, PETSC_COMM_WORLD);
312 double const surface_energy = _process_data.surface_energy;
313 MPI_Allreduce(&surface_energy, &_process_data.surface_energy, 1,
314 MPI_DOUBLE, MPI_SUM, PETSC_COMM_WORLD);
315 double const pressure_work = _process_data.pressure_work;
316 MPI_Allreduce(&pressure_work, &_process_data.pressure_work, 1,
317 MPI_DOUBLE, MPI_SUM, PETSC_COMM_WORLD);
321 "Elastic energy: {:g} Surface energy: {:g} Pressure work: {:g} at "
323 _process_data.elastic_energy, _process_data.surface_energy,
324 _process_data.pressure_work, t);
325 if (_process_data.propagating_pressurized_crack)
327 INFO(
"Pressure: {:g} at time: {:g} ", _process_data.pressure, t);
332template <
int DisplacementDim>
334 std::vector<GlobalVector*>
const& x,
335 std::vector<GlobalVector*>
const& ,
const double t,
336 double const ,
const int process_id)
338 _process_data.crack_volume = 0.0;
340 if (isPhaseFieldProcess(process_id))
342 if (_process_data.propagating_pressurized_crack)
346 1 / _process_data.pressure);
351 std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
353 dof_tables.emplace_back(_local_to_global_index_map.get());
354 dof_tables.emplace_back(_local_to_global_index_map_single_component.get());
356 DBUG(
"PostNonLinearSolver crack volume computation.");
360 &LocalAssemblerInterface::computeCrackIntegral, _local_assemblers,
364 double const crack_volume = _process_data.crack_volume;
365 MPI_Allreduce(&crack_volume, &_process_data.crack_volume, 1, MPI_DOUBLE,
366 MPI_SUM, PETSC_COMM_WORLD);
369 INFO(
"Integral of crack: {:g}", _process_data.crack_volume);
371 if (_process_data.propagating_pressurized_crack)
373 _process_data.pressure_old = _process_data.pressure;
374 _process_data.pressure =
375 _process_data.injected_volume / _process_data.crack_volume;
376 _process_data.pressure_error =
377 std::abs(_process_data.pressure_old - _process_data.pressure) /
378 _process_data.pressure;
379 INFO(
"Internal pressure: {:g} and Pressure error: {:.4e}",
380 _process_data.pressure, _process_data.pressure_error);
384 _process_data.pressure);
388template <
int DisplacementDim>
401 if ((*_x_previous_timestep)[i] > _process_data.irreversible_threshold)
408template <
int DisplacementDim>
410 int const process_id)
const
412 return process_id == 1;
GlobalMatrix::IndexType GlobalIndexType
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Global vector based on Eigen vector.
void copyValues(std::vector< double > &u) const
void set(IndexType rowId, double v)
set entry
bool isAxiallySymmetric() const
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
virtual void preTimestep(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, double const t, double const delta_t)
virtual void initialize(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)
bool isPhaseFieldProcess(int const process_id) const
NumLib::LocalToGlobalIndexMap const & getDOFTable(const int process_id) const override
void updateConstraints(GlobalVector &lower, GlobalVector &upper, int const process_id) override
MathLib::MatrixSpecifications getMatrixSpecifications(const int process_id) const override
MeshLib::PropertyVector< double > * _nodal_forces
void postNonLinearSolverConcreteProcess(std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, const double t, double const dt, int const process_id) override
void postTimestepConcreteProcess(std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, const double t, const double delta_t, int const process_id) override
void constructDofTable() override
void assembleConcreteProcess(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b) override
void preTimestepConcreteProcess(std::vector< GlobalVector * > const &x, double const t, double const dt, const int process_id) override
bool isLinear() const override
void initializeBoundaryConditions(std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media) override
void initializeConcreteProcess(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order) override
Process specific initialization called by initialize().
void assembleWithJacobianConcreteProcess(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) override
PhaseFieldProcess(std::string name, MeshLib::Mesh &mesh, std::unique_ptr< ProcessLib::AbstractJacobianAssembler > &&jacobian_assembler, std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const ¶meters, unsigned const integration_order, std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > &&process_variables, PhaseFieldProcessData< DisplacementDim > &&process_data, SecondaryVariableCollection &&secondary_variables, bool const use_monolithic_scheme)
std::vector< std::size_t > const & getActiveElementIDs() const
Handles configuration of several secondary variables from the project file.
void assembleWithJacobian(std::size_t const mesh_item_id, LocalAssemblerInterface &local_assembler, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac)
void assemble(std::size_t const mesh_item_id, LocalAssemblerInterface &local_assembler, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, double const t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b)
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
void copy(PETScVector const &x, PETScVector &y)
void setLocalAccessibleVector(PETScVector const &x)
void scale(PETScVector &x, PetscScalar const a)
@ BY_LOCATION
Ordering data by spatial location.
GlobalSparsityPattern computeSparsityPattern(LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh)
Computes a sparsity pattern for the given inputs.
SecondaryVariableFunctions makeExtrapolator(const unsigned num_components, NumLib::Extrapolator &extrapolator, LocalAssemblerCollection const &local_assemblers, typename NumLib::ExtrapolatableLocalAssemblerCollection< LocalAssemblerCollection >::IntegrationPointValuesMethod integration_point_values_method)
static void executeSelectedMemberOnDereferenced(Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)
static void executeSelectedMemberDereferenced(Object &object, Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)
static void executeMemberOnDereferenced(Method method, Container const &container, Args &&... args)