OGS 6.2.0-97-g4a610c866
HydroMechanicsProcess.cpp
Go to the documentation of this file.
1 
10 #include "HydroMechanicsProcess.h"
11 
13 #include "MeshLib/ElementStatus.h"
14 #include "MeshLib/Elements/Utils.h"
15 #include "MeshLib/Mesh.h"
17 #include "MeshLib/Properties.h"
18 
21 
26 
27 #include "LocalAssembler/CreateLocalAssemblers.h"
31 
32 namespace ProcessLib
33 {
34 namespace LIE
35 {
36 namespace HydroMechanics
37 {
38 template <int GlobalDim>
40  MeshLib::Mesh& mesh,
41  std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
42  std::vector<std::unique_ptr<ParameterLib::ParameterBase>> const& parameters,
43  unsigned const integration_order,
44  std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>&&
45  process_variables,
47  SecondaryVariableCollection&& secondary_variables,
48  NumLib::NamedFunctionCaller&& named_function_caller,
49  bool const use_monolithic_scheme)
50  : Process(mesh, std::move(jacobian_assembler), parameters,
51  integration_order, std::move(process_variables),
52  std::move(secondary_variables), std::move(named_function_caller),
53  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<int> vec_fracture_mat_IDs;
58  std::vector<std::vector<MeshLib::Element*>> vec_vec_fracture_elements;
59  std::vector<std::vector<MeshLib::Element*>>
60  vec_vec_fracture_matrix_elements;
61  std::vector<std::vector<MeshLib::Node*>> vec_vec_fracture_nodes;
62  std::vector<std::pair<std::size_t, std::vector<int>>>
63  vec_branch_nodeID_matIDs;
64  std::vector<std::pair<std::size_t, std::vector<int>>>
65  vec_junction_nodeID_matIDs;
67  mesh, _vec_matrix_elements, vec_fracture_mat_IDs,
68  vec_vec_fracture_elements, vec_vec_fracture_matrix_elements,
69  vec_vec_fracture_nodes, vec_branch_nodeID_matIDs,
70  vec_junction_nodeID_matIDs);
72  vec_vec_fracture_elements[0].begin(),
73  vec_vec_fracture_elements[0].end());
76  vec_vec_fracture_matrix_elements[0].begin(),
77  vec_vec_fracture_matrix_elements[0].end());
79  vec_vec_fracture_nodes[0].begin(),
80  vec_vec_fracture_nodes[0].end());
81 
82  if (!_vec_fracture_elements.empty())
83  {
84  // set fracture property assuming a fracture forms a straight line
85  setFractureProperty(GlobalDim,
87  *_process_data.fracture_property);
88  }
89 
90  //
91  // If Neumann BCs for the displacement_jump variable are required they need
92  // special treatment because of the levelset function. The implementation
93  // exists in the version 6.1.0 (e54815cc07ee89c81f953a4955b1c788595dd725)
94  // and was removed due to lack of applications.
95  //
96 
97  if (!_process_data.deactivate_matrix_in_flow)
98  {
99  _process_data.p_element_status =
100  std::make_unique<MeshLib::ElementStatus>(&mesh);
101  }
102  else
103  {
104  auto range =
105  MeshLib::MeshInformation::getValueBounds<int>(mesh, "MaterialIDs");
106  std::vector<int> vec_p_inactive_matIDs;
107  for (int matID = range.first; matID <= range.second; matID++)
108  {
109  if (std::find(vec_fracture_mat_IDs.begin(),
110  vec_fracture_mat_IDs.end(),
111  matID) == vec_fracture_mat_IDs.end())
112  {
113  vec_p_inactive_matIDs.push_back(matID);
114  }
115  }
116  _process_data.p_element_status =
117  std::make_unique<MeshLib::ElementStatus>(&mesh,
118  vec_p_inactive_matIDs);
119 
120  const int monolithic_process_id = 0;
121  ProcessVariable const& pv_p =
122  getProcessVariables(monolithic_process_id)[0];
123  _process_data.p0 = &pv_p.getInitialCondition();
124  }
125 }
126 
127 template <int GlobalDim>
129 {
130  //------------------------------------------------------------
131  // prepare mesh subsets to define DoFs
132  //------------------------------------------------------------
133  // for extrapolation
135  std::make_unique<MeshLib::MeshSubset>(_mesh, _mesh.getNodes());
136  // pressure
138  _process_data.p_element_status->getActiveElements());
140  std::make_unique<MeshLib::MeshSubset>(_mesh, _mesh_nodes_p);
141  // regular u
143  std::make_unique<MeshLib::MeshSubset>(_mesh, _mesh.getNodes());
144  if (!_vec_fracture_nodes.empty())
145  {
146  // u jump
148  std::make_unique<MeshLib::MeshSubset>(_mesh, _vec_fracture_nodes);
149  }
150 
151  // Collect the mesh subsets in a vector.
152  std::vector<MeshLib::MeshSubset> all_mesh_subsets;
153  std::vector<int> vec_n_components;
154  std::vector<std::vector<MeshLib::Element*> const*> vec_var_elements;
155  // pressure
156  vec_n_components.push_back(1);
157  all_mesh_subsets.emplace_back(*_mesh_subset_nodes_p);
158  if (!_process_data.deactivate_matrix_in_flow)
159  {
160  vec_var_elements.push_back(&_mesh.getElements());
161  }
162  else
163  {
164  // TODO set elements including active nodes for pressure.
165  // cannot use ElementStatus
166  vec_var_elements.push_back(&_vec_fracture_matrix_elements);
167  }
168  // regular displacement
169  vec_n_components.push_back(GlobalDim);
170  std::generate_n(std::back_inserter(all_mesh_subsets), GlobalDim,
171  [&]() { return *_mesh_subset_matrix_nodes; });
172  vec_var_elements.push_back(&_vec_matrix_elements);
173  if (!_vec_fracture_nodes.empty())
174  {
175  // displacement jump
176  vec_n_components.push_back(GlobalDim);
177  std::generate_n(std::back_inserter(all_mesh_subsets), GlobalDim,
178  [&]() { return *_mesh_subset_fracture_nodes; });
179  vec_var_elements.push_back(&_vec_fracture_matrix_elements);
180  }
181 
182  INFO("[LIE/HM] creating a DoF table");
184  std::make_unique<NumLib::LocalToGlobalIndexMap>(
185  std::move(all_mesh_subsets),
186  vec_n_components,
187  vec_var_elements,
189 
190  DBUG("created %d DoF", _local_to_global_index_map->size());
191 }
192 
193 template <int GlobalDim>
195  NumLib::LocalToGlobalIndexMap const& dof_table,
196  MeshLib::Mesh const& mesh,
197  unsigned const integration_order)
198 {
199  assert(mesh.getDimension() == GlobalDim);
200  INFO("[LIE/HM] creating local assemblers");
201  const int monolithic_process_id = 0;
206  mesh.getElements(), dof_table,
207  // use displacment process variable for shapefunction order
209  monolithic_process_id)[1].get().getShapeFunctionOrder(),
210  _local_assemblers, mesh.isAxiallySymmetric(), integration_order,
211  _process_data);
212 
213  auto mesh_prop_sigma_xx = MeshLib::getOrCreateMeshProperty<double>(
214  const_cast<MeshLib::Mesh&>(mesh), "stress_xx",
216  mesh_prop_sigma_xx->resize(mesh.getNumberOfElements());
217  _process_data.mesh_prop_stress_xx = mesh_prop_sigma_xx;
218 
219  auto mesh_prop_sigma_yy = MeshLib::getOrCreateMeshProperty<double>(
220  const_cast<MeshLib::Mesh&>(mesh), "stress_yy",
222  mesh_prop_sigma_yy->resize(mesh.getNumberOfElements());
223  _process_data.mesh_prop_stress_yy = mesh_prop_sigma_yy;
224 
225  auto mesh_prop_sigma_zz = MeshLib::getOrCreateMeshProperty<double>(
226  const_cast<MeshLib::Mesh&>(mesh), "stress_zz",
228  mesh_prop_sigma_zz->resize(mesh.getNumberOfElements());
229  _process_data.mesh_prop_stress_zz = mesh_prop_sigma_zz;
230 
231  auto mesh_prop_sigma_xy = MeshLib::getOrCreateMeshProperty<double>(
232  const_cast<MeshLib::Mesh&>(mesh), "stress_xy",
234  mesh_prop_sigma_xy->resize(mesh.getNumberOfElements());
235  _process_data.mesh_prop_stress_xy = mesh_prop_sigma_xy;
236 
237  if (GlobalDim == 3)
238  {
239  auto mesh_prop_sigma_xz = MeshLib::getOrCreateMeshProperty<double>(
240  const_cast<MeshLib::Mesh&>(mesh), "stress_xz",
242  mesh_prop_sigma_xz->resize(mesh.getNumberOfElements());
243  _process_data.mesh_prop_stress_xz = mesh_prop_sigma_xz;
244 
245  auto mesh_prop_sigma_yz = MeshLib::getOrCreateMeshProperty<double>(
246  const_cast<MeshLib::Mesh&>(mesh), "stress_yz",
248  mesh_prop_sigma_yz->resize(mesh.getNumberOfElements());
249  _process_data.mesh_prop_stress_yz = mesh_prop_sigma_yz;
250  }
251 
252  auto mesh_prop_epsilon_xx = MeshLib::getOrCreateMeshProperty<double>(
253  const_cast<MeshLib::Mesh&>(mesh), "strain_xx",
255  mesh_prop_epsilon_xx->resize(mesh.getNumberOfElements());
256  _process_data.mesh_prop_strain_xx = mesh_prop_epsilon_xx;
257 
258  auto mesh_prop_epsilon_yy = MeshLib::getOrCreateMeshProperty<double>(
259  const_cast<MeshLib::Mesh&>(mesh), "strain_yy",
261  mesh_prop_epsilon_yy->resize(mesh.getNumberOfElements());
262  _process_data.mesh_prop_strain_yy = mesh_prop_epsilon_yy;
263 
264  auto mesh_prop_epsilon_zz = MeshLib::getOrCreateMeshProperty<double>(
265  const_cast<MeshLib::Mesh&>(mesh), "strain_zz",
267  mesh_prop_epsilon_zz->resize(mesh.getNumberOfElements());
268  _process_data.mesh_prop_strain_zz = mesh_prop_epsilon_zz;
269 
270  auto mesh_prop_epsilon_xy = MeshLib::getOrCreateMeshProperty<double>(
271  const_cast<MeshLib::Mesh&>(mesh), "strain_xy",
273  mesh_prop_epsilon_xy->resize(mesh.getNumberOfElements());
274  _process_data.mesh_prop_strain_xy = mesh_prop_epsilon_xy;
275 
276  if (GlobalDim == 3)
277  {
278  auto mesh_prop_epsilon_xz = MeshLib::getOrCreateMeshProperty<double>(
279  const_cast<MeshLib::Mesh&>(mesh), "strain_xz",
281  mesh_prop_epsilon_xz->resize(mesh.getNumberOfElements());
282  _process_data.mesh_prop_strain_xz = mesh_prop_epsilon_xz;
283 
284  auto mesh_prop_epsilon_yz = MeshLib::getOrCreateMeshProperty<double>(
285  const_cast<MeshLib::Mesh&>(mesh), "strain_yz",
287  mesh_prop_epsilon_yz->resize(mesh.getNumberOfElements());
288  _process_data.mesh_prop_strain_yz = mesh_prop_epsilon_yz;
289  }
290 
291  auto mesh_prop_velocity = MeshLib::getOrCreateMeshProperty<double>(
292  const_cast<MeshLib::Mesh&>(mesh), "velocity",
294  mesh_prop_velocity->resize(mesh.getNumberOfElements() * 3);
295  _process_data.mesh_prop_velocity = mesh_prop_velocity;
296 
297  if (!_vec_fracture_elements.empty())
298  {
299  auto mesh_prop_levelset = MeshLib::getOrCreateMeshProperty<double>(
300  const_cast<MeshLib::Mesh&>(mesh), "levelset1",
302  mesh_prop_levelset->resize(mesh.getNumberOfElements());
303  for (MeshLib::Element const* e : _mesh.getElements())
304  {
305  if (e->getDimension() < GlobalDim)
306  {
307  continue;
308  }
309 
310  std::vector<FractureProperty*> fracture_props(
311  {_process_data.fracture_property.get()});
312  std::vector<JunctionProperty*> junction_props;
313  std::unordered_map<int, int> fracID_to_local({{0, 0}});
314  std::vector<double> levelsets = uGlobalEnrichments(
315  fracture_props, junction_props, fracID_to_local,
316  Eigen::Vector3d(e->getCenterOfGravity().getCoords()));
317  (*mesh_prop_levelset)[e->getID()] = levelsets[0];
318  }
319 
320  auto mesh_prop_w_n = MeshLib::getOrCreateMeshProperty<double>(
321  const_cast<MeshLib::Mesh&>(mesh), "w_n",
323  mesh_prop_w_n->resize(mesh.getNumberOfElements());
324  auto mesh_prop_w_s = MeshLib::getOrCreateMeshProperty<double>(
325  const_cast<MeshLib::Mesh&>(mesh), "w_s",
327  mesh_prop_w_s->resize(mesh.getNumberOfElements());
328  _process_data.mesh_prop_w_n = mesh_prop_w_n;
329  _process_data.mesh_prop_w_s = mesh_prop_w_s;
330 
331  auto mesh_prop_b = MeshLib::getOrCreateMeshProperty<double>(
332  const_cast<MeshLib::Mesh&>(mesh), "aperture",
334  mesh_prop_b->resize(mesh.getNumberOfElements());
335  auto const mesh_prop_matid = materialIDs(mesh);
336  if (!mesh_prop_matid)
337  {
338  OGS_FATAL("Could not access MaterialIDs property from mesh.");
339  }
340  auto const& frac = _process_data.fracture_property;
341  for (MeshLib::Element const* e : _mesh.getElements())
342  {
343  if (e->getDimension() == GlobalDim)
344  {
345  continue;
346  }
347  if ((*mesh_prop_matid)[e->getID()] != frac->mat_id)
348  {
349  continue;
350  }
352  x.setElementID(e->getID());
353  (*mesh_prop_b)[e->getID()] = frac->aperture0(0, x)[0];
354  }
355  _process_data.mesh_prop_b = mesh_prop_b;
356 
357  auto mesh_prop_k_f = MeshLib::getOrCreateMeshProperty<double>(
358  const_cast<MeshLib::Mesh&>(mesh), "k_f",
360  mesh_prop_k_f->resize(mesh.getNumberOfElements());
361  _process_data.mesh_prop_k_f = mesh_prop_k_f;
362 
363  auto mesh_prop_fracture_stress_shear =
364  MeshLib::getOrCreateMeshProperty<double>(
365  const_cast<MeshLib::Mesh&>(mesh), "f_stress_s",
367  mesh_prop_fracture_stress_shear->resize(mesh.getNumberOfElements());
368  _process_data.mesh_prop_fracture_stress_shear =
369  mesh_prop_fracture_stress_shear;
370 
371  auto mesh_prop_fracture_stress_normal =
372  MeshLib::getOrCreateMeshProperty<double>(
373  const_cast<MeshLib::Mesh&>(mesh), "f_stress_n",
375  mesh_prop_fracture_stress_normal->resize(mesh.getNumberOfElements());
376  _process_data.mesh_prop_fracture_stress_normal =
377  mesh_prop_fracture_stress_normal;
378 
379  auto mesh_prop_fracture_shear_failure =
380  MeshLib::getOrCreateMeshProperty<double>(
381  const_cast<MeshLib::Mesh&>(mesh), "f_shear_failure",
383  mesh_prop_fracture_shear_failure->resize(mesh.getNumberOfElements());
384  _process_data.mesh_prop_fracture_shear_failure =
385  mesh_prop_fracture_shear_failure;
386 
387  auto mesh_prop_nodal_w = MeshLib::getOrCreateMeshProperty<double>(
388  const_cast<MeshLib::Mesh&>(mesh), "nodal_w",
389  MeshLib::MeshItemType::Node, GlobalDim);
390  mesh_prop_nodal_w->resize(mesh.getNumberOfNodes() * GlobalDim);
391  _process_data.mesh_prop_nodal_w = mesh_prop_nodal_w;
392 
393  auto mesh_prop_nodal_b = MeshLib::getOrCreateMeshProperty<double>(
394  const_cast<MeshLib::Mesh&>(mesh), "nodal_aperture",
396  mesh_prop_nodal_b->resize(mesh.getNumberOfNodes());
397  _process_data.mesh_prop_nodal_b = mesh_prop_nodal_b;
398 
399  if (GlobalDim == 3)
400  {
401  auto mesh_prop_w_s2 = MeshLib::getOrCreateMeshProperty<double>(
402  const_cast<MeshLib::Mesh&>(mesh), "w_s2",
404  mesh_prop_w_s2->resize(mesh.getNumberOfElements());
405  _process_data.mesh_prop_w_s2 = mesh_prop_w_s2;
406 
407  auto mesh_prop_fracture_stress_shear2 =
408  MeshLib::getOrCreateMeshProperty<double>(
409  const_cast<MeshLib::Mesh&>(mesh), "f_stress_s2",
411  mesh_prop_fracture_stress_shear2->resize(
412  mesh.getNumberOfElements());
413  _process_data.mesh_prop_fracture_stress_shear2 =
414  mesh_prop_fracture_stress_shear2;
415  }
416 
417  auto mesh_prop_nodal_p = MeshLib::getOrCreateMeshProperty<double>(
418  const_cast<MeshLib::Mesh&>(mesh), "pressure_interpolated",
420  mesh_prop_nodal_p->resize(mesh.getNumberOfNodes());
421  _process_data.mesh_prop_nodal_p = mesh_prop_nodal_p;
422 
423  _process_data.mesh_prop_nodal_forces =
424  MeshLib::getOrCreateMeshProperty<double>(
425  const_cast<MeshLib::Mesh&>(mesh), "NodalForces",
426  MeshLib::MeshItemType::Node, GlobalDim);
427  assert(_process_data.mesh_prop_nodal_forces->size() ==
428  GlobalDim * mesh.getNumberOfNodes());
429 
430  _process_data.mesh_prop_nodal_forces_jump =
431  MeshLib::getOrCreateMeshProperty<double>(
432  const_cast<MeshLib::Mesh&>(mesh), "NodalForcesJump",
433  MeshLib::MeshItemType::Node, GlobalDim);
434  assert(_process_data.mesh_prop_nodal_forces_jump->size() ==
435  GlobalDim * mesh.getNumberOfNodes());
436 
437  _process_data.mesh_prop_hydraulic_flow =
438  MeshLib::getOrCreateMeshProperty<double>(
439  const_cast<MeshLib::Mesh&>(mesh), "HydraulicFlow",
441  assert(_process_data.mesh_prop_hydraulic_flow->size() ==
442  mesh.getNumberOfNodes());
443  }
444 }
445 
446 template <int GlobalDim>
448  const double t, GlobalVector const& x, int const process_id)
449 {
450  DBUG("Compute the secondary variables for HydroMechanicsProcess.");
451  const auto& dof_table = getDOFTable(process_id);
452 
453  {
454  ProcessLib::ProcessVariable const& pv =
455  getProcessVariables(process_id)[0];
456 
460  dof_table, t, x, _coupled_solutions);
461  }
462 
463  // Copy displacement jumps in a solution vector to mesh property
464  // Remark: the copy is required because mesh properties for primary
465  // variables are set during output and are not ready yet when this function
466  // is called.
467  int g_variable_id = 0;
468  {
469  const int monolithic_process_id = 0;
470  auto const& pvs = getProcessVariables(monolithic_process_id);
471  auto const it =
472  std::find_if(pvs.begin(), pvs.end(), [](ProcessVariable const& pv) {
473  return pv.getName() == "displacement_jump1";
474  });
475  if (it == pvs.end())
476  {
477  OGS_FATAL(
478  "Didn't find expected 'displacement_jump1' process "
479  "variable.");
480  }
481  g_variable_id = static_cast<int>(std::distance(pvs.begin(), it));
482  }
483 
484  MathLib::LinAlg::setLocalAccessibleVector(x);
485 
486  const int monolithic_process_id = 0;
487  ProcessVariable& pv_g =
488  this->getProcessVariables(monolithic_process_id)[g_variable_id];
489  auto const num_comp = pv_g.getNumberOfComponents();
490  auto& mesh_prop_g = *MeshLib::getOrCreateMeshProperty<double>(
491  _mesh, pv_g.getName(), MeshLib::MeshItemType::Node, num_comp);
492  for (int component_id = 0; component_id < num_comp; ++component_id)
493  {
494  auto const& mesh_subset = dof_table.getMeshSubset(
495  g_variable_id, component_id);
496  auto const mesh_id = mesh_subset.getMeshID();
497  for (auto const* node : mesh_subset.getNodes())
498  {
500  node->getID());
501 
502  auto const global_index =
503  dof_table.getGlobalIndex(l, g_variable_id, component_id);
504  mesh_prop_g[node->getID() * num_comp + component_id] =
505  x[global_index];
506  }
507  }
508 
509  // compute nodal w and aperture
510  auto const& R = _process_data.fracture_property->R;
511  auto const& b0 = _process_data.fracture_property->aperture0;
512  MeshLib::PropertyVector<double>& vec_w = *_process_data.mesh_prop_nodal_w;
513  MeshLib::PropertyVector<double>& vec_b = *_process_data.mesh_prop_nodal_b;
514 
515  auto compute_nodal_aperture = [&](std::size_t const node_id,
516  double const w_n) {
517  // skip aperture computation for element-wise defined b0 because there
518  // are jumps on the nodes between the element's values.
519  if (dynamic_cast<ParameterLib::MeshElementParameter<double> const*>(
520  &b0))
521  {
522  return std::numeric_limits<double>::quiet_NaN();
523  }
524 
526  x.setNodeID(node_id);
527  return w_n + b0(/*time independent*/ 0, x)[0];
528  };
529 
530  Eigen::VectorXd g(GlobalDim), w(GlobalDim);
531  for (MeshLib::Node const* node : _vec_fracture_nodes)
532  {
533  auto const node_id = node->getID();
534  g.setZero();
535  for (int k = 0; k < GlobalDim; k++)
536  {
537  g[k] = mesh_prop_g[node_id * GlobalDim + k];
538  }
539 
540  w.noalias() = R * g;
541  for (int k = 0; k < GlobalDim; k++)
542  {
543  vec_w[node_id * GlobalDim + k] = w[k];
544  }
545 
546  vec_b[node_id] = compute_nodal_aperture(node_id, w[GlobalDim - 1]);
547  }
548 }
549 
550 template <int GlobalDim>
552 {
553  return false;
554 }
555 
556 template <int GlobalDim>
558  const double t, GlobalVector const& x, GlobalMatrix& M, GlobalMatrix& K,
559  GlobalVector& b)
560 {
561  DBUG("Assemble HydroMechanicsProcess.");
562 
563  std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
564  dof_table = {std::ref(*_local_to_global_index_map)};
565  // Call global assembler for each local assembly item.
568  dof_table, t, x, M, K, b, _coupled_solutions);
569 }
570 
571 template <int GlobalDim>
573  const double t, GlobalVector const& x, GlobalVector const& xdot,
574  const double dxdot_dx, const double dx_dx, GlobalMatrix& M, GlobalMatrix& K,
575  GlobalVector& b, GlobalMatrix& Jac)
576 {
577  DBUG("AssembleWithJacobian HydroMechanicsProcess.");
578 
579  const int process_id =
581  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
582 
583  // Call global assembler for each local assembly item.
584  std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
585  dof_table = {std::ref(*_local_to_global_index_map)};
588  _local_assemblers, pv.getActiveElementIDs(), dof_table, t, x,
589  xdot, dxdot_dx, dx_dx, M, K, b, Jac, _coupled_solutions);
590 
591  auto copyRhs = [&](int const variable_id, auto& output_vector) {
592  transformVariableFromGlobalVector(b, variable_id,
594  output_vector, std::negate<double>());
595  };
596  copyRhs(0, *_process_data.mesh_prop_hydraulic_flow);
597  copyRhs(1, *_process_data.mesh_prop_nodal_forces);
598  copyRhs(2, *_process_data.mesh_prop_nodal_forces_jump);
599 }
600 
601 template <int GlobalDim>
603  GlobalVector const& x, double const t, double const dt,
604  const int process_id)
605 {
606  DBUG("PreTimestep HydroMechanicsProcess.");
607 
608  _process_data.dt = dt;
609  _process_data.t = t;
610 
611  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
612 
616  x, t, dt);
617 }
618 
619 // ------------------------------------------------------------------------------------
620 // template instantiation
621 // ------------------------------------------------------------------------------------
622 template class HydroMechanicsProcess<2>;
623 template class HydroMechanicsProcess<3>;
624 
625 } // namespace HydroMechanics
626 } // namespace LIE
627 } // namespace ProcessLib
ParameterLib::Parameter< double > const & getInitialCondition() const
std::size_t getNumberOfNodes() const
Get the number of nodes.
Definition: Mesh.h:99
void setFractureProperty(int const dim, MeshLib::Element const &e, FractureProperty &frac_prop)
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_nodes_p
MeshLib::Mesh & _mesh
Definition: Process.h:262
virtual NumLib::LocalToGlobalIndexMap const & getDOFTable(const int) const
Definition: Process.h:116
std::vector< Node * > getBaseNodes(std::vector< Element *> const &elements)
Definition: Utils.h:25
HydroMechanicsProcess(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< GlobalDim > &&process_data, SecondaryVariableCollection &&secondary_variables, NumLib::NamedFunctionCaller &&named_function_caller, bool const use_monolithic_scheme)
void setElementID(std::size_t element_id)
virtual void computeSecondaryVariable(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, const double t, GlobalVector const &x, CoupledSolutionsForStaggeredScheme const *coupled_xs)
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_all_nodes
Definition: Process.h:263
HydroMechanicsProcessData< GlobalDim > _process_data
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)
Definition: MeshUtils.cpp:213
static void executeSelectedMemberOnDereferenced(Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables(const int process_id) const
Definition: Process.h:124
std::vector< std::unique_ptr< LocalAssemblerInterface > > _local_assemblers
void setNodeID(std::size_t node_id)
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_matrix_nodes
A parameter represented by a mesh property vector.
void assembleConcreteProcess(const double t, GlobalVector const &x, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b) override
void assembleWithJacobianConcreteProcess(const double t, GlobalVector const &x, GlobalVector const &xdot, const double dxdot_dx, const double dx_dx, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) override
void createLocalAssemblers(std::vector< MeshLib::Element *> const &mesh_elements, NumLib::LocalToGlobalIndexMap const &dof_table, const unsigned shapefunction_order, std::vector< std::unique_ptr< LocalAssemblerInterface >> &local_assemblers, ExtraCtorArgs &&... extra_ctor_args)
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition: Mesh.h:105
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, GlobalVector const &x, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, CoupledSolutionsForStaggeredScheme const *const cpl_xs)
Definition of the Mesh class.
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, GlobalVector const &x, GlobalVector const &xdot, const double dxdot_dx, const double dx_dx, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac, CoupledSolutionsForStaggeredScheme const *const cpl_xs)
bool isAxiallySymmetric() const
Definition: Mesh.h:137
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition: Process.h:265
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)
Builds expression trees of named functions dynamically at runtime.
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)
static void executeMemberDereferenced(Object &object, Method method, Container const &container, Args &&... args)
static void executeSelectedMemberDereferenced(Object &object, Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)
void preTimestepConcreteProcess(GlobalVector const &x, double const t, double const dt, const int) override
std::string const & getName() const
int getNumberOfComponents() const
Returns the number of components of the process variable.
PropertyVector< int > const * materialIDs(Mesh const &mesh)
Definition: Mesh.cpp:403
void initializeConcreteProcess(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order) override
Process specific initialization called by initialize().
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition: Mesh.h:108
void computeSecondaryVariableConcrete(double const t, GlobalVector const &x, int const process_id) override
void transformVariableFromGlobalVector(GlobalVector const &input_vector, int const variable_id, NumLib::LocalToGlobalIndexMap const &local_to_global_index_map, MeshLib::PropertyVector< double > &output_vector, Functor mapFunction)
Definition: DOFTableUtil.h:59
CoupledSolutionsForStaggeredScheme * _coupled_solutions
Definition: Process.h:280
Handles configuration of several secondary variables from the project file.
Definition of the class Properties that implements a container of properties.
unsigned getDimension() const
Returns the dimension of the mesh (determined by the maximum dimension over all elements).
Definition: Mesh.h:81
Definition of the MeshInformation class.
std::size_t getNumberOfElements() const
Get the number of elements.
Definition: Mesh.h:96
#define OGS_FATAL(fmt,...)
Definition: Error.h:63
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_fracture_nodes
const bool _use_monolithic_scheme
Definition: Process.h:276
VectorMatrixAssembler _global_assembler
Definition: Process.h:274
Definition of the ElementStatus class.
Ordering data by component type.
std::vector< MeshLib::Element * > _vec_fracture_matrix_elements
std::vector< std::size_t > & getActiveElementIDs() const