20 namespace HeatTransportBHE
24 constexpr
int global_dim = 3;
25 constexpr
int element_dim = 1;
31 .template topLeftCorner<global_dim, element_dim>();
34 template <
typename ShapeFunction,
typename IntegrationMethod,
typename BHEType>
38 bool const is_axially_symmetric,
39 unsigned const integration_order,
41 : _process_data(process_data),
42 _integration_method(integration_order),
44 _element_id(e.getID()),
50 unsigned const n_integration_points =
53 _ip_data.reserve(n_integration_points);
56 auto const shape_matrices =
58 3 >(e, is_axially_symmetric,
62 for (
unsigned ip = 0; ip < n_integration_points; ip++)
64 auto const& sm = shape_matrices[ip];
69 sm.integralMeasure * sm.detJ});
85 static constexpr
int max_num_thermal_exchange_terms = 5;
87 for (
int idx_bhe_unknowns = 0; idx_bhe_unknowns <
bhe_unknowns;
90 typename ShapeMatricesType::template MatrixType<
92 matBHE_loc_R = ShapeMatricesType::template MatrixType<
97 for (
unsigned ip = 0; ip < n_integration_points; ip++)
100 auto const& w =
_ip_data[ip].integration_weight;
102 auto const& R =
_bhe.thermalResistance(idx_bhe_unknowns);
104 matBHE_loc_R += N.transpose() * N * (1 / R) * w;
116 if (idx_bhe_unknowns < max_num_thermal_exchange_terms)
118 _bhe.template assembleRMatrices<ShapeFunction::NPOINTS>(
136 template <
typename ShapeFunction,
typename IntegrationMethod,
typename BHEType>
140 double const ,
double const ,
141 std::vector<double>
const& ,
142 std::vector<double>
const& ,
143 std::vector<double>& local_M_data, std::vector<double>& local_K_data,
144 std::vector<double>& )
146 auto local_M = MathLib::createZeroedMatrix<BheLocalMatrixType>(
147 local_M_data, local_matrix_size, local_matrix_size);
148 auto local_K = MathLib::createZeroedMatrix<BheLocalMatrixType>(
149 local_K_data, local_matrix_size, local_matrix_size);
151 unsigned const n_integration_points =
152 _integration_method.getNumberOfPoints();
154 auto const& pipe_heat_capacities = _bhe.pipeHeatCapacities();
155 auto const& pipe_heat_conductions = _bhe.pipeHeatConductions();
156 auto const& pipe_advection_vectors =
157 _bhe.pipeAdvectionVectors(_element_direction);
158 auto const& cross_section_areas = _bhe.crossSectionAreas();
161 for (
unsigned ip = 0; ip < n_integration_points; ip++)
163 auto& ip_data = _ip_data[ip];
165 auto const& w = ip_data.integration_weight;
166 auto const& N = ip_data.N;
167 auto const& dNdx = ip_data.dNdx;
170 for (
int idx_bhe_unknowns = 0; idx_bhe_unknowns < bhe_unknowns;
174 auto const& mass_coeff = pipe_heat_capacities[idx_bhe_unknowns];
175 auto const& lambda = pipe_heat_conductions[idx_bhe_unknowns];
176 auto const& advection_vector =
177 pipe_advection_vectors[idx_bhe_unknowns];
180 double const advection_coefficient =
181 _rotation_matrix_1D_to_3D.dot(advection_vector);
183 auto const& A = cross_section_areas[idx_bhe_unknowns];
185 int const single_bhe_unknowns_index =
187 single_bhe_unknowns_size * idx_bhe_unknowns;
190 .template block<single_bhe_unknowns_size,
191 single_bhe_unknowns_size>(
192 single_bhe_unknowns_index, single_bhe_unknowns_index)
193 .noalias() += N.transpose() * N * mass_coeff * A * w;
198 .template block<single_bhe_unknowns_size,
199 single_bhe_unknowns_size>(
200 single_bhe_unknowns_index, single_bhe_unknowns_index)
201 .noalias() += dNdx.transpose() * dNdx * lambda * A * w;
203 int constexpr element_dim = 1;
205 .template block<single_bhe_unknowns_size,
206 single_bhe_unknowns_size>(
207 single_bhe_unknowns_index, single_bhe_unknowns_index)
209 advection_coefficient * N.transpose() *
210 dNdx.template topLeftCorner<element_dim,
211 ShapeFunction::NPOINTS>() *
217 local_K.template block<bhe_unknowns_size, bhe_unknowns_size>(
218 bhe_unknowns_index, bhe_unknowns_index) += _R_matrix;
222 .template block<bhe_unknowns_size, soil_temperature_size>(
223 bhe_unknowns_index, soil_temperature_index)
224 .noalias() += _R_pi_s_matrix;
226 .template block<soil_temperature_size, bhe_unknowns_size>(
227 soil_temperature_index, bhe_unknowns_index)
228 .noalias() += _R_pi_s_matrix.transpose();
232 .template block<soil_temperature_size, soil_temperature_size>(
233 soil_temperature_index, soil_temperature_index)
234 .noalias() += _bhe.number_of_grout_zones * _R_s_matrix;
const T * getCoords() const
const RotationMatrix & getRotationMatrixToGlobal() const
return a rotation matrix converting to global coordinates
virtual const Node * getNode(unsigned idx) const =0
virtual constexpr unsigned getDimension() const =0
Get dimension of the mesh element.
Eigen::Vector3d _element_direction
static constexpr int bhe_unknowns
ShapeMatrixPolicyType< ShapeFunction, 3 > ShapeMatricesType
IntegrationMethod _integration_method
ShapeMatricesType::template MatrixType< soil_temperature_size, soil_temperature_size > _R_s_matrix
static constexpr int single_bhe_unknowns_size
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data
ShapeMatricesType::template MatrixType< bhe_unknowns_size, soil_temperature_size > _R_pi_s_matrix
std::vector< IntegrationPointDataBHE< ShapeMatricesType >, Eigen::aligned_allocator< IntegrationPointDataBHE< ShapeMatricesType > > > _ip_data
ShapeMatricesType::template MatrixType< bhe_unknowns_size, bhe_unknowns_size > _R_matrix
static constexpr int soil_temperature_size
HeatTransportBHELocalAssemblerBHE(HeatTransportBHELocalAssemblerBHE const &)=delete
static constexpr int bhe_unknowns_size
std::vector< typename ShapeMatricesType::ShapeMatrices, Eigen::aligned_allocator< typename ShapeMatricesType::ShapeMatrices > > initShapeMatrices(MeshLib::Element const &e, bool const is_axially_symmetric, IntegrationMethod const &integration_method)
Eigen::Vector3d compute1Dto3DRotationMatrix(MeshLib::Element const &e)
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N