27namespace HydroMechanics
29template <
int DisplacementDim>
33 std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
34 std::vector<std::unique_ptr<ParameterLib::ParameterBase>>
const& parameters,
35 unsigned const integration_order,
36 std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>&&
40 bool const use_monolithic_scheme)
41 :
Process(std::move(name), mesh, std::move(jacobian_assembler), parameters,
42 integration_order, std::move(process_variables),
43 std::move(secondary_variables), use_monolithic_scheme),
44 _process_data(std::move(process_data))
53 std::make_unique<MeshLib::IntegrationPointWriter>(
59 std::make_unique<MeshLib::IntegrationPointWriter>(
68 std::make_unique<MeshLib::IntegrationPointWriter>(
69 "strain_rate_variable_ip", 1, integration_order,
74template <
int DisplacementDim>
80template <
int DisplacementDim>
83 const int process_id)
const
87 if (process_id == _process_data.mechanics_related_process_id)
89 auto const& l = *_local_to_global_index_map;
90 return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
91 &l.getGhostIndices(), &this->_sparsity_pattern};
95 auto const& l = *_local_to_global_index_map_with_base_nodes;
96 return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
97 &l.getGhostIndices(), &_sparsity_pattern_with_linear_element};
100template <
int DisplacementDim>
104 _mesh_subset_all_nodes = std::make_unique<MeshLib::MeshSubset>(
105 _mesh, _mesh.getNodes(), _process_data.use_taylor_hood_elements);
109 _mesh_subset_base_nodes = std::make_unique<MeshLib::MeshSubset>(
110 _mesh, _base_nodes, _process_data.use_taylor_hood_elements);
114 std::vector<MeshLib::MeshSubset> all_mesh_subsets_single_component{
115 *_mesh_subset_all_nodes};
116 _local_to_global_index_map_single_component =
117 std::make_unique<NumLib::LocalToGlobalIndexMap>(
118 std::move(all_mesh_subsets_single_component),
122 if (_process_data.isMonolithicSchemeUsed())
125 std::vector<MeshLib::MeshSubset> all_mesh_subsets{
126 *_mesh_subset_base_nodes};
129 const int monolithic_process_id = 0;
130 std::generate_n(std::back_inserter(all_mesh_subsets),
131 getProcessVariables(monolithic_process_id)[1]
133 .getNumberOfGlobalComponents(),
134 [&]() {
return *_mesh_subset_all_nodes; });
136 std::vector<int>
const vec_n_components{1, DisplacementDim};
137 _local_to_global_index_map =
138 std::make_unique<NumLib::LocalToGlobalIndexMap>(
139 std::move(all_mesh_subsets), vec_n_components,
141 assert(_local_to_global_index_map);
146 const int process_id = 1;
147 std::vector<MeshLib::MeshSubset> all_mesh_subsets;
148 std::generate_n(std::back_inserter(all_mesh_subsets),
149 getProcessVariables(process_id)[0]
151 .getNumberOfGlobalComponents(),
152 [&]() {
return *_mesh_subset_all_nodes; });
154 std::vector<int>
const vec_n_components{DisplacementDim};
155 _local_to_global_index_map =
156 std::make_unique<NumLib::LocalToGlobalIndexMap>(
157 std::move(all_mesh_subsets), vec_n_components,
162 std::vector<MeshLib::MeshSubset> all_mesh_subsets_base_nodes{
163 *_mesh_subset_base_nodes};
164 _local_to_global_index_map_with_base_nodes =
165 std::make_unique<NumLib::LocalToGlobalIndexMap>(
166 std::move(all_mesh_subsets_base_nodes),
171 *_local_to_global_index_map_with_base_nodes, _mesh);
173 assert(_local_to_global_index_map);
174 assert(_local_to_global_index_map_with_base_nodes);
178template <
int DisplacementDim>
182 unsigned const integration_order)
190 auto add_secondary_variable = [&](std::string
const& name,
191 int const num_components,
192 auto get_ip_values_function)
194 _secondary_variables.addSecondaryVariable(
198 std::move(get_ip_values_function)));
201 add_secondary_variable(
"sigma",
203 DisplacementDim>::RowsAtCompileTime,
204 &LocalAssemblerIF::getIntPtSigma);
206 add_secondary_variable(
"epsilon",
208 DisplacementDim>::RowsAtCompileTime,
209 &LocalAssemblerIF::getIntPtEpsilon);
211 add_secondary_variable(
"velocity", DisplacementDim,
212 &LocalAssemblerIF::getIntPtDarcyVelocity);
219 add_secondary_variable);
221 _process_data.pressure_interpolated =
222 MeshLib::getOrCreateMeshProperty<double>(
226 _process_data.principal_stress_vector[0] =
227 MeshLib::getOrCreateMeshProperty<double>(
228 const_cast<MeshLib::Mesh&
>(mesh),
"principal_stress_vector_1",
231 _process_data.principal_stress_vector[1] =
232 MeshLib::getOrCreateMeshProperty<double>(
233 const_cast<MeshLib::Mesh&
>(mesh),
"principal_stress_vector_2",
236 _process_data.principal_stress_vector[2] =
237 MeshLib::getOrCreateMeshProperty<double>(
238 const_cast<MeshLib::Mesh&
>(mesh),
"principal_stress_vector_3",
241 _process_data.principal_stress_values =
242 MeshLib::getOrCreateMeshProperty<double>(
246 _process_data.permeability = MeshLib::getOrCreateMeshProperty<double>(
257 *_local_to_global_index_map);
260template <
int DisplacementDim>
262 std::map<
int, std::shared_ptr<MaterialPropertyLib::Medium>>
const& media)
264 if (_process_data.isMonolithicSchemeUsed())
266 const int process_id_of_hydromechanics = 0;
267 initializeProcessBoundaryConditionsAndSourceTerms(
268 *_local_to_global_index_map, process_id_of_hydromechanics, media);
274 const int hydraulic_process_id = 0;
275 initializeProcessBoundaryConditionsAndSourceTerms(
276 *_local_to_global_index_map_with_base_nodes, hydraulic_process_id,
280 const int mechanical_process_id = 1;
281 initializeProcessBoundaryConditionsAndSourceTerms(
282 *_local_to_global_index_map, mechanical_process_id, media);
285template <
int DisplacementDim>
287 const double t,
double const dt, std::vector<GlobalVector*>
const& x,
288 std::vector<GlobalVector*>
const& x_prev,
int const process_id,
291 DBUG(
"Assemble the equations for HydroMechanics");
297 std::vector<NumLib::LocalToGlobalIndexMap const*> dof_table = {
298 _local_to_global_index_map.get()};
302 getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K,
306template <
int DisplacementDim>
309 const double t,
double const dt, std::vector<GlobalVector*>
const& x,
310 std::vector<GlobalVector*>
const& x_prev,
int const process_id,
314 bool const use_monolithic_scheme = _process_data.isMonolithicSchemeUsed();
315 if (use_monolithic_scheme)
318 "Assemble the Jacobian of HydroMechanics for the monolithic "
324 if (process_id == _process_data.hydraulic_process_id)
327 "Assemble the Jacobian equations of liquid fluid process in "
328 "HydroMechanics for the staggered scheme.");
333 "Assemble the Jacobian equations of mechanical process in "
334 "HydroMechanics for the staggered scheme.");
338 auto const dof_tables = getDOFTables(x.size());
341 _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev,
342 process_id, M, K, b, Jac);
344 auto copyRhs = [&](
int const variable_id,
auto& output_vector)
346 if (use_monolithic_scheme)
348 transformVariableFromGlobalVector(b, variable_id, *dof_tables[0],
350 std::negate<double>());
354 transformVariableFromGlobalVector(b, 0, *dof_tables[process_id],
356 std::negate<double>());
359 if (process_id == _process_data.hydraulic_process_id)
361 copyRhs(0, *_hydraulic_flow);
363 if (process_id == _process_data.mechanics_related_process_id)
365 copyRhs(1, *_nodal_forces);
369template <
int DisplacementDim>
371 std::vector<GlobalVector*>
const& x,
double const t,
double const dt,
372 const int process_id)
374 DBUG(
"PreTimestep HydroMechanicsProcess.");
376 if (hasMechanicalProcess(process_id))
379 &LocalAssemblerIF::preTimestep, _local_assemblers,
380 getActiveElementIDs(), *_local_to_global_index_map, *x[process_id],
385template <
int DisplacementDim>
387 std::vector<GlobalVector*>
const& x,
388 std::vector<GlobalVector*>
const& x_prev,
double const t,
double const dt,
389 const int process_id)
391 if (process_id != _process_data.hydraulic_process_id)
396 DBUG(
"PostTimestep HydroMechanicsProcess.");
399 &LocalAssemblerIF::postTimestep, _local_assemblers,
400 getActiveElementIDs(), getDOFTables(x.size()), x, x_prev, t, dt,
404template <
int DisplacementDim>
406 std::vector<GlobalVector*>
const& x,
407 std::vector<GlobalVector*>
const& x_prev,
const double t,
double const dt,
408 const int process_id)
410 DBUG(
"PostNonLinearSolver HydroMechanicsProcess.");
414 &LocalAssemblerIF::postNonLinearSolver, _local_assemblers,
415 getActiveElementIDs(), getDOFTables(x.size()), x, x_prev, t, dt,
419template <
int DisplacementDim>
423 int const process_id)
426 if (process_id != _process_data.mechanics_related_process_id)
431 DBUG(
"Set initial conditions of HydroMechanicsProcess.");
434 &LocalAssemblerIF::setInitialConditions, _local_assemblers,
435 getActiveElementIDs(), getDOFTables(x.size()), x, t, process_id);
438template <
int DisplacementDim>
440 double const t,
double const dt, std::vector<GlobalVector*>
const& x,
443 if (process_id != _process_data.hydraulic_process_id)
448 DBUG(
"Compute the secondary variables for HydroMechanicsProcess.");
451 &LocalAssemblerIF::computeSecondaryVariable, _local_assemblers,
452 getActiveElementIDs(), getDOFTables(x.size()), t, dt, x, x_prev,
456template <
int DisplacementDim>
457std::tuple<NumLib::LocalToGlobalIndexMap*, bool>
460 const bool manage_storage =
false;
461 return std::make_tuple(_local_to_global_index_map_single_component.get(),
465template <
int DisplacementDim>
469 if (hasMechanicalProcess(process_id))
471 return *_local_to_global_index_map;
475 return *_local_to_global_index_map_with_base_nodes;
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Global vector based on Eigen vector.
bool isAxiallySymmetric() const
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
unsigned getDimension() const
Returns the dimension of the mesh (determined by the maximum dimension over all elements).
Properties & getProperties()
void assembleConcreteProcess(const double t, double const, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b) override
std::vector< std::unique_ptr< LocalAssemblerIF > > _local_assemblers
HydroMechanicsProcess(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, HydroMechanicsProcessData< DisplacementDim > &&process_data, SecondaryVariableCollection &&secondary_variables, bool const use_monolithic_scheme)
void computeSecondaryVariableConcrete(double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_prev, const int process_id) override
void preTimestepConcreteProcess(std::vector< GlobalVector * > const &x, double const t, double const dt, const int process_id) override
void initializeBoundaryConditions(std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media) override
NumLib::LocalToGlobalIndexMap const & getDOFTable(const int process_id) const override
void assembleWithJacobianConcreteProcess(const double t, double const, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) override
void constructDofTable() override
bool isLinear() const override
void postNonLinearSolverConcreteProcess(std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, const double t, double const dt, int const process_id) override
MeshLib::PropertyVector< double > * _nodal_forces
void setInitialConditionsConcreteProcess(std::vector< GlobalVector * > &x, double const t, int const process_id) override
MeshLib::PropertyVector< double > * _hydraulic_flow
MathLib::MatrixSpecifications getMatrixSpecifications(const int process_id) const override
void initializeConcreteProcess(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order) override
Process specific initialization called by initialize().
std::tuple< NumLib::LocalToGlobalIndexMap *, bool > getDOFTableForExtrapolatorData() const override
void postTimestepConcreteProcess(std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, const double t, const double dt, int const process_id) override
std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > _integration_point_writer
const bool _use_monolithic_scheme
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)
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
std::vector< Node * > getBaseNodes(std::vector< Element * > const &elements)
@ 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)
void setIPDataInitialConditions(std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > const &_integration_point_writer, MeshLib::Properties const &mesh_properties, LocalAssemblersVector &local_assemblers)
void createLocalAssemblersHM(std::vector< MeshLib::Element * > const &mesh_elements, NumLib::LocalToGlobalIndexMap const &dof_table, std::vector< std::unique_ptr< LocalAssemblerInterface > > &local_assemblers, ProviderOrOrder const &provider_or_order, ExtraCtorArgs &&... extra_ctor_args)
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)
virtual std::vector< double > getSigma() const =0
virtual std::vector< double > getStrainRateVariable() const =0
virtual std::vector< double > getEpsilon() const =0