OGS
StokesFlowProcess.cpp
Go to the documentation of this file.
1 
11 #include "StokesFlowProcess.h"
12 
13 #include <cassert>
14 
15 #include "MeshLib/Elements/Utils.h"
17 
18 namespace ProcessLib
19 {
20 namespace StokesFlow
21 {
22 template <int GlobalDim>
24  std::string name,
25  MeshLib::Mesh& mesh,
26  std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
27  std::vector<std::unique_ptr<ParameterLib::ParameterBase>> const& parameters,
28  unsigned const integration_order,
29  std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>&&
30  process_variables,
31  StokesFlowProcessData&& process_data,
32  SecondaryVariableCollection&& secondary_variables,
33  bool const use_monolithic_scheme)
34  : Process(std::move(name), mesh, std::move(jacobian_assembler), parameters,
35  integration_order, std::move(process_variables),
36  std::move(secondary_variables), use_monolithic_scheme),
37  _process_data(std::move(process_data))
38 {
39 }
40 
41 template <int GlobalDim>
43 {
44  // Create single component dof in every of the mesh's nodes.
45  _mesh_subset_all_nodes =
46  std::make_unique<MeshLib::MeshSubset>(_mesh, _mesh.getNodes());
47  // Create single component dof in the mesh's base nodes.
48  _base_nodes = MeshLib::getBaseNodes(_mesh.getElements());
49  _mesh_subset_base_nodes =
50  std::make_unique<MeshLib::MeshSubset>(_mesh, _base_nodes);
51 
52  // TODO move the two data members somewhere else.
53  // for extrapolation of secondary variables
54  std::vector<MeshLib::MeshSubset> all_mesh_subsets_single_component{
55  *_mesh_subset_all_nodes};
56  _local_to_global_index_map_single_component =
57  std::make_unique<NumLib::LocalToGlobalIndexMap>(
58  std::move(all_mesh_subsets_single_component),
59  // by location order is needed for output
61 
62  assert(_use_monolithic_scheme);
63  {
64  // For vector variables, in this case liquid velocity.
65  int const process_id = 0;
66  std::vector<MeshLib::MeshSubset> all_mesh_subsets(
67  getProcessVariables(process_id)[0]
68  .get()
69  .getNumberOfGlobalComponents(),
70  *_mesh_subset_all_nodes);
71 
72  // For scalar variables, in this case pressure.
73  all_mesh_subsets.push_back(*_mesh_subset_base_nodes);
74 
75  std::vector<int> const vec_n_components{GlobalDim, 1};
76 
77  _local_to_global_index_map =
78  std::make_unique<NumLib::LocalToGlobalIndexMap>(
79  std::move(all_mesh_subsets), vec_n_components,
81  assert(_local_to_global_index_map);
82  }
83 }
84 
85 template <int GlobalDim>
87  NumLib::LocalToGlobalIndexMap const& dof_table,
88  MeshLib::Mesh const& mesh,
89  unsigned const integration_order)
90 {
91  ProcessLib::createLocalAssemblersStokes<GlobalDim, LocalAssemblerData>(
92  mesh.getElements(), dof_table, _local_assemblers,
93  NumLib::IntegrationOrder{integration_order}, mesh.isAxiallySymmetric(),
94  _process_data);
95 
96  _process_data.pressure_interpolated =
97  MeshLib::getOrCreateMeshProperty<double>(
98  const_cast<MeshLib::Mesh&>(mesh), "pressure_interpolated",
100 }
101 
102 template <int GlobalDim>
104 {
105  assert(_use_monolithic_scheme);
106  {
107  int const process_id = 0;
108  initializeProcessBoundaryConditionsAndSourceTerms(
109  *_local_to_global_index_map, process_id);
110  }
111 }
112 
113 template <int GlobalDim>
116  const int /*process_id*/) const
117 {
118  assert(_use_monolithic_scheme);
119  {
120  auto const& l = *_local_to_global_index_map;
121 
122  return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
123  &l.getGhostIndices(), &this->_sparsity_pattern};
124  }
125 }
126 
127 template <int GlobalDim>
129  const double t, double const dt, std::vector<GlobalVector*> const& x,
130  std::vector<GlobalVector*> const& xdot, int const process_id,
132 {
133  DBUG("Assemble StokesFlowProcess.");
134 
135  std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
136  dof_tables;
137  assert(_use_monolithic_scheme);
138  {
139  dof_tables.push_back(std::ref(*_local_to_global_index_map));
140  }
141 
142  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
143  // Call global assembler for each local assembly item.
145  _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers,
146  pv.getActiveElementIDs(), dof_tables, t, dt, x, xdot, process_id, M, K,
147  b);
148 }
149 
150 template <int GlobalDim>
152  const double /*t*/, double const /*dt*/,
153  std::vector<GlobalVector*> const& /*x*/,
154  std::vector<GlobalVector*> const& /*xdot*/, int const /*process_id*/,
155  GlobalMatrix& /*M*/, GlobalMatrix& /*K*/, GlobalVector& /*b*/,
156  GlobalMatrix& /*Jac*/)
157 {
158  OGS_FATAL(
159  "Assembly of Jacobian matrix has not yet been implemented for "
160  "StokesFlowProcess.");
161 }
162 
163 template <int GlobalDim>
165  double const t,
166  double const dt,
167  std::vector<GlobalVector*> const& x,
168  GlobalVector const& x_dot,
169  int const process_id)
170 {
171  if (process_id != 0)
172  {
173  return;
174  }
175 
176  std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
177  dof_tables.reserve(x.size());
178  assert(_use_monolithic_scheme);
179  {
180  dof_tables.push_back(_local_to_global_index_map.get());
181  }
182 
183  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
186  _local_assemblers, pv.getActiveElementIDs(), dof_tables, t, dt, x,
187  x_dot, process_id);
188 }
189 
190 template <int GlobalDim>
192  std::vector<GlobalVector*> const& x,
193  const double t,
194  const double dt,
195  int const process_id)
196 {
197  if (process_id != 0)
198  {
199  return;
200  }
201 
202  std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
203  dof_tables.reserve(x.size());
204  assert(_use_monolithic_scheme);
205  {
206  dof_tables.push_back(_local_to_global_index_map.get());
207  }
208 
209  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
212  pv.getActiveElementIDs(), dof_tables, x, t, dt);
213 }
214 
215 template <int GlobalDim>
217  const int /*process_id*/) const
218 {
219  assert(_use_monolithic_scheme);
220  {
221  return *_local_to_global_index_map;
222  }
223 }
224 
225 template class StokesFlowProcess<2>;
226 } // namespace StokesFlow
227 } // namespace ProcessLib
#define OGS_FATAL(...)
Definition: Error.h:26
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition: Logging.h:29
Global vector based on Eigen vector.
Definition: EigenVector.h:28
bool isAxiallySymmetric() const
Definition: Mesh.h:131
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition: Mesh.h:103
virtual void postTimestep(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, double const t, double const dt)
virtual void computeSecondaryVariable(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_dot, int const process_id)
std::vector< std::size_t > const & getActiveElementIDs() const
Handles configuration of several secondary variables from the project file.
StokesFlowProcess(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, StokesFlowProcessData &&process_data, SecondaryVariableCollection &&secondary_variables, bool const use_monolithic_scheme)
void initializeConcreteProcess(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order) override
Process specific initialization called by initialize().
void computeSecondaryVariableConcrete(double const, double const, std::vector< GlobalVector * > const &x, GlobalVector const &, int const) override
void assembleWithJacobianConcreteProcess(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, GlobalMatrix &Jac) 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
NumLib::LocalToGlobalIndexMap const & getDOFTable(const int) const override
void postTimestepConcreteProcess(std::vector< GlobalVector * > const &x, const double t, const double dt, int const process_id) override
MathLib::MatrixSpecifications getMatrixSpecifications(const int) const override
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)
static const double t
std::vector< Node * > getBaseNodes(std::vector< Element * > const &elements)
Definition: Utils.h:26
@ BY_LOCATION
Ordering data by spatial location.
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)