Off-Design Modeling and Operational Optimization of Trans-Critical Carbon Dioxide Heat Pumps

With a contribution of more than 51% of the total worldwide energy consumption, the heating and cooling sector plays a major role in the transition toward more sustainable energy systems and in the mitigation of climate change [1]. In this context, the installed capacity of heat pumps increased significantly in recent years and is projected to grow even more in the future, as shown in Fig. 1 for the two major heat pump applications in Europe, namely, district heating and industrial heating [2]. The former refers to heat pumps for space and water heating; the latter refers to heat pumps used for heat recovery and energy efficiency improvements in different industrial sectors.

Different heat pump configurations and refrigerants are currently used depending on the end applications. After the seminal work of Lorentzen [3] and Nekså [4], natural refrigerants such as carbon dioxide (CO₂) have been gaining traction in recent years due to their low global warming potential, toxicity, and flammability compared to traditional refrigerants, such as ammonia or R134a amongst many others [5]. However, having a considerably lower critical temperature (31 °C) than conventional refrigerants, CO₂ is subject to a number of challenges when applying it in a subcritical cycle: (i) limited operational temperature range at the sink side of the cycle; (ii) low enthalpy of vaporization; (iii) heat rejection temperature similar to the refrigerant's critical temperature; and (iv) ambient temperature in summer at similar level as the refrigerant's critical temperature. All these factors translate into limitations and poor performance under certain operating conditions [4]. These challenges are tackled by designing and operating CO₂ heat pumps through a trans-critical cycle, where the maximum pressure is significantly higher than the critical pressure of CO₂ (71.1 bar). This avoids being limited in the heat delivery temperatures, hence widening the operation range and improving the performance of the thermodynamic process.

Several studies dealt with the modeling and analysis of trans-critical CO₂ heat pumps in the past. For example, Sarkar et al. provided a thorough assessment of different system designs and carried out an extensive experimental campaign showing that trans-critical CO₂ heat pumps can be used for a wide range of applications ranging from residential to commercial heating to agriculture [6]. Dai et al. evaluated heat pumps integrated with mechanical subcooling in terms of energy, exergy, and economic performance [7], whereas Wang et al. investigated the feasibility and performance of heat pumps integrated with thermal energy storage for space heating [8]. Recently, He et al. presented a modified trans-critical CO₂ heat pump system with new water flow configuration for residential space heating [9], whereas Wang et al. proposed a novel concept for simultaneous space heating and cooling [10]. Cao et al. reviewed the developments of air source trans-critical CO₂ heat pumps using direct-heated type and recirculating-heated type, whereas Yao et al. performed a comparative study of upgraded heat pump systems with different heat sinks [11]. Lo Basso et al. studied the potential role of trans-critical CO₂ heat pumps within a solar cooling system for building services [12]. Despite the latest advancements, a gap exists concerning the modeling and optimization of trans-critical CO₂ heat pumps operating under different operating conditions throughout the year [13].

While heat pumps are typically designed for specific nominal operating points, their performance is affected by variable boundary conditions, such as the heat duty and the water supply temperature required by the final end user and the ambient temperature [14]. We address this challenge by presenting a modeling and optimization tool to predict and optimize the operation of trans-critical CO₂ heat pumps in off-design conditions. Whereas the model developed is general, it is applied for a custom-made product developed at MAN Energy Solutions Schweiz AG, which features higher performance with respect to conventional heat pumps [15]. We present the mathematical model of such a novel thermodynamic cycle, which describes the behavior of CO₂ across a wide range of thermodynamic conditions. A high level of accuracy is reached by the model thanks to the detailed modeling of the heat exchangers (HEXs) and of the compressor, which takes into account effects of motor-cooling, leakages and performance limits. Furthermore, the model maximizes the coefficient of performance (COP) of the heat pump via an accurate and computationally efficient heuristic optimization algorithm. This allows (i) to improve the overall heat pump performance during the year, hence reducing the operational costs, and (ii) to investigate the tradeoff between the optimal design under nominal and under actual operating conditions.

This paper is structured as follows. Section 2 describes the CO₂ heat pump configuration considered in this work. Section 3 presents the mathematical model and the optimization algorithm, and Sec. 4 presents the validation of the mathematical model. Section 5 discusses the results by means of two case studies. Finally, Sec. 6 summarizes the study and draws conclusions.

The trans-critical heat pump considered in this work consists of a compressor, a gas cooler HEX, an expansion manifold (i.e., expander stage followed by an expansion valve), and an evaporator HEX. The main difference with respect to a conventional subcritical cycle is the fact that the working fluid (also called refrigerant) is compressed to a pressure above the critical point, hence to its supercritical thermodynamic conditions. This implies that the conventional condensation at the sink side of the process is replaced by the sensible cooling of the refrigerant.

Figure 2 shows a process flow schematic of the considered CO₂ heat pump, which is also a custom-made product developed at MAN Energy Solutions Schweiz AG [15]. This system improves on conventional trans-critical CO₂ heat pumps as follows:

• The purpose of the recuperator HEX is to pre-heat the refrigerant after the evaporator, on the low-pressure side of the cycle, by using the excess heat carried after the gas cooler, on the high-pressure side of the cycle. In fact, the recuperator is only needed when the end user requires high-return temperatures on the gas cooler water return side. For low return temperatures, the recuperator does not improve the heat pump COP and simply increases the investment cost of the system [6]. Compared to the gas cooler, the recuperator has a much smaller heat exchange area, because of the much higher temperature difference required across the HEX, and the lower heat duty.

• An additional heat exchanger, the rejection HEX, may be required upstream the expander to meet the required volume flow conditions. Heat is rejected to a source flow prior to entering the evaporator, and the trans-critical CO₂ is cooled down to suitable conditions at the expander suction.

• An expander stage is deployed between the rejection HEX (if required) and the expansion valve to mechanically exploit the pressure difference between the high- and low-pressure sides of the thermodynamic cycle. The expander stage is mounted on the same shaft as the compressor, thus reducing the overall mechanical shaft power of the compressor. The expander stage rotational speed is therefore imposed by the compressor motor.

This section presents an overview of the mathematical models of all relevant heat pump components (see Sec. 2), as well as the optimization algorithm used to maximize the heat pump COP under varying boundary conditions.

The major model inputs can be grouped into three categories:

• Load boundary conditions, which characterize the end-use application of the heat pump, independently of its specific components. These include the overall heat load that has to be met during the year, the temperatures of the heat recovery fluid (i.e., the fluid on the sink side, water here) required at the inlet and the outlet of the heat pump, and the temperature and pressure of the cold source of the heat pump (water here).

• Turbomachinery operating maps, which include the nondimensional maps of both the compressor and the expander.

• HEX design, which describes the exchange area of the HEXs, as well as their operating conditions at the design point.

The overall modeling and optimization procedure is shown in Fig. 3. First, the aforementioned input data are defined, alongside an initial guess of the decision variables of the optimization problem, i.e., the operation variables of the heat pump. Based on this, the model determines the conditions of CO₂ along the heat pump cycle, the COP of the system, as well as the penalties associated with all individual components (see Sec. 3.3). The modeling steps are embedded within an optimization procedure, which determines the optimal combination of operation variables that maximizes the COP of the system while minimizing the penalties of the individual components. The optimization procedure is based on a genetic algorithm, which solves the model iteratively and stops when predefined convergence criteria are met (see Sec. 3.3).

The considered HEXs include the gas cooler, the recuperator, and the rejection HEX (see Fig. 2), which are all modeled through the same approach. Various approaches can be found in the literature for modeling HEXs in CO₂ heat pumps. For example, Sarkar et al. [16], Li and Wang [17], and Ye et al. [18] modeled tube HEX based on their detailed geometrical characteristics. However, HEXs currently used in CO₂ heat pumps are often characterized by complex geometries, which results in models with high computational complexity. To reduce the computational complexity, a novel model is developed that enables the calculation of the heat transfer coefficients in off-design conditions based on the behavior of the HEXs at the design point. Such a novel approach is detailed in the following for a generic HEX. It requires the following input data:

• the type of HEX and its exchange area, A

• the inlet temperatures, T, and pressures, p, of both the hot (H) and the cold (C) fluids within the HEX

• performance of the HEX at the design point, which includes the pressure drops, $Δp$⁠, and the inlet and outlet temperature, pressure, and mass flow rates, m, of both the hot and cold fluids, hence the total heat load of the unit, Q

Here, Q_i is the heat load, U_i is the overall heat transfer coefficient, A_i is the heat exchange area, and LMTD is the logarithmic mean temperature difference across the ith segment of the HEX, which is computed by using $ΔTi=TH,i−TC,i$⁠. $mH$ and $mC$ are the mass flow rates of the hot and cold fluids, respectively, which are the same through all segments; $ΔhH,i$ and $ΔhC,i$ are the enthalpy variations of the hot and cold fluids across the ith segment of the HEX, respectively. $αH,i$ and $αC,i$ are heat transfer coefficients of the hot and cold fluids, respectively, across the ith segment of the HEX; they are calculated by using the values of the fluid thermal conductivity, λ, viscosity, μ, and Prandtl number, Pr. The coefficient B is the base heat transfer coefficient, which only depends on the type of fluid and on the type and geometry of the HEX; a and b are constant parameters that depend on the type of fluid and on the type of HEX, and are reported in Table 1.

The following calculation steps are performed for all segments $i∈{1,…,N}$⁠.

• The heat load, $Qi=Q/N$⁠, is equal for all segments and is known from the total heat load of the system.

• The known values of mass flow rates at design conditions, $mH$ and $mC$⁠, are used to compute the design enthalpy changes for both fluids, $ΔhH,i$ and $ΔhC,i$⁠, via Eq. (2). The known inlet temperatures and pressures are used to determine the enthalpy changes, hence the enthalpy profiles across the HEX. Such enthalpy profiles are combined with the linear pressure profiles assumed between the inlet and outlet fluid pressures at design conditions to determine the temperature profiles, $ΔTi$⁠, via the NIST standard reference database [19]. The temperature profiles are then used to determine LMTD_i via Eq. (4).

• With known pressure and temperature (or enthalpy) profiles, the thermodynamic properties of the HEX fluids, namely, λ, μ, and Pr, are computed via the NIST standard reference database [19] for both fluids. Furthermore, a first guess for the value of the base heat transfer coefficients, $BH$ and $BC$⁠, is chosen based on the type of fluids and on the type and geometry of the HEX. These quantities are used in Eq. (5) to determine the heat transfer coefficients, $αH,i$ and $αC,i$⁠, of both the hot and the cold fluids. The fluid heat transfer coefficients are then used to compute the overall heat transfer coefficient of the HEX, U_i.

• The heat load, the LMTD, and the overall heat transfer coefficients are used to determine the area of the ith segment, A_i. The total area calculated from the sum of the segment areas is compared against the actual area of the HEX, A. The values of the base heat transfer coefficients of bot fluids are changed iteratively until the calculated HEX area matches the actual one.

• By using the values of the base heat transfer coefficients, of the off-design heat load, and of the HEX exchange area, Eqs. (1) to (5) are used to determine the mass flow rate, pressure, and temperature profiles along the HEX in off-design conditions. For both the cold and hot fluids, the pressure drops along the HEX are computed as

  $Δp=Δpdes(mmdes)2$

  (6)

  where the subscript “des” denotes design conditions.

This iterative procedure proves to model the heat transfer across all HEX at the same time with a high enough accuracy and a low computational complexity. The accuracy of the model is assessed by comparing the total heat load of the system and the HEX pressure drops against manufacturer data. The model validation is shown in Sec. 4.

Here, turbomachinery include the compressor and the expander of the trans-critical CO₂ heat pump (see Fig. 2), which use CO₂ as the working fluid. The turbomachinery performance is often modeled through nondimensional maps of work versus flow characteristics, which are typically provided by manufacturers. Here, the input data required by the turbomachinery models are:

• the type and geometry of the turbomachinery

• the nondimensional performance maps of the turbomachinery in terms of work coefficient, μ₀, and flow coefficient, $ϕ$

• the inlet conditions of the working fluid, namely temperature, pressure, and mass flow rate

• the rotational speed, ω, of the turbomachinery.

The following calculation steps are performed to estimate the outlet conditions and the performance for both compressor [20] and the expander units [21]:

• The rotational speed (in revolutions per minute, rpm) and the diameter of the turbomachinery, D, are used to determine the tip wheel velocity, u

  $u=πω60D$

  (7)
• The tip wheel velocity and the known volume flow rate of the working fluid at inlet conditions, V, are used to compute the flow coefficient

  $ϕ=VuD2$

  (8)
• The work coefficient μ₀ and the polytropic efficiency η_p are determined through the turbomachinery nondimensional map by knowing the flow coefficient and the inlet stage Mach number. The total enthalpy difference across the turbomachinery is then calculated as

  $Δhtot=μ0u2$

  (9)

• The remaining outlet conditions of the turbomachinery are calculated by following the polytropic head corresponding to η_p and illustrated in Fig. 4. The polytropic head is divided into N segments of equal enthalpy difference $Δhi$ defining the polytropic path. Here, N = 30 proves to be the optimal tradeoff between accuracy and computational complexity.

• Starting from the inlet conditions (h₁ and p₁ in Fig. 4) and from the outlet enthalpy (h₂ in Fig. 4), the polytropic head is calculated for all segments $i∈{1,…,N}$ by first calculating the entropy at State 3, $s3,i$⁠, as a function of the inlet pressure, $p1,i$⁠, and the outlet enthalpy, $h2,i$⁠. The outlet entropy $s2,i$ is calculated through Eqs. (10) and (11) for compressor [20] and expander [21], respectively

  $s2,i=s1,i+(1−ηp)(s3,i−s1,i)$

  (10)

  $s2,i=s3,i+s1,i−s3,iηp$

  (11)

  The inlet conditions of segment i are set equal to the outlet conditions of segment i – 1, and Eqs. (10) and (11) are applied for all segments $i∈{1,…,N}$⁠. The outlet conditions of segment N define the outlet conditions of the turbomachinery unit.

• Once the outlet entropy of segment i is known, the remaining properties, i.e. the temperature and pressure profiles, are calculated as functions of entropy and enthalpy.

Here, $x∈ℝX$ is the vector defining the decision variables of the optimization problem, with X being the dimensionality of x; $Pi∈ℝ$ are the penalties of the individual components, and M the number of penalties considered, namely: (1) the compressor feasibility, (2) the system target load, (3) the area of the rejection HEX, (4) the area of the evaporator, and (5) the motor torque feasibility (see Fig. 3). The inequality constraints, g_j, include the lower- and upper-bound constraints for all optimization variables, which can range between a minimum and a maximum value. The equality constraints, h_k, include the equations describing the behavior of HEXs and turbomachinery presented in Secs. 3.2.1 and 3.2.2, respectively.

The decision variables of the optimization problem are:

• the mass flow rate of CO₂ in the loop (kg/s)

• the speed of the compressor (rpm)

• the outlet temperature of the recuperator HEX (if present) (°C)

• the temperature difference across the rejection HEX (if present) (°C)

• the mass flow rate of cold source (kg/s)

• the evaporator HEX pressure (bar).

The penalty functions are defined by multiplying a penalty coefficient, which is scaled according to the magnitude of the objective function, by the following discrepancies:

• (P1) Compressor feasibility: Difference between the work coefficient of the compressor obtained by the model and its lower and upper values; it quantifies the violation of the operation range of the compressor.

• (P2) System target load: Difference between the heat load of the gas cooler obtained by the model with that of the required heat load demand from the system.

• (P3) Rejection HEX penalty: Difference between the area of rejection HEX obtained by the model and its minimum and maximum values.

• (P4) Evaporator heat load penalty: Difference between the heat load of the evaporator obtained by the model and its maximum and minimum values.

• (P5) Motor torque feasibility: Difference between the torque required by the compressor and maximum allowed torque that the motor is able to deliver.

The optimization solver adopted is a genetic algorithm based on a covariance matrix adaptation evolution strategy [22]. The tolerance of the solver is the difference between the values of the decision variables in two subsequent iterations and is set to $10−4$⁠. The computational time required to meet the convergence criteria highly depends on the complexity of the system, i.e., presence or not of the rejection and recuperator HEXs. This increases both the number of optimization variables, hence of iterations required, and the time required per iteration. For the cases considered here, which include the recuperator HEX, the average time per iteration is close to 3 s, and requires approximately 1,100 iterations to meet the convergence criteria. This results in a typical run time shorter than 1 h.

The validation of the model is presented in Fig. 5 for the system's most relevant components, namely, the gas cooler, the evaporator, and the compressor. No validation is shown for the expander due the difficulty in measuring its performance under real operating conditions. The model results are assessed in terms of relative errors with respect to manufacturer data for different technology designs and for a variety of operating conditions. The overall model performance is assessed via the mean absolute percentage error (MAPE) [23] across all considered designs and conditions (i.e., the average of the absolute values of the relative errors).

The models of the gas cooler and the recuperator HEXs are assessed by comparing the total heat provided, Q, against manufacturer data. The validation results are presented in Fig. 5(a), which shows the relative error for two technology designs, simply denoted as Design 1 and Design 2, and four load points (LP), which are typical for the different seasons. These load points are selected as they cover the majority of the HEXs operation, with the heat load ranging from about 15 MWth to about 30 MWth of heat across case studies and load points.

Smaller errors are generally obtained for the gas cooler than for the evaporator, resulting in a MAPE of 0.29% for the former and 1.31% for the latter. Both values are deemed satisfactory for the applications of interest and improve the performance of earlier work presented for the gas cooler HEX [24].

Higher relative errors are observed when considering the pressure drops along the HEXs. MAPEs of 8.8% and 13% are obtained for the water and the CO₂ side of the gas cooler, respectively, and MAPES of 34% and 2.1% are obtained for the water and the CO₂ side of the evaporator, respectively. On the one hand, such larger relative discrepancies are due to the simple model adopted for pressure losses, given by Eq. (6). On the other hand, they correspond to absolute discrepancies smaller than about 100 kPa, which are deemed acceptable for the applications of interest, especially when comparing them to the operating pressures of the trans-critical CO₂ heat pump (in the range of 115–140 bar).

The model of the compressor is assessed by comparing the discharge temperature, the discharge pressure, and the shaft power against manufacturer data. The validation results are presented in Fig. 5(b), which shows the relative error for two technology designs and the aforementioned load points. These load points are selected as they cover the majority of the compressor operation, with the discharge temperature ranging between 390 and 450 K, the discharge pressure between 115 and 145 bar, and the shaft power between 6 and 10 MW electrical power across all considered ambient temperatures and load points.

Mean absolute percentage errors of 0.02%, 0.23%, and 0.09% are obtained for the discharge temperature, the discharge pressure and the shaft power, respectively, highlighting the excellent predicting performance of the compressor model.

In this section, the off-design model presented in Sec. 3 is applied to two real-world case studies to determine the optimal operational strategy that maximizes the heat pump performance, as well as its sensitivity to the most relevant decision variables, hereafter denoted as control variables.

The most relevant control variables are the compressor rotational speed and the loop resistance, which define the pressure rise across the compressor and ultimately its operating point. Whereas the rotational speed determines the power consumption of the compressor motor, and the resistance created by the expander manifold (see Fig. 2) determines the mass flowrate and the pressure ratio of the loop. Since the pressure ratio of the expander is given by its design characteristics and the rotational speed, the only way to control the overall resistance of the system is by adjusting the position of the expansion valve (hereafter referred as “valve”) downstream the expander stage. For a given rotational speed, closing the valve results in increasing the loop resistance (i.e., higher pressure ratio) and lowering the mass flowrate; this can be done until the surge line of the compressor is reached. In contrast, opening the valve results in a lower pressure ratio and a higher mass flowrate in the loop; this can be done until the compressor reaches choke conditions.

The optimal value of the valve position is thus a compromise, which takes into account the performance characteristics of all the components of the system, namely, the turbomachinery efficiencies and the performance of the HEXs. Although the compressor performance can be easily predicted through its characteristic map, the variation of the heat transfer coefficients of the HEXs due to changes in the thermal properties of supercritical CO₂, mostly in the gas cooler, adds another degree of complexity when operating trans-critical CO₂ heat pumps compared to conventional subcritical heat pumps.

The considered case studies include (1) industrial heat generation using a large body of water as heat source and (2) high-temperature district heating using the waste heat from data center (DC) as heat source.

Large bodies of water such as rivers, lakes, or seawater has been used for many years in the power industry to condense steam in Rankine cycles. The know-how gathered from this industry can be tapped when dealing with the heat source of heat pumps. This allows to use large flow rates of water, which result in lower temperature differences in the evaporator HEX, thus in higher evaporation temperatures hence higher COP. Moreover, the small amplitudes and slow variation in source temperature throughout the year facilitates the system operation, as it allows to gradually change the compressor rotational speed and the valve position as load conditions change.

The first considered case study describes an industry located near a large body of water, which installs a trans-critical CO₂ heat pump to heat up water from 40 °C to 95 °C for process heating application. Due to regulatory constraints based on realistic enquiries, it is allowed to extract a maximum water flowrate of 1,700 kg/s, which is fed via a circulation pump to the evaporator HEX. Given the fluctuations in the heating demand of the industrial end user during the year, the capability of the heat pump to meet different heat loads is a key feature of the system installed. Moreover, it is necessary not only to estimate the expected performance and operation of the system but also to predict its operational limitations in terms of heat output, both on the upper and lower bounds (i.e., part-load), throughout the year. The operation range of the heat pump (also called envelope), as well as its performance and control variables, for all possible operating conditions is shown in Fig. 6; the normalized compressor map selected and optimized for this case study is shown in Fig. 7, whereas the parameters of the HEXs are reported in Table 2. Different operating conditions, namely, winter, summer, and minimum load operations, are explicitly reported.

The operation range of the heat pump is delimited on the lower heat load end by the compressor design, as instabilities arise when operating the compressor at low mass flow rates near surge conditions. On the upper end, the operation is limited by the maximum allowed torque of the compressor motor, which limits the electrical power input to the compressor hence the heat load of the heat pump. Similarly, on the high end of the operation range, the maximum attainable flowrate of the compressor may be limited by choke conditions. Therefore, the size and shape of the heat pump envelope depend on the design criteria of the compressor. The maximum heat produced by the heat pump can also be limited by the maximum heat that can be extracted from the water source, in this case a lake or river. As shown in Fig. 6, for source temperatures lower than 3 °C, the maximum heat output is limited by the potential risk of freezing of the source water. For warmer source temperatures, the compressor remains the limiting component.

Figure 6 shows the interdependence of expansion valve position and compressor speed. The higher the heat load, and the higher the required temperature lift, the higher the compressor speed to increase the necessary pressure ratio. This also leads to a wider valve opening to increase the CO₂ mass flowrate.

Two major trends are observed when looking at the heat pump operation:

• Significant increase in COP for warmer source temperature. For a given heat load, the lower temperature lift that the compressor carries out results in a lower electrical consumption hence higher COP. While this is expected, the results allow to quantify accurately the relative performance deviation depending on the seasonal conditions. This means that operation cost variation throughout the year for varying operating conditions can be derived from the model results.

• For different source temperatures, the performance of the system is optimal at approximately 80% of the maximum heat load (see Fig. 8).

Starting from the maximum heat load, the reduction in heat load results in an increase in COP until the 80% load capacity. This is due to lower exergetic losses in the HEXs while keeping the efficiency of the turbomachinery near their best efficiency at design point. When further decreasing the heat load, the difference between operating point and design conditions results in lower turbomachinery efficiency, which offsets the lower exergy losses in the HEXs, and reduces the heat pump COP. The effect is accentuated by the reduction in turbulence in the HEX as a result of the mass flow rates reduction, which leads to higher exergetic losses and further penalizes the COP.

The second application presented in this study is the use of waste heat available from a DC as the heat source for a local high-temperature district heating (DH) network [25]. This application enables sector-coupling and exploits both the high performance and operational flexibility of trans-critical CO₂ heat pumps by providing both heat and cold to two different consumers at the same time [26,27]. The optimal profiles of the supply and return temperatures and of the control variables during the year are presented in Fig. 9. The DC returns water at 15 °C to the CO₂ heat pump, delivering 17 MWth of chilled water at 7 °C; this is used by the computer room air handler to cool down the servers room (also called white space). A constant cooling demand of the DC is assumed over the year according to the base load of a typical medium size DC. Accordingly, a constant heat load is delivered to the local DH network of a middle size European city (+200,000 inhabitants) at temperatures between 85 °C and 120 °C. Since the heat demand of DH networks varies seasonally, the coupling with the DC cooling demand is mostly possible for the base load of the DH heating demand.

The normalized compressor map selected and optimized for this case study is shown in Fig. 10, whereas the parameters of the HEXs are reported in Table 3. Compared to case study 1, the narrower range of operation resulting from the constant demand allows to select a more efficient compressor with a more limited operation range (see Fig. 10). Concerning the HEXs parameters, the major difference with respect to case study 1 is given by the evaporator (see Table 3). The higher temperature difference of the source side in the evaporator (8 K compared to approx. 3 K in case study 1) allows to increase the LMTD and to reduce the size of the component for equal heat load (see Eq. (1)).

The direct correlation between supply temperature and rotational speed of the compressor is directly linked to the discharge pressure that the compressor delivers in order to meet both the supply temperature to the DH (i.e., 85–120 °C) and the cooling load to the DC (i.e., 17 MWth) at 7 °C. As shown in Fig. 9, the higher supply temperatures required during the winter months result in higher rotational speeds. Moreover, the relatively small variation in speed throughout the year due to constant cooling supply allows to design the turbomachinery with a more restricted operating range, hence optimizing the efficiency over a narrower operating envelope

The variation in valve position between summer and winter operation is explained by the lower loop resistance (i.e., lower compressor discharge pressure) required to satisfy the thermal demands during summer, and therefore the higher CO₂ mass flowrate is necessary to deliver the cooling demand, maximizing both heat pump performance and profitability in operation.

Since the heat pump system delivers a constant cooling load of 17 MWth (i.e., dictated by the DC demand), the power input and the heat pump COP solely depend on the operating temperatures, which vary on the DH network side in a seasonal fashion. As shown in Fig. 9, the lower required temperature lift in summer allows to increase the COP of the system, fluctuating between 2.87 and 3.06 depending on the season.

The sensitivity of heat pump performance on boundary conditions is key when design a heat pump, and especially the turbomachinery. The developed model helps the decision-making process at design phase with the goal of optimizing the end user profitability. As an example, the higher profit expected from selling heat in winter leads to designs of the compressor and the expander that allows for the highest possible COP in winter, possibly reducing the COP during other periods of the year. This implies that the power consumption during summer operation could be further reduced if both compressor and expander designs were optimized specifically for these boundary conditions.

This study presents the model of an optimization tool developed and used to predict the operation of trans-critical CO₂ heat pumps in off-design conditions. As the operation of these systems strongly depend on the variable end user heat demand, it is key to predict the system performance for varying boundary conditions and to assess optimal operation strategies to maximize the heat pump performance, hence its profitability, under realistic conditions. The tool allows at the same time high accuracy and fast computational performance and can in principle be used for any heat and cold industrial energy demand. Here, it is presented by referring to two case studies of interest, namely, (i) heat generation using a large body of water as heat source and (ii) district heating using waste heat from data centers as heat source.

The analysis of the optimal heat pump operation shows that the rotational speed of the compressor and the throttle valve opening are the most relevant control variables to maximize the heat pump performance over the year. The proposed model determines the set-point values of the main control variables to maximize the heat pump performance under different operating conditions, while complying with the operational limits imposed by the individual heat pump components. Such components include the compressor, the expander and the HEXs. As a result of the optimization procedure, the optimal values of the COP, of the compressor speed, and of the valve opening are shown for the complete operation envelopes of the system, i.e., for a wide range of source temperature and net heat load. It is worth noting that the maximum heat pump efficiency is predicted at about 80% of max heat load capacity, which suggests to maximize the operating time for these conditions whenever possible.

Overall, the model can be used to assess and modify the operational points set by different industrial end users during the year. This allows to (i) evaluate the effectiveness of a system design under realistic operating conditions and (ii) maximize the economic and environmental performance of an end user, while accounting for both design and off-design operations.