Effects of the Working Fluid Charge in Organic Rankine Cycle Power Systems: Numerical and Experimental Analyses

Abstract

It is well known that organic Rankine cycle (ORC) power systems often operate in conditions differing from the nominal design point due to variations of the heat source and heat sink conditions. Similar to a vapor compression cycle, the system operation (e.g., subcooling level, pump cavitation) and performance (e.g., heat exchanger effectiveness) of an ORC are affected by the working fluid charge. This chapter presents a discussion of the effects of the charge inventory in ORC systems. In particular, both numerical and experimental aspects are presented. The importance of properly predicting the total amount of working fluid charge for optimizing design and off-design conditions is highlighted. Furthermore, an overview on state-of-the-art modeling approaches is also presented.

1. Introduction

Organic Rankine cycle (ORC) systems are widely acknowledged as one of the most suitable technologies for harvesting medium- to low-grade heat sources (i.e., below 300°C) from both renewable sources (e.g., geothermal and solar) and waste heat [1]. Nowadays, the total power capacity installed worldwide is estimated to be above 2.7 GWe, as reported by Tartière et al. in the ORC World Map [2, 3]. Furthermore, there has been a continuous increase in research activity related to ORC-based power systems over the last decades that has demonstrated the value of this technology and its potential to improve energy sustainability [4].

Generally, ORCs are designed and optimized for steady-state operating conditions, unless mobile applications are considered. However, in the vast majority of stationary applications (e.g., solar thermal power, geothermal, waste heat recovery, and combined heat and power), the heat source and heat sink conditions are subjected to fluctuations and the ORC systems often operate at part-load or off-design conditions. To this end, numerical and experimental studies on off-design ORC performance have been published in the scientific literature, as outlined by Dickes et al. [5]. Similar to vapor compression systems, ORC performance is sensitive to working fluid charge and an incorrect charge level affects the system performance, especially under part-load conditions. Yet, very limited work can be found on the impact of working fluid charge on off-design operations. In fact, besides the research conducted by the Authors [6, 7, 8], only two additional studies could be found in the literature. Liu et al. [9] discussed the effect of different working fluid charge masses on the system performance. Experimental and numerical analyses were carried out to identify the optimal working fluid charge. Pan et al. [10] analyzed the impact of the working fluid charge on the length distribution of the different zones in the heat exchangers under several operating conditions. Furthermore, different void fraction models and heat transfer correlations were compared to assess their impact on the working fluid mass estimations in the heat exchangers.

The present chapter discusses the relationship between working fluid charge and performance in ORC systems under different operating conditions. Furthermore, the fundamentals of deterministic numerical methods to account for working fluid charge in ORC simulations are also described.

2. Importance of working fluid charge in ORC systems

A general sub-critical ORC with an internal heat exchanger (or regenerator), as shown in Figure 1, is considered as an example to understand the effect of working fluid charge. Although this chapter focuses on sub-critical ORC systems, the principles and the methodologies discussed hereafter can be extended to transcritical cycles and other cycle configurations [11]. By referring to Figure 1, the system consists of three heat exchangers (e.g., brazed plate heat exchangers), a pump (either centrifugal or volumetric type), a positive displacement expander (e.g., scroll type), a liquid receiver (or buffer tank), and line sets (or pipelines) between each component. Depending on the type of expander (e.g., single- and twin-screw, swash plate, or scroll), a dedicated oil-injection loop with an additional pump may be present in the system (see Figure 1). An oil separator is typically installed after the expander to avoid accumulation of the lubricant oil inside the heat exchangers, where hot spots may occur.

Figure 1.

General schematic of a regenerative organic Rankine cycle system. An independent lubricant oil loop is also included for completeness. Typical location of the sensors is also shown.

The operation of a sub-critical ORC encompasses several aspects that are common to refrigeration systems, such as the control of the subcooling at the condenser outlet and the degree of superheating at the evaporator outlet, single- and two-phase heat transfer in the heat exchangers, as well as working fluid leaks. The amount of working fluid charged into the system has a direct effect on the performance, the flexibility, as well as on the operational costs of the unit. The costs associated with the working fluid and possible working fluid losses due to leaks can be significant, especially in large-scale systems (up to 10% of the plant cost) [12].

With respect to the system behavior, when the operating conditions shift from design to off-design, the working fluid charge migrates within the system, altering the transition between single-to two-phase heat transfer regimes inside the heat exchangers. If the system is over-charged, typically the length (or, equivalently, the spatial fraction) of the subcooled zone in the evaporator increases, reducing both the capacity and degree of superheating. However, if the system is under-charged, one of the drawbacks is the decrease of subcooling at the condenser outlet, which may cause the pump to cavitate [13]. Additionally, a liquid receiver can be installed before the pump, as shown in Figure 1, in order to compensate for charge migration during off-design conditions, especially in large scale ORC systems. Nevertheless, a non-proper charge level may completely empty or flood the liquid receiver.

The estimation of the optimum charge level is particularly challenging when additional aspects of the normal ORC operation are considered. In particular, lubricant oil may be pre-mixed with the working fluid to ensure the correct functioning of positive displacement expanders. In other cases, oil-separators are utilized to separate the working fluid from the lubricant oil at the expander outlet. However, oil entrainment in the vapor and oil solubility affect the distribution of charge in the system. Furthermore, non-condensable gases can also be present inside the system, which can impact the operation of the ORCs by altering the condensing pressure [14].

Based on the aforementioned reasoning, appropriate numerical methodologies are required to predict the performance of the system and to account for the total working fluid charge in a deterministic way such that the system model resembles the actual behavior of the system. Experimental results are also necessary to validate such numerical methods and to conduct additional analyses. Therefore, in the following sections, the state-of-the-art numerical methods utilized to simulate ORC systems are described and their limitations are also highlighted.

3. Charge-sensitive ORC modeling approaches

Modeling of ORC systems is typically performed for steady-state conditions using purely thermodynamic considerations. Assumptions are made regarding the minimum temperature difference (i.e., the pinch point) between the working fluid and the heat source and heat sink fluids. The degrees of subcooling at the pump inlet and superheat at the expander inlet are set, as well as the condensing and evaporating temperatures. Constant values of isentropic efficiency for the pump and expander are often assumed. Then, the cycle is solved by computing the enthalpy of the working fluid at each state point and assuming the pump and expander to be adiabatic. This type of simplified model works well as an initial working fluid screening tool and to broadly assess cycle performance trends. However, in such a model, the physical characteristics of the system components are not taken into account. It follows that a more detailed model is necessary to study the behavior of a real ORC system in off-design conditions. Such a detailed model also typically estimates the steady-state response of the system, although dynamic models have been constructed from similar considerations [15]. For the scope of this chapter, the following discussions will be limited to steady-state modeling.

3.1. Definition of a mechanistic model for ORC system

Detailed off-design models of ORC systems are built by connecting together different sub-models for each component of the system and are implicitly solved by driving to zero a number of residuals to ensure that the solution is within a physical domain. Such detailed models are based upon the extensive literature available for vapor compression cooling and heat pumping systems [16, 17]. A recent overview of off-design steady-state performance studies and modeling applied to ORC systems can be found in [8].

Generally, an ORC system model can be regarded as mechanistic if the inputs and the known parameters are similar to those an ORC operator would know in practice. The model then uses physical principles and empirical or semi-empirical component models to simulate system performance given these parameters, inputs and boundary conditions. An ORC system model requires knowledge of the pump and expander displacements, the heat exchangers geometry, and the total system volume. Such a model receives as inputs the working fluid type, the inlet temperature and mass flow rate of the source and sink fluids, the rotational speeds of the pump and expander, and the condenser exit subcooling or the total working fluid charge in the system. A truly mechanistic model is charge-sensitive, meaning the total refrigerant charge is known, but the condenser exit subcooling is determined by the ORC operation [5]. In fact, if the condenser exit subcooling is fixed within the simulation, an assumption is made regarding the system state, which is not known in a real system. Outputs of the model include rotating equipment efficiency, system efficiency, heat transfer rates, condensing and evaporating temperatures, and net power production. In addition, a Second Law analysis can also be applied to estimate the system irreversibilities [18].

A detailed ORC model such as the one described above can give information about off-design performance, working fluid charge sensitivity, migration of charge between operating points, pressure drops, and the impact of working fluid transport properties on the heat transfer processes. Although, it only applies to the exact system being modeled, it provides a general framework (see for example ORCmKit library [19]) whereby other systems can be simulated by adjusting the geometric parameters of the model.

While constructing a charge-sensitive model, a number of challenges, both numerical and thermophysical, need to be considered:

• Since the entire ORC model is obtained by connecting several sub-models representing the system components and a number of residuals are necessary to ensure a physically-meaningful solution, a robust solution scheme is required to perform charge-sensitive simulations with reasonable computational efforts.

• Due to the presence of multiple heat transfer mechanisms, pressure drops, and flow regimes inside the heat exchangers, it is important to properly identify the convective heat transfer coefficients in multi-zone heat exchangers and understand how the working fluid charge can improve the reliability of estimating such coefficients.

• The charge estimation inside the heat exchangers relies on the proper knowledge of spatial fraction occupied by each zone (i.e., single-phase liquid or vapor, two-phase) as well as the estimation of the working fluid density. In the case of single-phase zones, such calculation is straightforward. However, in the case of two-phase zones, the density depends on the pressure, temperature, quality, and flow regime through the void fraction. Void fraction models directly impact the charge estimation.

• Presence of a liquid receiver in the system.

• Difficulties in accounting all the system volumes (e.g., valves, sight glasses, filters, sensors, and line sets).

• Estimation of working fluid dissolved in lubricant oil.

• Correctness of the model in predicting the total working fluid charge with different working fluids.

• Uncertainties associated with the exact amount of working fluid charge present in the actual system.

In order to address the majority of these modeling aspects, general guidelines for developing charge-sensitive models are discussed in the following section. Among these aspects, heat transfer correlations and the length of each zone inside the heat exchangers account for the majority of the inaccuracy in charge estimations. However, it will be shown in Section 4.3 that the solubility of working fluid in oil can also be significant. A charge tuning scheme can be applied to improve the accuracy of charge estimation, as outlined in Section 3.3. Experimental results from a test case will be used to address the last three bullet points in Section 4.

3.2. Overall charge-sensitive model description

As aforementioned, the individual components are properly arranged to form the overall cycle model. By referring to the ORC schematic of Figure 1, two general approaches can be identified to assemble the cycle models, which are shown as block diagrams in Figures 2 and 3, respectively. By looking at the two diagrams, it appears evident that the major cycle components, that is, heat exchangers, pump and expander, are arranged in the same fashion. The major difference is the approach adopted to account for the line sets connecting the active components and, therefore, the associated pressure drops and heat losses. The approach shown in Figure 2 retains the physical description of the cycle with the line set sub-models placed between two consecutive cycle components. It follows that an accurate description of each line set, including internal and external diameters as well as equivalent length, is required. Whereas, in the block diagram of Figure 3, the line set losses are lumped into two single artificial components representing high- and low-pressure sides of the systems. These lumped components are placed at the outlet of the evaporator and the condenser, respectively.

Figure 2.

Block diagram showing how the individual component models are arranged to form the overall ORC model. Variables are updated as they pass through each component (adapted from [7]).

Figure 3.

Block diagram showing the solver architecture of the charge-sensitive ORC model proposed by Dickes et al. [8]. The inputs are indicated in blue, the outputs in brown, the parameters in green, the iteration variables in red and the model residuals in violet.

Each block diagram also determines the required initial guess values and the residual to be driven to zero to solve the system. In particular, in the case of the block diagram of Figure 2, the five initial guesses are the pump inlet specific enthalpy ( $h_{\mathit{\text{pump}},\mathit{suc}}$ ), the pump inlet and outlet pressures ( $p_{\mathit{\text{pump}},\mathit{suc}}$ , $p_{\mathit{\text{pump}},\mathit{dis}}$ ), the expander inlet specific enthalpy ( $h_{\text{exp},\mathit{suc}}$ ), and the expander outlet pressure ( $p_{\text{exp},\mathit{suc}}$ ). Such guesses are chosen because of the presence of the regenerator in the cycle and to ensure stability of the code. Given the set of guess values, an iteration of the cycle model consists of solving the components in the following order: pump and pump discharge line set, expander and expander discharge line set, regenerator and regenerator exit line sets for liquid and vapor sides, evaporator and expander suction line set, and condenser and pump suction line set. The solution of the cycle model is enforced by five residuals. In particular, four residuals ensure the continuity of thermodynamic states between two consecutive components. These are:

$$\begin{array}{l}{\text{Δ}}_1=h_{pump,suc}-h_{lineset,pump,out};{\text{Δ}}_2=p_{pump,suc}-p_{lineset,pump,out}\hfill \\ \begin{array}{ll}{\text{Δ}}_3=h_{\exp ,suc}-h_{lineset,\exp ,out};\hfill & {\text{Δ}}_4=p_{\exp ,suc}-p_{lineset,\exp ,out}\hfill \end{array}\hfill \end{array}$$ E1

The fifth residual is associated with the total system charge, which is estimated by summing all the charge contributions:

$${\Delta }_{5,\mathit{ch}\mathit{arge}}=m_{\mathit{\text{evap}}}+m_{\mathit{\text{cond}}}+m_{\mathit{\text{regen}}}+m_{\mathit{\text{pump}}}+m_{\exp }+m_{\mathit{tank}}+\sum _i{m}_{\mathit{\text{lineset}},\mathit{i}}-m_{\mathit{wf}}$$ E2

where the mass of working fluid for each component is given by $m_{\mathit{wf},\mathit{\text{component}}}=V_{\mathit{\text{component}}}{\rho }_{\mathit{wf}}$ . Particular attention has to be given to computing the averaged density of the working fluid, ${\rho }_{\mathit{wf}}$ , to account for single-, two-phase conditions for both pure working fluids or mixtures, as will be discussed in more detail in Section 3.3. Furthermore, it can also be mentioned that the fifth residual associated with the total charge, ${\Delta }_{5,\mathit{\text{charge}}}$ , is usually replaced by the difference between the imposed and calculated condenser exit subcooling in many off-design simulation models proposed in literature.

In the case of the block diagram of Figure 3, the model iterates on the evaporator outlet specific enthalpy ( $h_{\mathit{\text{evap}},\mathit{ex}}$ ), pump inlet and outlet pressures ( $p_{\mathit{\text{pump}},\mathit{su}}$ , $p_{\mathit{\text{pump}},\mathit{ex}}$ ), and condenser outlet subcooling ( $\Delta T_{\mathit{sc}}$ ). The residuals to be minimized are the condenser outlet specific enthalpy, the evaporator outlet specific enthalpy, the expander rotational speed, and the total system charge:

$$\begin{array}{lll}{\Delta }_1=1-\frac{h_{\mathit{cd},\mathit{ex},2}}{h_{\mathit{cd},\mathit{ex},2}};& {\Delta }_2=1-\frac{h_{\mathit{ev},\mathit{ex},2}}{h_{\mathit{ev},\mathit{ex},2}};& \begin{array}{ll}{\Delta }_3=1-\frac{N_{\text{exp},2}}{N_{\exp }};& {\Delta }_4=1-\frac{m_{\mathit{wf},2}}{m_{\mathit{wf}}}\end{array}\end{array}$$ E3

where the total working fluid charge is calculated by employing an analogous expression to that one appearing in Eq. (2).

In both solution schemes, the resulting problem to be solved is multi-dimensional, leading to possible convergence issues when the charge is imposed. To increase the robustness of the model, a multi-stage solver can be employed to run the simulations. Dickes et al. [8] described these ad-hoc solvers in detail. The thermophysical properties of the working fluid can be retrieved either from CoolProp [20] or REFPROP [21].

In order to simulate the entire cycle, different sub-models are required to characterize each of the components. In the context of charge-sensitive modeling, the heat exchangers typically have the largest volumes in the system and are subjected to varying flow regimes, as outlined in Section 3.3. However, pump, expander, liquid receiver, and line sets also require proper modeling. Particular emphasis is given to estimate the working fluid charge in each of these components. For an extensive description and examples of these models, the reader is invited to referred to [5, 6, 7, 8, 19].

As a general overview, the rotating equipment is modeled by using one or more of the following approaches: (i) a performance map based on experimental data (black box models); (ii) physics-based models with empirically determined parameters (gray box models); (iii) entirely physics-based models (white box models). Semi-empirical models are usually preferred as a compromise between physical-characteristics (e.g., under- over-expansion/compression, pressure drops, heat transfer and mechanical losses) and computational cost. Regarding the charge estimation of rotating equipment, the mass of working fluid is computed by knowing the internal volume of the machine and by computing an average density. This simplified approach is usually reasonable given the relatively small internal volumes of pumps and expanders compared to heat exchangers and other volumes of the system. However, especially in the case of oil lubricated positive displacement expanders, high pressures and temperatures may lead to dissolved refrigerant in the oil that could be accounted for using solubility data of the refrigerant-oil mixture [16, 22, 23].

A liquid receiver (or buffer tank) can be installed in an ORC system at the condenser outlet. Under normal operating conditions, the liquid receiver ensures a saturated liquid at the condenser outlet and serves as a mass damping device during off-design conditions. In most cases, the liquid receiver is modeled by neglecting heat losses to the environment as well as potential hydrostatic effects due to the height of liquid (i.e., the pressure is considered to be uniform inside the liquid receiver). Nevertheless, it is not straightforward to predict the mass of liquid stored inside the tank at any operating condition with a steady-state model. In an experimental setup, a level sensor could be used to monitor the liquid level inside the tank. However, in commercial ORC systems, this may not be a viable option. Recently, Dickes et al. [8] proposed a liquid receiver model by introducing four hypotheses that resulted in the following constraints for calculating the mass of working fluid inside the liquid receiver:

$$m_{\text{receiver}}=\left\{\begin{array}{l}{V}_{\text{tank}}{\rho }_{\text{tank}}\\ V_{\text{tank}}\left[L_{\text{tank}}{\rho }_{\mathrm{sat},\mathrm{l}}+\left(1-L_{\text{tank}}\right){\rho }_{\mathrm{sat},\mathrm{l}}\right]\\ V_{\text{tank}}{\rho }_{\mathrm{sat},\mathrm{v}}\\ V_{\text{tank}}{\rho }_{\text{tank}}\end{array}\right.\begin{array}{l}\mathit{\text{if}} h_{\text{tank},\mathrm{su}}<h_{\mathrm{sat}}\left(p_{\text{tank},\mathrm{su}},x=0\right)\\ \mathit{\text{if}} x_{\text{tank},\mathrm{su}}=0\\ \mathit{\text{if}} x_{\text{tank},\mathrm{su}}\in \left[0,1\right]\\ \mathit{\text{if}} h_{\text{tank},\mathrm{su}}>h_{\mathrm{sat}}\left(p_{\text{tank},\mathrm{su}},x=1\right)\end{array}$$ E4

where $V_{\mathit{\text{tank}}}$ is the total volume of the liquid receiver, $L_{\mathit{\text{tank}}}$ is the liquid level, ${\rho }_{\mathit{sat},l}$ and ${\rho }_{\mathit{sat},v}$ are the saturated liquid and vapor densities at the supply pressure, respectively. Under normal operating conditions and optimal charge, the receiver should be partially filled with liquid. If the charge is not proper or the ORC system is operating at strong off-design conditions, the receiver can be either full or empty of liquid. Furthermore, such mathematical formulation may result in numerical issues due to a discontinuity in the mass estimation when the fluid reaches its saturated liquid state (i.e., $x_{\mathit{\text{tank}},\mathit{su}}=0$ ) [8].

The line sets (or pipelines) connect the different system components and carry the working fluid. The line sets consist of piping, valves, fittings, sight glasses, filters and other elements and are also associated with pressure drops and heat losses. As mentioned at the beginning of this section, two different modeling approaches can be employed to estimate the working fluid charge carried by the line sets as well as the thermodynamic states at inlet and outlet of each line set. By referring to Figure 2, a total of six line sets is considered. In particular, each line set can carry single-phase or two-phase working fluid and they are modeled as an equivalent tube having inner and outer diameters. An equivalent length that accounts for straight sections ( $L_{\mathit{\text{straight}}}$ ) and all the fittings ( $L_{\mathit{eq},\mathit{\text{fittings}}}$ ) is computed by using the method of loss coefficients proposed by Munson et al. [24]. Distinction is made in the calculation of pressure drops to account for single and two-phase flow conditions. In the case of single-phase flow, the pressure drop for a certain line set having a certain internal diameter is computed by introducing the Darcy friction factor:

$$\Delta p_{\mathit{\text{lineset}},1\phi ,\mathit{ID}=\mathit{\text{const}}}=f\frac{\sum L_{\mathit{\text{straight}}}+\sum L_{\mathit{eq},\mathit{\text{fittings}}}}{\mathit{ID}}\left(\frac{{\rho }_{\mathit{\text{lineset}},\mathit{su}}{v}^2}{2}\right)$$ E5

where $v$ is the working fluid velocity. Under two-phase conditions, the Lockhart-Martinelli method for two-phase frictional pressure drop in tubes [25] is used. That is:

$$\Delta p_{\mathit{\text{lineset}},2\phi ,\mathit{ID}=\mathit{\text{const}}}=-{\left(\frac{\mathit{dp}}{\mathit{dz}}\right)}_F\left(\sum L_{\mathit{\text{straight}}}+\sum L_{\mathit{eq},\mathit{\text{fittings}}}\right)$$ E6

The heat losses though the line sets to the surroundings can be calculated with the effectiveness-NTU method [17], as given in Eq. (7), where ${\mathit{NTU}}_{\mathit{\text{lineset}}}$ can be either specified or obtained by regression using experimental data.

$$T_{\mathit{\text{lineset}},\mathit{ex}}=T_{\mathit{amb}}-e^{\left(-{\mathit{NTU}}_{\mathit{\text{lineset}}}\right)}\left(T_{\mathit{amb}}-T_{\mathit{\text{lineset}},\mathit{su}}\right)$$ E7

A more simplified approach is proposed in Figure 3, where the line sets and associated losses are lumped in the high- and the low-pressure lines by using single fictitious components placed at the outlet of the evaporator and the condenser, respectively. The pressure losses are computed as a linear function of the working fluid kinetic energy with two coefficients to be calibrated using experimental data. That is:

$$\Delta p_{\mathit{\text{lineset}}}=K\left(\frac{m_{\mathit{wf}}^2}{{\rho }_{\mathit{\text{lineset}},\mathit{su}}}\right)+B$$ E8

The ambient losses are modeled by introducing an overall heat transfer coefficient of the line set ( ${\mathit{UA}}_{\mathit{\text{lineset}}}$ ), which is the sum of the thermal resistance of the pipe, the thermal resistance associated with insulation and the convective $\mathit{UA}$ values associated with the inside and outside of the line set:

$$Q_{\mathit{\text{lineset}}}={\mathit{UA}}_{\mathit{\text{lineset}}}\left(T_{\mathit{\text{lineset}},\mathit{su}}-T_{\mathit{amb}}\right)$$ E9

$${\mathit{UA}}_{\mathit{\text{lineset}}}={\left(\frac{1}{{\mathit{UA}}_{\mathit{ID}}}+\frac{1}{{\mathit{UA}}_{\mathit{OD}}}+R_{\mathit{\text{lineset}}}+R_{\mathit{\text{insulation}}}\right)}^{-1}$$ E10

The working fluid charge of each line set is calculated in a manner that is similar to the approach that will be presented for each zone of the heat exchangers (see Section 3.3). In the case of single-phase working fluid, the mass of working fluid is calculated as

$$m_{\mathit{\text{lineset}}}=V_{\mathit{\text{lineset}}}{\rho }_{\mathit{\text{lineset}}}$$ E11

where ${\rho }_{\mathit{\text{lineset}}}$ is the average between inlet and outlet densities of a line set. If a line set is considered adiabatic, then ${\rho }_{\mathit{\text{lineset}}}\equiv {\rho }_{\mathit{\text{lineset}},\mathit{su}}$ . Under two-phase conditions, the average density of the working fluid is obtained by introducing an average void fraction over the length of the line set ( ${\alpha }_{\mathit{\text{lineset}}}$ ). Mathematically, this can be expressed as:

$${\rho }_{\mathit{\text{lineset}}}={\alpha }_{\mathit{\text{lineset}}}{\rho }_{\mathit{sat},v}+\left(1-{\alpha }_{\mathit{\text{lineset}}}\right){\rho }_{\mathit{sat},l}$$ E12

3.3. Working fluid charge in heat exchangers

An ORC system operates between a heat source and heat sink by using heat exchangers (HEX) for the heat input and rejection processes, respectively. The working fluid undergoes phase changes in both evaporator and condenser and different phases co-exist during these processes. An internal heat exchanger (or regenerator) can be present in the system to improve the cycle thermal efficiency [11]. The heat exchangers account for the majority of the volume in an ORC system and therefore, it is particularly important to have reliable models to predict the heat transfer rate, zone lengths, and working fluid charge.

To predict the different zones inside an HEX, three different numerical methods can be used: moving-boundary [7, 8], finite-volume (or discretized) [16], and hybrid approach [26]. The main advantage of moving-boundary models is their fast computation time and good accuracy. Finite-volume models provide greater spatial resolution but can be computationally expensive. Hybrid models represent a compromise between the former methods. By considering the moving-boundary model as an example, the heat exchanger can be divided in multiple zones, as shown in Figure 4. Each zone is characterized by a heat transfer coefficient and heat transfer area through which a certain heat transfer process occurs. By knowing the inlet conditions of both streams, the effective heat transfer rate between the hot and cold streams is calculated by enforcing that the total surface area occupied by the different zones is equal to the geometrical surface area of the heat exchanger [8]. It follows that the total working fluid mass inside the heat exchanger ( $m_{\mathit{HEX}}$ ) can be calculated as the sum of the working fluid mass associated with each zone:

$$m_{\mathit{wf},\mathit{HEX}}=\sum _{i=1}^N{\rho }_i\left(w_i{V}_{\mathit{HEX}}\right)$$ E13

where $w_i$ is the volume fraction occupied by the i-th zone in the heat exchanger, and ${\rho }_i$ is the mean density of working fluid in the i-th zone. In order to ensure a correct estimation of the working fluid charge, the knowledge of the zone division (or in other words, correct estimation of the heat transfer rates) and the evaluation of the working fluid density under both single- and two-phase conditions are required. In particular, the latter implies the selection of the proper void fraction model that characterizes both boiling and condensing processes.

Figure 4.

Schematic of counterflow HEX model simulated by moving boundary method showing phase boundaries (solid lines) and zone boundaries (dashed lines).

Various heat exchanger types are employed in ORC systems (e.g., brazed plate heat exchangers, fin and tubes, shell and tubes, etc.) and a number of state-of-the-art heat transfer correlations are available in the literature that are used to predict the convective heat transfer coefficients in each zone. However, these correlations are often purely empirical and calibrated to fit experimental data and multiple correlations may be available for the same type of heat exchanger and working fluid. Furthermore, these correlations are typically derived for refrigeration systems, which present different operating conditions than ORCs. It follows that such correlations are likely to over-predict or under-predict the actual heat transfer coefficient within a certain zone, even though the overall predictions of the heat transfer rate may be reasonable. Nevertheless, errors in predictions of heat transfer within individual zones can lead to significant errors in charge estimations [8]. To improve the accuracy in predicting the heat transfer coefficients in the different zones, empirical parameters of these correlations should be fitted with experimental results by using identification methods [8].

The accuracy of the heat transfer predictions alone does not ensure the correctness of the charge estimation. In fact, as already indicated in Eq. (14), under two-phase conditions, the density is a function of the thermodynamic conditions as well as the flow pattern characterizing each phase. The void fraction, $\alpha$ , is an essential parameter associated with two-phase flow and it is related to the fluid quality, $x$ . Based on the void fraction, the working fluid mass of a two-phase zone having length $L$ and cross section $A_c$ is calculated as:

$$m_{\mathit{wf},2\phi }=m_{\mathit{wf},v}+m_{\mathit{wf},l}=A_c\left({\rho }_v\underset{0}{\overset{L}{\int }}\alpha \left(x\right)\mathit{dl}+{\rho }_l\underset{0}{\overset{L}{\int }}\left[1-\alpha \left(x\right)\right]\mathit{dl}\right)$$ E14

To be further noted is that the void fraction $\alpha \left(x\right)$ is integrated along the length of the zone as a function of the working fluid quality. A common simplification is to assume a uniform heat flux in the two-phase zone, that is, a linear evolution of the quality along the length of the zone. However, a more appropriate approach would require evaluation of the effective spatial evolution of the quality in the heat exchanger for example by discretizing the two-phase zone into multiple sub-cells (e.g., by employing a hybrid heat exchanger model) [8]. Similar to the convective heat transfer correlations, several void fraction models can be found in the scientific literature [8, 10] that lead to different values of the mean density under two-phase conditions. An example of investigation of the combined effect of different convective heat transfer correlations and void fraction models on the total charge estimation is shown in Figure 5 that is based on the work done by Dickes et al. [8].

Figure 5.

Working fluid charge predictions for different void fraction models as well as three adjusting methods to tune the convective heat transfer correlations [8].

The study considered a total of 360 charge-sensitive simulations including 40 experimental data points, 4 tuning methods, 5 void fraction correlations, and 72 heat exchanger models. It was concluded that the total working fluid charge estimation is mainly affected by the adjusting (or tuning) method applied to the convective heat transfer coefficient rather than the void fraction model. However, a trade-off between the two aspects exists since higher accuracy in the thermal performance predictions led to a large scatter of the charge inventory predictions.

4. Experimental considerations on working fluid charge with pure working fluids and zeotropic mixtures

In Section 3, the numerical aspects to be considered while developing a charge-sensitive ORC model have been discussed. However, in order to understand the usefulness of the model in predicting a real ORC behavior, it is necessary to consider experimental data. To this end, in the following sub-sections, several aspects related to analyzing different charge levels using both experimental and charge-sensitive model results are presented.

4.1. Data reduction and charge uncertainties

Steady-state operation of an ORC system can be detected by adopting for example the conditions proposed by Woodland et al. [15]. Alternatively, the methodology described by Dickes et al. [5] could also be employed. Due to the fact that experimental data is subject to different uncertainties, possible errors or sensor malfunction, a post-treatment data can be performed to identify possible outliers by using the open-source GPExp library [27]. Furthermore, a reconciliation method can also be applied to correct the measured data affected by measurement error propagation to satisfy physical system constraints (e.g., energy balances) [28].

Uncertainty propagation from the measured system properties to the calculated quantities is usually performed according to the method described by Figliola and Beasley [29]. In particular, it is important to discuss the uncertainties associated with the total system charges, which are often not thoroughly accounted for when reporting experimental results of ORC systems. In general, there are several sources of uncertainty in the measurement of the total system charge, such as inaccuracy of the scale used to weight the supply cylinder, leakage of working fluid from the system during operation, small portions of working fluid occupying charging hoses, small losses in system charge due to purging of charging hoses, stacked uncertainties due to the summing of several incremental charging steps. In experimental test rigs, the uncertainty due to purging of system charging hoses can be eliminated by leaving all charging hoses connected and full of working fluid during all tests and incremental charge conditions. Leakage of working fluid from the system during operation can be neglected after ensuring that repeated tests at the same conditions, but several days apart, yield approximately the same condenser subcooling. The main sources of uncertainty in the working fluid charge level are the inaccuracy of the scale and the stacked uncertainties due to several incremental steps in adding or subtracting charge. As an example, if a scale has a known resolution of ±0.005 kg, and a measurement uncertainty of ten times the resolution or ± 0.045 kg is conservatively assumed, the uncertainty of the total system charge ( $u_{\mathit{\text{charge}}}$ ) can be estimated as a function of the number of consecutive charging steps, $n_{\mathit{\text{steps}}}$ , as:

$$u_{\mathit{\text{charge}}}=\sqrt{{\left(0.045\right)}^2·n_{\mathit{\text{steps}}}}$$ E15

As reported by Woodland [6], given a system charge of 18.1 kg, the relative uncertainty in system charge assuming a maximum of 8 charge steps is still less than 1%. Therefore, a high degree of confidence can be placed on the measured system charge level. More challenging is determining the actual charge inside the heat exchangers, which would require active measurements of their weight during operation.

In the case of a zeotropic mixture, for example, R134a-R245fa (0.625–0.375), uncertainty in the charge level of each component of the binary mixture results in some uncertainties in the concentration of the mixture, which is an important parameter in computing the thermophysical properties of the mixture. However, even with the conservative estimate on the charge uncertainty previously mentioned, the resulting relative uncertainty in the concentration is ≤0.004, which results in a relative uncertainty of 1% or less, depending on which component concentration is considered. Depending on the mixture considered, the sensitivity of the mixture thermodynamic properties to the concentration can be much smaller for single-phase states than the sensitivity to temperature. For example, in the case of R134a-R245fa, the sensitivity of enthalpy to the concentration uncertainty is two orders of magnitude less than the sensitivity to the temperature uncertainty for single-phase states [6]. As a result, the measured concentration of the working fluid is not a dominant factor in measurement uncertainty in this case.

4.2. Effect of charge variations in an ORC system

To conduct a working fluid charge study, it is desired to explore a wide range of working conditions, different working fluids, and several charge levels. In the current available literature, only the study published by Woodland [6] meets all of these criteria, and, therefore, it is used hereafter as an example. In particular, three working fluids were chosen to investigate low grade waste heat recovery by means of an ORC: R134a, R245fa, and the zeotropic mixture R134a-R245fa (0.625–0.375). In addition, two heat source temperatures, and several charge levels of the working fluids were considered. The range of conditions tested experimentally for each working fluid is summarized in Table 1. A complete description of the ORC system can be found in [6] and its schematic is shown in Figure 1. The effect of the working fluid charge on the ORC system performance is shown in Figure 6, where the Second Law efficiency, defined as the ratio of the system net power to the total exergy rate available, is plotted as a function of the subcooling level for different operating conditions. It can be seen that for a given working fluid and heat source inlet temperature, there exists an optimum charge level that maximizes the Second Law efficiency (i.e., power output) of the system. Furthermore, for each working fluid, the efficiency tends to be higher at lower subcooling due to lower condensing pressures dictated by lower fluid charge.

Working fluid	Charge [kg]	Source temperature [°C]	Expander speed range [rpm]	Pump speed range [rpm]
R134a	18.1	80, 100	2800–3000	1000
	20.4	80, 100	1200–3000	1000
	22.7	100	1500–2300	1000
R245fa	27.4	80, 100	1000–2000	300–500
	29.5	80, 100	1000–2000	300–500
	31.8	80, 100	1000–2000	300–500
R134a-R245fa (0.625–0.375)	24.9	80, 100	1000–1500	500–1000
	27.2	80, 100	1000–2500	500–1000
	29.5	80, 100	1000–2200	500–1000

Table 1.

Summary of experimental test matrix to investigate the effect of working fluid charge [6].

Figure 6.

Measured second law efficiency versus condenser exit subcooling for each working fluid. The contour shows the working fluid charge and the size of the markers is proportional to the expander rotational speed (data from [6]).

4.3. Total working fluid charge model predictions

One of the main purposes of a charge-sensitive model is to predict the total charge of a certain ORC system for which it has been developed and its behavior (e.g., evaporating and condensing temperatures, subcooling and superheating). This needs to hold true also under the circumstances that the system is under-charged or over-charged as well as if the working fluid is changed, that is, drop-in replacement. To this end, the same ORC system considered in Section 4.2 is used as an example [6]. The charge-sensitive model is employed to simulate three different working fluids (i.e., R134a, R245fa, and R134a-R245fa (0.635–0.375)) at three different charge levels, without any tuning. The untuned model prediction of the system charge versus the measured charge for all the test points is shown in Figure 7. It can be seen that even without tuning, the model predictions are reasonably accurate due to the fact that the heat exchanger volumes as well as the liquid line sets have been accounted for. These liquid volumes remain constant across all experimental tests and therefore set the baseline charge level of the system. However, the model does not predict the variability in charge well for all cases. For R134a, the trend in system charge is captured well, but with increased underprediction of the charge at higher charge levels. For R245fa, the trend is not captured well. A general increase is observed, but the underprediction is significant. For the zeotropic mixture, no trend is visible at all in the charge prediction. Variability in the predicted charge is due to the changing liquid lengths in the heat exchangers, the void fraction of two-phase flow, and the solubility of the working fluid in the lubricant oil. In particular, if the boundary conditions of the ORC are fixed then an increase in the working fluid charge in the system leads to a decrease in the two-phase zone in the condenser (resp. subcooled zone), as shown in Figure 8(a). If the void fraction model and the heat transfer coefficient of the two-phase zone are not estimated correctly, the charge of the subcooled region can be largely under predicted due to the rapid increase in liquid level inside the condenser, as shown in Figure 8(b). Solubility of the working fluid in the oil was not considered in the model. However, it can have a significant impact on prediction of the total system charge. In fact, if an oil separator is present in the ORC system, the working fluid dissolved in the lubricant oil could be up to 10% of the total charge depending on the solubility of the working fluid-lubricant oil mixture [6].

Figure 7.

Parity plot showing the ability of the charge-sensitive model to predict the total working fluid charge without any charge tuning.

4.4. Charge tuning schemes

As shown in Figure 8, the estimation of working fluid charge in system modeling is usually biased due to inaccurate estimation of system volumes, ambiguous flow patterns under two-phase flow conditions, solubility of working fluid in the lubricant oil, etc. A charge-tuning scheme can be used to eliminate the bias. Although such methods have been widely used for vapor compression heat pumping cycles [16, 17], they have not been employed in charge-sensitive ORC system modeling. An empirical two-point charge tuning equation that incorporates a linear function of the length of the subcooled section in the condenser ( $L_{\mathit{\text{cond}},\mathit{sc}}$ ) can be employed:

$$m_{\mathit{wf},\mathit{\text{meas}}}-m_{\mathit{wf},\mathit{\text{calc}}}=C_{m,0}+C_{m,1}{L}_{\mathit{\text{cond}},\mathit{sc}}$$ E16

where $C_{m,0}$ and $C_{m,1}$ are the regression parameters to be estimated by linear regression. The estimated correction term from Eq. (16) is added to the total working fluid charge computed by the model to give an unbiased estimate [17].

Figure 8.

Heat exchanger length fraction (a) and breakdown of component charge (b) versus condenser exit subcooling when the total charge is imposed as model input. Pump speed = 1000 rpm, expander speed = 2000 rpm, source temperature = 100°C, source and sink fluid flow rates are 0.45 and 1 kg/s, respectively. The charge calculation is untuned.

5. Conclusions

In this chapter, the effects of the total working fluid charge on the performance and operation of an ORC system have been discussed. In particular, both numerical and experimental have been used to demonstrate the importance of considering the charge inventory when analyzing off-design conditions of the ORC system. Under-charging or over-charging the system has a direct effect on the spatial distribution of the different zones inside the heat exchangers. Additional aspects such as oil solubility can also have non-negligible impact on the system performance. A detailed mechanistic model of the ORC allows identification of optimal charge level. However, charge-tuning is always required to account for inaccuracies and uncertainties associated with the modeling assumptions. In fact, two-phase heat transfer correlations and void fraction models directly affect the charge inventory estimations. The analysis of the charge inventory of ORC systems is also particularly important when considering drop-in replacements such as HFOs and their blends. The size of the heat exchangers as well as the piping system require an optimal working fluid charge to ensure proper behavior of the ORC at both design and off-design conditions.

Conflict of interest

The Authors declare no conflict of interest.

Nomenclature

A	area, [m2]
c	tuning coefficient, [−]
C	regression coefficient, [−]
h	specific enthalpy, [kJ/kg]
m	mass, [kg]
N	rotational speed, [rpm]
Nu	Nusselt number, [−]
p	pressure, [Pa]
Q	heat transfer rate, [W]
Subscripts
calc	calculated
cond	condenser
evap	evaporator
ex	exit
exp.	expander
R	Thermal resistance, [K/W]
Re	Reynolds number, [−]
T	temperature, [K]
UA	heat conductance, [W/K]
V	volume, [m3]
x	quality, [−]
α	void fraction, [−]
ρ	density, [kg/m3]
meas	measured
sc	subcooled
su	supply
wf	working fluid