GSHP integrated with DH

[AlexanderKlimko]

I'm currently working on my master thesis in power generation through clean conversion technologies.

The project aim to provide economical benefits by improving the efficiency of residential ground source heat pumps (GSHP) using the return flow of the district heating network in Sweden. During the summer months the return flow of the district heating network is routed into the boreholes in order to increase the ground temperature. By increasing the temperature in the boreholes, the goal is to increase the coefficient of performance (COP) and amount of energy able to be retrieved of the GSHP. This is intended to decrease the compressor power during operation of the GSHP in the winter months, resulting in lower costs for the operator of the GSHP. Further, by enabling heat to be stored in the boreholes, economical benefits for the cogeneration plant operators are sought through sales of otherwise wasted heat during summer months.

I aim to simulate such a system using a combination of TESpy for the GSHP and pyg-function to simulate the thermal response in the boreholes. The model should be able to simulate the operation of such a proposed system for a time period of 40 years. Hence the part simulated using TESpy is done for each hour of simulation passing results into pyg-function and vice versa for the next time step.

However, I am running into an issue when trying to change the heat transfer fluid (HTF) for the borehole thermal energy storage (BTES) from 100% water to a mixture of 75% water and 25% ethanol. As ground temperatures go below 0°C in the winter, it's needed to add ethanol to the water in order to prevent freezing. I've tried using enthalpy instead of temperature for the HTF to no avail. The solver won't converge.

I'd be very happy for any comments, hints or tips on my code!

from tespy.networks import Network
from tespy.components import (
    CycleCloser, Compressor, Valve, SimpleHeatExchanger, Pump, Condenser, HeatExchanger
)
from tespy.connections import Connection, Ref
import pandas as pd
import CoolProp.CoolProp as CP
import simulation as sim
import btes as bt

def initialise_tes_network(input_data: pd.DataFrame, T_g: float, print_results: bool = False) -> tuple[Network, float, float]:
    """
    Initializes and creates a network for a thermal energy system, setting up components, connections,
    parameters, and solving the system for the given conditions.

    :param input_data: DataFrame containing system inputs, including heat demand, temperatures, 
                       and other relevant parameters for the thermal energy system.
    :param T_g: Ground temperature in °C, used as the BTES forward flow temperature.
    :param print_results: Boolean flag indicating whether to print the results of the solved network (default: False).
    
    :return: A tuple containing:
        - Solved Network object representing the thermal energy system.
        - Mass flow rate in the BTES loop (kg/s).
        - Heat exchanged by the BTES (Wh).
    """
    
    # Create a network object
    nw = Network(
        T_unit='C', p_unit='bar', h_unit='kJ / kg', m_unit='kg / s'
    )

    K = 273.15  # Kelvin conversion constant
    # OBSERVE: ALL HEAT EXCHANGERS ARE COUNTER CURRENT

    # -------------------------------------------------------------------------
    # Heat demand/supply
    # -------------------------------------------------------------------------
    Q_cons = input_data.iloc[0, 1]

    # -------------------------------------------------------------------------
    # Temperatures
    # -------------------------------------------------------------------------
    T_cons_ff = 60  # Consumer forward flow temperature
    T_cons_bf = 50  # Consumer backward flow temperature
    T_btes_ff = T_g  # BTES forward flow temperature
    T_dh_ff = input_data.iloc[0, 3]  # District heating forward flow temperature

    # -------------------------------------------------------------------------
    # Fluids
    # -------------------------------------------------------------------------
    wf_hp = 'R410A'  # Working fluid in the heat pump
    wf_btes = 'water', 'ethanol'  # Working fluid in BTES
    wf_cons = 'water'  # Working fluid for the consumer
    wf_dh = 'water'  # Working fluid for district heating

    # Validate properties
    rho_btes = CP.PropsSI('D', 'T', T_g + K, 'P', 1.2e5, f"{wf_btes[0]}[0.75]&{wf_btes[1]}[0.25]")   # Density
    cp_btes = CP.PropsSI('C', 'T', T_g + K, 'P', 1.2e5, f"{wf_btes[0]}[0.75]&{wf_btes[1]}[0.25]")    # Specific heat
    h_btes = CP.PropsSI('H', 'T', T_g + K, 'P', 1.2e5, f"{wf_btes[0]}[0.75]&{wf_btes[1]}[0.25]")/1000     # Enthalpy
    print(f"Density: {rho_btes} kg/m³, Specific heat: {cp_btes} J/kg·K, Enthalpy: {h_btes} kJ/kg")

    # -------------------------------------------------------------------------
    # Consumer components
    # -------------------------------------------------------------------------
    p_cons = Pump('Consumer recirculation pump')
    cons = SimpleHeatExchanger('Consumer')
    cc2 = CycleCloser('Consumer cycle closer')

    # -------------------------------------------------------------------------
    # Heat pump components
    # -------------------------------------------------------------------------
    co = Condenser('Condenser')
    cp = Compressor('Compressor')
    va = Valve('Valve')
    ev = HeatExchanger('Evaporator')
    cc1 = CycleCloser('Heat pump cycle closer')

    # -------------------------------------------------------------------------
    # BTES components
    # -------------------------------------------------------------------------
    btes = SimpleHeatExchanger('BTES')
    cc3 = CycleCloser('BTES cycle closer')
    p_btes = Pump('BTES recirculation pump')

    # -------------------------------------------------------------------------
    # Temporary components (District heating)
    # -------------------------------------------------------------------------
    he = HeatExchanger('Heat exchanger district heating')
    dh = SimpleHeatExchanger('District Heating')
    cc4 = CycleCloser('District heating cycle closer')
    p_dh = Pump('District heating pump')

    # -------------------------------------------------------------------------
    # Connections
    # -------------------------------------------------------------------------
    # Heat pump connections
    c0 = Connection(va, 'out1', cc1, 'in1', label='0')
    c1 = Connection(cc1, 'out1', ev, 'in2', label='1')
    c2 = Connection(ev, 'out2', cp, 'in1', label='2')
    c3 = Connection(cp, 'out1', co, 'in1', label='3')
    c4 = Connection(co, 'out1', va, 'in1', label='4')

    # Consumer connections
    c5 = Connection(p_cons, 'out1', co, 'in2', label='5')
    c6 = Connection(co, 'out2', cons, 'in1', label='6')
    c7 = Connection(cons, 'out1', cc2, 'in1', label='7')
    c8 = Connection(cc2, 'out1', p_cons, 'in1', label='8')

    # BTES connections
    c9 = Connection(ev, 'out1', he, 'in2', label='9')
    c10 = Connection(he, 'out2', cc3, 'in1', label='10')
    c11 = Connection(cc3, 'out1', btes, 'in1', label='11')
    c12 = Connection(btes, 'out1', p_btes, 'in1', label='12')
    c13 = Connection(p_btes, 'out1', ev, 'in1', label='13')

    # District heating connections
    c14 = Connection(he, 'out1', cc4, 'in1', label='14')
    c15 = Connection(cc4, 'out1', dh, 'in1', label='15')
    c16 = Connection(dh, 'out1', p_dh, 'in1', label='16')
    c17 = Connection(p_dh, 'out1', he, 'in1', label='17')

    # Add connections to the network
    nw.add_conns(c5, c6, c7, c8)
    nw.add_conns(c0, c1, c2, c3, c4)
    nw.add_conns(c9, c10, c11, c12, c13)
    nw.add_conns(c14, c15, c16, c17)

    # -------------------------------------------------------------------------
    # Parametrisation
    # -------------------------------------------------------------------------
    # Set fluid types
    c3.set_attr(fluid={wf_hp: 1})       # Sets the fluid to 100% R410A
    c6.set_attr(fluid={wf_cons: 1})     # Sets the fluid to 100% Water
    c12.set_attr(fluid={wf_btes[0]: 0.75, wf_btes[1]: 0.25})    # Fluid composition, index 0: water, index 1: ethanol.
    c16.set_attr(fluid={wf_dh: 1})      # Sets the fluid to 100% Water

    # Set pressure ratios and temperature differences
    co.set_attr(pr1=0.98, pr2=0.98)     # Pressure loss 2% on hot and cold side
    ev.set_attr(pr1=0.98, pr2=0.98, ttd_l=2.0)
    cons.set_attr(pr=0.98)              # Pressure loss 2% over the consumer
    btes.set_attr(pr=0.98)              # Pressure loss 2% over the BTES
    he.set_attr(pr1=0.98, pr2=0.98)
    dh.set_attr(pr=0.98)

    # Set pump and compressor efficiencies
    p_cons.set_attr(eta_s=0.5)
    p_btes.set_attr(eta_s=0.5)
    p_dh.set_attr(eta_s=0.5)
    cp.set_attr(eta_s=0.85)

    # Set temperatures
    c6.set_attr(T=T_cons_ff)    # Consumer feed temperature
    c7.set_attr(T=T_cons_bf)    # Consumer backflow temperature
    c12.set_attr(h=h_btes)      # Enthalpy at connection 12 after exiting BTES
    c16.set_attr(T=T_dh_ff)     # District heating feed temperature
    c14.set_attr(T=Ref(c17, 1, -5))     # District heating backflow temperature

    # Set pressures
    p_cond = CP.PropsSI('P', 'Q', 1, 'T', K + T_cons_ff + 3, wf_hp) * 1e-5
    c4.set_attr(p=p_cond)  # Pressure after condenser in the HP loop
    p_ev = CP.PropsSI('P', 'Q', 0, 'T', K + T_btes_ff - 2, wf_hp) * 1e-5
    c2.set_attr(p=p_ev)    # Pressure after evaporator in the HP loop
    c6.set_attr(p=2)       # Pressure in the consumer loop
    c12.set_attr(p=1.2)    # Pressure in the BTES loop
    c17.set_attr(p=1, m=0.01)

    # Set heat demand/supply
    cons.set_attr(Q=Q_cons)

    # Set fluid state
    c2.set_attr(x=1)

    # Solve the network and print results if enabled
    nw.solve('design', print_results=print_results)
    nw.save('data/network')
    m_flow_btes_loop = c12.m.val
    q_btes = btes.Q.val

    if print_results:
        nw.print_results()
    
    return nw, m_flow_btes_loop, q_btes

def main():
     # Read the input file
    m_input = sim.import_data('1_day_input_sample.csv')

    # Initialise the TESpy network
    network_object = initialise_tes_network(m_input, 10, print_results=True)

if __name__ == '__main__':
    main()

[fwitte]

Hi @AlexanderKlimko,

thanks for reaching out! TESPy has implementations for the following mixing rules

• ideal gas (everything is assumed gas phase)
• ideal gas with water beeing able to condense (default, everything is gas, the water condenses up until the partial pressure of the water in the gas is in saturation given the mixture temperature)
• incompressible (simple mixing rule, which is a good guess for incompressible liquid mixing but not considering any interaction between the mixed fluids).

Because you are trying to use pure water mixed with ethanol, you cannot use any of the present mixing rules however, because they will always call the pure fluid properties from CoolProp. Those limit water to ~0 °C. Instead, you can either use CoolProp's inbuilt mixture backend for your specification (fluid={"Water[0.75]&Enthanol[0.25]: 1}). This may not work however, I have never tested it, because the mixture backend of CoolProp does not work well (that is also the reason, why I implemented custom implementations for mixtures in tespy). Or, you can use the incompressibles backend, which should work: http://www.coolprop.org/fluid_properties/Incompressibles.html#id136. See my answer to this discussion on how to use it in tespy: #473.

On a side note: I would be interested in the results of your thesis, if you like to share them with me (I'll keep them confidential!)

Best

Francesco

[AlexanderKlimko]

Hi @fwitte,

Thank you for the reply. I've tried several methods in order for a water-ethanol mixture to provide a convergent solution. I've found that if you find a convergent solution with non-mixtures such as water the solver is able to find for temperatures > 0°C using a water-ethanol mixture. However, I've not found any way to use such a mixture for temperatures below 0°C that provides a convergent solution. I've instead resorted to using c12.set_attr(fluid={'INCOMP::MEG[0.25]': 1}) as heat transfer fluid for the borehole branch. Basically, the trick is to find a convergent solution and later adjust parameters.

I can share my repository for my thesis with you!

Kindly,

Alexander

[fwitte]

Hi Alexander,

thank you very much for the feedback, I have received the invitation.

Have a nice weekend

Francesco