template<typename ShapeFunction, typename IntegrationMethod, typename BHEType>
class ProcessLib::HeatTransportBHE::HeatTransportBHELocalAssemblerBHE< ShapeFunction, IntegrationMethod, BHEType >
Definition at line 28 of file HeatTransportBHELocalAssemblerBHE.h.
|
| HeatTransportBHELocalAssemblerBHE (HeatTransportBHELocalAssemblerBHE const &)=delete |
|
| HeatTransportBHELocalAssemblerBHE (HeatTransportBHELocalAssemblerBHE &&)=delete |
|
| HeatTransportBHELocalAssemblerBHE (MeshLib::Element const &e, BHEType const &bhe, bool const is_axially_symmetric, unsigned const integration_order, HeatTransportBHEProcessData &process_data) |
|
void | assemble (double const, double const, std::vector< double > const &, std::vector< double > const &, std::vector< double > &, std::vector< double > &, std::vector< double > &) override |
|
Eigen::Map< const Eigen::RowVectorXd > | getShapeMatrix (const unsigned integration_point) const override |
| Provides the shape matrix at the given integration point. More...
|
|
virtual | ~LocalAssemblerInterface ()=default |
|
void | setInitialConditions (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, double const t, bool const use_monolithic_scheme, int const process_id) |
|
virtual void | initialize (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table) |
|
virtual void | preAssemble (double const, double const, std::vector< double > const &) |
|
virtual void | assembleForStaggeredScheme (double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_xdot, int const process_id, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data) |
|
virtual void | assembleWithJacobian (double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_xdot, const double dxdot_dx, const double dx_dx, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data) |
|
virtual void | assembleWithJacobianForStaggeredScheme (double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_xdot, const double dxdot_dx, const double dx_dx, int const process_id, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data) |
|
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) |
|
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) |
|
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) |
|
void | postNonLinearSolver (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, GlobalVector const &xdot, double const t, double const dt, bool const use_monolithic_scheme, int const process_id) |
|
virtual std::vector< double > | interpolateNodalValuesToIntegrationPoints (std::vector< double > const &) |
|
virtual Eigen::Vector3d | getFlux (MathLib::Point3d const &, double const, std::vector< double > const &) const |
|
virtual Eigen::Vector3d | getFlux (MathLib::Point3d const &, double const, std::vector< std::vector< double >> const &) const |
| Fits to staggered scheme. More...
|
|
template<typename ShapeFunction , typename IntegrationMethod , typename BHEType >
Definition at line 35 of file HeatTransportBHELocalAssemblerBHE-impl.h.
48 assert(e.getDimension() == 1);
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});
76 Eigen::Map<Eigen::Vector3d const>(e.getNode(0)->getCoords(), 3);
78 Eigen::Map<Eigen::Vector3d const>(e.getNode(1)->getCoords(), 3);
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>(
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
Eigen::Vector3d const _rotation_matrix_1D_to_3D
static constexpr int single_bhe_unknowns_size
std::size_t const _element_id
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data
HeatTransportBHEProcessData & _process_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
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
References ProcessLib::HeatTransportBHE::HeatTransportBHELocalAssemblerBHE< ShapeFunction, IntegrationMethod, BHEType >::_bhe, ProcessLib::HeatTransportBHE::HeatTransportBHELocalAssemblerBHE< ShapeFunction, IntegrationMethod, BHEType >::_element_direction, ProcessLib::HeatTransportBHE::HeatTransportBHELocalAssemblerBHE< ShapeFunction, IntegrationMethod, BHEType >::_integration_method, ProcessLib::HeatTransportBHE::HeatTransportBHELocalAssemblerBHE< ShapeFunction, IntegrationMethod, BHEType >::_ip_data, ProcessLib::HeatTransportBHE::HeatTransportBHELocalAssemblerBHE< ShapeFunction, IntegrationMethod, BHEType >::_R_matrix, ProcessLib::HeatTransportBHE::HeatTransportBHELocalAssemblerBHE< ShapeFunction, IntegrationMethod, BHEType >::_R_pi_s_matrix, ProcessLib::HeatTransportBHE::HeatTransportBHELocalAssemblerBHE< ShapeFunction, IntegrationMethod, BHEType >::_R_s_matrix, ProcessLib::HeatTransportBHE::HeatTransportBHELocalAssemblerBHE< ShapeFunction, IntegrationMethod, BHEType >::_secondary_data, ProcessLib::HeatTransportBHE::HeatTransportBHELocalAssemblerBHE< ShapeFunction, IntegrationMethod, BHEType >::bhe_unknowns, ProcessLib::HeatTransportBHE::HeatTransportBHELocalAssemblerBHE< ShapeFunction, IntegrationMethod, BHEType >::bhe_unknowns_size, MathLib::TemplatePoint< T, DIM >::getCoords(), MeshLib::Element::getDimension(), MeshLib::Element::getNode(), NumLib::initShapeMatrices(), ProcessLib::HeatTransportBHE::SecondaryData< ShapeMatrixType >::N, ProcessLib::HeatTransportBHE::HeatTransportBHELocalAssemblerBHE< ShapeFunction, IntegrationMethod, BHEType >::single_bhe_unknowns_size, and ProcessLib::HeatTransportBHE::HeatTransportBHELocalAssemblerBHE< ShapeFunction, IntegrationMethod, BHEType >::soil_temperature_size.