OGS 6.2.1-97-g73d1aeda3
HeatTransportBHEProcess.cpp
Go to the documentation of this file.
1 
11 
12 #include <cassert>
15 
18 
21 
22 namespace ProcessLib
23 {
24 namespace HeatTransportBHE
25 {
27  std::string name,
28  MeshLib::Mesh& mesh,
29  std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
30  std::vector<std::unique_ptr<ParameterLib::ParameterBase>> const& parameters,
31  unsigned const integration_order,
32  std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>&&
33  process_variables,
34  HeatTransportBHEProcessData&& process_data,
35  SecondaryVariableCollection&& secondary_variables,
36  NumLib::NamedFunctionCaller&& named_function_caller)
37  : Process(std::move(name), mesh, std::move(jacobian_assembler), parameters,
38  integration_order, std::move(process_variables),
39  std::move(secondary_variables), std::move(named_function_caller)),
40  _process_data(std::move(process_data)),
41  _bheMeshData(getBHEDataInMesh(mesh))
42 {
43  if (_bheMeshData.BHE_mat_IDs.size() !=
45  {
46  OGS_FATAL(
47  "The number of the given BHE properties (%d) are not consistent "
48  "with the number of BHE groups in the mesh (%d).",
50  _bheMeshData.BHE_mat_IDs.size());
51  }
52 
53  auto material_ids = MeshLib::materialIDs(mesh);
54  if (material_ids == nullptr)
55  {
56  OGS_FATAL("Not able to get material IDs! ");
57  }
58 
60 
61  // create a map from a material ID to a BHE ID
62  for (int i = 0; i < static_cast<int>(_bheMeshData.BHE_mat_IDs.size()); i++)
63  {
64  // fill in the map structure
66  i;
67  }
68 }
69 
71 {
72  // Create single component dof in every of the mesh's nodes.
74  std::make_unique<MeshLib::MeshSubset>(_mesh, _mesh.getNodes());
75 
76  //
77  // Soil temperature variable defined on the whole mesh.
78  //
80  std::make_unique<MeshLib::MeshSubset>(_mesh, _mesh.getNodes());
81  std::vector<MeshLib::MeshSubset> all_mesh_subsets{*_mesh_subset_soil_nodes};
82 
83  std::vector<std::vector<MeshLib::Element*> const*> vec_var_elements;
84  vec_var_elements.push_back(&(_mesh.getElements()));
85 
86  std::vector<int> vec_n_components{
87  1}; // one component for the soil temperature variable.
88 
89  //
90  // BHE nodes with BHE type dependend number of variables.
91  //
92  int const n_BHEs = _process_data._vec_BHE_property.size();
93  assert(n_BHEs == static_cast<int>(_bheMeshData.BHE_mat_IDs.size()));
94  assert(n_BHEs == static_cast<int>(_bheMeshData.BHE_nodes.size()));
95  assert(n_BHEs == static_cast<int>(_bheMeshData.BHE_elements.size()));
96 
97  // the BHE nodes need to be cherry-picked from the vector
98  for (int i = 0; i < n_BHEs; i++)
99  {
100  auto const number_of_unknowns =
101  visit([](auto const& bhe) { return bhe.number_of_unknowns; },
103  auto const& bhe_nodes = _bheMeshData.BHE_nodes[i];
104  auto const& bhe_elements = _bheMeshData.BHE_elements[i];
105 
106  // All the BHE nodes have additional variables.
107  _mesh_subset_BHE_nodes.push_back(
108  std::make_unique<MeshLib::MeshSubset const>(_mesh, bhe_nodes));
109 
110  std::generate_n(
111  std::back_inserter(all_mesh_subsets),
112  // Here the number of components equals to the
113  // number of unknowns on the BHE
114  number_of_unknowns,
115  [& ms = _mesh_subset_BHE_nodes.back()]() { return *ms; });
116 
117  vec_n_components.push_back(number_of_unknowns);
118  vec_var_elements.push_back(&bhe_elements);
119  }
120 
122  std::make_unique<NumLib::LocalToGlobalIndexMap>(
123  std::move(all_mesh_subsets),
124  vec_n_components,
125  vec_var_elements,
127 
128  // in case of debugging the dof table, activate the following line
129  // std::cout << *_local_to_global_index_map << "\n";
130 }
131 
133  NumLib::LocalToGlobalIndexMap const& dof_table,
134  MeshLib::Mesh const& mesh,
135  unsigned const integration_order)
136 {
137  // Quick access map to BHE's through element ids.
138  std::unordered_map<std::size_t, BHE::BHETypes*> element_to_bhe_map;
139  int const n_BHEs = _process_data._vec_BHE_property.size();
140  for (int i = 0; i < n_BHEs; i++)
141  {
142  auto const& bhe_elements = _bheMeshData.BHE_elements[i];
143  for (auto const& e : bhe_elements)
144  {
145  element_to_bhe_map[e->getID()] =
147  }
148  }
149 
150  assert(mesh.getDimension() == 3);
153  mesh.getElements(), dof_table, _local_assemblers, element_to_bhe_map,
154  mesh.isAxiallySymmetric(), integration_order, _process_data);
155 
156  // Create BHE boundary conditions for each of the BHEs
158 }
159 
161  GlobalVector const& x,
162  GlobalMatrix& M,
163  GlobalMatrix& K,
164  GlobalVector& b)
165 {
166  DBUG("Assemble HeatTransportBHE process.");
167 
168  const int process_id = 0;
169  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
170 
171  std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
172  dof_table = {std::ref(*_local_to_global_index_map)};
173  // Call global assembler for each local assembly item.
176  pv.getActiveElementIDs(), dof_table, t, x, M, K, b,
178 }
179 
181  const double /*t*/, GlobalVector const& /*x*/, GlobalVector const& /*xdot*/,
182  const double /*dxdot_dx*/, const double /*dx_dx*/, GlobalMatrix& /*M*/,
183  GlobalMatrix& /*K*/, GlobalVector& /*b*/, GlobalMatrix& /*Jac*/)
184 {
185  OGS_FATAL(
186  "HeatTransportBHE: analytical Jacobian assembly is not implemented");
187 }
188 
190  const double t, GlobalVector const& x, int const process_id)
191 {
192  DBUG("Compute heat flux for HeatTransportBHE process.");
193 
194  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
198  getDOFTable(process_id), t, x, _coupled_solutions);
199 }
200 
202  std::vector<std::vector<MeshLib::Node*>> const& all_bhe_nodes)
203 {
204  const int process_id = 0;
205  auto& bcs = _boundary_conditions[process_id];
206 
207  int const n_BHEs = static_cast<int>(_process_data._vec_BHE_property.size());
208 
209  // for each BHE
210  for (int bhe_i = 0; bhe_i < n_BHEs; bhe_i++)
211  {
212  auto const& bhe_nodes = all_bhe_nodes[bhe_i];
213  // find the variable ID
214  // the soil temperature is 0-th variable
215  // the BHE temperature is therefore bhe_i + 1
216  const int variable_id = bhe_i + 1;
217 
218  // Bottom and top nodes w.r.t. the z coordinate.
219  auto const bottom_top_nodes = std::minmax_element(
220  begin(bhe_nodes), end(bhe_nodes),
221  [&](auto const& a, auto const& b) {
222  return a->getCoords()[2] < b->getCoords()[2];
223  });
224  auto const bc_bottom_node_id = (*bottom_top_nodes.first)->getID();
225  auto const bc_top_node_id = (*bottom_top_nodes.second)->getID();
226 
227  auto get_global_bhe_bc_indices =
228  [&](std::size_t const node_id,
229  std::pair<int, int> const& in_out_component_id) {
230  return std::make_pair(
231  _local_to_global_index_map->getGlobalIndex(
232  {_mesh.getID(), MeshLib::MeshItemType::Node, node_id},
233  variable_id, in_out_component_id.first),
234  _local_to_global_index_map->getGlobalIndex(
236  variable_id, in_out_component_id.second));
237  };
238 
239  auto createBCs = [&](auto& bhe) {
240  for (auto const& in_out_component_id :
241  bhe.inflow_outflow_bc_component_ids)
242  {
243  // Top, inflow.
244  bcs.addBoundaryCondition(
246  get_global_bhe_bc_indices(bc_top_node_id,
247  in_out_component_id),
248  [&bhe](double const T, double const t) {
249  return bhe.updateFlowRateAndTemperature(T, t);
250  }));
251 
252  // Bottom, outflow.
253  bcs.addBoundaryCondition(
255  get_global_bhe_bc_indices(bc_bottom_node_id,
256  in_out_component_id)));
257  }
258  };
259  visit(createBCs, _process_data._vec_BHE_property[bhe_i]);
260  }
261 }
262 } // namespace HeatTransportBHE
263 } // namespace ProcessLib
HeatTransportBHEProcess(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, HeatTransportBHEProcessData &&process_data, SecondaryVariableCollection &&secondary_variables, NumLib::NamedFunctionCaller &&named_function_caller)
MeshLib::Mesh & _mesh
Definition: Process.h:287
virtual NumLib::LocalToGlobalIndexMap const & getDOFTable(const int) const
Definition: Process.h:117
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:288
std::vector< std::vector< MeshLib::Element * > > BHE_elements
Definition: MeshUtils.h:36
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:125
std::size_t getID() const
Get id of the mesh.
Definition: Mesh.h:123
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)
std::vector< BoundaryConditionCollection > _boundary_conditions
Definition: Process.h:331
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
bool isAxiallySymmetric() const
Definition: Mesh.h:137
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition: Process.h:290
void createBHEBoundaryConditionTopBottom(std::vector< std::vector< MeshLib::Node *>> const &all_bhe_nodes)
Builds expression trees of named functions dynamically at runtime.
void computeSecondaryVariableConcrete(double const t, GlobalVector const &x, int const process_id) override
std::vector< std::unique_ptr< HeatTransportBHELocalAssemblerInterface > > _local_assemblers
static void executeSelectedMemberDereferenced(Object &object, Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)
void initializeConcreteProcess(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order) override
Process specific initialization called by initialize().
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_soil_nodes
void createLocalAssemblers(std::vector< MeshLib::Element *> const &mesh_elements, NumLib::LocalToGlobalIndexMap const &dof_table, std::vector< std::unique_ptr< LocalAssemblerInterface >> &local_assemblers, ExtraCtorArgs &&... extra_ctor_args)
std::unique_ptr< BHEBottomDirichletBoundaryCondition > createBHEBottomDirichletBoundaryCondition(std::pair< GlobalIndexType, GlobalIndexType > &&in_out_global_indices)
std::vector< std::vector< MeshLib::Node * > > BHE_nodes
Definition: MeshUtils.h:37
std::vector< std::unique_ptr< MeshLib::MeshSubset const > > _mesh_subset_BHE_nodes
PropertyVector< int > const * materialIDs(Mesh const &mesh)
Definition: Mesh.cpp:403
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition: Mesh.h:108
CoupledSolutionsForStaggeredScheme * _coupled_solutions
Definition: Process.h:305
Handles configuration of several secondary variables from the project file.
void assembleConcreteProcess(const double t, GlobalVector const &x, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b) override
unsigned getDimension() const
Returns the dimension of the mesh (determined by the maximum dimension over all elements).
Definition: Mesh.h:81
#define OGS_FATAL(fmt,...)
Definition: Error.h:63
BHEMeshData getBHEDataInMesh(MeshLib::Mesh const &mesh)
Definition: MeshUtils.cpp:49
VectorMatrixAssembler _global_assembler
Definition: Process.h:299
std::unique_ptr< BHEInflowDirichletBoundaryCondition< BHEUpdateCallback > > createBHEInflowDirichletBoundaryCondition(std::pair< GlobalIndexType, GlobalIndexType > &&in_out_global_indices, BHEUpdateCallback bhe_update_callback)
Ordering data by component type.
std::vector< std::size_t > & getActiveElementIDs() const