Merge pull request #2818 from ChaofanChen/BHE_1P

Added single pipe flow process based on the HeatTransportBHE

Documentation/ProjectFile/prj/processes/process/HEAT_TRANSPORT_BHE/borehole_heat_exchangers/borehole_heat_exchanger/flow_and_temperature_control/BuildingPowerCurveConstantFlow/c_BuildingPowerCurveConstantFlow.md

+It is the BHE boundary condition of various heat load from buildings. Meanwhile, the COP of the heat pump is intergrated and used for calculating the transit heat load on the BHE.

Documentation/ProjectFile/prj/processes/process/HEAT_TRANSPORT_BHE/borehole_heat_exchangers/borehole_heat_exchanger/flow_and_temperature_control/BuildingPowerCurveConstantFlow/t_cop_heating_curve.md

+It is the COP curve of the heat pump at different outflow temperature from the BHE.

Documentation/ProjectFile/prj/processes/process/HEAT_TRANSPORT_BHE/borehole_heat_exchangers/borehole_heat_exchanger/flow_and_temperature_control/BuildingPowerCurveConstantFlow/t_flow_rate.md

+It is the flow rate of the BHE and it is a constant value in this BC.

Documentation/ProjectFile/prj/processes/process/HEAT_TRANSPORT_BHE/borehole_heat_exchangers/borehole_heat_exchanger/flow_and_temperature_control/BuildingPowerCurveConstantFlow/t_power_curve.md

+It is the building heat load curve.

ProcessLib/HeatTransportBHE/BHE/BHECommonCoaxial.cpp

@@ -87,7 +87,8 @@ BHECommonCoaxial::pipeHeatConductions() const
 }
 
 std::array<Eigen::Vector3d, BHECommonCoaxial::number_of_unknowns>
-BHECommonCoaxial::pipeAdvectionVectors() const
+BHECommonCoaxial::pipeAdvectionVectors(
+    Eigen::Vector3d const& /*elem_direction*/) const
 {
     double const rho_r = refrigerant.density;
     double const Cp_r = refrigerant.specific_heat_capacity;
@@ -130,6 +131,26 @@ BHECommonCoaxial::calcThermalResistances(double const Nu_inner_pipe,
     return getThermalResistances(R_gs, R_ff, R_fg);
 }
 
+std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>
+BHECommonCoaxial::getBHEInflowDirichletBCNodesAndComponents(
+    std::size_t const top_node_id,
+    std::size_t const /*bottom_node_id*/,
+    int const in_component_id) const
+{
+    return {std::make_pair(top_node_id, in_component_id),
+            std::make_pair(top_node_id, in_component_id + 1)};
+}
+
+std::optional<
+    std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>>
+BHECommonCoaxial::getBHEBottomDirichletBCNodesAndComponents(
+    std::size_t const bottom_node_id, int const in_component_id,
+    int const out_component_id) const
+{
+    return {{std::make_pair(bottom_node_id, in_component_id),
+             std::make_pair(bottom_node_id, out_component_id)}};
+}
+
 void BHECommonCoaxial::updateHeatTransferCoefficients(double const flow_rate)
 {
     auto const tm_flow_properties_annulus =

ProcessLib/HeatTransportBHE/BHE/BHECommonCoaxial.h

@@ -10,6 +10,7 @@
 #pragma once
 
 #include <Eigen/Eigen>
+#include <optional>
 #include "BHECommon.h"
 #include "FlowAndTemperatureControl.h"
 #include "PipeConfigurationCoaxial.h"
@@ -48,10 +49,22 @@ public:
     static constexpr std::pair<int, int> inflow_outflow_bc_component_ids[] = {
         {0, 1}};
 
+    std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>
+    getBHEInflowDirichletBCNodesAndComponents(
+        std::size_t const top_node_id,
+        std::size_t const /*bottom_node_id*/,
+        int const in_component_id) const;
+
+    std::optional<
+        std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>>
+    getBHEBottomDirichletBCNodesAndComponents(std::size_t const bottom_node_id,
+                                              int const in_component_id,
+                                              int const out_component_id) const;
+
     std::array<double, number_of_unknowns> pipeHeatConductions() const;
 
-    std::array<Eigen::Vector3d, number_of_unknowns> pipeAdvectionVectors()
-        const;
+    std::array<Eigen::Vector3d, number_of_unknowns> pipeAdvectionVectors(
+        Eigen::Vector3d const& /*elem_direction*/) const;
 
     double cross_section_area_inner_pipe, cross_section_area_annulus,
         cross_section_area_grout;

ProcessLib/HeatTransportBHE/BHE/BHETypes.h

@@ -11,6 +11,7 @@
 #pragma once
 
 #include <variant>
+#include "BHE_1P.h"
 #include "BHE_1U.h"
 #include "BHE_2U.h"
 #include "BHE_CXA.h"
@@ -22,7 +23,7 @@ namespace HeatTransportBHE
 {
 namespace BHE
 {
-using BHETypes = std::variant<BHE_1U, BHE_CXA, BHE_CXC, BHE_2U>;
+using BHETypes = std::variant<BHE_1U, BHE_CXA, BHE_CXC, BHE_2U, BHE_1P>;
 }  // end of namespace BHE
 }  // end of namespace HeatTransportBHE
 }  // end of namespace ProcessLib

ProcessLib/HeatTransportBHE/BHE/BHE_1P.cpp

+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2020, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ */
+
+#include "BHE_1P.h"
+
+#include <boost/math/constants/constants.hpp>
+#include "FlowAndTemperatureControl.h"
+#include "Physics.h"
+#include "ThermoMechanicalFlowProperties.h"
+
+namespace ProcessLib
+{
+namespace HeatTransportBHE
+{
+namespace BHE
+{
+BHE_1P::BHE_1P(BoreholeGeometry const& borehole,
+               RefrigerantProperties const& refrigerant,
+               GroutParameters const& grout,
+               FlowAndTemperatureControl const& flowAndTemperatureControl,
+               PipeConfiguration1PType const& pipes,
+               bool const use_python_bcs)
+    : BHECommon{borehole, refrigerant, grout, flowAndTemperatureControl,
+                use_python_bcs},
+      _pipe(pipes)
+{
+    _thermal_resistances.fill(std::numeric_limits<double>::quiet_NaN());
+
+    // Initialize thermal resistances.
+    auto values = visit(
+        [&](auto const& control) {
+            return control(refrigerant.reference_temperature,
+                           0. /* initial time */);
+        },
+        flowAndTemperatureControl);
+    updateHeatTransferCoefficients(values.flow_rate);
+}
+
+std::array<double, BHE_1P::number_of_unknowns> BHE_1P::pipeHeatCapacities()
+    const
+{
+    double const rho_r = refrigerant.density;
+    double const specific_heat_capacity = refrigerant.specific_heat_capacity;
+    double const rho_g = grout.rho_g;
+    double const porosity_g = grout.porosity_g;
+    double const heat_cap_g = grout.heat_cap_g;
+
+    return {{
+        /*pipe*/ rho_r * specific_heat_capacity,
+        /*grout*/ (1.0 - porosity_g) * rho_g * heat_cap_g,
+    }};
+}
+
+std::array<double, BHE_1P::number_of_unknowns> BHE_1P::pipeHeatConductions()
+    const
+{
+    double const lambda_r = refrigerant.thermal_conductivity;
+    double const rho_r = refrigerant.density;
+    double const Cp_r = refrigerant.specific_heat_capacity;
+    double const alpha_L = _pipe.longitudinal_dispersion_length;
+    double const porosity_g = grout.porosity_g;
+    double const lambda_g = grout.lambda_g;
+
+    // Here we calculate the laplace coefficients in the governing
+    // equations of the BHE.
+    return {{
+        // pipe, Eq. 19
+        (lambda_r + rho_r * Cp_r * alpha_L * _flow_velocity),
+        // grout, Eq. 21
+        (1.0 - porosity_g) * lambda_g,
+    }};
+}
+
+std::array<Eigen::Vector3d, BHE_1P::number_of_unknowns>
+BHE_1P::pipeAdvectionVectors(Eigen::Vector3d const& elem_direction) const
+{
+    double const& rho_r = refrigerant.density;
+    double const& Cp_r = refrigerant.specific_heat_capacity;
+    Eigen::Vector3d adv_vector = rho_r * Cp_r * _flow_velocity * elem_direction;
+
+    return {// pipe, Eq. 19
+            adv_vector,
+            // grout, Eq. 21
+            {0, 0, 0}};
+}
+
+double BHE_1P::compute_R_gs(double const chi, double const R_g)
+{
+    return (1 - chi) * R_g;
+}
+
+void BHE_1P::updateHeatTransferCoefficients(double const flow_rate)
+
+{
+    auto const tm_flow_properties = calculateThermoMechanicalFlowPropertiesPipe(
+        _pipe.single_pipe, borehole_geometry.length, refrigerant, flow_rate);
+
+    _flow_velocity = tm_flow_properties.velocity;
+    _thermal_resistances =
+        calcThermalResistances(tm_flow_properties.nusselt_number);
+}
+
+// Nu is the Nusselt number.
+std::array<double, BHE_1P::number_of_unknowns> BHE_1P::calcThermalResistances(
+    double const Nu)
+{
+    constexpr double pi = boost::math::constants::pi<double>();
+
+    double const lambda_r = refrigerant.thermal_conductivity;
+    double const lambda_g = grout.lambda_g;
+    double const lambda_p = _pipe.single_pipe.wall_thermal_conductivity;
+
+    // thermal resistances due to advective flow of refrigerant in the pipe
+    double const R_adv_i1 = 1.0 / (Nu * lambda_r * pi);
+
+    // thermal resistance due to thermal conductivity of the pipe wall material
+    double const R_con_a = std::log(_pipe.single_pipe.outsideDiameter() /
+                                    _pipe.single_pipe.diameter) /
+                           (2.0 * pi * lambda_p);
+
+    // thermal resistances of the grout
+    double const D = borehole_geometry.diameter;
+    double const pipe_outside_diameter = _pipe.single_pipe.outsideDiameter();
+
+    double const chi = std::log(std::sqrt(D * D + pipe_outside_diameter *
+                                                      pipe_outside_diameter) /
+                                std::sqrt(2) / pipe_outside_diameter) /
+                       std::log(D / pipe_outside_diameter);
+    double const R_g =
+        std::log(D / pipe_outside_diameter) / 2 / (pi * lambda_g);
+
+    double const R_con_b = chi * R_g;
+
+    // thermal resistances due to grout-soil exchange
+    double const R_gs = compute_R_gs(chi, R_g);
+
+    // Eq. 29 and 30
+    double const R_fg = R_adv_i1 + R_con_a + R_con_b;
+
+    return {{R_fg, R_gs}};
+}
+
+std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>
+BHE_1P::getBHEInflowDirichletBCNodesAndComponents(
+    std::size_t const top_node_id,
+    std::size_t const bottom_node_id,
+    int const in_component_id) const
+{
+    return {std::make_pair(top_node_id, in_component_id),
+            std::make_pair(bottom_node_id, in_component_id)};
+}
+
+std::optional<
+    std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>>
+BHE_1P::getBHEBottomDirichletBCNodesAndComponents(
+    std::size_t const /*bottom_node_id*/,
+    int const /*in_component_id*/,
+    int const /*out_component_id*/) const
+{
+    return {};
+}
+
+std::array<double, BHE_1P::number_of_unknowns> BHE_1P::crossSectionAreas() const
+{
+    return {{_pipe.single_pipe.area(),
+             borehole_geometry.area() - _pipe.single_pipe.outsideArea()}};
+}
+
+double BHE_1P::updateFlowRateAndTemperature(double const T_out,
+                                            double const current_time)
+{
+    auto values =
+        visit([&](auto const& control) { return control(T_out, current_time); },
+              flowAndTemperatureControl);
+    updateHeatTransferCoefficients(values.flow_rate);
+    return values.temperature;
+}
+}  // namespace BHE
+}  // namespace HeatTransportBHE
+}  // namespace ProcessLib

ProcessLib/HeatTransportBHE/BHE/BHE_1P.h

+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2020, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ */
+
+#pragma once
+
+#include <Eigen/Eigen>
+#include <optional>
+
+#include "BaseLib/Error.h"
+
+#include "BHECommon.h"
+#include "FlowAndTemperatureControl.h"
+#include "PipeConfiguration1PType.h"
+
+namespace ProcessLib
+{
+namespace HeatTransportBHE
+{
+namespace BHE
+{
+/**
+ * The BHE_1P class is the realization of single-pipe type of Borehole Heate
+ * Exchanger. In this class, the pipe heat capacity, pipe heat conductiion, pie
+ * advection vectors are intialized according to the geometry of the single-pipe
+ * type of BHE. For this type of BHE, 2 primary unknowns are assigned on the 1D
+ * BHE elements. They are the temperature in the pipe T_p, and temperature of
+ * the grout zone sorrounding the single pipe T_g. These two primary varaibles
+ * are solved according to heat convection and conduction equations on the pipes
+ * and also in the grout zones. The interaction of the 1P type of BHE and the
+ * sorrounding soil is regulated through the thermal resistance values, which
+ * are calculated specifically during the initialization of the class.
+ */
+class BHE_1P final : public BHECommon
+{
+public:
+    BHE_1P(BoreholeGeometry const& borehole,
+           RefrigerantProperties const& refrigerant,
+           GroutParameters const& grout,
+           FlowAndTemperatureControl const& flowAndTemperatureControl,
+           PipeConfiguration1PType const& pipes,
+           bool const use_python_bcs);
+
+    static constexpr int number_of_unknowns = 2;
+    static constexpr int number_of_grout_zones = 1;
+
+    std::array<double, number_of_unknowns> pipeHeatCapacities() const;
+
+    std::array<double, number_of_unknowns> pipeHeatConductions() const;
+
+    std::array<Eigen::Vector3d, number_of_unknowns> pipeAdvectionVectors(
+        Eigen::Vector3d const& elem_direction) const;
+
+    template <int NPoints,
+              typename SingleUnknownMatrixType,
+              typename RMatrixType,
+              typename RPiSMatrixType,
+              typename RSMatrixType>
+    void assembleRMatrices(
+        int const idx_bhe_unknowns,
+        Eigen::MatrixBase<SingleUnknownMatrixType> const& matBHE_loc_R,
+        Eigen::MatrixBase<RMatrixType>& R_matrix,
+        Eigen::MatrixBase<RPiSMatrixType>& R_pi_s_matrix,
+        Eigen::MatrixBase<RSMatrixType>& R_s_matrix) const
+    {
+        // Here we are looping over two resistance terms
+        // First PHI_fg is the resistance between pipe and grout
+        // Second PHI_gs is the resistance between grout and soil
+        switch (idx_bhe_unknowns)
+        {
+            case 0:  // PHI_fg
+                R_matrix.block(0, NPoints, NPoints, NPoints) +=
+                    -1.0 * matBHE_loc_R;
+                R_matrix.block(NPoints, 0, NPoints, NPoints) +=
+                    -1.0 * matBHE_loc_R;
+
+                R_matrix.block(0, 0, NPoints, NPoints) +=
+                    matBHE_loc_R;  // K_i/o
+                R_matrix.block(NPoints, NPoints, NPoints, NPoints) +=
+                    matBHE_loc_R;  // K_fg
+                return;
+            case 1:  // PHI_gs
+                R_s_matrix += matBHE_loc_R;
+
+                R_pi_s_matrix.block(NPoints, 0, NPoints, NPoints) +=
+                    -1.0 * matBHE_loc_R;
+
+                R_matrix.block(NPoints, NPoints, NPoints, NPoints) +=
+                    matBHE_loc_R;  // K_fg
+                return;
+            default:
+                OGS_FATAL(
+                    "Error!!! In the function BHE_1P::assembleRMatrices, "
+                    "the index of bhe resistance term is out of range! ");
+        }
+    }
+
+    /// Return the inflow temperature for the boundary condition.
+    double updateFlowRateAndTemperature(double T_out, double current_time);
+
+    double thermalResistance(int const unknown_index) const
+    {
+        return _thermal_resistances[unknown_index];
+    }
+
+    static constexpr std::pair<int, int> inflow_outflow_bc_component_ids[] = {
+        {0, 1}};
+
+    std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>
+    getBHEInflowDirichletBCNodesAndComponents(std::size_t const top_node_id,
+                                              std::size_t const bottom_node_id,
+                                              int const in_component_id) const;
+
+    std::optional<
+        std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>>
+    getBHEBottomDirichletBCNodesAndComponents(
+        std::size_t const /*bottom_node_id*/,
+        int const /*in_component_id*/,
+        int const /*out_component_id*/) const;
+
+public:
+    std::array<double, number_of_unknowns> crossSectionAreas() const;
+
+    void updateHeatTransferCoefficients(double const flow_rate);
+
+protected:
+    PipeConfiguration1PType const _pipe;
+
+    /// Flow velocity inside the pipes. Depends on the flow_rate.
+    double _flow_velocity = std::numeric_limits<double>::quiet_NaN();
+
+private:
+    std::array<double, number_of_unknowns> calcThermalResistances(
+        double const Nu);
+
+    double compute_R_gs(double const chi, double const R_g);
+
+private:
+    /// PHI_fg, PHI_gs;
+    /// Here we store the thermal resistances needed for computation of the heat
+    /// exchange coefficients in the governing equations of BHE.
+    std::array<double, number_of_unknowns> _thermal_resistances;
+};
+}  // namespace BHE
+}  // namespace HeatTransportBHE
+}  // namespace ProcessLib

ProcessLib/HeatTransportBHE/BHE/BHE_1U.cpp

@@ -84,7 +84,7 @@ std::array<double, BHE_1U::number_of_unknowns> BHE_1U::pipeHeatConductions()
 }
 
 std::array<Eigen::Vector3d, BHE_1U::number_of_unknowns>
-BHE_1U::pipeAdvectionVectors() const
+BHE_1U::pipeAdvectionVectors(Eigen::Vector3d const& /*elem_direction*/) const
 {
     double const& rho_r = refrigerant.density;
     double const& Cp_r = refrigerant.specific_heat_capacity;
@@ -228,6 +228,27 @@ std::array<double, BHE_1U::number_of_unknowns> BHE_1U::calcThermalResistances(
     // -------------------------------------------------------------------------
 }
 
+std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>
+BHE_1U::getBHEInflowDirichletBCNodesAndComponents(
+    std::size_t const top_node_id,
+    std::size_t const /*bottom_node_id*/,
+    int const in_component_id) const
+{
+    return {std::make_pair(top_node_id, in_component_id),
+            std::make_pair(top_node_id, in_component_id + 1)};
+}
+
+std::optional<
+    std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>>
+BHE_1U::getBHEBottomDirichletBCNodesAndComponents(
+    std::size_t const bottom_node_id,
+    int const in_component_id,
+    int const out_component_id) const
+{
+    return {{std::make_pair(bottom_node_id, in_component_id),
+             std::make_pair(bottom_node_id, out_component_id)}};
+}
+
 std::array<double, BHE_1U::number_of_unknowns> BHE_1U::crossSectionAreas() const
 {
     return {{_pipes.inlet.area(), _pipes.outlet.area(),

ProcessLib/HeatTransportBHE/BHE/BHE_1U.h

@@ -11,6 +11,7 @@
 #pragma once
 
 #include <Eigen/Eigen>
+#include <optional>
 
 #include "BaseLib/Error.h"
 
@@ -55,8 +56,8 @@ public:
 
     std::array<double, number_of_unknowns> pipeHeatConductions() const;
 
-    std::array<Eigen::Vector3d, number_of_unknowns> pipeAdvectionVectors()
-        const;
+    std::array<Eigen::Vector3d, number_of_unknowns> pipeAdvectionVectors(
+        Eigen::Vector3d const& /*elem_direction*/) const;
 
     template <int NPoints,
               typename SingleUnknownMatrixType,
@@ -138,6 +139,17 @@ public:
     static constexpr std::pair<int, int> inflow_outflow_bc_component_ids[] = {
         {0, 1}};
 
+    std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>
+    getBHEInflowDirichletBCNodesAndComponents(
+        std::size_t const top_node_id,
+        std::size_t const /*bottom_node_id*/,
+        int const in_component_id) const;
+    std::optional<
+        std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>>
+    getBHEBottomDirichletBCNodesAndComponents(std::size_t const bottom_node_id,
+                                              int const in_component_id,
+                                              int const out_component_id) const;
+
 public:
     std::array<double, number_of_unknowns> crossSectionAreas() const;
 

ProcessLib/HeatTransportBHE/BHE/BHE_2U.cpp

@@ -94,7 +94,7 @@ std::array<double, BHE_2U::number_of_unknowns> BHE_2U::pipeHeatConductions()
 }
 
 std::array<Eigen::Vector3d, BHE_2U::number_of_unknowns>
-BHE_2U::pipeAdvectionVectors() const
+BHE_2U::pipeAdvectionVectors(Eigen::Vector3d const& /*elem_direction*/) const
 {
     double const rho_r = refrigerant.density;
     double const Cp_r = refrigerant.specific_heat_capacity;
@@ -258,6 +258,27 @@ std::array<double, BHE_2U::number_of_unknowns> BHE_2U::calcThermalResistances(
     // -------------------------------------------------------------------------
 }
 
+std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>
+BHE_2U::getBHEInflowDirichletBCNodesAndComponents(
+    std::size_t const top_node_id,
+    std::size_t const /*bottom_node_id*/,
+    int const in_component_id) const
+{
+    return {std::make_pair(top_node_id, in_component_id),
+            std::make_pair(top_node_id, in_component_id + 2)};
+}
+
+std::optional<
+    std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>>
+BHE_2U::getBHEBottomDirichletBCNodesAndComponents(
+    std::size_t const bottom_node_id,
+    int const in_component_id,
+    int const out_component_id) const
+{
+    return {{std::make_pair(bottom_node_id, in_component_id),
+             std::make_pair(bottom_node_id, out_component_id)}};
+}
+
 std::array<double, BHE_2U::number_of_unknowns> BHE_2U::crossSectionAreas() const
 {
     return {{

ProcessLib/HeatTransportBHE/BHE/BHE_2U.h

@@ -12,6 +12,7 @@
 
 #include <Eigen/Eigen>
 
+#include <optional>
 #include "BaseLib/Error.h"
 
 #include "BHECommon.h"
@@ -55,8 +56,8 @@ public:
 
     std::array<double, number_of_unknowns> pipeHeatConductions() const;
 
-    std::array<Eigen::Vector3d, number_of_unknowns> pipeAdvectionVectors()
-        const;
+    std::array<Eigen::Vector3d, number_of_unknowns> pipeAdvectionVectors(
+        Eigen::Vector3d const& /*elem_direction*/) const;
 
     template <int NPoints,
               typename SingleUnknownMatrixType,
@@ -200,6 +201,18 @@ public:
     static constexpr std::pair<int, int> inflow_outflow_bc_component_ids[] = {
         {0, 2}, {1, 3}};
 
+    std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>
+    getBHEInflowDirichletBCNodesAndComponents(
+        std::size_t const top_node_id,
+        std::size_t const /*bottom_node_id*/,
+        int const in_component_id) const;
+
+    std::optional<
+        std::array<std::pair<std::size_t /*node_id*/, int /*component*/>, 2>>
+    getBHEBottomDirichletBCNodesAndComponents(std::size_t const bottom_node_id,
+                                              int const in_component_id,
+                                              int const out_component_id) const;
+
 public:
     std::array<double, number_of_unknowns> crossSectionAreas() const;
 

ProcessLib/HeatTransportBHE/BHE/BHE_CXC.h

@@ -12,8 +12,6 @@
 
 #include <Eigen/Eigen>
 
-#include "BaseLib/Error.h"
-
 #include "BHECommonCoaxial.h"
 
 namespace ProcessLib

ProcessLib/HeatTransportBHE/BHE/BuildingPowerCurves.h

+/**
+ * \file
+ *
+ * \copyright
+ * Copyright (c) 2012-2020, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ */
+
+#pragma once
+
+namespace MathLib
+{
+class PiecewiseLinearInterpolation;
+}
+
+namespace ProcessLib
+{
+namespace HeatTransportBHE
+{
+namespace BHE
+{
+struct BuildingPowerCurves
+{
+    MathLib::PiecewiseLinearInterpolation const& power_curve;
+    MathLib::PiecewiseLinearInterpolation const& cop_heating_curve;
+};
+}  // namespace BHE
+}  // namespace HeatTransportBHE
+}  // namespace ProcessLib

ProcessLib/HeatTransportBHE/BHE/CreateBHE1PType.cpp

+/**
+ * \file
+ *
+ * \copyright
+ * Copyright (c) 2012-2020, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ */
+
+#include "BaseLib/ConfigTree.h"
+#include "CreateBHEUType.h"
+
+#include "BHE_1P.h"
+#include "CreateFlowAndTemperatureControl.h"
+
+namespace ProcessLib
+{
+namespace HeatTransportBHE
+{
+namespace BHE
+{
+static std::tuple<BoreholeGeometry,
+                  RefrigerantProperties,
+                  GroutParameters,
+                  FlowAndTemperatureControl,
+                  PipeConfiguration1PType,
+                  bool>
+parseBHE1PTypeConfig(
+    BaseLib::ConfigTree const& config,
+    std::map<std::string,
+             std::unique_ptr<MathLib::PiecewiseLinearInterpolation>> const&
+        curves)
+{
+    // if the BHE is using python boundary condition
+    auto const bhe_if_use_python_bc_conf =
+        //! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__borehole_heat_exchangers__borehole_heat_exchanger__use_bhe_pipe_network}
+        config.getConfigParameter<bool>("use_bhe_pipe_network", false);
+
+    if (bhe_if_use_python_bc_conf)
+    {
+        DBUG("BHE 1P using python boundary conditions.");
+    }
+
+    auto const borehole_geometry =
+        //! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__borehole_heat_exchangers__borehole_heat_exchanger__borehole}
+        createBoreholeGeometry(config.getConfigSubtree("borehole"));
+
+    //! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__borehole_heat_exchangers__borehole_heat_exchanger__pipes}
+    auto const& pipes_config = config.getConfigSubtree("pipes");
+    //! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__borehole_heat_exchangers__borehole_heat_exchanger__pipes__inlet}
+    Pipe const inlet_pipe = createPipe(pipes_config.getConfigSubtree("inlet"));
+
+    const auto pipe_longitudinal_dispersion_length =
+        //! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__borehole_heat_exchangers__borehole_heat_exchanger__pipes__longitudinal_dispersion_length}
+        pipes_config.getConfigParameter<double>(
+            "longitudinal_dispersion_length");
+    PipeConfiguration1PType const pipes{inlet_pipe,
+                                        pipe_longitudinal_dispersion_length};
+
+    //! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__borehole_heat_exchangers__borehole_heat_exchanger__grout}
+    auto const grout = createGroutParameters(config.getConfigSubtree("grout"));
+
+    auto const refrigerant =
+        //! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__borehole_heat_exchangers__borehole_heat_exchanger__refrigerant}
+        createRefrigerantProperties(config.getConfigSubtree("refrigerant"));
+
+    auto const flowAndTemperatureControl = createFlowAndTemperatureControl(
+        //! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__borehole_heat_exchangers__borehole_heat_exchanger__flow_and_temperature_control}
+        config.getConfigSubtree("flow_and_temperature_control"),
+        curves,
+        refrigerant);
+
+    return {borehole_geometry,         refrigerant, grout,
+            flowAndTemperatureControl, pipes,       bhe_if_use_python_bc_conf};
+}
+
+template <typename T_BHE>
+T_BHE createBHE1PType(
+    BaseLib::ConfigTree const& config,
+    std::map<std::string,
+             std::unique_ptr<MathLib::PiecewiseLinearInterpolation>> const&
+        curves)
+{
+    auto SinglePipeType = parseBHE1PTypeConfig(config, curves);
+    return {std::get<0>(SinglePipeType), std::get<1>(SinglePipeType),
+            std::get<2>(SinglePipeType), std::get<3>(SinglePipeType),
+            std::get<4>(SinglePipeType), std::get<5>(SinglePipeType)};
+}
+
+template BHE_1P createBHE1PType<BHE_1P>(
+    BaseLib::ConfigTree const& config,
+    std::map<std::string,
+             std::unique_ptr<MathLib::PiecewiseLinearInterpolation>> const&
+        curves);
+}  // namespace BHE
+}  // namespace HeatTransportBHE
+}  // namespace ProcessLib

ProcessLib/HeatTransportBHE/BHE/CreateBHE1PType.h

+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2020, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ */
+
+#pragma once
+
+#include <map>
+#include <memory>
+#include <string>
+
+namespace BaseLib
+{
+class ConfigTree;
+}
+namespace ProcessLib
+{
+namespace HeatTransportBHE
+{
+namespace BHE
+{
+template <typename T_BHE>
+T_BHE createBHE1PType(
+    BaseLib::ConfigTree const& config,
+    std::map<std::string,
+             std::unique_ptr<MathLib::PiecewiseLinearInterpolation>> const&
+        curves);
+}  // namespace BHE
+}  // namespace HeatTransportBHE
+}  // namespace ProcessLib

ProcessLib/HeatTransportBHE/BHE/CreateBHEUType.cpp

@@ -11,6 +11,7 @@
 #include "CreateBHEUType.h"
 #include "BaseLib/ConfigTree.h"
 
+#include "BHE_1P.h"
 #include "BHE_1U.h"
 #include "BHE_2U.h"
 #include "CreateFlowAndTemperatureControl.h"

ProcessLib/HeatTransportBHE/BHE/CreateFlowAndTemperatureControl.cpp

@@ -12,6 +12,7 @@
 #include "BaseLib/ConfigTree.h"
 #include "MathLib/InterpolationAlgorithms/PiecewiseLinearInterpolation.h"
 
+#include "BuildingPowerCurves.h"
 #include "CreateFlowAndTemperatureControl.h"
 #include "RefrigerantProperties.h"
 
@@ -86,6 +87,31 @@ FlowAndTemperatureControl createFlowAndTemperatureControl(
                                       refrigerant.specific_heat_capacity,
                                       refrigerant.density};
     }
+
+    if (type == "BuildingPowerCurveConstantFlow")
+    {
+        auto const& power_curve = *BaseLib::getOrError(
+            curves,
+            //! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__borehole_heat_exchangers__borehole_heat_exchanger__flow_and_temperature_control__BuildingPowerCurveConstantFlow__power_curve}
+            config.getConfigParameter<std::string>("power_curve"),
+            "Required power curve not found.");
+
+        auto const& cop_heating_curve = *BaseLib::getOrError(
+            curves,
+            //! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__borehole_heat_exchangers__borehole_heat_exchanger__flow_and_temperature_control__BuildingPowerCurveConstantFlow__cop_heating_curve}
+            config.getConfigParameter<std::string>("cop_heating_curve"),
+            "Required power curve not found.");
+
+        BuildingPowerCurves const building_power_curves{power_curve,
+                                                        cop_heating_curve};
+
+        //! \ogs_file_param{prj__processes__process__HEAT_TRANSPORT_BHE__borehole_heat_exchangers__borehole_heat_exchanger__flow_and_temperature_control__BuildingPowerCurveConstantFlow__flow_rate}
+        auto const flow_rate = config.getConfigParameter<double>("flow_rate");
+
+        return BuildingPowerCurveConstantFlow{
+            building_power_curves, flow_rate,
+            refrigerant.specific_heat_capacity, refrigerant.density};
+    }
     OGS_FATAL("FlowAndTemperatureControl type '%s' is not implemented.",
               type.c_str());
 }

ProcessLib/HeatTransportBHE/BHE/FlowAndTemperatureControl.h

@@ -11,6 +11,7 @@
 #pragma once
 
 #include <variant>
+#include "BuildingPowerCurves.h"
 #include "MathLib/InterpolationAlgorithms/PiecewiseLinearInterpolation.h"
 
 namespace ProcessLib
@@ -81,10 +82,34 @@ struct PowerCurveConstantFlow
     double density;
 };
 
+struct BuildingPowerCurveConstantFlow
+{
+    FlowAndTemperature operator()(double const T_out, double const time) const
+    {
+        double const power = building_power_curves.power_curve.getValue(time);
+        double const cop =
+            building_power_curves.cop_heating_curve.getValue(T_out);
+
+        if (power == 0)
+        {
+            return {0.0, T_out};
+        }
+        return {flow_rate,
+                power * (cop - 1) / cop / flow_rate / heat_capacity / density +
+                    T_out};
+    }
+    BuildingPowerCurves const building_power_curves;
+
+    double flow_rate;
+    double heat_capacity;
+    double density;
+};
+
 using FlowAndTemperatureControl = std::variant<TemperatureCurveConstantFlow,
                                                FixedPowerConstantFlow,
                                                FixedPowerFlowCurve,
-                                               PowerCurveConstantFlow>;
+                                               PowerCurveConstantFlow,
+                                               BuildingPowerCurveConstantFlow>;
 }  // namespace BHE
 }  // namespace HeatTransportBHE
 }  // namespace ProcessLib