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
new file mode 100644
index 0000000000000000000000000000000000000000..c75385a76ff80457576de35b754ec7e101b37a07
--- /dev/null
+++ b/Documentation/ProjectFile/prj/processes/process/HEAT_TRANSPORT_BHE/t_linear.md
@@ -0,0 +1 @@
+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.
diff --git a/Documentation/ProjectFile/prj/processes/process/HEAT_TRANSPORT_BHE/t_use_algebraic_bc.md b/Documentation/ProjectFile/prj/processes/process/HEAT_TRANSPORT_BHE/t_use_algebraic_bc.md
new file mode 100644
index 0000000000000000000000000000000000000000..b735f6a1fa8a5f198a743639514fd1bf9fa675e9
--- /dev/null
+++ b/Documentation/ProjectFile/prj/processes/process/HEAT_TRANSPORT_BHE/t_use_algebraic_bc.md
@@ -0,0 +1 @@
+The default setting is false. When algebraic boundary condition is activated, the inflow temperatures at the top and bottom of the BHE are no longer imposed as a Dirichlet type boundary condition. Instead, the linear equation system will be enlarged to accommodate the numerical constrains. The effect of this feature is to accelerate BHE simulations.
diff --git a/Documentation/ProjectFile/prj/processes/process/HEAT_TRANSPORT_BHE/t_weighting_factor.md b/Documentation/ProjectFile/prj/processes/process/HEAT_TRANSPORT_BHE/t_weighting_factor.md
new file mode 100644
index 0000000000000000000000000000000000000000..b557a842c7f1a4426ee3fab4a0614cfdfd2d11dd
--- /dev/null
+++ b/Documentation/ProjectFile/prj/processes/process/HEAT_TRANSPORT_BHE/t_weighting_factor.md
@@ -0,0 +1 @@
+This value is the weighting factor used to impose an algebraic boundary condition in the **HeatTransportBHE** process. The default value is 100. This value is only effective, when the option **\<use_algebraic_bc\>** is set to true. This value is multiplied with the max value of the column wise inner product and then imposed in the additional rows at the end of the linear equation system, so that the unknown temperatures at the top and bottom of the BHE are imposed with certain desired values.
diff --git a/MathLib/LinAlg/LinAlg.cpp b/MathLib/LinAlg/LinAlg.cpp
index 5e688a3fd1ee252324d11c4d3c765d5a77852d80..8be69e838047222034cd61b5f856fcffc3027d5e 100644
--- a/MathLib/LinAlg/LinAlg.cpp
+++ b/MathLib/LinAlg/LinAlg.cpp
@@ -167,6 +167,34 @@ void matMultAdd(PETScMatrix const& A, PETScVector const& v1,
v3.getRawVector());
}
+void linearSysNormalize(PETScMatrix const& /*A*/, PETScMatrix& /*new_A*/,
+ PETScVector const& /*b*/, PETScVector& /*new_b*/)
+{
+ // The following block is deactivated, because there is no tests yet for the
+ // normalization operation in PETSc. This will be a task for later.
+ /*
+ assert(&A != &new_A);
+ assert(&b != &new_b);
+
+ PetscInt n_rows(0);
+ PetscInt n_cols(0);
+ MatGetSize(A.getRawMatrix(), &n_rows, &n_cols);
+ // only when A matrix is not square
+ if (n_rows != n_cols)
+ {
+ // new_b = A^T * b
+ MatMultTranspose(A.getRawMatrix(), b.getRawVector(),
+ new_b.getRawVector());
+ // new_A = A^T * A
+ MatTranspose(A.getRawMatrix(), MAT_INITIAL_MATRIX,
+ &(new_A.getRawMatrix()));
+ }
+ */
+ OGS_FATAL(
+ "Normalization operation is not implemented yet for PETSc library! "
+ "Program terminated.");
+}
+
void finalizeAssembly(PETScMatrix& A)
{
A.finalizeAssembly(MAT_FINAL_ASSEMBLY);
@@ -313,6 +341,25 @@ void matMultAdd(EigenMatrix const& A, EigenVector const& v1,
v2.getRawVector() + A.getRawMatrix() * v1.getRawVector();
}
+void linearSysNormalize(EigenMatrix const& A, EigenMatrix& new_A,
+ EigenVector const& b, EigenVector& new_b)
+{
+ // make sure that new_A and new_b are not the same memory
+ assert(&A != &new_A);
+ assert(&b != &new_b);
+
+ if (A.getRawMatrix().rows() == A.getRawMatrix().cols())
+ {
+ WARN(
+ "The number of rows and columns are the same for the LHS matrix."
+ "Are you sure you still need to normalize the LHS matrix and RHS "
+ "vector? ");
+ }
+
+ new_b.getRawVector() = A.getRawMatrix().transpose() * b.getRawVector();
+ new_A.getRawMatrix() = A.getRawMatrix().transpose() * A.getRawMatrix();
+}
+
void finalizeAssembly(EigenMatrix& x)
{
x.getRawMatrix().makeCompressed();
diff --git a/MathLib/LinAlg/LinAlg.h b/MathLib/LinAlg/LinAlg.h
index f982306db82122b41c3224cb2348b13257dd8750..206395e8f907e2e72d3d30b5388030312cfdd7e9 100644
--- a/MathLib/LinAlg/LinAlg.h
+++ b/MathLib/LinAlg/LinAlg.h
@@ -196,6 +196,11 @@ void matMult(PETScMatrix const& A, PETScVector const& x, PETScVector& y);
void matMultAdd(PETScMatrix const& A, PETScVector const& v1,
PETScVector const& v2, PETScVector& v3);
+// new_A = A^T * A
+// new_b = A^T * b
+void linearSysNormalize(PETScMatrix const& A, PETScMatrix& new_A,
+ PETScVector const& b, PETScVector& new_b);
+
void finalizeAssembly(PETScMatrix& A);
void finalizeAssembly(PETScVector& x);
@@ -264,7 +269,10 @@ void matMult(EigenMatrix const& A, EigenVector const& x, EigenVector& y);
// v3 = A*v1 + v2
void matMultAdd(EigenMatrix const& A, EigenVector const& v1,
EigenVector const& v2, EigenVector& v3);
-
+// new_A = A^T * A
+// new_b = A^T * b
+void linearSysNormalize(EigenMatrix const& A, EigenMatrix& new_A,
+ EigenVector const& b, EigenVector& new_b);
void finalizeAssembly(EigenMatrix& x);
void finalizeAssembly(EigenVector& A);
diff --git a/NumLib/ODESolver/EquationSystem.h b/NumLib/ODESolver/EquationSystem.h
index 643d1bdb57a501c7fccd75818ec9f35b1864f57e..c851ed8f92cdc09796c63498dd8142ff2251f04f 100644
--- a/NumLib/ODESolver/EquationSystem.h
+++ b/NumLib/ODESolver/EquationSystem.h
@@ -39,6 +39,15 @@ public:
*/
virtual bool isLinear() const = 0;
+ /*! Check whether normalization of A and rhs is required.
+ *
+ * \remark
+ * In some processes, a normalization operation is required, to calculate
+ * A^T * A, and overwrite A; also calculate A^T * rhs and overwrite rhs.
+ * This parameter reflect whether such operation is required.
+ */
+ virtual bool requiresNormalization() const = 0;
+
/*! Prepares a new iteration in the solution process of this equation.
*
* \param iter the current iteration number, starting from 1.
diff --git a/NumLib/ODESolver/MatrixTranslator.cpp b/NumLib/ODESolver/MatrixTranslator.cpp
index 5be630904e119bcafc50ba743010a1a5f6ea84d3..c86edbba043b5e45aa24bd391ead4e28c9f574e5 100644
--- a/NumLib/ODESolver/MatrixTranslator.cpp
+++ b/NumLib/ODESolver/MatrixTranslator.cpp
@@ -43,6 +43,25 @@ void MatrixTranslatorGeneral<ODESystemTag::FirstOrderImplicitQuasilinear>::
NumLib::GlobalVectorProvider::provider.releaseVector(tmp);
}
+void MatrixTranslatorGeneral<ODESystemTag::FirstOrderImplicitQuasilinear>::
+ normalizeAandRhs(GlobalMatrix& A, GlobalVector& b) const
+{
+ namespace LinAlg = MathLib::LinAlg;
+
+ // check whether A is square?
+
+ GlobalMatrix new_A(A);
+ GlobalVector new_b(b);
+ LinAlg::copy(A, new_A);
+ LinAlg::copy(b, new_b);
+ // rhs = A^T * rhs
+ // A = A^T * A
+ LinAlg::linearSysNormalize(A, new_A, b, new_b);
+
+ LinAlg::copy(new_A, A);
+ LinAlg::copy(new_b, b);
+}
+
void MatrixTranslatorGeneral<ODESystemTag::FirstOrderImplicitQuasilinear>::
computeResidual(GlobalMatrix const& M, GlobalMatrix const& K,
GlobalVector const& b, double const dt,
diff --git a/NumLib/ODESolver/MatrixTranslator.h b/NumLib/ODESolver/MatrixTranslator.h
index b3752586c0bde9e2bf7b72ef26d9d717e78b5dee..60cc11dfcf5ee11acc2426ba1fb11ff78b8dc653 100644
--- a/NumLib/ODESolver/MatrixTranslator.h
+++ b/NumLib/ODESolver/MatrixTranslator.h
@@ -49,6 +49,10 @@ public:
const GlobalVector& b, const GlobalVector& x_prev,
GlobalVector& rhs) const = 0;
+ //! Computes \f$ A = A^T \cdot A \f$, and
+ //! also \f$ rhs = A^T \cdot rhs \f$.
+ virtual void normalizeAandRhs(GlobalMatrix& A, GlobalVector& b) const = 0;
+
/*! Computes \c res from \c M, \c K, \c b, \f$ \hat x \f$ and \f$ x_N \f$.
* You might also want read the remarks on
* \ref concept_time_discretization "time discretization".
@@ -105,6 +109,10 @@ public:
const GlobalVector& b, const GlobalVector& x_prev,
GlobalVector& rhs) const override;
+ //! Computes \f$ A = A^T \cdot A \f$, and
+ //! also \f$ rhs = A^T \cdot rhs \f$.
+ void normalizeAandRhs(GlobalMatrix& A, GlobalVector& b) const override;
+
//! Computes \f$ r = M \cdot \hat x + K \cdot x_C - b \f$.
void computeResidual(GlobalMatrix const& M, GlobalMatrix const& K,
GlobalVector const& b, double dt,
diff --git a/NumLib/ODESolver/NonlinearSolver.cpp b/NumLib/ODESolver/NonlinearSolver.cpp
index 324c11905c646a6aa7159e2e70795abec62a2e2a..59c0886b55bfc7038cf25188999cb25475cc852d 100644
--- a/NumLib/ODESolver/NonlinearSolver.cpp
+++ b/NumLib/ODESolver/NonlinearSolver.cpp
@@ -167,6 +167,11 @@ NonlinearSolverStatus NonlinearSolver<NonlinearSolverTag::Picard>::solve(
sys.assemble(x_new, x_prev, process_id);
sys.getA(A);
sys.getRhs(*x_prev[process_id], rhs);
+
+ // Normalize the linear equation system, if required
+ if (sys.requiresNormalization())
+ sys.getAandRhsNormalized(A, rhs);
+
INFO("[time] Assembly took {:g} s.", time_assembly.elapsed());
// Subtract non-equilibrium initial residuum if set
diff --git a/NumLib/ODESolver/NonlinearSystem.h b/NumLib/ODESolver/NonlinearSystem.h
index e3b20f6e12d2a41aea49260afcbf3a0b26a94e16..82d9681525eb94e9f4dba8ddaf9dd3155b53456f 100644
--- a/NumLib/ODESolver/NonlinearSystem.h
+++ b/NumLib/ODESolver/NonlinearSystem.h
@@ -121,6 +121,12 @@ public:
virtual void getRhs(GlobalVector const& x_prev,
GlobalVector& rhs) const = 0;
+ //! Writes the A_transposed times A into \c A
+ //! and also writes A_transposed times rhs into \c rhs
+ //! \pre getA() and getRhs must have been called before.
+ virtual void getAandRhsNormalized(GlobalMatrix& A,
+ GlobalVector& rhs) const = 0;
+
//! Pre-compute known solutions and possibly store them internally.
virtual void computeKnownSolutions(GlobalVector const& x,
int const process_id) = 0;
diff --git a/NumLib/ODESolver/TimeDiscretizedODESystem.h b/NumLib/ODESolver/TimeDiscretizedODESystem.h
index 4fffd07cbcbed200093087a3eaff2adcd2dfdc51..e93751dd9d8c5fd32a8319d6d5fd3093fffa9a7a 100644
--- a/NumLib/ODESolver/TimeDiscretizedODESystem.h
+++ b/NumLib/ODESolver/TimeDiscretizedODESystem.h
@@ -120,6 +120,11 @@ public:
bool isLinear() const override { return _ode.isLinear(); }
+ bool requiresNormalization() const override
+ {
+ return _ode.requiresNormalization();
+ }
+
void preIteration(const unsigned iter, GlobalVector const& x) override
{
_ode.preIteration(iter, x);
@@ -207,6 +212,11 @@ public:
_mat_trans->computeRhs(*_M, *_K, *_b, x_prev, rhs);
}
+ void getAandRhsNormalized(GlobalMatrix& A, GlobalVector& rhs) const override
+ {
+ _mat_trans->normalizeAandRhs(A, rhs);
+ }
+
void computeKnownSolutions(GlobalVector const& x,
int const process_id) override;
@@ -217,6 +227,11 @@ public:
bool isLinear() const override { return _ode.isLinear(); }
+ bool requiresNormalization() const override
+ {
+ return _ode.requiresNormalization();
+ }
+
void preIteration(const unsigned iter, GlobalVector const& x) override
{
_ode.preIteration(iter, x);
diff --git a/ProcessLib/HeatTransportBHE/BHE/BHECommon.h b/ProcessLib/HeatTransportBHE/BHE/BHECommon.h
index 85acb09d18e7b264148eaa39a2af6ae93f90b20e..a084529ae249b26b9091e1ce47b1e93d8ee8ba21 100644
--- a/ProcessLib/HeatTransportBHE/BHE/BHECommon.h
+++ b/ProcessLib/HeatTransportBHE/BHE/BHECommon.h
@@ -43,6 +43,11 @@ struct BHECommon
GroutParameters const grout;
FlowAndTemperatureControl const flowAndTemperatureControl;
bool const use_python_bcs;
+ constexpr bool isPowerBC() const
+ {
+ return std::visit([](auto const& ftc) { return ftc.is_power_bc; },
+ flowAndTemperatureControl);
+ }
};
} // end of namespace BHE
} // end of namespace HeatTransportBHE
diff --git a/ProcessLib/HeatTransportBHE/BHE/CreateFlowAndTemperatureControl.cpp b/ProcessLib/HeatTransportBHE/BHE/CreateFlowAndTemperatureControl.cpp
index c290cfa6eaa1c2588ec237d92a956ee6d49f25e5..ced0293262dccd2aefb5932b7856882b76e9ae72 100644
--- a/ProcessLib/HeatTransportBHE/BHE/CreateFlowAndTemperatureControl.cpp
+++ b/ProcessLib/HeatTransportBHE/BHE/CreateFlowAndTemperatureControl.cpp
@@ -67,6 +67,7 @@ FlowAndTemperatureControl createFlowAndTemperatureControl(
//! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__borehole_heat_exchangers__borehole_heat_exchanger__flow_and_temperature_control__FixedPowerConstantFlow__flow_rate}
auto const flow_rate = config.getConfigParameter<double>("flow_rate");
+
return FixedPowerConstantFlow{flow_rate, power,
refrigerant.specific_heat_capacity,
refrigerant.density};
diff --git a/ProcessLib/HeatTransportBHE/BHE/FlowAndTemperatureControl.h b/ProcessLib/HeatTransportBHE/BHE/FlowAndTemperatureControl.h
index b0b815a57c9793baadab29d6b36452c46ad69be1..06f477527120518ceaa84c82c52ea60570cf9483 100644
--- a/ProcessLib/HeatTransportBHE/BHE/FlowAndTemperatureControl.h
+++ b/ProcessLib/HeatTransportBHE/BHE/FlowAndTemperatureControl.h
@@ -35,6 +35,7 @@ struct TemperatureCurveConstantFlow
}
double flow_rate;
MathLib::PiecewiseLinearInterpolation const& temperature_curve;
+ static constexpr bool is_power_bc = false;
};
struct TemperatureCurveFlowCurve
@@ -47,6 +48,7 @@ struct TemperatureCurveFlowCurve
}
MathLib::PiecewiseLinearInterpolation const& flow_rate_curve;
MathLib::PiecewiseLinearInterpolation const& temperature_curve;
+ static constexpr bool is_power_bc = false;
};
struct FixedPowerConstantFlow
@@ -60,6 +62,7 @@ struct FixedPowerConstantFlow
double power; // Value is expected to be in Watt.
double heat_capacity;
double density;
+ static constexpr bool is_power_bc = true;
};
struct FixedPowerFlowCurve
@@ -74,6 +77,7 @@ struct FixedPowerFlowCurve
double power; // Value is expected to be in Watt.
double heat_capacity;
double density;
+ static constexpr bool is_power_bc = true;
};
struct PowerCurveConstantFlow
@@ -92,6 +96,7 @@ struct PowerCurveConstantFlow
double flow_rate;
double heat_capacity;
double density;
+ static constexpr bool is_power_bc = true;
};
struct PowerCurveFlowCurve
@@ -111,6 +116,7 @@ struct PowerCurveFlowCurve
double heat_capacity;
double density;
+ static constexpr bool is_power_bc = true;
};
struct BuildingPowerCurveConstantFlow
@@ -134,6 +140,7 @@ struct BuildingPowerCurveConstantFlow
double flow_rate;
double heat_capacity;
double density;
+ static constexpr bool is_power_bc = true;
};
using FlowAndTemperatureControl = std::variant<TemperatureCurveConstantFlow,
diff --git a/ProcessLib/HeatTransportBHE/CreateHeatTransportBHEProcess.cpp b/ProcessLib/HeatTransportBHE/CreateHeatTransportBHEProcess.cpp
index 2e215b8c31eb6228858650bddf727f81230ce2ec..2a83a931a8de0713c58965e18a15115635077ed7 100644
--- a/ProcessLib/HeatTransportBHE/CreateHeatTransportBHEProcess.cpp
+++ b/ProcessLib/HeatTransportBHE/CreateHeatTransportBHEProcess.cpp
@@ -104,6 +104,39 @@ std::unique_ptr<Process> createHeatTransportBHEProcess(
//! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__use_server_communication}
config.getConfigParameter<bool>("use_server_communication", false);
+ auto const using_algebraic_bc =
+ //! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__use_algebraic_bc}
+ config.getConfigParameter<bool>("use_algebraic_bc", false);
+
+ auto const weighting_factor =
+ //! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__weighting_factor}
+ config.getConfigParameter<float>("weighting_factor", 100.0);
+
+ auto const is_linear =
+ //! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__linear}
+ config.getConfigParameter<bool>("linear", false);
+ if (is_linear)
+ {
+ if (!using_algebraic_bc)
+ {
+ OGS_FATAL(
+ "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.")
+ }
+ else
+ {
+ WARN(
+ "You specified that the process simulated by OGS is linear. "
+ "With that optimization the process will be assembled only "
+ "once and the non-linear solver will do only one iteration per "
+ "time step. No 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!");
+ }
+ }
+
for (
auto const& bhe_config :
//! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__borehole_heat_exchangers__borehole_heat_exchanger}
@@ -218,7 +251,8 @@ std::unique_ptr<Process> createHeatTransportBHEProcess(
HeatTransportBHEProcessData process_data(
std::move(media_map), std::move(bhes), py_object, using_tespy,
- using_server_communication);
+ using_server_communication,
+ {using_algebraic_bc, weighting_factor, is_linear});
SecondaryVariableCollection secondary_variables;
diff --git a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp
index 5b2a54b1110f5fd3d1fa6c6ca23693c812efce4a..4e7add016613aaf0e48c351facef95dded91637b 100644
--- a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp
+++ b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp
@@ -172,6 +172,13 @@ void HeatTransportBHEProcess::assembleConcreteProcess(
_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)
+ {
+ algebraicBcConcreteProcess(t, dt, x, x_prev, process_id, M, K, b);
+ }
+
+ //_global_output(t, process_id, M, K, b);
}
void HeatTransportBHEProcess::assembleWithJacobianConcreteProcess(
@@ -340,16 +347,120 @@ void HeatTransportBHEProcess::postTimestepConcreteProcess(
}
}
+void HeatTransportBHEProcess::algebraicBcConcreteProcess(
+ [[maybe_unused]] const double t, double const /*dt*/,
+ [[maybe_unused]] std::vector<GlobalVector*> const& x,
+ std::vector<GlobalVector*> const& /*xprev*/, int const /*process_id*/,
+ [[maybe_unused]] GlobalMatrix& M, [[maybe_unused]] GlobalMatrix& K,
+ [[maybe_unused]] GlobalVector& b)
+{
+#ifndef USE_PETSC
+ auto M_normal = M.getRawMatrix();
+ auto K_normal = K.getRawMatrix();
+ auto n_original_rows = K_normal.rows();
+ auto const n_BHE_bottom_pairs = _vec_bottom_BHE_node_indices.size();
+ auto const n_BHE_top_pairs = _vec_top_BHE_node_indices.size();
+
+ // apply weighting factor based on the max value from column wise inner
+ // product and scale it with user defined value
+ const double w_val =
+ _process_data._algebraic_BC_Setting._weighting_factor *
+ (Eigen::RowVectorXd::Ones(K_normal.rows()) * K_normal.cwiseAbs())
+ .maxCoeff();
+
+ M_normal.conservativeResize(
+ M_normal.rows() + n_BHE_bottom_pairs + n_BHE_top_pairs,
+ M_normal.cols());
+ K_normal.conservativeResize(
+ K_normal.rows() + n_BHE_bottom_pairs + n_BHE_top_pairs,
+ K_normal.cols());
+
+ for (std::size_t i = 0; i < n_BHE_bottom_pairs; i++)
+ {
+ Eigen::SparseVector<double> M_Plus(M_normal.cols());
+ M_Plus.setZero();
+ M_normal.row(n_original_rows + i) = M_Plus;
+
+ Eigen::SparseVector<double> K_Plus(K_normal.cols());
+ K_Plus.setZero();
+
+ auto const [bhe_idx, first_BHE_bottom_index, second_BHE_bottom_index] =
+ _vec_bottom_BHE_node_indices[i];
+
+ K_Plus.insert(first_BHE_bottom_index) = w_val;
+ K_Plus.insert(second_BHE_bottom_index) = -w_val;
+
+ K_normal.row(n_original_rows + i) = K_Plus;
+ }
+
+ auto b_normal = b.getRawVector();
+ Eigen::SparseVector<double> b_Plus(b_normal.rows() + n_BHE_bottom_pairs +
+ n_BHE_top_pairs);
+ b_Plus.setZero();
+
+ // Copy values from the original column vector to the modified one
+ for (int i = 0; i < b_normal.innerSize(); ++i)
+ {
+ b_Plus.insert(i) = b_normal.coeff(i);
+ }
+
+ for (std::size_t i = 0; i < n_BHE_top_pairs; i++)
+ {
+ Eigen::SparseVector<double> M_Plus(M_normal.cols());
+ M_Plus.setZero();
+ M_normal.row(n_original_rows + n_BHE_bottom_pairs + i) = M_Plus;
+
+ Eigen::SparseVector<double> K_Plus(K_normal.cols());
+ K_Plus.setZero();
+
+ auto const [bhe_idx, first_BHE_top_index, second_BHE_top_index] =
+ _vec_top_BHE_node_indices[i];
+
+ auto first_BHE_top_index_pair = first_BHE_top_index;
+ auto second_BHE_top_index_pair = second_BHE_top_index;
+
+ K_Plus.insert(first_BHE_top_index_pair) =
+ w_val; // for power BC, the inflow node must be positive
+ K_Plus.insert(second_BHE_top_index_pair) =
+ -w_val; // for power BC, the outflow node must be negative
+
+ K_normal.row(n_original_rows + n_BHE_bottom_pairs + i) = K_Plus;
+
+ // get the delta_T value here
+ double const T_out = (*x[0])[second_BHE_top_index_pair];
+
+ auto calculate_delta_T = [&](auto& bhe)
+ {
+ auto const T_in = bhe.updateFlowRateAndTemperature(T_out, t);
+ return T_in - T_out;
+ };
+ auto delta_T = std::visit(calculate_delta_T,
+ _process_data._vec_BHE_property[bhe_idx]);
+
+ b_Plus.insert(n_original_rows + n_BHE_bottom_pairs + i) =
+ delta_T * w_val;
+ }
+
+ M.getRawMatrix() = M_normal;
+ K.getRawMatrix() = K_normal;
+ b.getRawVector() = b_Plus;
+#else
+ OGS_FATAL(
+ "The Algebraic Boundary Condition is not implemented for use with "
+ "PETsc Library! Simulation will be terminated.");
+#endif
+}
+
void HeatTransportBHEProcess::createBHEBoundaryConditionTopBottom(
std::vector<std::vector<MeshLib::Node*>> const& all_bhe_nodes)
{
const int process_id = 0;
auto& bcs = _boundary_conditions[process_id];
- int const n_BHEs = static_cast<int>(_process_data._vec_BHE_property.size());
+ std::size_t const n_BHEs = _process_data._vec_BHE_property.size();
// for each BHE
- for (int bhe_i = 0; bhe_i < n_BHEs; bhe_i++)
+ for (std::size_t bhe_i = 0; bhe_i < n_BHEs; bhe_i++)
{
auto const& bhe_nodes = all_bhe_nodes[bhe_i];
// find the variable ID
@@ -434,6 +545,20 @@ void HeatTransportBHEProcess::createBHEBoundaryConditionTopBottom(
nodes_and_components[1].second));
};
+ auto get_global_bhe_bc_indices_with_bhe_idx =
+ [&](std::size_t bhe_idx,
+ std::array<
+ std::pair<std::size_t /*node_id*/, int /*component*/>, 2>
+ nodes_and_components)
+ {
+ return std::make_tuple(
+ bhe_idx,
+ get_global_index(nodes_and_components[0].first,
+ nodes_and_components[0].second),
+ get_global_index(nodes_and_components[1].first,
+ nodes_and_components[1].second));
+ };
+
auto createBCs =
[&, bc_top_node_id = bhe_boundary_nodes[0]->getID(),
bc_bottom_node_id = bhe_boundary_nodes[1]->getID()](auto& bhe)
@@ -442,10 +567,11 @@ void HeatTransportBHEProcess::createBHEBoundaryConditionTopBottom(
bhe.inflow_outflow_bc_component_ids)
{
if (bhe.use_python_bcs ||
- _process_data._use_server_communication)
+ this->_process_data._use_server_communication)
// call BHEPythonBoundarycondition
{
- if (_process_data.py_bc_object) // the bc object exist
+ if (this->_process_data
+ .py_bc_object) // the bc object exist
{
// apply the customized top, inflow BC.
bcs.addBoundaryCondition(
@@ -466,16 +592,33 @@ void HeatTransportBHEProcess::createBHEBoundaryConditionTopBottom(
}
else
{
- // Top, inflow, normal case
- bcs.addBoundaryCondition(
- createBHEInflowDirichletBoundaryCondition(
- get_global_bhe_bc_indices(
- bhe.getBHEInflowDirichletBCNodesAndComponents(
- bc_top_node_id, bc_bottom_node_id,
- in_out_component_id.first)),
- [&bhe](double const T, double const t) {
- return bhe.updateFlowRateAndTemperature(T, t);
- }));
+ if (this->_process_data._algebraic_BC_Setting
+ ._use_algebraic_bc &&
+ bhe.isPowerBC())
+ {
+ // for algebraic_bc method, record the pair of indices
+ // in a separate vector
+ _vec_top_BHE_node_indices.push_back(
+ get_global_bhe_bc_indices_with_bhe_idx(
+ bhe_i,
+ {{{bc_top_node_id, in_out_component_id.first},
+ {bc_top_node_id,
+ in_out_component_id.second}}}));
+ }
+ else
+ {
+ // Top, inflow, normal case
+ bcs.addBoundaryCondition(
+ createBHEInflowDirichletBoundaryCondition(
+ get_global_bhe_bc_indices(
+ bhe.getBHEInflowDirichletBCNodesAndComponents(
+ bc_top_node_id, bc_bottom_node_id,
+ in_out_component_id.first)),
+ [&bhe](double const T, double const t) {
+ return bhe.updateFlowRateAndTemperature(T,
+ t);
+ }));
+ }
}
auto const bottom_nodes_and_components =
@@ -484,9 +627,12 @@ void HeatTransportBHEProcess::createBHEBoundaryConditionTopBottom(
in_out_component_id.first,
in_out_component_id.second);
- if (bottom_nodes_and_components)
+ if (bottom_nodes_and_components &&
+ !this->_process_data._algebraic_BC_Setting
+ ._use_algebraic_bc)
{
- // Bottom, outflow, all cases
+ // Bottom, outflow, all cases | not needed for algebraic_bc
+ // method
bcs.addBoundaryCondition(
createBHEBottomDirichletBoundaryCondition(
get_global_bhe_bc_indices(
@@ -495,6 +641,19 @@ void HeatTransportBHEProcess::createBHEBoundaryConditionTopBottom(
{bc_bottom_node_id,
in_out_component_id.second}}})));
}
+ else if (bottom_nodes_and_components &&
+ this->_process_data._algebraic_BC_Setting
+ ._use_algebraic_bc)
+ {
+ // for algebraic_bc method, record the pair of indices in a
+ // separate vector
+ _vec_bottom_BHE_node_indices.push_back(
+ get_global_bhe_bc_indices_with_bhe_idx(
+ bhe_i,
+ {{{bc_bottom_node_id, in_out_component_id.first},
+ {bc_bottom_node_id,
+ in_out_component_id.second}}}));
+ }
}
};
visit(createBCs, _process_data._vec_BHE_property[bhe_i]);
diff --git a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h
index e59d3af65fcc6ecc1cc3cc5d45273c5ef0ea6865..a50d8f66e2ba1f2687791887629617dc9ec2d85d 100644
--- a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h
+++ b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h
@@ -39,7 +39,18 @@ public:
//! \name ODESystem interface
//! @{
- bool isLinear() const override { return false; }
+ bool isLinear() const override
+ {
+ return _process_data._algebraic_BC_Setting._is_linear;
+ }
+
+ bool requiresNormalization() const override
+ {
+ // In the current setup, when using algebraic bc,
+ // then normalization is always required
+ return _process_data._algebraic_BC_Setting._use_algebraic_bc;
+ }
+ //! @}
void computeSecondaryVariableConcrete(double const t, double const dt,
std::vector<GlobalVector*> const& x,
@@ -76,6 +87,12 @@ private:
const double t, const double dt,
int const process_id) override;
+ void algebraicBcConcreteProcess(const double t, double const dt,
+ std::vector<GlobalVector*> const& x,
+ std::vector<GlobalVector*> const& xdot,
+ int const process_id, GlobalMatrix& M,
+ GlobalMatrix& K, GlobalVector& b);
+
NumLib::IterationResult postIterationConcreteProcess(
GlobalVector const& x) override;
@@ -89,6 +106,18 @@ private:
std::vector<std::unique_ptr<MeshLib::MeshSubset const>>
_mesh_subset_BHE_soil_nodes;
+ // a vector of tuple structure containing the indices of BHE top nodes,
+ // used only for algebraic boundary conditions
+ // first object is the index of BHE
+ // second and third object is the global indices of a pair of unknowns,
+ // pointing to the inflow and outflow temperature
+ std::vector<std::tuple<std::size_t, GlobalIndexType, GlobalIndexType>>
+ _vec_top_BHE_node_indices;
+ // a vector of tuple structure containing the indices of BHE bottom nodes,
+ // used only for algebraic boundary conditions
+ // same structure as the top node vector
+ std::vector<std::tuple<std::size_t, GlobalIndexType, GlobalIndexType>>
+ _vec_bottom_BHE_node_indices;
std::unique_ptr<MeshLib::MeshSubset const> _mesh_subset_soil_nodes;
diff --git a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcessData.h b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcessData.h
index 134b42c5a5415d7fccc53d4727e5e906e4972322..5aa319a139c2f2b28bde86555a5d49e62d7437ce 100644
--- a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcessData.h
+++ b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcessData.h
@@ -23,6 +23,15 @@ class Element;
namespace ProcessLib::HeatTransportBHE
{
+struct AlgebraicBCSetting
+{
+ const bool _use_algebraic_bc;
+
+ const double _weighting_factor;
+
+ const bool _is_linear;
+};
+
struct HeatTransportBHEProcessData final
{
HeatTransportBHEProcessData(
@@ -31,12 +40,14 @@ struct HeatTransportBHEProcessData final
BHEInflowPythonBoundaryConditionPythonSideInterface* py_bc_object_ =
nullptr,
const bool use_tespy = false,
- const bool use_server_communication = false)
+ const bool use_server_communication = false,
+ AlgebraicBCSetting algebraicBCSetting = {false, 100.0, 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)
+ _use_server_communication(use_server_communication),
+ _algebraic_BC_Setting(algebraicBCSetting)
{
}
MaterialPropertyLib::MaterialSpatialDistributionMap media_map;
@@ -52,5 +63,7 @@ struct HeatTransportBHEProcessData final
const bool _use_tespy;
const bool _use_server_communication;
+
+ AlgebraicBCSetting const _algebraic_BC_Setting;
};
} // namespace ProcessLib::HeatTransportBHE
diff --git a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE-impl.h b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE-impl.h
index e4fb2d485fe66eb8490b213f7d3355f195ea4337..f4c44748c0d294e4f50c626adeadab1489df81d1 100644
--- a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE-impl.h
+++ b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE-impl.h
@@ -234,8 +234,8 @@ void HeatTransportBHELocalAssemblerBHE<ShapeFunction, BHEType>::
auto local_rhs = MathLib::createZeroedVector<BheLocalVectorType>(
local_rhs_data, local_matrix_size);
- std::vector<double> local_M_data(local_Jac_data.size());
- std::vector<double> local_K_data(local_Jac_data.size());
+ std::vector<double> local_M_data;
+ std::vector<double> local_K_data;
assemble(t, dt, local_x, local_x_prev, local_M_data, local_K_data,
local_rhs_data /*not going to be used*/);
diff --git a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h
index fc1443e76f99470bc426e1ed73dbb11e93b0359c..6b9b5b42e58b705c5318e7d101c595a0a7d7d9af 100644
--- a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h
+++ b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h
@@ -211,8 +211,8 @@ void HeatTransportBHELocalAssemblerSoil<ShapeFunction>::assembleWithJacobian(
auto local_rhs = MathLib::createZeroedVector<NodalVectorType>(
local_rhs_data, local_matrix_size);
- std::vector<double> local_M_data(local_Jac_data.size());
- std::vector<double> local_K_data(local_Jac_data.size());
+ std::vector<double> local_M_data;
+ std::vector<double> local_K_data;
assemble(t, dt, local_x, local_x_prev, local_M_data, local_K_data,
local_rhs_data /*not going to be used*/);
diff --git a/ProcessLib/HeatTransportBHE/Tests.cmake b/ProcessLib/HeatTransportBHE/Tests.cmake
index 7334aa699df821e36ac68813d652e6affa2e09b1..19f4fc36d7b1b5453c27a516f7eed11d685e7d85 100644
--- a/ProcessLib/HeatTransportBHE/Tests.cmake
+++ b/ProcessLib/HeatTransportBHE/Tests.cmake
@@ -1,3 +1,73 @@
+AddTest(
+ NAME HeatTransportBHE_1U_3D_bhe_sandwich
+ PATH Parabolic/T/3D_BHE_Sandwich
+ EXECUTABLE ogs
+ EXECUTABLE_ARGS sandwich.prj
+ WRAPPER time
+ TESTER vtkdiff
+ REQUIREMENTS NOT OGS_USE_MPI
+ RUNTIME 10
+ DIFF_DATA
+ sandwich_ts_10_t_600.000000.vtu sandwich_ts_10_t_600.000000.vtu temperature_BHE1 temperature_BHE1 0 5e-15
+ 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_Newton
+ PATH Parabolic/T/3D_BHE_Sandwich
+ EXECUTABLE ogs
+ EXECUTABLE_ARGS sandwich_newton.xml
+ WRAPPER time
+ TESTER vtkdiff
+ REQUIREMENTS NOT OGS_USE_MPI
+ RUNTIME 10
+ DIFF_DATA
+ sandwich_ts_10_t_600.000000.vtu sandwich_newton_ts_10_t_600.000000.vtu temperature_BHE1 temperature_BHE1 0 1e-14
+ sandwich_ts_10_t_600.000000.vtu sandwich_newton_ts_10_t_600.000000.vtu temperature_soil temperature_soil 0 1e-13
+)
+
+AddTest(
+ NAME HeatTransportBHE_1U_3D_bhe_sandwich_algebraicBC
+ PATH Parabolic/T/3D_BHE_Sandwich
+ EXECUTABLE ogs
+ EXECUTABLE_ARGS sandwich_algebraicBC.xml
+ WRAPPER time
+ TESTER vtkdiff
+ REQUIREMENTS NOT OGS_USE_MPI
+ RUNTIME 20
+ DIFF_DATA
+ sandwich_ts_10_t_600.000000.vtu sandwich_algebraic_bc_ts_10_t_600.000000.vtu temperature_BHE1 temperature_BHE1 0 1e-6
+ sandwich_ts_10_t_600.000000.vtu sandwich_algebraic_bc_ts_10_t_600.000000.vtu temperature_soil temperature_soil 0 5e-9
+)
+
+AddTest(
+ NAME HeatTransportBHE_1U_3D_bhe_sandwich_fixed_power
+ PATH Parabolic/T/3D_BHE_Sandwich
+ RUNTIME 100
+ EXECUTABLE ogs
+ EXECUTABLE_ARGS sandwich_fixed_power.prj
+ WRAPPER time
+ TESTER vtkdiff
+ REQUIREMENTS NOT OGS_USE_MPI
+ DIFF_DATA
+ sandwich_fixed_power_ts_10_t_600.000000.vtu sandwich_fixed_power_ts_10_t_600.000000.vtu temperature_BHE1 temperature_BHE1 0 5e-15
+ sandwich_fixed_power_ts_10_t_600.000000.vtu sandwich_fixed_power_ts_10_t_600.000000.vtu temperature_soil temperature_soil 0 1e-13
+)
+
+AddTest(
+ NAME HeatTransportBHE_1U_3D_bhe_sandwich_fixed_power_algebraicBC
+ PATH Parabolic/T/3D_BHE_Sandwich
+ RUNTIME 100
+ EXECUTABLE ogs
+ EXECUTABLE_ARGS sandwich_fixed_power_algebraicBC.xml
+ WRAPPER time
+ TESTER vtkdiff
+ REQUIREMENTS NOT OGS_USE_MPI
+ DIFF_DATA
+ sandwich_fixed_power_ts_10_t_600.000000.vtu sandwich_fixed_power_algebraic_bc_ts_10_t_600.000000.vtu temperature_BHE1 temperature_BHE1 0 5e-3
+ sandwich_fixed_power_ts_10_t_600.000000.vtu sandwich_fixed_power_algebraic_bc_ts_10_t_600.000000.vtu temperature_soil temperature_soil 0 1e-6
+)
+
AddTest(
NAME HeatTransportBHE_1U_3D_beier_sandbox
PATH Parabolic/T/3D_Beier_sandbox
@@ -26,6 +96,20 @@ AddTest(
beier_sandbox_ts_10_t_600.000000.vtu beier_sandbox_newton_ts_10_t_600.000000.vtu temperature_soil temperature_soil 0 1e-13
)
+AddTest(
+ NAME HeatTransportBHE_1U_3D_beier_sandbox_algebraicBC
+ PATH Parabolic/T/3D_Beier_sandbox
+ EXECUTABLE ogs
+ EXECUTABLE_ARGS beier_sandbox_algebraicBC.xml
+ WRAPPER time
+ TESTER vtkdiff
+ REQUIREMENTS NOT OGS_USE_MPI
+ RUNTIME 20
+ DIFF_DATA
+ beier_sandbox_ts_10_t_600.000000.vtu beier_sandbox_algebraic_bc_ts_10_t_600.000000.vtu temperature_BHE1 temperature_BHE1 0 5e-7
+ beier_sandbox_ts_10_t_600.000000.vtu beier_sandbox_algebraic_bc_ts_10_t_600.000000.vtu temperature_soil temperature_soil 0 5e-10
+)
+
AddTest(
NAME HeatTransportBHE_1U_beier_sandbox_fixed_power_constant_flow
PATH Parabolic/T/3D_Beier_sandbox
@@ -40,6 +124,20 @@ AddTest(
fixed_power_constant_flow_ts_10_t_600.000000.vtu fixed_power_constant_flow_ts_10_t_600.000000.vtu temperature_soil temperature_soil 0 1e-13
)
+AddTest(
+ NAME HeatTransportBHE_1U_beier_sandbox_fixed_power_constant_flow_algebraicBC
+ PATH Parabolic/T/3D_Beier_sandbox
+ RUNTIME 220
+ EXECUTABLE ogs
+ EXECUTABLE_ARGS fixed_power_constant_flow_algebraicBC.xml
+ WRAPPER time
+ TESTER vtkdiff
+ REQUIREMENTS NOT OGS_USE_MPI
+ DIFF_DATA
+ fixed_power_constant_flow_ts_10_t_600.000000.vtu fixed_power_constant_flow_algebraic_bc_ts_10_t_600.000000.vtu temperature_BHE1 temperature_BHE1 0 1e-4
+ fixed_power_constant_flow_ts_10_t_600.000000.vtu fixed_power_constant_flow_algebraic_bc_ts_10_t_600.000000.vtu temperature_soil temperature_soil 0 5e-9
+)
+
AddTest(
NAME HeatTransportBHE_coaxial_pipe_3D_deep_BHE_CXA
PATH Parabolic/T/3D_deep_BHE
@@ -85,7 +183,7 @@ AddTest(
AddTest(
NAME HeatTransportBHE_3D_BHE_groundwater_advection
PATH Parabolic/T/3D_BHE_GW_advection
- RUNTIME 8
+ RUNTIME 4
EXECUTABLE ogs
EXECUTABLE_ARGS BHE_GW_advection.prj
WRAPPER time
@@ -176,5 +274,5 @@ AddTest(
)
if(NOT OGS_USE_PETSC)
- NotebookTest(NOTEBOOKFILE Parabolic/T/BHE_1P/pipe_flow_ebhe.md RUNTIME 600)
+ NotebookTest(NOTEBOOKFILE Parabolic/T/BHE_1P/pipe_flow_ebhe.md RUNTIME 200)
endif()
diff --git a/ProcessLib/Process.h b/ProcessLib/Process.h
index 3b2f076105dfe583ea709b33677caf4bfaaf2c8b..b8c4fd9ed64b77eb76548a8296a3b2392bd1f436 100644
--- a/ProcessLib/Process.h
+++ b/ProcessLib/Process.h
@@ -192,6 +192,8 @@ public:
{
}
+ bool requiresNormalization() const override { return false; }
+
protected:
std::vector<NumLib::LocalToGlobalIndexMap const*> getDOFTables(
int const number_of_processes) const;
diff --git a/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich.prj b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich.prj
new file mode 100644
index 0000000000000000000000000000000000000000..d07938d11aefda2deeb191a6606a59b15ff689f3
--- /dev/null
+++ b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich.prj
@@ -0,0 +1,242 @@
+<?xml version="1.0" encoding="ISO-8859-1"?>
+<OpenGeoSysProject>
+ <meshes>
+ <mesh>sandwich.vtu</mesh>
+ </meshes>
+ <processes>
+ <process>
+ <name>HeatTransportBHE</name>
+ <type>HEAT_TRANSPORT_BHE</type>
+ <integration_order>2</integration_order>
+ <process_variables>
+ <process_variable>temperature_soil</process_variable>
+ <process_variable>temperature_BHE1</process_variable>
+ </process_variables>
+ <borehole_heat_exchangers>
+ <borehole_heat_exchanger>
+ <type>1U</type>
+ <flow_and_temperature_control>
+ <type>TemperatureCurveConstantFlow</type>
+ <flow_rate>1e-3</flow_rate>
+ <temperature_curve>fluid_temperature</temperature_curve>
+ </flow_and_temperature_control>
+ <borehole>
+ <length>3.0</length>
+ <diameter>0.2</diameter>
+ </borehole>
+ <grout>
+ <density>2190.0</density>
+ <porosity>0.0</porosity>
+ <specific_heat_capacity>1781</specific_heat_capacity>
+ <thermal_conductivity>1.4</thermal_conductivity>
+ </grout>
+ <pipes>
+ <inlet>
+ <diameter> 0.02</diameter>
+ <wall_thickness>0.003</wall_thickness>
+ <wall_thermal_conductivity>0.43</wall_thermal_conductivity>
+ </inlet>
+ <outlet>
+ <diameter>0.02</diameter>
+ <wall_thickness>0.003</wall_thickness>
+ <wall_thermal_conductivity>0.43</wall_thermal_conductivity>
+ </outlet>
+ <distance_between_pipes>0.075</distance_between_pipes>
+ <longitudinal_dispersion_length>0.001</longitudinal_dispersion_length>
+ </pipes>
+ <refrigerant>
+ <density>1000</density>
+ <viscosity>0.005</viscosity>
+ <specific_heat_capacity>4000</specific_heat_capacity>
+ <thermal_conductivity>0.5</thermal_conductivity>
+ <reference_temperature>293.15</reference_temperature>
+ </refrigerant>
+ </borehole_heat_exchanger>
+ </borehole_heat_exchangers>
+ </process>
+ </processes>
+ <media>
+ <medium id="0">
+ <phases>
+ <phase>
+ <type>AqueousLiquid</type>
+ <properties>
+ <property>
+ <name>phase_velocity</name>
+ <type>Constant</type>
+ <value>0 0 0</value>
+ </property>
+ <property>
+ <name>specific_heat_capacity</name>
+ <type>Constant</type>
+ <value>3795</value>
+ </property>
+ <property>
+ <name>density</name>
+ <type>Constant</type>
+ <value>1000</value>
+ </property>
+ </properties>
+ </phase>
+ <phase>
+ <type>Solid</type>
+ <properties>
+ <property>
+ <name>specific_heat_capacity</name>
+ <type>Constant</type>
+ <value>2000</value>
+ </property>
+ <property>
+ <name>density</name>
+ <type>Constant</type>
+ <value>1000</value>
+ </property>
+ </properties>
+ </phase>
+ <phase>
+ <type>Gas</type>
+ <properties>
+ <property>
+ <name>specific_heat_capacity</name>
+ <type>Constant</type>
+ <value>1000</value>
+ </property>
+ <property>
+ <name>density</name>
+ <type>Constant</type>
+ <value>2500</value>
+ </property>
+ </properties>
+ </phase>
+ </phases>
+ <properties>
+ <property>
+ <name>porosity</name>
+ <type>Constant</type>
+ <value>0</value>
+ </property>
+ <property>
+ <name>thermal_conductivity</name>
+ <type>Constant</type>
+ <value>2.0</value>
+ </property>
+ <property>
+ <name>thermal_longitudinal_dispersivity</name>
+ <type>Constant</type>
+ <value>5</value>
+ </property>
+ <property>
+ <name>thermal_transversal_dispersivity</name>
+ <type>Constant</type>
+ <value>0.5</value>
+ </property>
+ </properties>
+ </medium>
+ </media>
+ <time_loop>
+ <processes>
+ <process ref="HeatTransportBHE">
+ <nonlinear_solver>basic_picard</nonlinear_solver>
+ <convergence_criterion>
+ <type>DeltaX</type>
+ <norm_type>NORM2</norm_type>
+ <reltol>1.e-10</reltol>
+ </convergence_criterion>
+ <time_discretization><type>BackwardEuler</type></time_discretization>
+ <time_stepping>
+ <type>FixedTimeStepping</type>
+ <t_initial> 0.0 </t_initial>
+ <t_end> 600 </t_end>
+ <timesteps>
+ <pair>
+ <repeat>1</repeat>
+ <delta_t>60</delta_t>
+ </pair>
+ </timesteps>
+ </time_stepping>
+ </process>
+ </processes>
+ <output>
+ <type>VTK</type>
+ <prefix>sandwich</prefix>
+ <timesteps>
+ <pair>
+ <repeat> 1 </repeat>
+ <each_steps> 1 </each_steps>
+ </pair>
+ </timesteps>
+ <variables>
+ <variable>temperature_soil</variable>
+ <variable>temperature_BHE1</variable>
+ </variables>
+ </output>
+ </time_loop>
+ <parameters>
+ <parameter>
+ <name>T0</name>
+ <type>Constant</type>
+ <value>293.15</value>
+ </parameter>
+ <parameter>
+ <name>T_soil_imposed</name>
+ <type>Constant</type>
+ <value>273.15</value>
+ </parameter>
+ <parameter>
+ <name>T0_BHE</name>
+ <type>Constant</type>
+ <values>293.15 293.15 293.15 293.15</values>
+ </parameter>
+ </parameters>
+ <process_variables>
+ <process_variable>
+ <name>temperature_soil</name>
+ <components>1</components>
+ <order>1</order>
+ <initial_condition>T0</initial_condition>
+ <boundary_conditions>
+ </boundary_conditions>
+ </process_variable>
+ <process_variable>
+ <name>temperature_BHE1</name>
+ <components>4</components>
+ <order>1</order>
+ <initial_condition>T0_BHE</initial_condition>
+ </process_variable>
+ </process_variables>
+ <nonlinear_solvers>
+ <nonlinear_solver>
+ <name>basic_picard</name>
+ <type>Picard</type>
+ <max_iter>100</max_iter>
+ <linear_solver>general_linear_solver</linear_solver>
+ </nonlinear_solver>
+ </nonlinear_solvers>
+ <linear_solvers>
+ <linear_solver>
+ <name>general_linear_solver</name>
+ <lis>-i cg -p jacobi -tol 1e-16 -maxiter 10000</lis>
+ <eigen>
+ <solver_type>BiCGSTAB</solver_type>
+ <precon_type>ILUT</precon_type>
+ <max_iteration_step>1000</max_iteration_step>
+ <error_tolerance>1e-16</error_tolerance>
+ </eigen>
+ <petsc>
+ <prefix>gw</prefix>
+ <parameters>-gw_ksp_type cg -gw_pc_type bjacobi -gw_ksp_rtol 1e-16 -gw_ksp_max_it 10000</parameters>
+ </petsc>
+ </linear_solver>
+ </linear_solvers>
+ <curves>
+ <curve>
+ <name>fluid_temperature</name>
+ <coords>
+ 0 864000
+ </coords>
+ <values>
+ 273.15 273.15
+ </values>
+ </curve>
+ </curves>
+</OpenGeoSysProject>
diff --git a/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich.vtu b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..8939db6e61d2a2934620fda039bbc07113c66719
--- /dev/null
+++ b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich.vtu
@@ -0,0 +1,38 @@
+<?xml version="1.0"?>
+<VTKFile type="UnstructuredGrid" version="0.1" byte_order="LittleEndian" header_type="UInt32">
+ <UnstructuredGrid>
+ <Piece NumberOfPoints="12" NumberOfCells="6">
+ <PointData>
+ </PointData>
+ <CellData>
+ <DataArray type="Int32" Name="MaterialIDs" format="ascii" RangeMin="0" RangeMax="1">
+ 0 0 0 1 1 1
+ </DataArray>
+ </CellData>
+ <Points>
+ <DataArray type="Float64" Name="Points" NumberOfComponents="3" format="ascii" RangeMin="0" RangeMax="3.1622776602">
+ 0 0 0 0 1 0
+ 1 0 0 0 0 1
+ 0 1 1 1 0 1
+ 0 0 2 0 1 2
+ 1 0 2 0 0 3
+ 0 1 3 1 0 3
+ </DataArray>
+ </Points>
+ <Cells>
+ <DataArray type="Int64" Name="connectivity" format="ascii" RangeMin="0" RangeMax="11">
+ 0 1 2 3 4 5
+ 3 4 5 6 7 8
+ 6 7 8 9 10 11
+ 9 6 6 3 3 0
+ </DataArray>
+ <DataArray type="Int64" Name="offsets" format="ascii" RangeMin="6" RangeMax="24">
+ 6 12 18 20 22 24
+ </DataArray>
+ <DataArray type="UInt8" Name="types" format="ascii" RangeMin="3" RangeMax="13">
+ 13 13 13 3 3 3
+ </DataArray>
+ </Cells>
+ </Piece>
+ </UnstructuredGrid>
+</VTKFile>
diff --git a/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_algebraicBC.xml b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_algebraicBC.xml
new file mode 100644
index 0000000000000000000000000000000000000000..5df447f6532f4c96135aa206d214b6a912a13bf1
--- /dev/null
+++ b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_algebraicBC.xml
@@ -0,0 +1,15 @@
+<?xml version="1.0" encoding="ISO-8859-1"?>
+<OpenGeoSysProjectDiff base_file="sandwich.prj">
+ <add sel="/*/processes/process">
+ <use_algebraic_bc>true</use_algebraic_bc>
+ </add>
+ <add sel="/*/processes/process">
+ <weighting_factor>100</weighting_factor>
+ </add>
+ <add sel="/*/processes/process">
+ <linear>true</linear>
+ </add>
+ <replace sel="/*/time_loop/output/prefix/text()">sandwich_algebraic_bc</replace>
+ <replace sel="/*/linear_solvers/linear_solver/eigen/solver_type/text()">CG</replace>
+ <replace sel="/*/linear_solvers/linear_solver/eigen/precon_type/text()">DIAGONAL</replace>
+</OpenGeoSysProjectDiff>
diff --git a/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_fixed_power.prj b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_fixed_power.prj
new file mode 100644
index 0000000000000000000000000000000000000000..99a3fc2be90cb1fc8f88b2fa88724f714b4923a7
--- /dev/null
+++ b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_fixed_power.prj
@@ -0,0 +1,231 @@
+<?xml version="1.0" encoding="ISO-8859-1"?>
+<OpenGeoSysProject>
+ <meshes>
+ <mesh>sandwich.vtu</mesh>
+ </meshes>
+ <processes>
+ <process>
+ <name>HeatTransportBHE</name>
+ <type>HEAT_TRANSPORT_BHE</type>
+ <integration_order>2</integration_order>
+ <process_variables>
+ <process_variable>temperature_soil</process_variable>
+ <process_variable>temperature_BHE1</process_variable>
+ </process_variables>
+ <borehole_heat_exchangers>
+ <borehole_heat_exchanger>
+ <type>1U</type>
+ <flow_and_temperature_control>
+ <type>FixedPowerConstantFlow</type>
+ <flow_rate>1e-3</flow_rate>
+ <power>10</power>
+ </flow_and_temperature_control>
+ <borehole>
+ <length>3.0</length>
+ <diameter>0.2</diameter>
+ </borehole>
+ <grout>
+ <density>2190.0</density>
+ <porosity>0.0</porosity>
+ <specific_heat_capacity>1781</specific_heat_capacity>
+ <thermal_conductivity>1.4</thermal_conductivity>
+ </grout>
+ <pipes>
+ <inlet>
+ <diameter> 0.02</diameter>
+ <wall_thickness>0.003</wall_thickness>
+ <wall_thermal_conductivity>0.43</wall_thermal_conductivity>
+ </inlet>
+ <outlet>
+ <diameter>0.02</diameter>
+ <wall_thickness>0.003</wall_thickness>
+ <wall_thermal_conductivity>0.43</wall_thermal_conductivity>
+ </outlet>
+ <distance_between_pipes>0.075</distance_between_pipes>
+ <longitudinal_dispersion_length>0.001</longitudinal_dispersion_length>
+ </pipes>
+ <refrigerant>
+ <density>1000</density>
+ <viscosity>0.005</viscosity>
+ <specific_heat_capacity>4000</specific_heat_capacity>
+ <thermal_conductivity>0.5</thermal_conductivity>
+ <reference_temperature>293.15</reference_temperature>
+ </refrigerant>
+ </borehole_heat_exchanger>
+ </borehole_heat_exchangers>
+ </process>
+ </processes>
+ <media>
+ <medium id="0">
+ <phases>
+ <phase>
+ <type>AqueousLiquid</type>
+ <properties>
+ <property>
+ <name>phase_velocity</name>
+ <type>Constant</type>
+ <value>0 0 0</value>
+ </property>
+ <property>
+ <name>specific_heat_capacity</name>
+ <type>Constant</type>
+ <value>3795</value>
+ </property>
+ <property>
+ <name>density</name>
+ <type>Constant</type>
+ <value>1000</value>
+ </property>
+ </properties>
+ </phase>
+ <phase>
+ <type>Solid</type>
+ <properties>
+ <property>
+ <name>specific_heat_capacity</name>
+ <type>Constant</type>
+ <value>2000</value>
+ </property>
+ <property>
+ <name>density</name>
+ <type>Constant</type>
+ <value>1000</value>
+ </property>
+ </properties>
+ </phase>
+ <phase>
+ <type>Gas</type>
+ <properties>
+ <property>
+ <name>specific_heat_capacity</name>
+ <type>Constant</type>
+ <value>1000</value>
+ </property>
+ <property>
+ <name>density</name>
+ <type>Constant</type>
+ <value>2500</value>
+ </property>
+ </properties>
+ </phase>
+ </phases>
+ <properties>
+ <property>
+ <name>porosity</name>
+ <type>Constant</type>
+ <value>0</value>
+ </property>
+ <property>
+ <name>thermal_conductivity</name>
+ <type>Constant</type>
+ <value>2.0</value>
+ </property>
+ <property>
+ <name>thermal_longitudinal_dispersivity</name>
+ <type>Constant</type>
+ <value>5</value>
+ </property>
+ <property>
+ <name>thermal_transversal_dispersivity</name>
+ <type>Constant</type>
+ <value>0.5</value>
+ </property>
+ </properties>
+ </medium>
+ </media>
+ <time_loop>
+ <processes>
+ <process ref="HeatTransportBHE">
+ <nonlinear_solver>basic_picard</nonlinear_solver>
+ <convergence_criterion>
+ <type>DeltaX</type>
+ <norm_type>NORM2</norm_type>
+ <reltol>1.e-6</reltol>
+ </convergence_criterion>
+ <time_discretization><type>BackwardEuler</type></time_discretization>
+ <time_stepping>
+ <type>FixedTimeStepping</type>
+ <t_initial> 0.0 </t_initial>
+ <t_end> 600 </t_end>
+ <timesteps>
+ <pair>
+ <repeat>10</repeat>
+ <delta_t>60</delta_t>
+ </pair>
+ </timesteps>
+ </time_stepping>
+ </process>
+ </processes>
+ <output>
+ <type>VTK</type>
+ <prefix>sandwich_fixed_power</prefix>
+ <timesteps>
+ <pair>
+ <repeat> 1 </repeat>
+ <each_steps> 1 </each_steps>
+ </pair>
+ </timesteps>
+ <variables>
+ <variable>temperature_soil</variable>
+ <variable>temperature_BHE1</variable>
+ </variables>
+ </output>
+ </time_loop>
+ <parameters>
+ <parameter>
+ <name>T0</name>
+ <type>Constant</type>
+ <value>293.15</value>
+ </parameter>
+ <parameter>
+ <name>T_soil_imposed</name>
+ <type>Constant</type>
+ <value>273.15</value>
+ </parameter>
+ <parameter>
+ <name>T0_BHE</name>
+ <type>Constant</type>
+ <values>293.15 293.15 293.15 293.15</values>
+ </parameter>
+ </parameters>
+ <process_variables>
+ <process_variable>
+ <name>temperature_soil</name>
+ <components>1</components>
+ <order>1</order>
+ <initial_condition>T0</initial_condition>
+ <boundary_conditions>
+ </boundary_conditions>
+ </process_variable>
+ <process_variable>
+ <name>temperature_BHE1</name>
+ <components>4</components>
+ <order>1</order>
+ <initial_condition>T0_BHE</initial_condition>
+ </process_variable>
+ </process_variables>
+ <nonlinear_solvers>
+ <nonlinear_solver>
+ <name>basic_picard</name>
+ <type>Picard</type>
+ <max_iter>100</max_iter>
+ <linear_solver>general_linear_solver</linear_solver>
+ </nonlinear_solver>
+ </nonlinear_solvers>
+ <linear_solvers>
+ <linear_solver>
+ <name>general_linear_solver</name>
+ <lis>-i cg -p jacobi -tol 1e-16 -maxiter 10000</lis>
+ <eigen>
+ <solver_type>BiCGSTAB</solver_type>
+ <precon_type>ILUT</precon_type>
+ <max_iteration_step>1000</max_iteration_step>
+ <error_tolerance>1e-16</error_tolerance>
+ </eigen>
+ <petsc>
+ <prefix>gw</prefix>
+ <parameters>-gw_ksp_type cg -gw_pc_type bjacobi -gw_ksp_rtol 1e-16 -gw_ksp_max_it 10000</parameters>
+ </petsc>
+ </linear_solver>
+ </linear_solvers>
+</OpenGeoSysProject>
diff --git a/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_fixed_power_algebraicBC.xml b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_fixed_power_algebraicBC.xml
new file mode 100644
index 0000000000000000000000000000000000000000..a39dedc70acce2f377f07535b160a2d09cf2153b
--- /dev/null
+++ b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_fixed_power_algebraicBC.xml
@@ -0,0 +1,15 @@
+<?xml version="1.0" encoding="ISO-8859-1"?>
+<OpenGeoSysProjectDiff base_file="sandwich_fixed_power.prj">
+ <add sel="/*/processes/process">
+ <use_algebraic_bc>true</use_algebraic_bc>
+ </add>
+ <add sel="/*/processes/process">
+ <weighting_factor>100</weighting_factor>
+ </add>
+ <add sel="/*/processes/process">
+ <linear>true</linear>
+ </add>
+ <replace sel="/*/time_loop/output/prefix/text()">sandwich_fixed_power_algebraic_bc</replace>
+ <replace sel="/*/linear_solvers/linear_solver/eigen/solver_type/text()">CG</replace>
+ <replace sel="/*/linear_solvers/linear_solver/eigen/precon_type/text()">DIAGONAL</replace>
+</OpenGeoSysProjectDiff>
diff --git a/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_fixed_power_ts_10_t_600.000000.vtu b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_fixed_power_ts_10_t_600.000000.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..58c13e2dc7c7d914e314bf50b1fc2f86c822602a
--- /dev/null
+++ b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_fixed_power_ts_10_t_600.000000.vtu
@@ -0,0 +1,27 @@
+<?xml version="1.0"?>
+<VTKFile type="UnstructuredGrid" version="1.0" byte_order="LittleEndian" header_type="UInt64" compressor="vtkZLibDataCompressor">
+ <UnstructuredGrid>
+ <FieldData>
+ <DataArray type="Int8" Name="OGS_VERSION" NumberOfTuples="5" format="appended" RangeMin="46" RangeMax="54" offset="0" />
+ </FieldData>
+ <Piece NumberOfPoints="12" NumberOfCells="6" >
+ <PointData>
+ <DataArray type="Float64" Name="temperature_BHE1" NumberOfComponents="4" format="appended" RangeMin="0" RangeMax="586.50506585" offset="64" />
+ <DataArray type="Float64" Name="temperature_soil" format="appended" RangeMin="293.14993086" RangeMax="293.15020812" offset="276" />
+ </PointData>
+ <CellData>
+ </CellData>
+ <Points>
+ <DataArray type="Float64" Name="Points" NumberOfComponents="3" format="appended" RangeMin="0" RangeMax="3.1622776602" offset="400" />
+ </Points>
+ <Cells>
+ <DataArray type="Int64" Name="connectivity" format="appended" RangeMin="" RangeMax="" offset="512" />
+ <DataArray type="Int64" Name="offsets" format="appended" RangeMin="" RangeMax="" offset="620" />
+ <DataArray type="UInt8" Name="types" format="appended" RangeMin="" RangeMax="" offset="696" />
+ </Cells>
+ </Piece>
+ </UnstructuredGrid>
+ <AppendedData encoding="base64">
+ _AQAAAAAAAAAAgAAAAAAAAAUAAAAAAAAADQAAAAAAAAA=eF4z0zPVMwIAAvgA+g==AQAAAAAAAAAAgAAAAAAAAIABAAAAAAAAfAAAAAAAAAA=eF6bIJF7e0ZokUPxa/nF04C0/+rfT2uCihxcKlZfBtEMFIL1zf+vzASa634pe9VUIK302PI9yNzzqytOUcP85ugn52cBzW3682TDFCD9KLbjB9hcW9uD1DDfuPLrydlAc0/nrNsxGUjvVDzHWAs09+UHwR3UMB8Aj3g/IQ==AQAAAAAAAAAAgAAAAAAAAGAAAAAAAAAAOwAAAAAAAAA=eF5Li1ngkB5U5HD//3PZNCD9DEp/WAsRv/cEwr8IpZutFoLFNbZB+EZQOo9hEVi8oQTCr4XSAJ2dKw8=AQAAAAAAAAAAgAAAAAAAACABAAAAAAAAMQAAAAAAAAA=eF5jYCAGfLDHz8cFcOlDp3GphwMHVC5cnQMaHwYcGLACDjRxmD6YOLo5HA4ArTYOng==AQAAAAAAAAAAgAAAAAAAAMAAAAAAAAAAMAAAAAAAAAA=eF5jYIAARijNBKWZoTQLlGYlIM4GpdmhNAcBcU4ozQWludHE2dBoZjQaBgAsKAB/AQAAAAAAAAAAgAAAAAAAADAAAAAAAAAAGAAAAAAAAAA=eF5jY4AAHigtBKVFoLQYlJaA0gAJcABnAQAAAAAAAAAAgAAAAAAAAAYAAAAAAAAADgAAAAAAAAA=eF7j5eVlZmYGAADbADE=
+ </AppendedData>
+</VTKFile>
diff --git a/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_newton.xml b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_newton.xml
new file mode 100644
index 0000000000000000000000000000000000000000..073327ee73b4d5bb328c30f241ae0f52f3b3ab37
--- /dev/null
+++ b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_newton.xml
@@ -0,0 +1,9 @@
+<?xml version="1.0" encoding="ISO-8859-1"?>
+<OpenGeoSysProjectDiff base_file="sandwich.prj">
+ <replace sel="/*/time_loop/processes/process/nonlinear_solver/text()">newton_nonlinear_solver</replace>
+ <replace sel="/*/time_loop/processes/process/convergence_criterion/type/text()">Residual</replace>
+ <replace sel="/*/time_loop/output/prefix/text()">sandwich_newton</replace>
+ <replace sel="/*/nonlinear_solvers/nonlinear_solver/name/text()">newton_nonlinear_solver</replace>
+ <replace sel="/*/nonlinear_solvers/nonlinear_solver/type/text()">Newton</replace>
+ <replace sel="/*/nonlinear_solvers/nonlinear_solver/max_iter/text()">10</replace>
+</OpenGeoSysProjectDiff>
diff --git a/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_ts_10_t_600.000000.vtu b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_ts_10_t_600.000000.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..cbe4bc0f70b43383e96ad4336a98490d7eedf3fd
--- /dev/null
+++ b/Tests/Data/Parabolic/T/3D_BHE_Sandwich/sandwich_ts_10_t_600.000000.vtu
@@ -0,0 +1,27 @@
+<?xml version="1.0"?>
+<VTKFile type="UnstructuredGrid" version="1.0" byte_order="LittleEndian" header_type="UInt64" compressor="vtkZLibDataCompressor">
+ <UnstructuredGrid>
+ <FieldData>
+ <DataArray type="Int8" Name="OGS_VERSION" NumberOfTuples="5" format="appended" RangeMin="46" RangeMax="54" offset="0" />
+ </FieldData>
+ <Piece NumberOfPoints="12" NumberOfCells="6" >
+ <PointData>
+ <DataArray type="Float64" Name="temperature_BHE1" NumberOfComponents="4" format="appended" RangeMin="0" RangeMax="565.91439316" offset="64" />
+ <DataArray type="Float64" Name="temperature_soil" format="appended" RangeMin="293.11463784" RangeMax="293.16174145" offset="276" />
+ </PointData>
+ <CellData>
+ </CellData>
+ <Points>
+ <DataArray type="Float64" Name="Points" NumberOfComponents="3" format="appended" RangeMin="0" RangeMax="3.1622776602" offset="408" />
+ </Points>
+ <Cells>
+ <DataArray type="Int64" Name="connectivity" format="appended" RangeMin="" RangeMax="" offset="520" />
+ <DataArray type="Int64" Name="offsets" format="appended" RangeMin="" RangeMax="" offset="628" />
+ <DataArray type="UInt8" Name="types" format="appended" RangeMin="" RangeMax="" offset="704" />
+ </Cells>
+ </Piece>
+ </UnstructuredGrid>
+ <AppendedData encoding="base64">
+ _AQAAAAAAAAAAgAAAAAAAAAUAAAAAAAAADQAAAAAAAAA=eF4z0zPVMwIAAvgA+g==AQAAAAAAAAAAgAAAAAAAAIABAAAAAAAAfAAAAAAAAAA=eF47P8/gU6hwocN5KO1zw/+ghFuRw2GPokMgmoFC8KCyx5MNaG5NpmrHUiBtL3+6WARo7jKfIj5ZKphfYvrGdJtQoYOvFHPSV6D53//tlBMAmisify1OkQrmp4EA0PyyfwYXXUQKHVj8hc5xA82duWjrclUqmA8AKOsxSw==AQAAAAAAAAAAgAAAAAAAAGAAAAAAAAAAQQAAAAAAAAA=eF5b583RfzWwyGH6e9a6aUFFDj1QmkfxZx9IPO8bF5hfAKX3V58Ci5vmSIH5elBaXqQWLP5S2ATMvw6lAf4DKNs=AQAAAAAAAAAAgAAAAAAAACABAAAAAAAAMQAAAAAAAAA=eF5jYCAGfLDHz8cFcOlDp3GphwMHVC5cnQMaHwYcGLACDjRxmD6YOLo5HA4ArTYOng==AQAAAAAAAAAAgAAAAAAAAMAAAAAAAAAAMAAAAAAAAAA=eF5jYIAARijNBKWZoTQLlGYlIM4GpdmhNAcBcU4ozQWludHE2dBoZjQaBgAsKAB/AQAAAAAAAAAAgAAAAAAAADAAAAAAAAAAGAAAAAAAAAA=eF5jY4AAHigtBKVFoLQYlJaA0gAJcABnAQAAAAAAAAAAgAAAAAAAAAYAAAAAAAAADgAAAAAAAAA=eF7j5eVlZmYGAADbADE=
+ </AppendedData>
+</VTKFile>
diff --git a/Tests/Data/Parabolic/T/3D_Beier_sandbox/beier_sandbox_algebraicBC.xml b/Tests/Data/Parabolic/T/3D_Beier_sandbox/beier_sandbox_algebraicBC.xml
new file mode 100644
index 0000000000000000000000000000000000000000..b0ee68e555c1a276d07fd8defa128662dfd36d87
--- /dev/null
+++ b/Tests/Data/Parabolic/T/3D_Beier_sandbox/beier_sandbox_algebraicBC.xml
@@ -0,0 +1,15 @@
+<?xml version="1.0" encoding="ISO-8859-1"?>
+<OpenGeoSysProjectDiff base_file="beier_sandbox.prj">
+ <add sel="/*/processes/process">
+ <use_algebraic_bc>true</use_algebraic_bc>
+ </add>
+ <add sel="/*/processes/process">
+ <weighting_factor>100</weighting_factor>
+ </add>
+ <add sel="/*/processes/process">
+ <linear>true</linear>
+ </add>
+ <replace sel="/*/time_loop/output/prefix/text()">beier_sandbox_algebraic_bc</replace>
+ <replace sel="/*/linear_solvers/linear_solver/eigen/solver_type/text()">CG</replace>
+ <replace sel="/*/linear_solvers/linear_solver/eigen/precon_type/text()">DIAGONAL</replace>
+</OpenGeoSysProjectDiff>
diff --git a/Tests/Data/Parabolic/T/3D_Beier_sandbox/fixed_power_constant_flow_algebraicBC.xml b/Tests/Data/Parabolic/T/3D_Beier_sandbox/fixed_power_constant_flow_algebraicBC.xml
new file mode 100644
index 0000000000000000000000000000000000000000..963f9fa17a1dae41c122ba7803a2b90a7e7b305a
--- /dev/null
+++ b/Tests/Data/Parabolic/T/3D_Beier_sandbox/fixed_power_constant_flow_algebraicBC.xml
@@ -0,0 +1,15 @@
+<?xml version="1.0" encoding="ISO-8859-1"?>
+<OpenGeoSysProjectDiff base_file="fixed_power_constant_flow.prj">
+ <add sel="/*/processes/process">
+ <use_algebraic_bc>true</use_algebraic_bc>
+ </add>
+ <add sel="/*/processes/process">
+ <weighting_factor>100</weighting_factor>
+ </add>
+ <add sel="/*/processes/process">
+ <linear>true</linear>
+ </add>
+ <replace sel="/*/time_loop/output/prefix/text()">fixed_power_constant_flow_algebraic_bc</replace>
+ <replace sel="/*/linear_solvers/linear_solver/eigen/solver_type/text()">CG</replace>
+ <replace sel="/*/linear_solvers/linear_solver/eigen/precon_type/text()">DIAGONAL</replace>
+</OpenGeoSysProjectDiff>
diff --git a/Tests/NumLib/ODEs.h b/Tests/NumLib/ODEs.h
index 98d7d61b869e772dc5bf241b436ad5f028c462f5..2f3d2e6623482cc492290a945c7030123a5908db 100644
--- a/Tests/NumLib/ODEs.h
+++ b/Tests/NumLib/ODEs.h
@@ -115,6 +115,8 @@ public:
return true;
}
+ bool requiresNormalization() const override { return false; }
+
std::size_t const N = 2;
};
@@ -205,6 +207,8 @@ public:
return false;
}
+ bool requiresNormalization() const override { return false; }
+
std::size_t const N = 1;
};
@@ -333,6 +337,8 @@ public:
return false;
}
+ bool requiresNormalization() const override { return false; }
+
std::size_t const N = 2;
};