OGS
HydroMechanicsProcess.cpp
Go to the documentation of this file.
1
12
13#include <range/v3/view/join.hpp>
14
22#include "MeshLib/Mesh.h"
23#include "MeshLib/Properties.h"
32
33namespace ProcessLib
34{
35namespace LIE
36{
37namespace HydroMechanics
38{
39template <int DisplacementDim>
41 std::string name,
42 MeshLib::Mesh& mesh,
43 std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
44 std::vector<std::unique_ptr<ParameterLib::ParameterBase>> const& parameters,
45 unsigned const integration_order,
46 std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>&&
47 process_variables,
49 SecondaryVariableCollection&& secondary_variables,
50 bool const use_monolithic_scheme)
51 : Process(std::move(name), mesh, std::move(jacobian_assembler), parameters,
52 integration_order, std::move(process_variables),
53 std::move(secondary_variables), use_monolithic_scheme),
54 _process_data(std::move(process_data))
55{
56 INFO("[LIE/HM] looking for fracture elements in the given mesh");
57 std::vector<std::pair<std::size_t, std::vector<int>>>
58 vec_branch_nodeID_matIDs;
59 std::vector<std::pair<std::size_t, std::vector<int>>>
60 vec_junction_nodeID_matIDs;
64 _vec_fracture_nodes, vec_branch_nodeID_matIDs,
65 vec_junction_nodeID_matIDs);
66
67 if (_vec_fracture_mat_IDs.size() !=
68 _process_data.fracture_properties.size())
69 {
71 "The number of the given fracture properties ({:d}) are not "
72 "consistent with the number of fracture groups in a mesh ({:d}).",
73 _process_data.fracture_properties.size(),
75 }
76
77 // create a map from a material ID to a fracture ID
78 auto max_frac_mat_id = std::max_element(_vec_fracture_mat_IDs.begin(),
80 _process_data.map_materialID_to_fractureID.resize(*max_frac_mat_id + 1);
81 for (unsigned i = 0; i < _vec_fracture_mat_IDs.size(); i++)
82 {
83 _process_data.map_materialID_to_fractureID[_vec_fracture_mat_IDs[i]] =
84 i;
85 }
86
87 // create a table of connected fracture IDs for each element
88 _process_data.vec_ele_connected_fractureIDs.resize(
89 mesh.getNumberOfElements());
90 for (unsigned i = 0; i < _vec_fracture_matrix_elements.size(); i++)
91 {
92 for (auto e : _vec_fracture_matrix_elements[i])
93 {
94 _process_data.vec_ele_connected_fractureIDs[e->getID()].push_back(
95 i);
96 }
97 }
98
99 // set fracture property
100 for (auto& fracture_prop : _process_data.fracture_properties)
101 {
102 // based on the 1st element assuming a fracture forms a straight line
104 DisplacementDim,
105 *_vec_fracture_elements[fracture_prop.fracture_id][0],
106 fracture_prop);
107 }
108
109 // set branches
110 for (auto const& [vec_branch_nodeID, matID] : vec_branch_nodeID_matIDs)
111 {
112 auto const master_matId = matID[0];
113 auto const slave_matId = matID[1];
114 auto& master_frac =
115 _process_data.fracture_properties
116 [_process_data.map_materialID_to_fractureID[master_matId]];
117 auto& slave_frac =
118 _process_data.fracture_properties
119 [_process_data.map_materialID_to_fractureID[slave_matId]];
120
121 master_frac.branches_master.push_back(createBranchProperty(
122 *mesh.getNode(vec_branch_nodeID), master_frac, slave_frac));
123
124 slave_frac.branches_slave.push_back(createBranchProperty(
125 *mesh.getNode(vec_branch_nodeID), master_frac, slave_frac));
126 }
127
128 // set junctions
129 transform(cbegin(vec_junction_nodeID_matIDs),
130 cend(vec_junction_nodeID_matIDs),
131 back_inserter(_vec_junction_nodes),
132 [&](auto& vec_junction_nodeID_matID)
133 {
134 return const_cast<MeshLib::Node*>(
135 _mesh.getNode(vec_junction_nodeID_matID.first));
136 });
137
138 for (std::size_t i = 0; i < vec_junction_nodeID_matIDs.size(); i++)
139 {
140 auto const& material_ids = vec_junction_nodeID_matIDs[i].second;
141 assert(material_ids.size() == 2);
142 std::array<int, 2> fracture_ids{
143 {_process_data.map_materialID_to_fractureID[material_ids[0]],
144 _process_data.map_materialID_to_fractureID[material_ids[1]]}};
145
146 _process_data.junction_properties.emplace_back(
147 i, *mesh.getNode(vec_junction_nodeID_matIDs[i].first),
148 fracture_ids);
149 }
150
151 // create a table of connected junction IDs for each element
152 _process_data.vec_ele_connected_junctionIDs.resize(
153 mesh.getNumberOfElements());
154 for (unsigned i = 0; i < vec_junction_nodeID_matIDs.size(); i++)
155 {
156 auto node = mesh.getNode(vec_junction_nodeID_matIDs[i].first);
157 for (auto id :
159 {
160 _process_data.vec_ele_connected_junctionIDs[id].push_back(i);
161 }
162 }
163
164 // create a table of junction node and connected elements
166 vec_junction_nodeID_matIDs.size());
167 for (unsigned i = 0; i < vec_junction_nodeID_matIDs.size(); i++)
168 {
169 auto node = mesh.getNode(vec_junction_nodeID_matIDs[i].first);
170 for (auto e : mesh.getElementsConnectedToNode(*node))
171 {
173 const_cast<MeshLib::Element*>(e));
174 }
175 }
176
177 //
178 // If Neumann BCs for the displacement_jump variable are required they need
179 // special treatment because of the levelset function. The implementation
180 // exists in the version 6.1.0 (e54815cc07ee89c81f953a4955b1c788595dd725)
181 // and was removed due to lack of applications.
182 //
183
184 if (!_process_data.deactivate_matrix_in_flow)
185 {
186 _process_data.p_element_status =
187 std::make_unique<MeshLib::ElementStatus>(&mesh);
188 }
189 else
190 {
191 auto const range =
193 if (!range)
194 {
195 OGS_FATAL(
196 "Could not get minimum/maximum ranges values for the "
197 "MaterialIDs property in the mesh '{:s}'.",
198 mesh.getName());
199 }
200
201 std::vector<int> vec_p_inactive_matIDs;
202 for (int matID = range->first; matID <= range->second; matID++)
203 {
204 if (std::find(_vec_fracture_mat_IDs.begin(),
206 matID) == _vec_fracture_mat_IDs.end())
207 {
208 vec_p_inactive_matIDs.push_back(matID);
209 }
210 }
211 _process_data.p_element_status =
212 std::make_unique<MeshLib::ElementStatus>(&mesh,
213 vec_p_inactive_matIDs);
214
215 const int monolithic_process_id = 0;
216 ProcessVariable const& pv_p =
217 getProcessVariables(monolithic_process_id)[0];
219 }
220 MeshLib::PropertyVector<int> const* material_ids(
221 mesh.getProperties().getPropertyVector<int>("MaterialIDs"));
222 _process_data.mesh_prop_materialIDs = material_ids;
223}
224
225template <int DisplacementDim>
227{
228 //------------------------------------------------------------
229 // prepare mesh subsets to define DoFs
230 //------------------------------------------------------------
231 // for extrapolation
232 _mesh_subset_all_nodes =
233 std::make_unique<MeshLib::MeshSubset>(_mesh, _mesh.getNodes());
234 // pressure
235 _mesh_nodes_p = MeshLib::getBaseNodes(
236 _process_data.p_element_status->getActiveElements());
237 _mesh_subset_nodes_p =
238 std::make_unique<MeshLib::MeshSubset>(_mesh, _mesh_nodes_p);
239 // regular u
240 _mesh_subset_matrix_nodes =
241 std::make_unique<MeshLib::MeshSubset>(_mesh, _mesh.getNodes());
242 // u jump
243 for (unsigned i = 0; i < _vec_fracture_nodes.size(); i++)
244 {
245 _mesh_subset_fracture_nodes.push_back(
246 std::make_unique<MeshLib::MeshSubset const>(
247 _mesh, _vec_fracture_nodes[i]));
248 }
249 // enrichment for junctions
250 _mesh_subset_junction_nodes =
251 std::make_unique<MeshLib::MeshSubset>(_mesh, _vec_junction_nodes);
252
253 // Collect the mesh subsets in a vector. (pressure, displacement,
254 // displacement_jump_fracture, and displacement_jump_junction)
255 std::vector<MeshLib::MeshSubset> all_mesh_subsets;
256 all_mesh_subsets.emplace_back(*_mesh_subset_nodes_p);
257 std::generate_n(std::back_inserter(all_mesh_subsets), DisplacementDim,
258 [&]() { return *_mesh_subset_matrix_nodes; });
259 for (auto const& ms : _mesh_subset_fracture_nodes)
260 {
261 std::generate_n(std::back_inserter(all_mesh_subsets),
262 DisplacementDim,
263 [&]() { return *ms; });
264 }
265 std::generate_n(std::back_inserter(all_mesh_subsets),
266 DisplacementDim,
267 [&]() { return *_mesh_subset_junction_nodes; });
268
269 // The corresponding number of components for (pressure, displacement,
270 // displacement_jump_fracture, and displacement_jump_junction).
271 std::vector<int> vec_n_components;
272 vec_n_components.push_back(1); // pressure
273 vec_n_components.insert(
274 vec_n_components.end(),
275 1 + _vec_fracture_mat_IDs.size() + _vec_junction_nodes.size(),
276 DisplacementDim); // all displacements
277
278 auto const all_fracture_matrix_elements = _vec_fracture_matrix_elements |
279 ranges::views::join |
280 ranges::to<std::vector>();
281 std::vector<std::vector<MeshLib::Element*> const*> vec_var_elements;
282 if (!_process_data.deactivate_matrix_in_flow)
283 {
284 vec_var_elements.push_back(&_mesh.getElements());
285 }
286 else
287 {
288 // TODO set elements including active nodes for pressure.
289 // cannot use ElementStatus
290 vec_var_elements.push_back(&all_fracture_matrix_elements);
291 }
292 vec_var_elements.push_back(&_vec_matrix_elements);
293 for (unsigned i = 0; i < _vec_fracture_matrix_elements.size(); i++)
294 {
295 vec_var_elements.push_back(&_vec_fracture_matrix_elements[i]);
296 }
297 for (unsigned i = 0; i < _vec_junction_fracture_matrix_elements.size(); i++)
298 {
299 vec_var_elements.push_back(&_vec_junction_fracture_matrix_elements[i]);
300 }
301
302 INFO("[LIE/HM] creating a DoF table");
303 _local_to_global_index_map =
304 std::make_unique<NumLib::LocalToGlobalIndexMap>(
305 std::move(all_mesh_subsets),
306 vec_n_components,
307 vec_var_elements,
309
310 DBUG("[LIE/HM] created {:d} DoF", _local_to_global_index_map->size());
311}
312
313template <int DisplacementDim>
315 MeshLib::Element const& e, MeshLib::Mesh& mesh)
316{
317 Eigen::Vector3d const pt(getCenterOfGravity(e).asEigenVector3d());
318 std::vector<FractureProperty*> e_fracture_props;
319 std::unordered_map<int, int> e_fracID_to_local;
320 unsigned tmpi = 0;
321 for (auto fid : _process_data.vec_ele_connected_fractureIDs[e.getID()])
322 {
323 e_fracture_props.push_back(&_process_data.fracture_properties[fid]);
324 e_fracID_to_local.insert({fid, tmpi++});
325 }
326 std::vector<JunctionProperty*> e_junction_props;
327 std::unordered_map<int, int> e_juncID_to_local;
328 tmpi = 0;
329 for (auto fid : _process_data.vec_ele_connected_junctionIDs[e.getID()])
330 {
331 e_junction_props.push_back(&_process_data.junction_properties[fid]);
332 e_juncID_to_local.insert({fid, tmpi++});
333 }
334 std::vector<double> const levelsets(uGlobalEnrichments(
335 e_fracture_props, e_junction_props, e_fracID_to_local, pt));
336
337 auto update_levelset_property = [&](unsigned const i, int const id,
338 unsigned const levelset_idx_offset,
339 unsigned const name_offset)
340 {
341 auto levelset_property = MeshLib::getOrCreateMeshProperty<double>(
342 const_cast<MeshLib::Mesh&>(mesh),
343 "levelset" + std::to_string(id + 1 + name_offset),
345 levelset_property->resize(mesh.getNumberOfElements());
346 (*levelset_property)[e.getID()] = levelsets[i + levelset_idx_offset];
347 };
348
349 for (unsigned i = 0; i < e_fracture_props.size(); i++)
350 {
351 update_levelset_property(i, e_fracture_props[i]->fracture_id, 0, 0);
352 }
353 for (unsigned i = 0; i < e_junction_props.size(); i++)
354 {
355 update_levelset_property(i, e_junction_props[i]->junction_id,
356 e_fracture_props.size(),
357 _process_data.fracture_properties.size());
358 }
359}
360
361template <int DisplacementDim>
363 NumLib::LocalToGlobalIndexMap const& dof_table,
364 MeshLib::Mesh const& mesh,
365 unsigned const integration_order)
366{
367 assert(mesh.getDimension() == DisplacementDim);
368 INFO("[LIE/HM] creating local assemblers");
370 DisplacementDim, HydroMechanicsLocalAssemblerMatrix,
373 mesh.getElements(), dof_table, _local_assemblers,
374 NumLib::IntegrationOrder{integration_order}, mesh.isAxiallySymmetric(),
375 _process_data);
376
377 auto add_secondary_variable = [&](std::string const& name,
378 int const num_components,
379 auto get_ip_values_function)
380 {
381 _secondary_variables.addSecondaryVariable(
382 name,
383 makeExtrapolator(num_components, getExtrapolator(),
384 _local_assemblers,
385 std::move(get_ip_values_function)));
386 };
387
388 add_secondary_variable("sigma",
390 DisplacementDim>::RowsAtCompileTime,
391 &LocalAssemblerInterface::getIntPtSigma);
392
393 add_secondary_variable("epsilon",
395 DisplacementDim>::RowsAtCompileTime,
396 &LocalAssemblerInterface::getIntPtEpsilon);
397
398 add_secondary_variable("velocity", DisplacementDim,
399 &LocalAssemblerInterface::getIntPtDarcyVelocity);
400
401 add_secondary_variable("fracture_velocity", DisplacementDim,
402 &LocalAssemblerInterface::getIntPtFractureVelocity);
403
404 add_secondary_variable("fracture_stress", DisplacementDim,
405 &LocalAssemblerInterface::getIntPtFractureStress);
406
407 add_secondary_variable("fracture_aperture", 1,
408 &LocalAssemblerInterface::getIntPtFractureAperture);
409
410 add_secondary_variable(
411 "fracture_permeability", 1,
412 &LocalAssemblerInterface::getIntPtFracturePermeability);
413
414 _process_data.element_stresses = MeshLib::getOrCreateMeshProperty<double>(
415 const_cast<MeshLib::Mesh&>(mesh), "sigma_avg",
418 DisplacementDim>::RowsAtCompileTime);
419
420 _process_data.element_velocities = MeshLib::getOrCreateMeshProperty<double>(
421 const_cast<MeshLib::Mesh&>(mesh), "velocity_avg",
422 MeshLib::MeshItemType::Cell, DisplacementDim);
423
424 for (MeshLib::Element const* e : _mesh.getElements())
425 {
426 if (e->getDimension() < DisplacementDim)
427 {
428 continue;
429 }
430
431 updateElementLevelSets(*e, const_cast<MeshLib::Mesh&>(mesh));
432 }
433
434 _process_data.element_local_jumps =
436 const_cast<MeshLib::Mesh&>(mesh), "local_jump_w_avg",
437 MeshLib::MeshItemType::Cell, DisplacementDim);
438
439 _process_data.element_fracture_stresses =
441 const_cast<MeshLib::Mesh&>(mesh), "fracture_stress_avg",
442 MeshLib::MeshItemType::Cell, DisplacementDim);
443
444 _process_data.element_fracture_velocities =
446 const_cast<MeshLib::Mesh&>(mesh), "fracture_velocity_avg",
447 MeshLib::MeshItemType::Cell, DisplacementDim);
448
450 const_cast<MeshLib::Mesh&>(mesh), "fracture_aperture_avg",
452
453 mesh_prop_b->resize(mesh.getNumberOfElements());
454 auto const& mesh_prop_matid = *_process_data.mesh_prop_materialIDs;
455 for (auto const& fracture_prop : _process_data.fracture_properties)
456 {
457 for (MeshLib::Element const* e : _mesh.getElements())
458 {
459 if (e->getDimension() == DisplacementDim)
460 {
461 continue;
462 }
463 if (mesh_prop_matid[e->getID()] != fracture_prop.mat_id)
464 {
465 continue;
466 }
467 // Mean value for the element. This allows usage of node based
468 // properties for aperture.
469 (*mesh_prop_b)[e->getID()] =
470 fracture_prop.aperture0
471 .getNodalValuesOnElement(*e, /*time independent*/ 0)
472 .mean();
473 }
474 }
475
476 auto mesh_prop_k_f = MeshLib::getOrCreateMeshProperty<double>(
477 const_cast<MeshLib::Mesh&>(mesh), "fracture_permeability_avg",
479 mesh_prop_k_f->resize(mesh.getNumberOfElements());
480 _process_data.mesh_prop_k_f = mesh_prop_k_f;
481
482 auto mesh_prop_fracture_shear_failure =
484 const_cast<MeshLib::Mesh&>(mesh), "f_shear_failure",
486 mesh_prop_fracture_shear_failure->resize(mesh.getNumberOfElements());
487 _process_data.mesh_prop_fracture_shear_failure =
488 mesh_prop_fracture_shear_failure;
489
490 auto mesh_prop_nodal_p = MeshLib::getOrCreateMeshProperty<double>(
491 const_cast<MeshLib::Mesh&>(mesh), "pressure_interpolated",
493 mesh_prop_nodal_p->resize(mesh.getNumberOfNodes());
494 _process_data.mesh_prop_nodal_p = mesh_prop_nodal_p;
495
496 _process_data.mesh_prop_nodal_forces =
498 const_cast<MeshLib::Mesh&>(mesh), "NodalForces",
499 MeshLib::MeshItemType::Node, DisplacementDim);
500 assert(_process_data.mesh_prop_nodal_forces->size() ==
501 DisplacementDim * mesh.getNumberOfNodes());
502
503 _process_data.mesh_prop_nodal_forces_jump =
505 const_cast<MeshLib::Mesh&>(mesh), "NodalForcesJump",
506 MeshLib::MeshItemType::Node, DisplacementDim);
507 assert(_process_data.mesh_prop_nodal_forces_jump->size() ==
508 DisplacementDim * mesh.getNumberOfNodes());
509
510 _process_data.mesh_prop_hydraulic_flow =
512 const_cast<MeshLib::Mesh&>(mesh), "MassFlowRate",
514 assert(_process_data.mesh_prop_hydraulic_flow->size() ==
515 mesh.getNumberOfNodes());
516 _process_data.mesh_prop_b = mesh_prop_b;
517}
518
519template <int DisplacementDim>
521 std::vector<GlobalVector*> const& x,
522 std::vector<GlobalVector*> const& x_prev, const double t, double const dt,
523 int const process_id)
524{
525 if (process_id == 0)
526 {
527 DBUG("PostTimestep HydroMechanicsProcess.");
528
531 _local_assemblers, getActiveElementIDs(), getDOFTables(x.size()), x,
532 x_prev, t, dt, process_id);
533 }
534
535 DBUG("Compute the secondary variables for HydroMechanicsProcess.");
536
537 const auto& dof_table = getDOFTable(process_id);
538
539 // Copy displacement jumps in a solution vector to mesh property
540 // Remark: the copy is required because mesh properties for primary
541 // variables are set during output and are not ready yet when this function
542 // is called.
543 int g_variable_id = 0;
544 {
545 const int monolithic_process_id = 0;
546 auto const& pvs = getProcessVariables(monolithic_process_id);
547 auto const it =
548 std::find_if(pvs.begin(), pvs.end(),
549 [](ProcessVariable const& pv)
550 { return pv.getName() == "displacement_jump1"; });
551 if (it == pvs.end())
552 {
553 OGS_FATAL(
554 "Didn't find expected 'displacement_jump1' process variable.");
555 }
556 g_variable_id = static_cast<int>(std::distance(pvs.begin(), it));
557 }
558
560
561 const int monolithic_process_id = 0;
562 ProcessVariable& pv_g =
563 this->getProcessVariables(monolithic_process_id)[g_variable_id];
564 auto const num_comp = pv_g.getNumberOfGlobalComponents();
565 auto& mesh_prop_g = *MeshLib::getOrCreateMeshProperty<double>(
566 _mesh, pv_g.getName(), MeshLib::MeshItemType::Node, num_comp);
567 for (int component_id = 0; component_id < num_comp; ++component_id)
568 {
569 auto const& mesh_subset =
570 dof_table.getMeshSubset(g_variable_id, component_id);
571 for (auto const& l : MeshLib::views::meshLocations(
572 mesh_subset, MeshLib::MeshItemType::Node))
573 {
574 auto const global_index =
575 dof_table.getGlobalIndex(l, g_variable_id, component_id);
576 auto const node_id = l.item_id;
577 mesh_prop_g[node_id * num_comp + component_id] =
578 (*x[process_id])[global_index];
579 }
580 }
581}
582
583template <int DisplacementDim>
585{
586 return false;
587}
588
589template <int DisplacementDim>
591 const double t, double const dt, std::vector<GlobalVector*> const& x,
592 std::vector<GlobalVector*> const& x_prev, int const process_id,
594{
595 DBUG("Assemble HydroMechanicsProcess.");
596
597 std::vector<NumLib::LocalToGlobalIndexMap const*> dof_table = {
598 _local_to_global_index_map.get()};
599 // Call global assembler for each local assembly item.
601 _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers,
602 dof_table, t, dt, x, x_prev, process_id, &M, &K, &b);
603}
604
605template <int DisplacementDim>
608 const double t, double const dt, std::vector<GlobalVector*> const& x,
609 std::vector<GlobalVector*> const& x_prev, int const process_id,
610 GlobalVector& b, GlobalMatrix& Jac)
611{
612 DBUG("AssembleWithJacobian HydroMechanicsProcess.");
613
614 // Call global assembler for each local assembly item.
615 std::vector<NumLib::LocalToGlobalIndexMap const*> dof_table = {
616 _local_to_global_index_map.get()};
619 _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
620 process_id, &b, &Jac);
621
622 auto copyRhs = [&](int const variable_id, auto& output_vector)
623 {
624 transformVariableFromGlobalVector(b, variable_id,
625 *_local_to_global_index_map,
626 output_vector, std::negate<double>());
627 };
628 copyRhs(0, *_process_data.mesh_prop_hydraulic_flow);
629 copyRhs(1, *_process_data.mesh_prop_nodal_forces);
630 copyRhs(2, *_process_data.mesh_prop_nodal_forces_jump);
631}
632
633template <int DisplacementDim>
635 std::vector<GlobalVector*> const& x, double const t, double const dt,
636 int const process_id)
637{
638 DBUG("PreTimestep HydroMechanicsProcess.");
639
642 getActiveElementIDs(), *_local_to_global_index_map, *x[process_id], t,
643 dt);
644}
645
646// ------------------------------------------------------------------------------------
647// template instantiation
648// ------------------------------------------------------------------------------------
649template class HydroMechanicsProcess<2>;
650template class HydroMechanicsProcess<3>;
651
652} // namespace HydroMechanics
653} // namespace LIE
654} // namespace ProcessLib
Definition of the ElementStatus class.
#define OGS_FATAL(...)
Definition Error.h:26
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:35
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:30
Definition of the MeshInformation class.
Definition of the class Properties that implements a container of properties.
Definition of the Mesh class.
Global vector based on Eigen vector.
Definition EigenVector.h:26
std::size_t getID() const
Returns the ID of the element.
Definition Element.h:91
bool isAxiallySymmetric() const
Definition Mesh.h:139
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition Mesh.h:111
unsigned getDimension() const
Returns the dimension of the mesh (determined by the maximum dimension over all elements).
Definition Mesh.h:90
const Node * getNode(std::size_t idx) const
Get the node with the given index.
Definition Mesh.h:93
Properties & getProperties()
Definition Mesh.h:136
const std::string getName() const
Get name of the mesh.
Definition Mesh.h:105
std::size_t getNumberOfNodes() const
Get the number of nodes.
Definition Mesh.h:102
std::vector< Element const * > const & getElementsConnectedToNode(std::size_t node_id) const
Definition Mesh.cpp:257
std::size_t getNumberOfElements() const
Get the number of elements.
Definition Mesh.h:99
PropertyVector< T > const * getPropertyVector(std::string_view name) const
static std::optional< std::pair< T, T > > const getValueBounds(MeshLib::PropertyVector< T > const &property)
HydroMechanicsProcess(std::string name, MeshLib::Mesh &mesh, std::unique_ptr< ProcessLib::AbstractJacobianAssembler > &&jacobian_assembler, std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, 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 postTimestepConcreteProcess(std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, int const process_id) override
void initializeConcreteProcess(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order) override
Process specific initialization called by initialize().
std::vector< std::vector< MeshLib::Node * > > _vec_fracture_nodes
std::vector< std::vector< MeshLib::Element * > > _vec_junction_fracture_matrix_elements
void preTimestepConcreteProcess(std::vector< GlobalVector * > const &x, double const t, double const dt, int const process_id) 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 updateElementLevelSets(MeshLib::Element const &e, MeshLib::Mesh &mesh)
std::vector< std::vector< MeshLib::Element * > > _vec_fracture_elements
HydroMechanicsProcessData< DisplacementDim > _process_data
void assembleWithJacobianConcreteProcess(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalVector &b, GlobalMatrix &Jac) override
std::vector< std::vector< MeshLib::Element * > > _vec_fracture_matrix_elements
virtual void postTimestep(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, int const process_id)
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)
ParameterLib::Parameter< double > const & getInitialCondition() const
std::string const & getName() const
int getNumberOfGlobalComponents() const
Returns the number of components of the process variable.
MeshLib::Mesh & _mesh
Definition Process.h:365
std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > const & getProcessVariables() const
Definition Process.h:156
Handles configuration of several secondary variables from the project file.
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)
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, GlobalVector *b, GlobalMatrix *Jac)
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
void setLocalAccessibleVector(PETScVector const &x)
Definition LinAlg.cpp:27
constexpr ranges::views::view_closure ids
For an element of a range view return its id.
Definition Mesh.h:227
auto meshLocations(Mesh const &mesh, MeshItemType const item_type)
Definition Mesh.h:238
PropertyVector< T > * getOrCreateMeshProperty(Mesh &mesh, std::string const &property_name, MeshItemType const item_type, int const number_of_components)
std::vector< Node * > getBaseNodes(std::vector< Element * > const &elements)
Definition Utils.h:26
@ BY_COMPONENT
Ordering data by component type.
void createLocalAssemblers(std::vector< MeshLib::Element * > const &mesh_elements, NumLib::LocalToGlobalIndexMap const &dof_table, std::vector< std::unique_ptr< LocalAssemblerInterface > > &local_assemblers, NumLib::IntegrationOrder const integration_order, ExtraCtorArgs &&... extra_ctor_args)
void setFractureProperty(int const dim, MeshLib::Element const &e, FractureProperty &frac_prop)
void getFractureMatrixDataInMesh(MeshLib::Mesh const &mesh, std::vector< MeshLib::Element * > &vec_matrix_elements, std::vector< int > &vec_fracture_mat_IDs, std::vector< std::vector< MeshLib::Element * > > &vec_fracture_elements, std::vector< std::vector< MeshLib::Element * > > &vec_fracture_matrix_elements, std::vector< std::vector< MeshLib::Node * > > &vec_fracture_nodes, std::vector< std::pair< std::size_t, std::vector< int > > > &vec_branch_nodeID_matIDs, std::vector< std::pair< std::size_t, std::vector< int > > > &vec_junction_nodeID_matIDs)
std::vector< double > uGlobalEnrichments(std::vector< FractureProperty * > const &frac_props, std::vector< JunctionProperty * > const &junction_props, std::unordered_map< int, int > const &fracID_to_local, Eigen::Vector3d const &x)
BranchProperty createBranchProperty(MeshLib::Node const &branchNode, FractureProperty const &master_frac, FractureProperty const &slave_frac)
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 executeMemberDereferenced(Object &object, Method method, Container const &container, Args &&... args)