Parametrization of a two stage heat pump with flashtank

[JanKochanski]

Hello,

I am fairly new to tespy and currently trying to model a heat pump with two stages and a flash tank/intermediate vessel. I've been going through the tutorials I could find in the documentation, which cover the parametrization as well as debugging and the heat pump modelling, but did not find a solution so far, that would solve the issues.

The model is sketched in the System Scheme and should be driven by the temperatures at the seawater inlet and outlet (in 4°C, out 2°C) as well as the given temperatures at the warm side (40°C in, 80°C out). The pressure ratios for both compressors should be equal, I am however not sure how to implement this best. I did not manage to make the approach as shown in the tutorial Build a Heat Pump Stepwise work in my settings.

My goal is to create a fairly flexible model of the heat pump, with the goal to test different refrigerants with the same setup later (which means using hardware-specific parameters in the next step).

I would be really grateful if I could receive some help with the parametrization, and I am always happy to get feedback regarding code and methodology!
Best regards, Jan

heatpumpV4.py

[fwitte]

Hi @JanKochanski,

thank you for reaching out! Here is some general advice:

• build your model incrementally until you have it completed, then you can go back an remove previous steps if necessary
• make the model in the class solve once trivially, then change the parameters to the desired set of parameters to solve the final model. If that does not work in a single step, do multiple steps in between
• there were some issues with your parameter setup:
  • the drum and droplet seperator both have specification of saturated steam/liquid implemented in them, no need to specify those
  • the pressure value at the end of the heat source and heat sink side was missing, the pump does not know to which pressure to increase to
  • I changed the condenser to a moving boundary model, because then we can impose the (true) pinch
  • I changed the subcooler to a monving boundary model to also impose the pinch, as it might be either at inlet or at outlet (theoretically at least)
  • for the first solve run I specify compressor outlet states instead of compressor efficiency and condensation/evaporation pressure/temperature levels instead of heat exchanger pinches as these are more simple equations
• you can implement arbitrary equations with the UserDefinedEquation class, I did that to impose pressure ratio equality of the compressors
• I hacked my changes into your model structure. I actually want to move to make a "TESPyModelTemplate" class or something like that, to enable easy plotting, parameter specification, compatibility with the optimization tools etc.. See also: https://tespy.readthedocs.io/en/main/integration/workflows.html#integration-model-class-template-label. This is one of the top points on the agenda for our tespy user meeting coming up in April in Brunswick: https://github.com/oemof/oemof/wiki/Meeting-2026.04%3A-tespy-user-meeting. You are welcome to participate!

If you need more frequent support you can also contact me directly :)

Best

Francesco

from tespy.networks import Network
from tespy.components import CycleCloser, Compressor, Valve, HeatExchanger, Merge, DropletSeparator, Condenser, MovingBoundaryHeatExchanger, Source, Sink, Pump, Drum
from tespy.connections import Connection
from tespy.tools import UserDefinedEquation
import math

class heatpump():
    def __init__(self, refrigerant):
        #init pressures
        self.p_high = 15
        self.p_low = 2

        #district heating side
        self.m_DH = 120.0  # in kg/s
        self.T_DH_return = 40.0  # in °C
        self.T_DH_supply = 80.0  # in °C

        #sea water side
        self.T_SW_warm = 4.0  # in °C
        self.T_SW_cold = 2.0   # in °C

        # initialize the network
        self.my_plant = Network()
        self.refrigerant = refrigerant
        # set the unitsystem for temperatures to °C and for pressure to bar
        self.my_plant.units.set_defaults(
            temperature="degC", pressure="bar", enthalpy="kJ/kg", heat="MW", power="MW")

        # Build the heat pump model from it's components
        self.CC = CycleCloser('cycle closer')
        self.Cond = MovingBoundaryHeatExchanger('condenser')
        self.Evp = HeatExchanger('evaporator')
        self.Drum = Drum('drum')
        self.SC = MovingBoundaryHeatExchanger('subcooler')
        self.ExVHP = Valve('Expansion Valve High Pressure')
        self.ExVLP = Valve('Expansion Valve Low Pressure')
        self.CompLP = Compressor('Compressor Low Pressure')
        self.CompHP = Compressor('Compressor High Pressure')
        self.IGC = Merge('Intermediate Gas Cooler')
        self.IV = DropletSeparator("Intermediate Vessel")

        # secundary side components
        self.DHP = Pump('District heating pump')
        self.DHWC = Source('District heating return side')
        self.DHWW = Sink("District heating supply side")

        self.SWP = Pump('Sea Water Pump')
        self.SWW = Source('Sea water warm side')
        self.SWC = Sink('Sea water cold side')


        # connections of sea water side
        self.C30 = Connection(self.SWW, 'out1',
                              self.Evp, 'in1',
                              label='Connection 30, Sea Water Warm Side to Evaporator')
        self.C31 = Connection(self.Evp, 'out1',
                              self.SWP, 'in1',
                              label='Connection 31, Evaporator to Sea Water Pump')
        self.C32 = Connection(self.SWP, 'out1',
                              self.SWC, 'in1',
                              label='Connection 32, Sea Water Pump to Sea Water Cold Side')
        self.my_plant.add_conns(self.C30,
                                self.C31,
                                self.C32)

        # connections of district heating side
        self.C20 = Connection(self.DHWC, 'out1',
                              self.SC, 'in2',
                              label='Connection 20, District Heating Return Side to Subcooler')
        self.C21 = Connection(self.SC, 'out2',
                              self.Cond, 'in2',
                              label='Connection 21, Subcooler to Condenser')
        self.C22 = Connection(self.Cond, 'out2',
                              self.DHP, 'in1',
                              label='Connection 22, Condenser to District Heating Pump')

        self.C23 = Connection(self.DHP, 'out1',
                              self.DHWW, 'in1',
                              label='Connection 23, District Heating Pump to District Heating Supply Side')

         # add district heating connections to the network
        self.my_plant.add_conns(self.C20,
                                self.C21,
                                self.C22,
                                self.C23)

        # connections of heat pump
        self.C1 = Connection(self.CC, 'out1',
                             self.Cond, 'in1',
                             label='Connection 1, Cyclecloser to Condenser')

        self.C2 = Connection(self.Cond, 'out1',
                             self.SC, 'in1',
                             label='Connection 2, Condenser to Subcooler')

        self.C3 = Connection(self.SC, 'out1',
                             self.ExVHP, 'in1',
                             label='Connection 3, Subcooler to Expansion Valve High Pressure')

        self.C4 = Connection(self.ExVHP, 'out1',
                             self.IV, 'in1',
                             label='Connection 4, Expansion Valve High Pressure to Intermediate Vessel')

        self.C5 = Connection(self.IV, 'out1', # saturated liquid
                             self.ExVLP, 'in1',
                             label='Connection 5, Intermediate Vessel to Expansion Valve Low Pressure')

        self.C6 = Connection(self.ExVLP, 'out1',
                             self.Drum, 'in1', # mixed state
                             label='Connection 6, Expansion Valve Low Pressure to Drum')

        self.C7 = Connection(self.Drum, 'out1', # saturated liquid
                             self.Evp, 'in2',
                             label='Connection 7, Drum to Evaporator')

        self.C8 = Connection(self.Evp, 'out2',
                             self.Drum, 'in2', # saturated gas
                             label='Connection 8, Evaporator to Drum')

        self.C9 = Connection(self.Drum, 'out2', # saturated gas
                             self.CompLP, 'in1',
                             label='Connection 9, Drum to Compressor Low Pressure')

        self.C10 = Connection(self.CompLP, 'out1',
                             self.IGC, 'in1',
                             label='Connection 10, Compressor Low Pressure to Intermediate Gas Cooler')

        self.C11 = Connection(self.IV, 'out2', # saturated gas
                             self.IGC, 'in2',
                             label='Connection 11, Intermediate Vessel to Intermediate Gas Cooler')

        self.C12 = Connection(self.IGC, 'out1',
                              self.CompHP, 'in1',
                              label='Connection 12, Intermediate Gas Cooler to Compressor High Pressure')

        self.C13 = Connection(self.CompHP, 'out1',
                              self.CC, 'in1',
                              label='Connection 13, Compressor High Pressure to Cyclecloser')

        # add heatpump connections to the network
        self.my_plant.add_conns(self.C1,
                                self.C2,
                                self.C3,
                                self.C4,
                                self.C5,
                                self.C6,
                                self.C7,
                                self.C8,
                                self.C9,
                                self.C10,
                                self.C11,
                                self.C12,
                                self.C13)


        # define component parameters
        self.CompHP.set_attr()#eta_s=0.8)
        self.CompLP.set_attr()#eta_s=0.8)

        self.Cond.set_attr(ttd_u=None, pr1 = 0.99, pr2 = 0.99)
        self.Evp.set_attr(ttd_l=None, pr1 = 0.99)  # needs to be hot side, cold side cannot have pressure loss
        self.SC.set_attr(ttd_l=None, pr1 = 0.99, pr2 = 0.99)        # not currently used

        self.SWP.set_attr(eta_s=0.75)
        self.DHP.set_attr(eta_s=0.75)

        # define connection parameters in the district heating side
        self.C20.set_attr(T=self.T_DH_return,
                          p=6,                 # Assumed pressure at district heating inlet
                          fluid={"water": 1},
                          m=self.m_DH)
        self.C21.set_attr(T=self.T_DH_return + 5)
        self.C22.set_attr(T=self.T_DH_supply)
        self.C23.set_attr(p=7)  # DH supply temperature

        # define connection parameters in the sea water side
        self.C30.set_attr(T=self.T_SW_warm,
                          fluid={"water": 1},  # Sea water warm side temperature
                          p=1.5) # Assumed pressure at sea water inlet

        self.C31.set_attr(T=self.T_SW_cold)
        self.C32.set_attr(p=1.5)

        # define connection parameters in the heatpump
        self.C1.set_attr(fluid={self.refrigerant: 1})
        self.C2.set_attr(x=0, T=90)
        self.C3.set_attr() # Temeprature should be defined by ttd_u of condenser
        self.C5.set_attr(T_bubble=50) # saturated liquid -> not required
        self.C7.set_attr(T=0) # Temperature should be defined by ttd_l of evaporator
        self.C8.set_attr(x=0.9) # not fully evaporated, otherwise the drum would not make sense
        self.C10.set_attr(td_dew=30)
        self.C13.set_attr(td_dew=30)#,
                          #p=self.calc_p_mid()) # initialize with "good guess" and calculate it after first solution
        self.my_plant.solve(mode='design')#, init_only=True)

    def set_final_parameters(self):
        self.C10.set_attr(td_dew=None)
        self.C13.set_attr(td_dew=None)
        self.CompHP.set_attr(eta_s=0.8)
        self.CompLP.set_attr(eta_s=0.8)

        def same_pressure_ratio_ude(ude):
            c1, c2, c3 = ude.conns
            return c2.p.val_SI ** 2 - c3.p.val_SI * c1.p.val_SI  # same as c2 / c1 - c3 / c2

        def same_pressure_ratio_ude_dependents(ude):
            return [c.p for c in ude.conns]

        ude = UserDefinedEquation(
            "same pressure ude",
            same_pressure_ratio_ude,
            same_pressure_ratio_ude_dependents,
            conns=[self.C9, self.C10, self.C13]
        )

        self.C5.set_attr(T_bubble=None)
        self.my_plant.add_ude(ude)

        self.C7.set_attr(T=None)
        self.Evp.set_attr(ttd_l=5)

        self.C2.set_attr(T=None)
        self.Cond.set_attr(td_pinch=5)

        self.C22.set_attr(T=None)
        self.SC.set_attr(td_pinch=5)

    def calc_p_mid(self):
        def is_valid_number(x):
            return isinstance(x, (int, float)) and not math.isnan(x)
         # recalculate intermediate pressure based on previous solution
        if is_valid_number(self.C1.p.val) and is_valid_number(self.C9.p.val):
            print(f"Using calculated high pressure: {self.C1.p.val:.2f} bar")
            print(f"Using calculated low pressure: {self.C9.p.val:.2f} bar")
            p = (self.C9.p.val*(self.C9.p.val/self.C1.p.val)**0.5)
        # for initialization with fixed pressures
        else:
            p = (self.p_low*(self.p_high/self.p_low)**0.5)
        print(f"Calculated intermediate pressure: {p:.2f} bar")
        return p

    def solve(self):
        self.my_plant.solve(mode='design')#, init_only=True)
        self.my_plant.print_results()
        # calculate and print performance indicators
        p_el = self.CompLP.P.val + self.CompHP.P.val
        q_heat = abs(self.Cond.Q.val)
        cop = q_heat / p_el
        print(f"Heating Power: {abs(self.Cond.Q.val):.2f} MW")
        print(f"COP: {cop:.2f}")
        print(f"Temperature after HP Compressor: {self.C1.T.val:.2f} °C")
        print(f"Massflow gaseous from IV (C11): {self.C11.m.val:.4f} kg/s")
        print(f"Massflow liquid from IV (C5): {self.C5.m.val:.4f} kg/s")
        print(f"Gas after 1st valve (C4): {self.C4.x.val:.4f}")

if __name__ == "__main__":
    hp = heatpump("R134a")
    hp.set_final_parameters()
    hp.solve()

[JanKochanski]

Hi @fwitte,

Thank you so much for your help! I saw the same approach of initializing and then solving once with trivial inputs in streamlit dashboard heatpumps codes, but was not quite sure if this step is always needed. Is this step done to make the solver start in an area of the solution space, where it can converge faster / is not going to diverge as easily?

Thank you again for your support, I highly appreciate it!
Best regards,
Jan

[fwitte]

Is this step done to make the solver start in an area of the solution space, where it can converge faster / is not going to diverge as easily?

Yes, that is exactly the reason for it. Some equations are harder (depend on more variables and require a good initial guess) than others. If there is a good initial guess (after the first solution since tespy always starts its simulation with the variables initialized from the previous solution), then these are no problem usually.