OGS
ProcessLib::HeatTransportBHE::HeatTransportBHEProcess Class Referencefinal

Detailed Description

Definition at line 24 of file HeatTransportBHEProcess.h.

#include <HeatTransportBHEProcess.h>

Inheritance diagram for ProcessLib::HeatTransportBHE::HeatTransportBHEProcess:
[legend]
Collaboration diagram for ProcessLib::HeatTransportBHE::HeatTransportBHEProcess:
[legend]

Public Member Functions

 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, BHEMeshData &&bhe_mesh_data)
 
void computeSecondaryVariableConcrete (double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_prev, int const process_id) override
 
ODESystem interface
bool isLinear () const override
 
bool requiresNormalization () const override
 
- Public Member Functions inherited from ProcessLib::Process
 Process (std::string name_, MeshLib::Mesh &mesh, std::unique_ptr< 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, SecondaryVariableCollection &&secondary_variables, const bool use_monolithic_scheme=true)
 
void preTimestep (std::vector< GlobalVector * > const &x, const double t, const double delta_t, const int process_id)
 Preprocessing before starting assembly for new timestep.
 
void postTimestep (std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, const double t, const double delta_t, int const process_id)
 Postprocessing after a complete timestep.
 
void postNonLinearSolver (std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, const double t, double const dt, int const process_id)
 
void preIteration (const unsigned iter, GlobalVector const &x) final
 
void computeSecondaryVariable (double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_prev, int const process_id)
 compute secondary variables for the coupled equations or for output.
 
NumLib::IterationResult postIteration (GlobalVector const &x) final
 
void initialize (std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
 
void setInitialConditions (std::vector< GlobalVector * > &process_solutions, std::vector< GlobalVector * > const &process_solutions_prev, double const t, int const process_id)
 
MathLib::MatrixSpecifications getMatrixSpecifications (const int process_id) const override
 
void updateDeactivatedSubdomains (double const time, const int process_id)
 
virtual bool isMonolithicSchemeUsed () const
 
virtual void extrapolateIntegrationPointValuesToNodes (const double, std::vector< GlobalVector * > const &, std::vector< GlobalVector * > &)
 
void preAssemble (const double t, double const dt, GlobalVector const &x) final
 
void assemble (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) final
 
void assembleWithJacobian (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) final
 
void preOutput (const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id)
 
std::vector< NumLib::IndexValueVector< GlobalIndexType > > const * getKnownSolutions (double const t, GlobalVector const &x, int const process_id) const final
 
virtual NumLib::LocalToGlobalIndexMap const & getDOFTable (const int) const
 
MeshLib::MeshgetMesh () const
 
std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > const & getProcessVariables () const
 
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables (const int process_id) const
 
std::vector< std::size_t > const & getActiveElementIDs () const
 
SecondaryVariableCollection const & getSecondaryVariables () const
 
std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > const & getIntegrationPointWriters () const
 
virtual Eigen::Vector3d getFlux (std::size_t, MathLib::Point3d const &, double const, std::vector< GlobalVector * > const &) const
 
virtual void solveReactionEquation (std::vector< GlobalVector * > &, std::vector< GlobalVector * > const &, double const, double const, NumLib::EquationSystem &, int const)
 
bool requiresNormalization () const override
 
- Public Member Functions inherited from ProcessLib::SubmeshAssemblySupport
virtual std::vector< std::vector< std::string > > initializeAssemblyOnSubmeshes (std::vector< std::reference_wrapper< MeshLib::Mesh > > const &meshes)
 
virtual ~SubmeshAssemblySupport ()=default
 

Private Member Functions

void constructDofTable () override
 
void initializeConcreteProcess (NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order) override
 Process specific initialization called by initialize().
 
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 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
 
void createBHEBoundaryConditionTopBottom (std::vector< std::vector< MeshLib::Node * > > const &all_bhe_nodes)
 
void preTimestepConcreteProcess (std::vector< GlobalVector * > const &x, const double t, const double dt, int const process_id) 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
 
void algebraicBcConcreteProcess (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)
 
NumLib::IterationResult postIterationConcreteProcess (GlobalVector const &x) override
 

Private Attributes

HeatTransportBHEProcessData _process_data
 
std::vector< std::unique_ptr< HeatTransportBHELocalAssemblerInterface > > _local_assemblers
 
std::vector< std::unique_ptr< MeshLib::MeshSubset const > > _mesh_subset_BHE_nodes
 
std::vector< std::unique_ptr< MeshLib::MeshSubset const > > _mesh_subset_BHE_soil_nodes
 
std::vector< std::tuple< std::size_t, GlobalIndexType, GlobalIndexType > > _vec_top_BHE_node_indices
 
std::vector< std::tuple< std::size_t, GlobalIndexType, GlobalIndexType > > _vec_bottom_BHE_node_indices
 
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_soil_nodes
 
const BHEMeshData _bheMeshData
 

Additional Inherited Members

- Public Attributes inherited from ProcessLib::Process
std::string const name
 
- Static Public Attributes inherited from ProcessLib::Process
static PROCESSLIB_EXPORT const std::string constant_one_parameter_name = "constant_one"
 
- Protected Member Functions inherited from ProcessLib::Process
std::vector< NumLib::LocalToGlobalIndexMap const * > getDOFTables (int const number_of_processes) const
 
NumLib::ExtrapolatorgetExtrapolator () const
 
NumLib::LocalToGlobalIndexMap const & getSingleComponentDOFTable () const
 
void initializeProcessBoundaryConditionsAndSourceTerms (const NumLib::LocalToGlobalIndexMap &dof_table, const int process_id, std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
 
void constructMonolithicProcessDofTable ()
 
void constructDofTableOfSpecifiedProcessStaggeredScheme (const int specified_process_id)
 
virtual std::tuple< NumLib::LocalToGlobalIndexMap *, bool > getDOFTableForExtrapolatorData () const
 
std::vector< GlobalIndexTypegetIndicesOfResiduumWithoutInitialCompensation () const override
 
- Protected Attributes inherited from ProcessLib::Process
MeshLib::Mesh_mesh
 
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_all_nodes
 
std::unique_ptr< NumLib::LocalToGlobalIndexMap_local_to_global_index_map
 
SecondaryVariableCollection _secondary_variables
 
CellAverageData cell_average_data_
 
std::unique_ptr< ProcessLib::AbstractJacobianAssembler_jacobian_assembler
 
VectorMatrixAssembler _global_assembler
 
const bool _use_monolithic_scheme
 
unsigned const _integration_order
 
std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > _integration_point_writer
 
GlobalSparsityPattern _sparsity_pattern
 
std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > _process_variables
 
std::vector< BoundaryConditionCollection_boundary_conditions
 

Constructor & Destructor Documentation

◆ HeatTransportBHEProcess()

ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::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,
BHEMeshData && bhe_mesh_data )

Definition at line 28 of file HeatTransportBHEProcess.cpp.

39 : Process(std::move(name), mesh, std::move(jacobian_assembler), parameters,
40 integration_order, std::move(process_variables),
41 std::move(secondary_variables)),
42 _process_data(std::move(process_data)),
43 _bheMeshData(std::move(bhe_mesh_data))
44{
45 if (_bheMeshData.BHE_mat_IDs.size() !=
47 {
49 "The number of the given BHE properties ({:d}) are not consistent "
50 "with the number of BHE groups in the mesh ({:d}).",
53 }
54
55 auto material_ids = MeshLib::materialIDs(mesh);
56 if (material_ids == nullptr)
57 {
58 OGS_FATAL("Not able to get material IDs! ");
59 }
60
62
63 // create a map from a material ID to a BHE ID
64 for (int i = 0; i < static_cast<int>(_bheMeshData.BHE_mat_IDs.size()); i++)
65 {
66 // fill in the map structure
68 i;
69 }
70}
#define OGS_FATAL(...)
Definition Error.h:26
std::string const name
Definition Process.h:362
Process(std::string name_, MeshLib::Mesh &mesh, std::unique_ptr< 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, SecondaryVariableCollection &&secondary_variables, const bool use_monolithic_scheme=true)
Definition Process.cpp:44
PropertyVector< int > const * materialIDs(Mesh const &mesh)
Definition Mesh.cpp:268

References _bheMeshData, ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::_map_materialID_to_BHE_ID, ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::_mesh_prop_materialIDs, _process_data, ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::_vec_BHE_property, ProcessLib::HeatTransportBHE::BHEMeshData::BHE_mat_IDs, MeshLib::materialIDs(), and OGS_FATAL.

Member Function Documentation

◆ algebraicBcConcreteProcess()

void ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::algebraicBcConcreteProcess ( 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 )
private

Definition at line 373 of file HeatTransportBHEProcess.cpp.

379{
380#ifndef USE_PETSC
381 auto M_normal = M.getRawMatrix();
382 auto K_normal = K.getRawMatrix();
383 auto n_original_rows = K_normal.rows();
384 auto const n_BHE_bottom_pairs = _vec_bottom_BHE_node_indices.size();
385 auto const n_BHE_top_pairs = _vec_top_BHE_node_indices.size();
386
387 // apply weighting factor based on the max value from column wise inner
388 // product and scale it with user defined value
389 const double w_val =
391 (Eigen::RowVectorXd::Ones(K_normal.rows()) * K_normal.cwiseAbs())
392 .maxCoeff();
393
394 M_normal.conservativeResize(
395 M_normal.rows() + n_BHE_bottom_pairs + n_BHE_top_pairs,
396 M_normal.cols());
397 K_normal.conservativeResize(
398 K_normal.rows() + n_BHE_bottom_pairs + n_BHE_top_pairs,
399 K_normal.cols());
400
401 for (std::size_t i = 0; i < n_BHE_bottom_pairs; i++)
402 {
403 Eigen::SparseVector<double> M_Plus(M_normal.cols());
404 M_Plus.setZero();
405 M_normal.row(n_original_rows + i) = M_Plus;
406
407 Eigen::SparseVector<double> K_Plus(K_normal.cols());
408 K_Plus.setZero();
409
410 auto const [bhe_idx, first_BHE_bottom_index, second_BHE_bottom_index] =
412
413 K_Plus.insert(first_BHE_bottom_index) = w_val;
414 K_Plus.insert(second_BHE_bottom_index) = -w_val;
415
416 K_normal.row(n_original_rows + i) = K_Plus;
417 }
418
419 auto b_normal = b.getRawVector();
420 Eigen::SparseVector<double> b_Plus(b_normal.rows() + n_BHE_bottom_pairs +
421 n_BHE_top_pairs);
422 b_Plus.setZero();
423
424 // Copy values from the original column vector to the modified one
425 for (int i = 0; i < b_normal.innerSize(); ++i)
426 {
427 b_Plus.insert(i) = b_normal.coeff(i);
428 }
429
430 for (std::size_t i = 0; i < n_BHE_top_pairs; i++)
431 {
432 Eigen::SparseVector<double> M_Plus(M_normal.cols());
433 M_Plus.setZero();
434 M_normal.row(n_original_rows + n_BHE_bottom_pairs + i) = M_Plus;
435
436 Eigen::SparseVector<double> K_Plus(K_normal.cols());
437 K_Plus.setZero();
438
439 auto const [bhe_idx, first_BHE_top_index, second_BHE_top_index] =
441
442 auto first_BHE_top_index_pair = first_BHE_top_index;
443 auto second_BHE_top_index_pair = second_BHE_top_index;
444
445 K_Plus.insert(first_BHE_top_index_pair) =
446 w_val; // for power BC, the inflow node must be positive
447 K_Plus.insert(second_BHE_top_index_pair) =
448 -w_val; // for power BC, the outflow node must be negative
449
450 K_normal.row(n_original_rows + n_BHE_bottom_pairs + i) = K_Plus;
451
452 // get the delta_T value here
453 double const T_out = (*x[0])[second_BHE_top_index_pair];
454
455 auto calculate_delta_T = [&](auto& bhe)
456 {
457 auto const T_in = bhe.updateFlowRateAndTemperature(T_out, t);
458 return T_in - T_out;
459 };
460 auto delta_T = std::visit(calculate_delta_T,
462
463 b_Plus.insert(n_original_rows + n_BHE_bottom_pairs + i) =
464 delta_T * w_val;
465 }
466
467 M.getRawMatrix() = M_normal;
468 K.getRawMatrix() = K_normal;
469 b.getRawVector() = b_Plus;
470#else
471 OGS_FATAL(
472 "The Algebraic Boundary Condition is not implemented for use with "
473 "PETsc Library! Simulation will be terminated.");
474#endif
475}
RawMatrixType & getRawMatrix()
RawVectorType & getRawVector()
return a raw Eigen vector object
std::vector< std::tuple< std::size_t, GlobalIndexType, GlobalIndexType > > _vec_bottom_BHE_node_indices
std::vector< std::tuple< std::size_t, GlobalIndexType, GlobalIndexType > > _vec_top_BHE_node_indices

References ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::_algebraic_BC_Setting, _process_data, ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::_vec_BHE_property, _vec_bottom_BHE_node_indices, _vec_top_BHE_node_indices, ProcessLib::HeatTransportBHE::AlgebraicBCSetting::_weighting_factor, and OGS_FATAL.

Referenced by assembleConcreteProcess().

◆ assembleConcreteProcess()

void ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::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 )
overrideprivatevirtual

Implements ProcessLib::Process.

Definition at line 184 of file HeatTransportBHEProcess.cpp.

188{
189 DBUG("Assemble HeatTransportBHE process.");
190
191 std::vector<NumLib::LocalToGlobalIndexMap const*> dof_table = {
193 // Call global assembler for each local assembly item.
196 getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K,
197 &b);
198 // Algebraic BC procedure.
200 {
201 algebraicBcConcreteProcess(t, dt, x, x_prev, process_id, M, K, b);
202 }
203
204 //_global_output(t, process_id, M, K, b);
205}
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:30
std::vector< std::unique_ptr< HeatTransportBHELocalAssemblerInterface > > _local_assemblers
void algebraicBcConcreteProcess(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)
std::vector< std::size_t > const & getActiveElementIDs() const
Definition Process.h:167
VectorMatrixAssembler _global_assembler
Definition Process.h:377
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition Process.h:368
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)
static void executeSelectedMemberDereferenced(Object &object, Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)

References ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::_algebraic_BC_Setting, ProcessLib::Process::_global_assembler, _local_assemblers, ProcessLib::Process::_local_to_global_index_map, _process_data, ProcessLib::HeatTransportBHE::AlgebraicBCSetting::_use_algebraic_bc, algebraicBcConcreteProcess(), ProcessLib::VectorMatrixAssembler::assemble(), DBUG(), NumLib::SerialExecutor::executeSelectedMemberDereferenced(), and ProcessLib::Process::getActiveElementIDs().

◆ assembleWithJacobianConcreteProcess()

void ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::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 )
overrideprivatevirtual

Implements ProcessLib::Process.

Definition at line 207 of file HeatTransportBHEProcess.cpp.

211{
212 DBUG("AssembleWithJacobian HeatTransportBHE process.");
213
214 std::vector<NumLib::LocalToGlobalIndexMap const*> dof_table = {
216
217 // Call global assembler for each local assembly item.
220 _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
221 process_id, &b, &Jac);
222}
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)

References ProcessLib::Process::_global_assembler, _local_assemblers, ProcessLib::Process::_local_to_global_index_map, ProcessLib::VectorMatrixAssembler::assembleWithJacobian(), DBUG(), NumLib::SerialExecutor::executeSelectedMemberDereferenced(), and ProcessLib::Process::getActiveElementIDs().

◆ computeSecondaryVariableConcrete()

void ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::computeSecondaryVariableConcrete ( double const t,
double const dt,
std::vector< GlobalVector * > const & x,
GlobalVector const & x_prev,
int const process_id )
overridevirtual

Reimplemented from ProcessLib::Process.

Definition at line 224 of file HeatTransportBHEProcess.cpp.

227{
228 DBUG("Compute heat flux for HeatTransportBHE process.");
229
230 std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
231 dof_tables.reserve(x.size());
232 std::generate_n(std::back_inserter(dof_tables), x.size(),
233 [&]() { return _local_to_global_index_map.get(); });
234
237 _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev,
238 process_id);
239}
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_prev, int const process_id)
static void executeSelectedMemberOnDereferenced(Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)

References _local_assemblers, ProcessLib::LocalAssemblerInterface::computeSecondaryVariable(), DBUG(), NumLib::SerialExecutor::executeSelectedMemberOnDereferenced(), and ProcessLib::Process::getActiveElementIDs().

◆ constructDofTable()

void ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::constructDofTable ( )
overrideprivatevirtual

This function is for general cases, in which all equations of the coupled processes have the same number of unknowns. For the general cases with the staggered scheme, all equations of the coupled processes share one DOF table hold by _local_to_global_index_map. Other cases can be considered by overloading this member function in the derived class.

Reimplemented from ProcessLib::Process.

Definition at line 72 of file HeatTransportBHEProcess.cpp.

73{
74 // Create single component dof in every of the mesh's nodes.
76 std::make_unique<MeshLib::MeshSubset>(_mesh, _mesh.getNodes());
77
78 //
79 // Soil temperature variable defined on the whole mesh.
80 //
82 std::make_unique<MeshLib::MeshSubset>(_mesh, _mesh.getNodes());
83 std::vector<MeshLib::MeshSubset> all_mesh_subsets{*_mesh_subset_soil_nodes};
84
85 std::vector<std::vector<MeshLib::Element*> const*> vec_var_elements;
86 vec_var_elements.push_back(&(_mesh.getElements()));
87
88 std::vector<int> vec_n_components{
89 1}; // one component for the soil temperature variable.
90
91 //
92 // BHE nodes with BHE type dependent number of variables.
93 //
94 int const n_BHEs = _process_data._vec_BHE_property.size();
95 assert(n_BHEs == static_cast<int>(_bheMeshData.BHE_mat_IDs.size()));
96 assert(n_BHEs == static_cast<int>(_bheMeshData.BHE_nodes.size()));
97 assert(n_BHEs == static_cast<int>(_bheMeshData.BHE_elements.size()));
98
99 // the BHE nodes need to be cherry-picked from the vector
100 for (int i = 0; i < n_BHEs; i++)
101 {
102 auto const number_of_unknowns =
103 visit([](auto const& bhe) { return bhe.number_of_unknowns; },
105 auto const& bhe_nodes = _bheMeshData.BHE_nodes[i];
106 auto const& bhe_elements = _bheMeshData.BHE_elements[i];
107
108 // All the BHE nodes have additional variables.
109 _mesh_subset_BHE_nodes.push_back(
110 std::make_unique<MeshLib::MeshSubset const>(_mesh, bhe_nodes));
111
112 std::generate_n(std::back_inserter(all_mesh_subsets),
113 // Here the number of components equals to the
114 // number of unknowns on the BHE
115 number_of_unknowns,
116 [&ms = _mesh_subset_BHE_nodes.back()]()
117 { return *ms; });
118
119 vec_n_components.push_back(number_of_unknowns);
120 vec_var_elements.push_back(&bhe_elements);
121 }
122
124 std::make_unique<NumLib::LocalToGlobalIndexMap>(
125 std::move(all_mesh_subsets),
126 vec_n_components,
127 vec_var_elements,
129
130 // in case of debugging the dof table, activate the following line
131 // std::cout << *_local_to_global_index_map << "\n";
132}
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition Mesh.h:106
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition Mesh.h:109
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_soil_nodes
std::vector< std::unique_ptr< MeshLib::MeshSubset const > > _mesh_subset_BHE_nodes
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_all_nodes
Definition Process.h:366
MeshLib::Mesh & _mesh
Definition Process.h:365
@ BY_COMPONENT
Ordering data by component type.
std::vector< std::vector< MeshLib::Node * > > BHE_nodes
Definition MeshUtils.h:38
std::vector< std::vector< MeshLib::Element * > > BHE_elements
Definition MeshUtils.h:37

References _bheMeshData, ProcessLib::Process::_local_to_global_index_map, ProcessLib::Process::_mesh, ProcessLib::Process::_mesh_subset_all_nodes, _mesh_subset_BHE_nodes, _mesh_subset_soil_nodes, _process_data, ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::_vec_BHE_property, ProcessLib::HeatTransportBHE::BHEMeshData::BHE_elements, ProcessLib::HeatTransportBHE::BHEMeshData::BHE_mat_IDs, ProcessLib::HeatTransportBHE::BHEMeshData::BHE_nodes, NumLib::BY_COMPONENT, MeshLib::Mesh::getElements(), and MeshLib::Mesh::getNodes().

◆ createBHEBoundaryConditionTopBottom()

void ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::createBHEBoundaryConditionTopBottom ( std::vector< std::vector< MeshLib::Node * > > const & all_bhe_nodes)
private

Definition at line 477 of file HeatTransportBHEProcess.cpp.

479{
480 const int process_id = 0;
481 auto& bcs = _boundary_conditions[process_id];
482
483 std::size_t const n_BHEs = _process_data._vec_BHE_property.size();
484
485 // for each BHE
486 for (std::size_t bhe_i = 0; bhe_i < n_BHEs; bhe_i++)
487 {
488 auto const& bhe_nodes = all_bhe_nodes[bhe_i];
489 // find the variable ID
490 // the soil temperature is 0-th variable
491 // the BHE temperature is therefore bhe_i + 1
492 const int variable_id = bhe_i + 1;
493
494 std::vector<MeshLib::Node*> bhe_boundary_nodes;
495
496 // cherry-pick the boundary nodes according to
497 // the number of connected line elements.
498 for (auto const& bhe_node : bhe_nodes)
499 {
500 // Count number of 1d elements connected with every BHE node.
501 auto const& connected_elements =
503 const std::size_t n_line_elements = std::count_if(
504 connected_elements.begin(), connected_elements.end(),
505 [](MeshLib::Element const* elem)
506 { return (elem->getDimension() == 1); });
507
508 if (n_line_elements == 1)
509 {
510 bhe_boundary_nodes.push_back(bhe_node);
511 }
512 }
513
514 if (bhe_boundary_nodes.size() != 2)
515 {
516 OGS_FATAL(
517 "Error!!! The BHE boundary nodes are not correctly found, "
518 "for every single BHE, there should be 2 boundary nodes.");
519 }
520
521 // For 1U, 2U, CXC, CXA type BHE, the node order in the boundary nodes
522 // vector should be rearranged according to its z coordinate in
523 // descending order. In these BHE types, the z coordinate on the top and
524 // bottom node is different. The BHE top node with a higher z coordinate
525 // should be placed at the first, while the BHE bottom node with a lower
526 // z coordinate should be placed at the second. For other horizontal BHE
527 // types e.g. 1P-type BHE, the z coordinate on the top and bottom node
528 // is identical. Thus the node order in the boundary nodes vector can
529 // not be rearranged according to its z coordinate. For these BHE types,
530 // the boundary node order is according to the default node id order in
531 // the model mesh.
532 // for 1P-type BHE
533 if ((*bhe_boundary_nodes[0])[2] == (*bhe_boundary_nodes[1])[2])
534 {
535 INFO(
536 "For 1P-type BHE, the BHE inflow and outflow "
537 "nodes are identified according to their mesh node id in "
538 "ascending order");
539 }
540 // for 1U, 2U, CXC, CXA type BHE
541 else
542 {
543 // swap the boundary nodes if the z coordinate of the
544 // first node is lower than it on the second node
545 if ((*bhe_boundary_nodes[0])[2] < (*bhe_boundary_nodes[1])[2])
546 {
547 std::swap(bhe_boundary_nodes[0], bhe_boundary_nodes[1]);
548 }
549 }
550
551 auto get_global_index =
552 [&](std::size_t const node_id, int const component)
553 {
554 return _local_to_global_index_map->getGlobalIndex(
556 variable_id, component);
557 };
558
559 auto get_global_bhe_bc_indices =
560 [&](std::array<
561 std::pair<std::size_t /*node_id*/, int /*component*/>, 2>
562 nodes_and_components)
563 {
564 return std::make_pair(
565 get_global_index(nodes_and_components[0].first,
566 nodes_and_components[0].second),
567 get_global_index(nodes_and_components[1].first,
568 nodes_and_components[1].second));
569 };
570
571 auto get_global_bhe_bc_indices_with_bhe_idx =
572 [&](std::size_t bhe_idx,
573 std::array<
574 std::pair<std::size_t /*node_id*/, int /*component*/>, 2>
575 nodes_and_components)
576 {
577 return std::make_tuple(
578 bhe_idx,
579 get_global_index(nodes_and_components[0].first,
580 nodes_and_components[0].second),
581 get_global_index(nodes_and_components[1].first,
582 nodes_and_components[1].second));
583 };
584
585 auto createBCs =
586 [&, bc_top_node_id = bhe_boundary_nodes[0]->getID(),
587 bc_bottom_node_id = bhe_boundary_nodes[1]->getID()](auto& bhe)
588 {
589 for (auto const& in_out_component_id :
590 bhe.inflow_outflow_bc_component_ids)
591 {
592 if (bhe.use_python_bcs ||
593 this->_process_data._use_server_communication)
594 // call BHEPythonBoundarycondition
595 {
596 if (this->_process_data
597 .py_bc_object) // the bc object exist
598 {
599 // apply the customized top, inflow BC.
600 bcs.addBoundaryCondition(
602 get_global_bhe_bc_indices(
603 bhe.getBHEInflowDirichletBCNodesAndComponents(
604 bc_top_node_id, bc_bottom_node_id,
605 in_out_component_id.first)),
606 bhe,
608 }
609 else
610 {
611 OGS_FATAL(
612 "The Python Boundary Condition was switched on, "
613 "but the data object does not exist! ");
614 }
615 }
616 else
617 {
620 bhe.isPowerBC())
621 {
622 // for algebraic_bc method, record the pair of indices
623 // in a separate vector
625 get_global_bhe_bc_indices_with_bhe_idx(
626 bhe_i,
627 {{{bc_top_node_id, in_out_component_id.first},
628 {bc_top_node_id,
629 in_out_component_id.second}}}));
630 }
631 else
632 {
633 // Top, inflow, normal case
634 bcs.addBoundaryCondition(
636 get_global_bhe_bc_indices(
637 bhe.getBHEInflowDirichletBCNodesAndComponents(
638 bc_top_node_id, bc_bottom_node_id,
639 in_out_component_id.first)),
640 [&bhe](double const T, double const t) {
641 return bhe.updateFlowRateAndTemperature(T,
642 t);
643 }));
644 }
645 }
646
647 auto const bottom_nodes_and_components =
648 bhe.getBHEBottomDirichletBCNodesAndComponents(
649 bc_bottom_node_id,
650 in_out_component_id.first,
651 in_out_component_id.second);
652
653 if (bottom_nodes_and_components &&
656 {
657 // Bottom, outflow, all cases | not needed for algebraic_bc
658 // method
659 bcs.addBoundaryCondition(
661 get_global_bhe_bc_indices(
662 {{{bc_bottom_node_id,
663 in_out_component_id.first},
664 {bc_bottom_node_id,
665 in_out_component_id.second}}})));
666 }
667 else if (bottom_nodes_and_components &&
670 {
671 // for algebraic_bc method, record the pair of indices in a
672 // separate vector
674 get_global_bhe_bc_indices_with_bhe_idx(
675 bhe_i,
676 {{{bc_bottom_node_id, in_out_component_id.first},
677 {bc_bottom_node_id,
678 in_out_component_id.second}}}));
679 }
680 }
681 };
682 visit(createBCs, _process_data._vec_BHE_property[bhe_i]);
683 }
684}
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:35
std::size_t getID() const
Get id of the mesh.
Definition Mesh.h:121
std::vector< Element const * > const & getElementsConnectedToNode(std::size_t node_id) const
Definition Mesh.cpp:256
std::vector< BoundaryConditionCollection > _boundary_conditions
Definition Process.h:405
std::unique_ptr< BHEInflowDirichletBoundaryCondition< BHEUpdateCallback > > createBHEInflowDirichletBoundaryCondition(std::pair< GlobalIndexType, GlobalIndexType > &&in_out_global_indices, BHEUpdateCallback bhe_update_callback)
std::unique_ptr< BHEBottomDirichletBoundaryCondition > createBHEBottomDirichletBoundaryCondition(std::pair< GlobalIndexType, GlobalIndexType > &&in_out_global_indices)
std::unique_ptr< BHEInflowPythonBoundaryCondition< BHEType > > createBHEInflowPythonBoundaryCondition(std::pair< GlobalIndexType, GlobalIndexType > &&in_out_global_indices, BHEType &bhe, BHEInflowPythonBoundaryConditionPythonSideInterface &py_bc_object)
BHEInflowPythonBoundaryConditionPythonSideInterface * py_bc_object
Python object computing BC values.

References ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::_algebraic_BC_Setting, ProcessLib::Process::_boundary_conditions, ProcessLib::Process::_local_to_global_index_map, ProcessLib::Process::_mesh, _process_data, ProcessLib::HeatTransportBHE::AlgebraicBCSetting::_use_algebraic_bc, ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::_vec_BHE_property, _vec_bottom_BHE_node_indices, _vec_top_BHE_node_indices, ProcessLib::HeatTransportBHE::createBHEBottomDirichletBoundaryCondition(), ProcessLib::HeatTransportBHE::createBHEInflowDirichletBoundaryCondition(), ProcessLib::createBHEInflowPythonBoundaryCondition(), MeshLib::Mesh::getElementsConnectedToNode(), MeshLib::Mesh::getID(), INFO(), MeshLib::Node, OGS_FATAL, and ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::py_bc_object.

Referenced by initializeConcreteProcess().

◆ initializeConcreteProcess()

void ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::initializeConcreteProcess ( NumLib::LocalToGlobalIndexMap const & dof_table,
MeshLib::Mesh const & mesh,
unsigned const integration_order )
overrideprivatevirtual

Process specific initialization called by initialize().

Implements ProcessLib::Process.

Definition at line 134 of file HeatTransportBHEProcess.cpp.

138{
139 // Quick access map to BHE's through element ids.
140 std::unordered_map<std::size_t, BHE::BHETypes*> element_to_bhe_map;
141 int const n_BHEs = _process_data._vec_BHE_property.size();
142 for (int i = 0; i < n_BHEs; i++)
143 {
144 auto const& bhe_elements = _bheMeshData.BHE_elements[i];
145 for (auto const& e : bhe_elements)
146 {
147 element_to_bhe_map[e->getID()] =
149 }
150 }
151
152 assert(mesh.getDimension() == 3);
154 HeatTransportBHELocalAssemblerSoil, HeatTransportBHELocalAssemblerBHE>(
155 mesh.getElements(), dof_table, _local_assemblers,
156 NumLib::IntegrationOrder{integration_order}, element_to_bhe_map,
157 mesh.isAxiallySymmetric(), _process_data);
158
159 // Create BHE boundary conditions for each of the BHEs
161
163 {
164 std::vector<std::size_t> const bhes_node_ids =
165 _bheMeshData.BHE_nodes | ranges::views::join |
166 ranges::views::transform([](auto const* const node)
167 { return node->getID(); }) |
168 ranges::to<std::vector>;
169
170 // all connected soil elements and also the BHE elements.
171 MeshLib::ElementSearch es{mesh};
172 es.searchByNodeIDs(bhes_node_ids);
173
176 mesh.getNumberOfElements(), false);
177 for (auto const id : es.getSearchedElementIDs())
178 {
180 }
181 }
182}
Element search class.
std::size_t searchByNodeIDs(const std::vector< std::size_t > &nodes)
Marks all elements connecting to any of the given nodes.
void createBHEBoundaryConditionTopBottom(std::vector< std::vector< MeshLib::Node * > > const &all_bhe_nodes)
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)

References _bheMeshData, _local_assemblers, ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::_mass_lumping, _process_data, ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::_vec_BHE_property, ProcessLib::HeatTransportBHE::BHEMeshData::BHE_elements, ProcessLib::HeatTransportBHE::BHEMeshData::BHE_nodes, createBHEBoundaryConditionTopBottom(), ProcessLib::HeatTransportBHE::createLocalAssemblers(), MeshLib::Mesh::getDimension(), MeshLib::Mesh::getElements(), MeshLib::Mesh::getNumberOfElements(), MeshLib::Mesh::isAxiallySymmetric(), ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::mass_lumping_soil_elements, and MeshLib::ElementSearch::searchByNodeIDs().

◆ isLinear()

◆ postIterationConcreteProcess()

NumLib::IterationResult ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::postIterationConcreteProcess ( GlobalVector const & x)
overrideprivatevirtual

Reimplemented from ProcessLib::Process.

Definition at line 241 of file HeatTransportBHEProcess.cpp.

243{
244 // if the process use python boundary condition
247
248 // Here the task is to get current time flowrate and flow temperature from
249 // TESPy and determine whether it converges.
250 auto const Tout_nodes_id =
252 const std::size_t n_bc_nodes = Tout_nodes_id.size();
253
254 for (std::size_t i = 0; i < n_bc_nodes; i++)
255 {
256 // read the T_out and store them in dataframe
258 x[Tout_nodes_id[i]];
259 }
260 // Transfer Tin and Tout to TESPy and return the results
261 auto const tespy_result = _process_data.py_bc_object->tespySolver(
263 std::get<1>(_process_data.py_bc_object->dataframe_network), // T_in
264 std::get<2>(_process_data.py_bc_object->dataframe_network)); // T_out
266 {
267 DBUG("Method `tespySolver' not overridden in Python script.");
268 }
269
270 // update the Tin and flow rate
271 for (std::size_t i = 0; i < n_bc_nodes; i++)
272 {
274 std::get<2>(tespy_result)[i];
276 std::get<3>(tespy_result)[i];
277 }
278 auto const tespy_has_converged = std::get<1>(tespy_result);
279 if (tespy_has_converged == true)
281
283}
std::tuple< double, std::vector< double >, std::vector< double >, std::vector< int >, std::vector< double > > dataframe_network
virtual std::tuple< bool, bool, std::vector< double >, std::vector< double > > tespySolver(double, std::vector< double > const &, std::vector< double > const &) const

References _process_data, ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::_use_tespy, ProcessLib::BHEInflowPythonBoundaryConditionPythonSideInterface::dataframe_network, DBUG(), ProcessLib::BHEInflowPythonBoundaryConditionPythonSideInterface::isOverriddenTespy(), ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::py_bc_object, NumLib::REPEAT_ITERATION, NumLib::SUCCESS, and ProcessLib::BHEInflowPythonBoundaryConditionPythonSideInterface::tespySolver().

◆ postTimestepConcreteProcess()

void ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::postTimestepConcreteProcess ( std::vector< GlobalVector * > const & x,
std::vector< GlobalVector * > const & x_prev,
const double t,
const double dt,
int const process_id )
overrideprivatevirtual

Reimplemented from ProcessLib::Process.

Definition at line 334 of file HeatTransportBHEProcess.cpp.

338{
339 if (_process_data.py_bc_object == nullptr ||
341 {
342 return;
343 }
344
345 auto& [time, Tin_value, Tout_value, Tout_nodes_ids, flowrate] =
347
348 // We found the problem that time != t, but it always equals the last
349 // step. The following line is to correct this, although we do not use
350 // it for server communication.
351 time = t;
352
353 auto const& solution = *x[process_id];
354
355 // Iterate through each BHE
356 const std::size_t n_bc_nodes = Tout_nodes_ids.size();
357 for (std::size_t i = 0; i < n_bc_nodes; i++)
358 {
359 // read the T_out and store them in dataframe
360 Tout_value[i] = solution[Tout_nodes_ids[i]];
361 }
362
363 // Transfer T_out to server_Communication
365 t, dt, Tin_value, Tout_value, flowrate);
368 {
369 DBUG("Method `serverCommunication' not overridden in Python script.");
370 }
371}
virtual void serverCommunicationPostTimestep(double, double, std::vector< double > const &, std::vector< double > const &, std::vector< double > const &) const

References _process_data, ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::_use_server_communication, ProcessLib::BHEInflowPythonBoundaryConditionPythonSideInterface::dataframe_network, DBUG(), ProcessLib::BHEInflowPythonBoundaryConditionPythonSideInterface::isOverriddenServerCommunicationPostTimestep(), ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::py_bc_object, and ProcessLib::BHEInflowPythonBoundaryConditionPythonSideInterface::serverCommunicationPostTimestep().

◆ preTimestepConcreteProcess()

void ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::preTimestepConcreteProcess ( std::vector< GlobalVector * > const & x,
const double t,
const double dt,
int const process_id )
overrideprivatevirtual

Reimplemented from ProcessLib::Process.

Definition at line 285 of file HeatTransportBHEProcess.cpp.

288{
289 if (_process_data.py_bc_object == nullptr ||
291 {
292 return;
293 }
294
295 auto& [time, Tin_value, Tout_value, Tout_nodes_ids, flowrate] =
297
298 // We found the problem that time != t, but it always equals the last
299 // step. The following line is to correct this, although we do not use
300 // it for server communication.
301 time = t;
302
303 auto const& solution = *x[process_id];
304
305 // Iterate through each BHE
306 const std::size_t n_bc_nodes = Tout_nodes_ids.size();
307 for (std::size_t i = 0; i < n_bc_nodes; i++)
308 {
309 // read the T_out and store them in dataframe
310 Tout_value[i] = solution[Tout_nodes_ids[i]];
311 }
312
313 // Transfer T_out to server_Communication and get back T_in and flowrate
314 auto const server_communication_result =
316 t, dt, Tin_value, Tout_value, flowrate);
319 {
320 DBUG("Method `serverCommunication' not overridden in Python script.");
321 }
322
323 auto const& [server_communication_Tin_value,
324 server_communication_flowrate] = server_communication_result;
325
326 std::copy(begin(server_communication_Tin_value),
327 end(server_communication_Tin_value),
328 begin(Tin_value));
329 std::copy(begin(server_communication_flowrate),
330 end(server_communication_flowrate),
331 begin(flowrate));
332}
virtual std::tuple< std::vector< double >, std::vector< double > > serverCommunicationPreTimestep(double, double, std::vector< double > const &, std::vector< double > const &, std::vector< double > const &) const

References _process_data, ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::_use_server_communication, ProcessLib::BHEInflowPythonBoundaryConditionPythonSideInterface::dataframe_network, DBUG(), ProcessLib::BHEInflowPythonBoundaryConditionPythonSideInterface::isOverriddenServerCommunicationPreTimestep(), ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::py_bc_object, and ProcessLib::BHEInflowPythonBoundaryConditionPythonSideInterface::serverCommunicationPreTimestep().

◆ requiresNormalization()

bool ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::requiresNormalization ( ) const
inlineoverride

Definition at line 48 of file HeatTransportBHEProcess.h.

49 {
50 // In the current setup, when using algebraic bc,
51 // then normalization is always required
53 }

References ProcessLib::HeatTransportBHE::HeatTransportBHEProcessData::_algebraic_BC_Setting, _process_data, and ProcessLib::HeatTransportBHE::AlgebraicBCSetting::_use_algebraic_bc.

Member Data Documentation

◆ _bheMeshData

const BHEMeshData ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::_bheMeshData
private

◆ _local_assemblers

std::vector<std::unique_ptr<HeatTransportBHELocalAssemblerInterface> > ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::_local_assemblers
private

◆ _mesh_subset_BHE_nodes

std::vector<std::unique_ptr<MeshLib::MeshSubset const> > ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::_mesh_subset_BHE_nodes
private

Definition at line 106 of file HeatTransportBHEProcess.h.

Referenced by constructDofTable().

◆ _mesh_subset_BHE_soil_nodes

std::vector<std::unique_ptr<MeshLib::MeshSubset const> > ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::_mesh_subset_BHE_soil_nodes
private

Definition at line 109 of file HeatTransportBHEProcess.h.

◆ _mesh_subset_soil_nodes

std::unique_ptr<MeshLib::MeshSubset const> ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::_mesh_subset_soil_nodes
private

Definition at line 123 of file HeatTransportBHEProcess.h.

Referenced by constructDofTable().

◆ _process_data

◆ _vec_bottom_BHE_node_indices

std::vector<std::tuple<std::size_t, GlobalIndexType, GlobalIndexType> > ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::_vec_bottom_BHE_node_indices
private

◆ _vec_top_BHE_node_indices

std::vector<std::tuple<std::size_t, GlobalIndexType, GlobalIndexType> > ProcessLib::HeatTransportBHE::HeatTransportBHEProcess::_vec_top_BHE_node_indices
private

The documentation for this class was generated from the following files: