Skip to content
Snippets Groups Projects
Commit dc212739 authored by wenqing's avatar wenqing
Browse files

[BHE] Rotate the velocity instead of rotating dNdx in the advection term

parent f73976db
No related branches found
No related tags found
No related merge requests found
...@@ -12,12 +12,23 @@ ...@@ -12,12 +12,23 @@
#include "HeatTransportBHELocalAssemblerBHE.h" #include "HeatTransportBHELocalAssemblerBHE.h"
#include "MathLib/LinAlg/Eigen/EigenMapTools.h" #include "MathLib/LinAlg/Eigen/EigenMapTools.h"
#include "MeshLib/ElementCoordinatesMappingLocal.h"
#include "NumLib/Fem/InitShapeMatrices.h" #include "NumLib/Fem/InitShapeMatrices.h"
namespace ProcessLib namespace ProcessLib
{ {
namespace HeatTransportBHE namespace HeatTransportBHE
{ {
Eigen::Vector3d compute1Dto3DRotationMatrix(MeshLib::Element const& e)
{
assert(e.getDimension() == 1);
const int global_dim = 3;
const MeshLib::ElementCoordinatesMappingLocal ele_local_coord(e,
global_dim);
return ele_local_coord.getRotationMatrixToGlobal().topLeftCorner(
global_dim, e.getDimension());
}
template <typename ShapeFunction, typename IntegrationMethod, typename BHEType> template <typename ShapeFunction, typename IntegrationMethod, typename BHEType>
HeatTransportBHELocalAssemblerBHE<ShapeFunction, IntegrationMethod, BHEType>:: HeatTransportBHELocalAssemblerBHE<ShapeFunction, IntegrationMethod, BHEType>::
HeatTransportBHELocalAssemblerBHE(MeshLib::Element const& e, HeatTransportBHELocalAssemblerBHE(MeshLib::Element const& e,
...@@ -28,7 +39,8 @@ HeatTransportBHELocalAssemblerBHE<ShapeFunction, IntegrationMethod, BHEType>:: ...@@ -28,7 +39,8 @@ HeatTransportBHELocalAssemblerBHE<ShapeFunction, IntegrationMethod, BHEType>::
: _process_data(process_data), : _process_data(process_data),
_integration_method(integration_order), _integration_method(integration_order),
_bhe(bhe), _bhe(bhe),
_element_id(e.getID()) _element_id(e.getID()),
_rotation_matrix_1D_to_3D(compute1Dto3DRotationMatrix(e))
{ {
// need to make sure that the BHE elements are one-dimensional // need to make sure that the BHE elements are one-dimensional
assert(e.getDimension() == 1); assert(e.getDimension() == 1);
...@@ -161,6 +173,11 @@ void HeatTransportBHELocalAssemblerBHE<ShapeFunction, IntegrationMethod, ...@@ -161,6 +173,11 @@ void HeatTransportBHELocalAssemblerBHE<ShapeFunction, IntegrationMethod,
auto const& lambda = pipe_heat_conductions[idx_bhe_unknowns]; auto const& lambda = pipe_heat_conductions[idx_bhe_unknowns];
auto const& advection_vector = auto const& advection_vector =
pipe_advection_vectors[idx_bhe_unknowns]; pipe_advection_vectors[idx_bhe_unknowns];
// R^T * v, 3D to 1D:
double const advection_coefficient =
_rotation_matrix_1D_to_3D.dot(advection_vector);
auto const& A = cross_section_areas[idx_bhe_unknowns]; auto const& A = cross_section_areas[idx_bhe_unknowns];
int const single_bhe_unknowns_index = int const single_bhe_unknowns_index =
...@@ -181,12 +198,13 @@ void HeatTransportBHELocalAssemblerBHE<ShapeFunction, IntegrationMethod, ...@@ -181,12 +198,13 @@ void HeatTransportBHELocalAssemblerBHE<ShapeFunction, IntegrationMethod,
single_bhe_unknowns_index, single_bhe_unknowns_index) single_bhe_unknowns_index, single_bhe_unknowns_index)
.noalias() += dNdx.transpose() * dNdx * lambda * A * w; .noalias() += dNdx.transpose() * dNdx * lambda * A * w;
// advection part // advection part
int constexpr element_dim = 1;
local_K local_K
.template block<single_bhe_unknowns_size, .template block<single_bhe_unknowns_size,
single_bhe_unknowns_size>( single_bhe_unknowns_size>(
single_bhe_unknowns_index, single_bhe_unknowns_index) single_bhe_unknowns_index, single_bhe_unknowns_index)
.noalias() += .noalias() += advection_coefficient * N.transpose() *
N.transpose() * advection_vector.transpose() * dNdx * A * w; dNdx.topLeftCorner(element_dim, N.size()) * A * w;
} }
} }
......
...@@ -92,6 +92,7 @@ private: ...@@ -92,6 +92,7 @@ private:
SecondaryData<typename ShapeMatrices::ShapeType> _secondary_data; SecondaryData<typename ShapeMatrices::ShapeType> _secondary_data;
Eigen::Vector3d _element_direction; Eigen::Vector3d _element_direction;
Eigen::Vector3d const _rotation_matrix_1D_to_3D;
typename ShapeMatricesType::template MatrixType<bhe_unknowns_size, typename ShapeMatricesType::template MatrixType<bhe_unknowns_size,
bhe_unknowns_size> bhe_unknowns_size>
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment