diff --git a/Documentation/ProjectFile/prj/processes/process/HEAT_TRANSPORT_BHE/t_linear.md b/Documentation/ProjectFile/prj/processes/process/HEAT_TRANSPORT_BHE/t_linear.md
index c75385a76ff80457576de35b754ec7e101b37a07..0f48474050f1cd6acf89a1e5db41ecc3acdf488c 100644
--- a/Documentation/ProjectFile/prj/processes/process/HEAT_TRANSPORT_BHE/t_linear.md
+++ b/Documentation/ProjectFile/prj/processes/process/HEAT_TRANSPORT_BHE/t_linear.md
@@ -1 +1,5 @@
-The default setting is false. When activated, non-linear iterations will be turned off. The simulation will just go through a single iteration within a time step. If it's set to true, the **\<use_algebraic_bc\>** option must be also set to true. Only when fluid properties does not change much with temperature, this feature can be activated. The effect of skipping the non-linear iteration is much faster simulation speed.
+The default setting is false. When activated, non-linear iterations will be turned off and the equation system for the soil part is only assembled once.
+If the `<use_algebraic_bc>` option is also set to true, the simulation will just go through a single iteration within a time step.
+Otherwise iterations are needed to fulfill the boundary conditions of the BHE.
+Only when fluid properties do not change much with temperature, this feature can be activated.
+The effect of skipping the non-linear iteration is much faster simulation speed.
diff --git a/ProcessLib/HeatTransportBHE/CreateHeatTransportBHEProcess.cpp b/ProcessLib/HeatTransportBHE/CreateHeatTransportBHEProcess.cpp
index bd724ac17b92656bba0ac3a22cbca16d5edbf1b3..eb2f140287495b8e337b46b3baff2b5c47bac2b0 100644
--- a/ProcessLib/HeatTransportBHE/CreateHeatTransportBHEProcess.cpp
+++ b/ProcessLib/HeatTransportBHE/CreateHeatTransportBHEProcess.cpp
@@ -195,10 +195,15 @@ std::unique_ptr<Process> createHeatTransportBHEProcess(
{
if (!using_algebraic_bc)
{
- OGS_FATAL(
+ WARN(
"You specified that the process simulated by OGS is linear. "
- "For the Heat-Transport-BHE process this can only be done "
- "together with setting the use_algebraic_bc option to true.")
+ "With that optimization the process will be assembled only "
+ "once and the non-linear solver will only do iterations per "
+ "time step to fulfill the BHE boundary conditions. No other "
+ "non-linearities will be resolved and OGS will not detect if "
+ "there are any non-linearities. It is your responsibility to "
+ "ensure that the assembled equation systems are linear, "
+ "indeed! There is no safety net!");
}
else
{
@@ -336,7 +341,7 @@ std::unique_ptr<Process> createHeatTransportBHEProcess(
HeatTransportBHEProcessData process_data(
std::move(media_map), std::move(bhes), py_object, using_tespy,
using_server_communication, mass_lumping,
- {using_algebraic_bc, weighting_factor, is_linear});
+ {using_algebraic_bc, weighting_factor}, is_linear);
SecondaryVariableCollection secondary_variables;
diff --git a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp
index 8528005f6f38feb56bb55c01295c9619fee6265e..576aad5db10d7ab6d73bd7105d3edff9ea93bc61 100644
--- a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp
+++ b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp
@@ -14,6 +14,7 @@
#include "BoundaryConditions/BHEBottomDirichletBoundaryCondition.h"
#include "BoundaryConditions/BHEInflowDirichletBoundaryCondition.h"
+#include "MathLib/LinAlg/SetMatrixSparsity.h"
#include "MeshLib/MeshSearch/ElementSearch.h"
#include "ProcessLib/BoundaryConditionAndSourceTerm/Python/BHEInflowPythonBoundaryCondition.h"
#include "ProcessLib/HeatTransportBHE/BHE/MeshUtils.h"
@@ -40,7 +41,8 @@ HeatTransportBHEProcess::HeatTransportBHEProcess(
integration_order, std::move(process_variables),
std::move(secondary_variables)),
_process_data(std::move(process_data)),
- _bheMeshData(std::move(bhe_mesh_data))
+ _bheMeshData(std::move(bhe_mesh_data)),
+ _asm_mat_cache{_process_data._is_linear, true /*use_monolithic_scheme*/}
{
if (_bheMeshData.BHE_mat_IDs.size() !=
_process_data._vec_BHE_property.size())
@@ -159,6 +161,33 @@ void HeatTransportBHEProcess::initializeConcreteProcess(
// Create BHE boundary conditions for each of the BHEs
createBHEBoundaryConditionTopBottom(_bheMeshData.BHE_nodes);
+ // Store BHE and soil elements to split the assembly and use the matrix
+ // cache in the linear case only for soil elements
+ if (_process_data._is_linear)
+ {
+ _bhes_element_ids = _bheMeshData.BHE_elements | ranges::views::join |
+ MeshLib::views::ids | ranges::to<std::vector>;
+
+ // sort bhe elements if needed
+ if (!std::is_sorted(_bhes_element_ids.begin(), _bhes_element_ids.end()))
+ {
+ std::sort(_bhes_element_ids.begin(), _bhes_element_ids.end());
+ }
+
+ _soil_element_ids = mesh.getElements() | MeshLib::views::ids |
+ ranges::to<std::vector>();
+
+ // sort soil elements if needed
+ if (!std::is_sorted(_soil_element_ids.begin(), _soil_element_ids.end()))
+ {
+ std::sort(_soil_element_ids.begin(), _soil_element_ids.end());
+ }
+
+ _soil_element_ids = ranges::views::set_difference(_soil_element_ids,
+ _bhes_element_ids) |
+ ranges::to<std::vector>();
+ }
+
if (_process_data._mass_lumping)
{
std::vector<std::size_t> const bhes_node_ids =
@@ -190,11 +219,42 @@ void HeatTransportBHEProcess::assembleConcreteProcess(
std::vector<NumLib::LocalToGlobalIndexMap const*> dof_table = {
_local_to_global_index_map.get()};
- // Call global assembler for each local assembly item.
- GlobalExecutor::executeSelectedMemberDereferenced(
- _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers,
- getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K,
- &b);
+
+ if (_process_data._is_linear)
+ {
+ if (!getActiveElementIDs().empty())
+ {
+ OGS_FATAL(
+ "Domain Deactivation is currently not implemnted with "
+ "linear optimization.");
+ }
+
+ auto const& spec = getMatrixSpecifications(process_id);
+
+ // use matrix cache for soil elements
+ _asm_mat_cache.assemble(t, dt, x, x_prev, process_id, &M, &K, &b,
+ dof_table, _global_assembler, _local_assemblers,
+ _soil_element_ids);
+
+ // reset the sparsity pattern for better performance in the BHE assembly
+ MathLib::setMatrixSparsity(M, *spec.sparsity_pattern);
+ MathLib::setMatrixSparsity(K, *spec.sparsity_pattern);
+
+ // Call global assembler for each local BHE assembly item.
+ GlobalExecutor::executeSelectedMemberDereferenced(
+ _global_assembler, &VectorMatrixAssembler::assemble,
+ _local_assemblers, _bhes_element_ids, dof_table, t, dt, x, x_prev,
+ process_id, &M, &K, &b);
+ }
+ else
+ {
+ // Call global assembler for each local assembly item.
+ GlobalExecutor::executeSelectedMemberDereferenced(
+ _global_assembler, &VectorMatrixAssembler::assemble,
+ _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x,
+ x_prev, process_id, &M, &K, &b);
+ }
+
// Algebraic BC procedure.
if (_process_data._algebraic_BC_Setting._use_algebraic_bc)
{
diff --git a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h
index 5ee4be7f1a3fe312366faf6f38498e8d1c553c57..e3090e51e72d9f3eae9037f8cd60d4712d177fa1 100644
--- a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h
+++ b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h
@@ -11,6 +11,7 @@
#pragma once
#include "HeatTransportBHEProcessData.h"
+#include "ProcessLib/Assembly/AssembledMatrixCache.h"
#include "ProcessLib/HeatTransportBHE/BHE/MeshUtils.h"
#include "ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHEProcessAssemblerInterface.h"
#include "ProcessLib/Process.h"
@@ -42,7 +43,10 @@ public:
//! @{
bool isLinear() const override
{
- return _process_data._algebraic_BC_Setting._is_linear;
+ // Linear Solver is called only once - valid only with
+ // _use_algebraic_bc and _is_linear
+ return _process_data._algebraic_BC_Setting._use_algebraic_bc &&
+ _process_data._is_linear;
}
bool requiresNormalization() const override
@@ -123,6 +127,12 @@ private:
std::unique_ptr<MeshLib::MeshSubset const> _mesh_subset_soil_nodes;
const BHEMeshData _bheMeshData;
+
+ AssembledMatrixCache _asm_mat_cache;
+
+ std::vector<std::size_t> _bhes_element_ids;
+
+ std::vector<std::size_t> _soil_element_ids;
};
} // namespace HeatTransportBHE
} // namespace ProcessLib
diff --git a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcessData.h b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcessData.h
index bb612187da7e71cab7b1106e324dce3223258d3b..f733bd945fabd6b3046644293842728d85c37c86 100644
--- a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcessData.h
+++ b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcessData.h
@@ -28,8 +28,6 @@ struct AlgebraicBCSetting
const bool _use_algebraic_bc;
const double _weighting_factor;
-
- const bool _is_linear;
};
struct HeatTransportBHEProcessData final
@@ -42,14 +40,16 @@ struct HeatTransportBHEProcessData final
const bool use_tespy = false,
const bool use_server_communication = false,
const bool mass_lumping = false,
- AlgebraicBCSetting algebraicBCSetting = {false, 100.0, false})
+ AlgebraicBCSetting algebraicBCSetting = {false, 100.0},
+ const bool is_linear = false)
: media_map(media_map_),
_vec_BHE_property(std::move(vec_BHEs_)),
py_bc_object(py_bc_object_),
_use_tespy(use_tespy),
_use_server_communication(use_server_communication),
_mass_lumping(mass_lumping),
- _algebraic_BC_Setting(algebraicBCSetting)
+ _algebraic_BC_Setting(algebraicBCSetting),
+ _is_linear(is_linear)
{
}
MaterialPropertyLib::MaterialSpatialDistributionMap media_map;
@@ -71,5 +71,7 @@ struct HeatTransportBHEProcessData final
std::vector<bool> mass_lumping_soil_elements;
AlgebraicBCSetting const _algebraic_BC_Setting;
+
+ const bool _is_linear;
};
} // namespace ProcessLib::HeatTransportBHE
diff --git a/ProcessLib/HeatTransportBHE/Tests.cmake b/ProcessLib/HeatTransportBHE/Tests.cmake
index 028d97295a5d76783ba2507fbd74daed07696b38..e14b64133f27cc99c59f4408ddfa7d2051a0c63d 100644
--- a/ProcessLib/HeatTransportBHE/Tests.cmake
+++ b/ProcessLib/HeatTransportBHE/Tests.cmake
@@ -12,6 +12,20 @@ AddTest(
sandwich_ts_10_t_600.000000.vtu sandwich_ts_10_t_600.000000.vtu temperature_soil temperature_soil 0 1e-13
)
+AddTest(
+ NAME HeatTransportBHE_1U_3D_bhe_sandwich_linear
+ PATH Parabolic/T/3D_BHE_Sandwich
+ EXECUTABLE ogs
+ EXECUTABLE_ARGS sandwich_linear.xml
+ WRAPPER time
+ TESTER vtkdiff
+ REQUIREMENTS NOT OGS_USE_MPI
+ RUNTIME 10
+ DIFF_DATA
+ sandwich_ts_10_t_600.000000.vtu sandwich_linear_ts_10_t_600.000000.vtu temperature_BHE1 temperature_BHE1 0 5e-15
+ sandwich_ts_10_t_600.000000.vtu sandwich_linear_ts_10_t_600.000000.vtu temperature_soil temperature_soil 0 1e-13
+)
+
AddTest(
NAME HeatTransportBHE_1U_3D_bhe_sandwich_Newton
PATH Parabolic/T/3D_BHE_Sandwich
@@ -110,6 +124,20 @@ AddTest(
beier_sandbox_ts_10_t_600.000000.vtu beier_sandbox_ts_10_t_600.000000.vtu temperature_soil temperature_soil 0 1e-13
)
+AddTest(
+ NAME HeatTransportBHE_1U_3D_beier_sandbox_linear
+ PATH Parabolic/T/3D_Beier_sandbox
+ EXECUTABLE ogs
+ EXECUTABLE_ARGS beier_sandbox_linear.xml
+ WRAPPER time
+ TESTER vtkdiff
+ REQUIREMENTS NOT OGS_USE_MPI
+ RUNTIME 20
+ DIFF_DATA
+ beier_sandbox_ts_10_t_600.000000.vtu beier_sandbox_linear_ts_10_t_600.000000.vtu temperature_BHE1 temperature_BHE1 0 5e-15
+ beier_sandbox_ts_10_t_600.000000.vtu beier_sandbox_linear_ts_10_t_600.000000.vtu temperature_soil temperature_soil 0 1e-13
+)
+
AddTest(
NAME HeatTransportBHE_1U_3D_beier_sandbox_Newton
PATH Parabolic/T/3D_Beier_sandbox
diff --git a/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_linear.xml b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_linear.xml
new file mode 100644
index 0000000000000000000000000000000000000000..ea5dbd1dd6934a3398361fcd83303a0f8bc6040e
--- /dev/null
+++ b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_linear.xml
@@ -0,0 +1,7 @@
+<?xml version="1.0" encoding="ISO-8859-1"?>
+<OpenGeoSysProjectDiff base_file="sandwich.prj">
+ <add sel="/*/processes/process">
+ <linear>true</linear>
+ </add>
+ <replace sel="/*/time_loop/output/prefix/text()">sandwich_linear</replace>
+</OpenGeoSysProjectDiff>
diff --git a/Tests/Data/Parabolic/T/3D_Beier_sandbox/beier_sandbox_linear.xml b/Tests/Data/Parabolic/T/3D_Beier_sandbox/beier_sandbox_linear.xml
new file mode 100644
index 0000000000000000000000000000000000000000..6f2c741bb14d60019b584afd4d589a0d98910b29
--- /dev/null
+++ b/Tests/Data/Parabolic/T/3D_Beier_sandbox/beier_sandbox_linear.xml
@@ -0,0 +1,7 @@
+<?xml version="1.0" encoding="ISO-8859-1"?>
+<OpenGeoSysProjectDiff base_file="beier_sandbox.prj">
+ <add sel="/*/processes/process">
+ <linear>true</linear>
+ </add>
+ <replace sel="/*/time_loop/output/prefix/text()">beier_sandbox_linear</replace>
+</OpenGeoSysProjectDiff>