From a875087aeddb998467fcefa997bae8fc2a59bd75 Mon Sep 17 00:00:00 2001 From: Boyan Meng <meng.boyan@ufz.de> Date: Sat, 12 Oct 2019 14:37:18 +0200 Subject: [PATCH] add the groundwater flow velocity and modify the assembly process --- .../HeatTransportBHELocalAssemblerSoil-impl.h | 29 ++++++++++++++----- .../HeatTransportBHELocalAssemblerSoil.h | 10 +++++-- .../IntegrationPointDataSoil.h | 16 ++++++++-- 3 files changed, 42 insertions(+), 13 deletions(-) diff --git a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h index 875430b76b3..1ae416e3dae 100644 --- a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h +++ b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h @@ -14,6 +14,7 @@ #include <vector> #include "MaterialLib/MPL/Medium.h" +#include "MaterialLib/MPL/Utils/FormEigenTensor.h" #include "MathLib/LinAlg/Eigen/EigenMapTools.h" #include "NumLib/Fem/ShapeMatrixPolicy.h" #include "ProcessLib/HeatTransportBHE/HeatTransportBHEProcessData.h" @@ -61,8 +62,10 @@ HeatTransportBHELocalAssemblerSoil<ShapeFunction, IntegrationMethod>:: auto const& sm = _shape_matrices[ip]; double const w = _integration_method.getWeightedPoint(ip).getWeight() * sm.integralMeasure * sm.detJ; - _ip_data.push_back( - {sm.N.transpose() * sm.N * w, sm.dNdx.transpose() * sm.dNdx * w}); + _ip_data.emplace_back( + sm.N, sm.dNdx, + _integration_method.getWeightedPoint(ip).getWeight() * + sm.integralMeasure * sm.detJ); _secondary_data.N[ip] = sm.N; } @@ -100,8 +103,9 @@ void HeatTransportBHELocalAssemblerSoil<ShapeFunction, IntegrationMethod>:: { pos.setIntegrationPoint(ip); auto& ip_data = _ip_data[ip]; - auto const& M_ip = ip_data.NTN_product_times_w; - auto const& K_ip = ip_data.dNdxTdNdx_product_times_w; + auto const& N = ip_data.N; + auto const& dNdx = ip_data.dNdx; + auto const& w = ip_data.integration_weight; auto const k_s = solid_phase @@ -130,16 +134,25 @@ void HeatTransportBHELocalAssemblerSoil<ShapeFunction, IntegrationMethod>:: .template value<double>(vars, pos, t); // for now only using the solid phase parameters - auto const velocity = + /*auto const velocity = + liquid_phase + .property(MaterialPropertyLib::PropertyType::phase_velocity) + .template value<GlobalDimVectorType>(vars, pos, t);*/ + auto const velocity_arr = liquid_phase .property(MaterialPropertyLib::PropertyType::phase_velocity) .value(vars, pos, t); - + auto velocity = Eigen::Map<Eigen::Matrix<double, 1, 3> const>( + std::get<MaterialPropertyLib::Vector>(velocity_arr).data(), 1, 3); // assemble Conductance matrix - local_K.noalias() += K_ip * k_s; + local_K.noalias() += + (dNdx.transpose() * dNdx * k_s + + N.transpose() * velocity * dNdx * density_f * heat_capacity_f) * + w; // assemble Mass matrix - local_M.noalias() += M_ip * density_s * heat_capacity_s; + local_M.noalias() += + N.transpose() * N * w * density_s * heat_capacity_s; } // debugging diff --git a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil.h b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil.h index 3018cae10a5..55c19392f8e 100644 --- a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil.h +++ b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil.h @@ -36,7 +36,12 @@ public: ShapeMatrixPolicyType<ShapeFunction, 3 /* GlobalDim */>; using NodalMatrixType = typename ShapeMatricesType::NodalMatrixType; using NodalVectorType = typename ShapeMatricesType::NodalVectorType; + using NodalRowVectorType = typename ShapeMatricesType::NodalRowVectorType; + using ShapeMatrices = typename ShapeMatricesType::ShapeMatrices; + using GlobalDimVectorType = typename ShapeMatricesType::GlobalDimVectorType; + using GlobalDimNodalMatrixType = + typename ShapeMatricesType::GlobalDimNodalMatrixType; HeatTransportBHELocalAssemblerSoil( HeatTransportBHELocalAssemblerSoil const&) = delete; @@ -68,8 +73,9 @@ private: HeatTransportBHEProcessData& _process_data; std::vector< - IntegrationPointDataSoil<NodalMatrixType>, - Eigen::aligned_allocator<IntegrationPointDataSoil<NodalMatrixType>>> + IntegrationPointDataSoil<NodalRowVectorType, GlobalDimNodalMatrixType>, + Eigen::aligned_allocator< + IntegrationPointDataSoil<NodalRowVectorType, GlobalDimNodalMatrixType>>> _ip_data; IntegrationMethod const _integration_method; diff --git a/ProcessLib/HeatTransportBHE/LocalAssemblers/IntegrationPointDataSoil.h b/ProcessLib/HeatTransportBHE/LocalAssemblers/IntegrationPointDataSoil.h index 1711addb634..a0b706978d7 100644 --- a/ProcessLib/HeatTransportBHE/LocalAssemblers/IntegrationPointDataSoil.h +++ b/ProcessLib/HeatTransportBHE/LocalAssemblers/IntegrationPointDataSoil.h @@ -14,11 +14,21 @@ namespace ProcessLib { namespace HeatTransportBHE { -template <typename NodalMatrixType> +template <typename NodalRowVectorType, typename GlobalDimNodalMatrixType> struct IntegrationPointDataSoil final { - NodalMatrixType const NTN_product_times_w; - NodalMatrixType const dNdxTdNdx_product_times_w; + IntegrationPointDataSoil(NodalRowVectorType N_, + GlobalDimNodalMatrixType dNdx_, + double const& integration_weight_) + : N(std::move(N_)), + dNdx(std::move(dNdx_)), + integration_weight(integration_weight_) + { + } + + NodalRowVectorType const N; + GlobalDimNodalMatrixType const dNdx; + double const integration_weight; EIGEN_MAKE_ALIGNED_OPERATOR_NEW; }; -- GitLab