24 template <
int DisplacementDim>
28 std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
29 std::vector<std::unique_ptr<ParameterLib::ParameterBase>>
const& parameters,
30 unsigned const integration_order,
31 std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>&&
35 bool const use_monolithic_scheme)
36 :
Process(std::move(
name), mesh, std::move(jacobian_assembler), parameters,
37 integration_order, std::move(process_variables),
38 std::move(secondary_variables), use_monolithic_scheme),
39 _process_data(std::move(process_data))
41 if (use_monolithic_scheme)
44 "Monolithic scheme is not implemented for the PhaseField process.");
51 template <
int DisplacementDim>
57 template <
int DisplacementDim>
60 const int process_id)
const
65 auto const& l = *_local_to_global_index_map;
66 return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
67 &l.getGhostIndices(), &this->_sparsity_pattern};
71 auto const& l = *_local_to_global_index_map_single_component;
72 return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
73 &l.getGhostIndices(), &_sparsity_pattern_with_single_component};
76 template <
int DisplacementDim>
83 return *_local_to_global_index_map;
87 return *_local_to_global_index_map_single_component;
90 template <
int DisplacementDim>
94 const int mechanics_process_id = 0;
95 constructDofTableOfSpecifiedProcessStaggeredScheme(mechanics_process_id);
99 std::vector<MeshLib::MeshSubset> all_mesh_subsets_single_component{
100 *_mesh_subset_all_nodes};
101 _local_to_global_index_map_single_component =
102 std::make_unique<NumLib::LocalToGlobalIndexMap>(
103 std::move(all_mesh_subsets_single_component),
107 assert(_local_to_global_index_map_single_component);
111 *_local_to_global_index_map_single_component, _mesh);
114 template <
int DisplacementDim>
118 unsigned const integration_order)
125 _secondary_variables.addSecondaryVariable(
128 DisplacementDim>::RowsAtCompileTime,
129 getExtrapolator(), _local_assemblers,
130 &LocalAssemblerInterface::getIntPtSigma));
132 _secondary_variables.addSecondaryVariable(
135 DisplacementDim>::RowsAtCompileTime,
136 getExtrapolator(), _local_assemblers,
137 &LocalAssemblerInterface::getIntPtEpsilon));
142 *_local_to_global_index_map);
145 template <
int DisplacementDim>
150 const int mechanical_process_id = 0;
151 initializeProcessBoundaryConditionsAndSourceTerms(
152 *_local_to_global_index_map, mechanical_process_id);
154 const int phasefield_process_id = 1;
155 initializeProcessBoundaryConditionsAndSourceTerms(
156 *_local_to_global_index_map_single_component, phasefield_process_id);
159 template <
int DisplacementDim>
161 const double t,
double const dt, std::vector<GlobalVector*>
const& x,
162 std::vector<GlobalVector*>
const& xdot,
int const process_id,
165 DBUG(
"Assemble PhaseFieldProcess.");
167 std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
174 "Assemble the equations of phase field in "
175 "PhaseFieldProcess for the staggered scheme.");
180 "Assemble the equations of deformation in "
181 "PhaseFieldProcess for the staggered scheme.");
183 dof_tables.emplace_back(*_local_to_global_index_map_single_component);
184 dof_tables.emplace_back(*_local_to_global_index_map);
195 template <
int DisplacementDim>
197 const double t,
double const dt, std::vector<GlobalVector*>
const& x,
198 std::vector<GlobalVector*>
const& xdot,
const double dxdot_dx,
202 std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
209 "Assemble the Jacobian equations of phase field in "
210 "PhaseFieldProcess for the staggered scheme.");
215 "Assemble the Jacobian equations of deformation in "
216 "PhaseFieldProcess for the staggered scheme.");
218 dof_tables.emplace_back(*_local_to_global_index_map);
219 dof_tables.emplace_back(*_local_to_global_index_map_single_component);
227 dxdot_dx, dx_dx, process_id, M, K, b, Jac);
232 std::transform(_nodal_forces->begin(), _nodal_forces->end(),
233 _nodal_forces->begin(), [](
double val) { return -val; });
237 template <
int DisplacementDim>
239 std::vector<GlobalVector*>
const& x,
double const t,
double const dt,
240 const int process_id)
242 DBUG(
"PreTimestep PhaseFieldProcess {:d}.", process_id);
244 _process_data.injected_volume = t;
246 _x_previous_timestep =
257 template <
int DisplacementDim>
259 std::vector<GlobalVector*>
const& x,
const double t,
260 const double ,
int const process_id)
262 if (isPhaseFieldProcess(process_id))
264 DBUG(
"PostTimestep PhaseFieldProcess.");
266 _process_data.elastic_energy = 0.0;
267 _process_data.surface_energy = 0.0;
268 _process_data.pressure_work = 0.0;
270 std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
273 dof_tables.emplace_back(*_local_to_global_index_map);
274 dof_tables.emplace_back(*_local_to_global_index_map_single_component);
277 getProcessVariables(process_id)[0];
280 &LocalAssemblerInterface::computeEnergy, _local_assemblers,
282 _process_data.elastic_energy, _process_data.surface_energy,
283 _process_data.pressure_work, _coupled_solutions);
285 INFO(
"Elastic energy: {:g} Surface energy: {:g} Pressure work: {:g} ",
286 _process_data.elastic_energy, _process_data.surface_energy,
287 _process_data.pressure_work);
291 template <
int DisplacementDim>
294 double const ,
const int process_id)
296 _process_data.crack_volume = 0.0;
298 if (!isPhaseFieldProcess(process_id))
300 std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
303 dof_tables.emplace_back(*_local_to_global_index_map);
304 dof_tables.emplace_back(*_local_to_global_index_map_single_component);
306 DBUG(
"PostNonLinearSolver crack volume computation.");
309 getProcessVariables(process_id)[0];
311 &LocalAssemblerInterface::computeCrackIntegral, _local_assemblers,
313 _process_data.crack_volume, _coupled_solutions);
315 INFO(
"Integral of crack: {:g}", _process_data.crack_volume);
317 if (_process_data.hydro_crack)
319 _process_data.pressure_old = _process_data.pressure;
320 _process_data.pressure =
321 _process_data.injected_volume / _process_data.crack_volume;
322 _process_data.pressure_error =
323 std::abs(_process_data.pressure_old - _process_data.pressure) /
324 _process_data.pressure;
325 INFO(
"Internal pressure: {:g} and Pressure error: {:.4e}",
326 _process_data.pressure, _process_data.pressure_error);
327 auto& u = *_coupled_solutions->coupled_xs[0];
329 _process_data.pressure);
334 if (_process_data.hydro_crack)
336 auto& u = *_coupled_solutions->coupled_xs[0];
338 1 / _process_data.pressure);
343 template <
int DisplacementDim>
356 if ((*_x_previous_timestep)[i] > _process_data.irreversible_threshold)
363 template <
int DisplacementDim>
365 int const process_id)
const
367 return process_id == 1;
GlobalMatrix::IndexType GlobalIndexType
void INFO(char const *fmt, Args const &... args)
void DBUG(char const *fmt, Args const &... args)
Global vector based on Eigen vector.
void copyValues(std::vector< double > &u) const
Copy vector values.
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 constructDofTable() override
void postTimestepConcreteProcess(std::vector< GlobalVector * > const &x, const double t, const double delta_t, int const process_id) override
void assembleConcreteProcess(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, 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 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 &xdot, const double dxdot_dx, const double dx_dx, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) override
void initializeBoundaryConditions() override
void postNonLinearSolverConcreteProcess(GlobalVector const &x, GlobalVector const &xdot, const double t, double const dt, int const process_id) 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< std::reference_wrapper< NumLib::LocalToGlobalIndexMap >> const &dof_tables, const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, const double dxdot_dx, const double dx_dx, 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< std::reference_wrapper< NumLib::LocalToGlobalIndexMap >> const &dof_tables, double const t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b)
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
void scale(PETScVector &x, double const a)
void copy(PETScVector const &x, PETScVector &y)
void setLocalAccessibleVector(PETScVector const &x)
@ 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)