From 54358a65ff84e7f78241fc94c385243be279fe11 Mon Sep 17 00:00:00 2001 From: Tom Fischer Date: Fri, 22 Mar 2024 13:32:04 +0100 Subject: [PATCH 1/3] [PL/LF] Add shape matrix cache to LFProcess --- ProcessLib/LiquidFlow/LiquidFlowProcess.cpp | 1 + ProcessLib/LiquidFlow/LiquidFlowProcess.h | 5 ++++- 2 files changed, 5 insertions(+), 1 deletion(-) diff --git a/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp b/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp index 23f7891cc4d..a7bf2ba9415 100644 --- a/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp +++ b/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp @@ -41,6 +41,7 @@ LiquidFlowProcess::LiquidFlowProcess( integration_order, std::move(process_variables), std::move(secondary_variables)), _process_data(std::move(process_data)), + _shape_matrix_cache{integration_order}, _surfaceflux(std::move(surfaceflux)), _is_linear(is_linear) { diff --git a/ProcessLib/LiquidFlow/LiquidFlowProcess.h b/ProcessLib/LiquidFlow/LiquidFlowProcess.h index 86b84e57e42..aac257a1e6f 100644 --- a/ProcessLib/LiquidFlow/LiquidFlowProcess.h +++ b/ProcessLib/LiquidFlow/LiquidFlowProcess.h @@ -14,9 +14,10 @@ #include -#include "LiquidFlowLocalAssembler.h" #include "LiquidFlowData.h" +#include "LiquidFlowLocalAssembler.h" #include "NumLib/DOF/LocalToGlobalIndexMap.h" +#include "NumLib/Fem/ShapeMatrixCache.h" #include "ProcessLib/Process.h" #include "ProcessLib/SurfaceFlux/SurfaceFluxData.h" @@ -108,6 +109,8 @@ class LiquidFlowProcess final : public Process std::vector> _local_assemblers; + NumLib::ShapeMatrixCache _shape_matrix_cache; + std::unique_ptr _surfaceflux; MeshLib::PropertyVector* _hydraulic_flow = nullptr; bool _is_linear = false; From 5d344ea3b990bb014e4c28f437555c30ebe0c9a2 Mon Sep 17 00:00:00 2001 From: Tom Fischer Date: Fri, 22 Mar 2024 14:54:49 +0100 Subject: [PATCH 2/3] [PL/LF] Pass shape matrix cache from LFProcess to LocalAssembler --- ProcessLib/LiquidFlow/LiquidFlowLocalAssembler.h | 8 ++++++-- ProcessLib/LiquidFlow/LiquidFlowProcess.cpp | 2 +- 2 files changed, 7 insertions(+), 3 deletions(-) diff --git a/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler.h b/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler.h index 71771a1b121..87650cd538b 100644 --- a/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler.h +++ b/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler.h @@ -24,6 +24,7 @@ #include "NumLib/Fem/FiniteElement/TemplateIsoparametric.h" #include "NumLib/Fem/InitShapeMatrices.h" #include "NumLib/Fem/Integration/GenericIntegrationMethod.h" +#include "NumLib/Fem/ShapeMatrixCache.h" #include "NumLib/Fem/ShapeMatrixPolicy.h" #include "ParameterLib/Parameter.h" #include "ProcessLib/LocalAssemblerInterface.h" @@ -87,10 +88,12 @@ class LiquidFlowLocalAssembler : public LiquidFlowLocalAssemblerInterface std::size_t const /*local_matrix_size*/, NumLib::GenericIntegrationMethod const& integration_method, bool const is_axially_symmetric, - LiquidFlowData const& process_data) + LiquidFlowData const& process_data, + NumLib::ShapeMatrixCache const& shape_matrix_cache) : _element(element), _integration_method(integration_method), - _process_data(process_data) + _process_data(process_data), + _shape_matrix_cache(shape_matrix_cache) { unsigned const n_integration_points = _integration_method.getNumberOfPoints(); @@ -227,6 +230,7 @@ class LiquidFlowLocalAssembler : public LiquidFlowLocalAssemblerInterface VelocityCacheType& darcy_velocity_at_ips) const; const LiquidFlowData& _process_data; + NumLib::ShapeMatrixCache const& _shape_matrix_cache; }; } // namespace LiquidFlow diff --git a/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp b/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp index a7bf2ba9415..3a5d2af3e8b 100644 --- a/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp +++ b/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp @@ -60,7 +60,7 @@ void LiquidFlowProcess::initializeConcreteProcess( ProcessLib::createLocalAssemblers( mesh_space_dimension, mesh.getElements(), dof_table, _local_assemblers, NumLib::IntegrationOrder{integration_order}, mesh.isAxiallySymmetric(), - _process_data); + _process_data, _shape_matrix_cache); _secondary_variables.addSecondaryVariable( "darcy_velocity", From 2b5594e1fe68fa2990e75dcbb792c19dc1aa639e Mon Sep 17 00:00:00 2001 From: Tom Fischer Date: Fri, 22 Mar 2024 15:16:58 +0100 Subject: [PATCH 3/3] [PL/LF] Use shape function from matrix cache --- .../LiquidFlow/LiquidFlowLocalAssembler-impl.h | 12 +++++++++--- ProcessLib/LiquidFlow/LiquidFlowLocalAssembler.h | 15 ++++++++------- 2 files changed, 17 insertions(+), 10 deletions(-) diff --git a/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler-impl.h b/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler-impl.h index 16d2116f6f4..4b3e41fffcc 100644 --- a/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler-impl.h +++ b/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler-impl.h @@ -143,12 +143,15 @@ void LiquidFlowLocalAssembler:: _process_data.element_rotation_matrices[_element.getID()].transpose() * _process_data.specific_body_force; + auto const& N = _shape_matrix_cache + .NsHigherOrder(); + for (unsigned ip = 0; ip < n_integration_points; ip++) { auto const& ip_data = _ip_data[ip]; double p = 0.; - NumLib::shapeFunctionInterpolate(local_x, ip_data.N, p); + NumLib::shapeFunctionInterpolate(local_x, N[ip], p); vars.liquid_phase_pressure = p; // Compute density: @@ -173,7 +176,7 @@ void LiquidFlowLocalAssembler:: // Assemble mass matrix, M local_M.noalias() += (porosity * ddensity_dpressure / fluid_density + storage) * - ip_data.N.transpose() * ip_data.N * ip_data.integration_weight; + N[ip].transpose() * N[ip] * ip_data.integration_weight; // Compute viscosity: auto const viscosity = @@ -292,11 +295,14 @@ void LiquidFlowLocalAssembler:: _process_data.element_rotation_matrices[_element.getID()].transpose() * _process_data.specific_body_force; + auto const& N = _shape_matrix_cache + .NsHigherOrder(); + for (unsigned ip = 0; ip < n_integration_points; ip++) { auto const& ip_data = _ip_data[ip]; double p = 0.; - NumLib::shapeFunctionInterpolate(local_x, ip_data.N, p); + NumLib::shapeFunctionInterpolate(local_x, N[ip], p); vars.liquid_phase_pressure = p; // Compute density: diff --git a/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler.h b/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler.h index 87650cd538b..86ae1838af3 100644 --- a/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler.h +++ b/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler.h @@ -37,14 +37,12 @@ namespace LiquidFlow template struct IntegrationPointData final { - explicit IntegrationPointData(NodalRowVectorType const& N_, - GlobalDimNodalMatrixType const& dNdx_, + explicit IntegrationPointData(GlobalDimNodalMatrixType const& dNdx_, double const& integration_weight_) - : N(N_), dNdx(dNdx_), integration_weight(integration_weight_) + : dNdx(dNdx_), integration_weight(integration_weight_) { } - NodalRowVectorType const N; GlobalDimNodalMatrixType const dNdx; double const integration_weight; @@ -115,7 +113,7 @@ class LiquidFlowLocalAssembler : public LiquidFlowLocalAssemblerInterface for (unsigned ip = 0; ip < n_integration_points; ip++) { _ip_data.emplace_back( - shape_matrices[ip].N, shape_matrices[ip].dNdx, + shape_matrices[ip].dNdx, _integration_method.getWeightedPoint(ip).getWeight() * shape_matrices[ip].integralMeasure * shape_matrices[ip].detJ * aperture_size); @@ -138,10 +136,13 @@ class LiquidFlowLocalAssembler : public LiquidFlowLocalAssemblerInterface Eigen::Map getShapeMatrix( const unsigned integration_point) const override { - auto const& N = _ip_data[integration_point].N; + auto const& N = + _shape_matrix_cache + .NsHigherOrder(); // assumes N is stored contiguously in memory - return Eigen::Map(N.data(), N.size()); + return Eigen::Map( + N[integration_point].data(), N[integration_point].size()); } std::vector const& getIntPtDarcyVelocity(