**Thermal Systems**

Editors

**Ivan CK Tam Brian Agnew**

MDPI • Basel • Beijing • Wuhan • Barcelona • Belgrade • Manchester • Tokyo • Cluj • Tianjin

*Editors* Ivan CK Tam Associate Professor in Marine Engineering Design & Technology Newcastle Research & Innovation Institute (NewRIIS) Newcastle University Singapore

Brian Agnew Professor of Energy and the Environment Newcastle Centre for Railway Research (NewRail) Newcastle University UK

*Editorial Office* MDPI St. Alban-Anlage 66 4052 Basel, Switzerland

This is a reprint of articles from the Special Issue published online in the open access journal *Energies* (ISSN 1996-1073) (available at: https://www.mdpi.com/journal/energies/special issues/TPS).

For citation purposes, cite each article independently as indicated on the article page online and as indicated below:

LastName, A.A.; LastName, B.B.; LastName, C.C. Article Title. *Journal Name* **Year**, *Volume Number*, Page Range.

**ISBN 978-3-03943-841-9 (Hbk) ISBN 978-3-03943-842-6 (PDF)**

© 2021 by the authors. Articles in this book are Open Access and distributed under the Creative Commons Attribution (CC BY) license, which allows users to download, copy and build upon published articles, as long as the author and publisher are properly credited, which ensures maximum dissemination and a wider impact of our publications.

The book as a whole is distributed by MDPI under the terms and conditions of the Creative Commons license CC BY-NC-ND.

## **Contents**


### **Younghyeon Kim, Seokyeon Im and Jaeyoung Han**


### **About the Editors**

**Ivan CK Tam** (Associate Professor, Dr.) is an Associate Professor at the Newcastle University with a strong track record of leading innovative and design projects. He has research interest in the combustion process, exhaust emission control, energy management and renewable energy technology. His recent research interest is the application of cryogenic technology in the use of liquefied natural gas and the organic Rankine cycles.

**Brian Agnew** (Prof. Dr.) joined the Department of Mechanical Engineering at Newcastle University in 1984 with research interest in heat transfer, internal combustion engines and thermal systems. Subsequently, he was appointed as Professor of Energy and the Environment at the School of the Built Environment, Northumbria University. He will continue as a Guest Member of Staff at Newcastle University until September 2022.

### *Editorial* **Thermal Systems—An Overview**

**Ivan CK Tam 1,\* and Brian Agnew <sup>2</sup>**


Received: 3 December 2020; Accepted: 30 December 2020; Published: 31 December 2020

We live in interesting times in which life as we know it is being threatened by human-made changes to the atmosphere we live. On the global scale, concern is focused on climate change due to greenhouse gas emissions and atmospheric pollution produced by combustion processes. The increase in global warming, added to the scarcity of fossil fuels, has motivated the development of new technologies to improve the efficiency of existing processes in power plants. To meet the dual challenges presented by these factors, consideration needs to be given to energy efficiency and pollution reduction in transport and energy conversion processes. A possible approach is through development of new ideas and innovative processes to current practices. Among the available options, multi-generation processes such as trigeneration cycles, battery storage systems, solar power plants and heat pumps have been widely studied as they potentially allow for greater efficiency, lower costs, and reduced emissions. On the other hand, some researchers have been working to increase the potential of energy generation processes through heat recovery with steam generators, organic Rankine cycles, and absorption chillers. This Special Issue is a collection of fundamental or applied and numerical or experimental investigations. Many new concepts in thermal systems and energy utilization are explored, discussed, and published as original research papers in the "Thermal Systems".

The first paper, presented by Ochoa et al. [1], offers an extensive thermo-economic analysis of a heat recovery steam generation system integrated with an absorption refrigeration chiller and a gas micro-turbine. The effect of compressor inlet air temperature on thermo-economic performance of trigeneration systems was studied and analyzed in detail based on a validated model. They found some operational conditions where exergy was highly destroyed due to the exergy inefficiencies of the equipment such as combustion chamber, microturbine, compressor, evaporator, heat exchanger and generator which are found to be important as exergo-economic factors.

In another investigation, Ochoa et al. [2] present an analysis of a waste heat recovery system based on the organic Rankine cycle from the exhaust gases of an internal combustion engine. They studied the exergy destroyed values and the rate of fuel exergy, product exergy, and loss exergy. They found exergo-economic analysis was a powerful method to identify the correct allocation of the irreversibility and the real potential for improvement between components.

Zhang et al. [3] perform an experimental investigation to enhance the working performance and temperature control of electric vehicle batteries through a thermal management system with a heat pipe and thermoelectric cooler. Heat pipes with high thermal conductivity were used to accelerate dissipating heat on the surface of the battery with an additional thermoelectric cooler to increase discharge rate. The findings support the results generated from engineering simulation and show that the combined system can effectively reduce the surface temperature of a battery.

Qian et al. [4] propose the application of oscillating heat pipes to reduce thermal damage in an abrasive milling tool. Heat pipes are passive heat transfer devices with excellent heat transport capacity and they are applied to the machining process to enhance heat transfer. The experimental investigation studied the effects of centrifugal acceleration, heat flux, and working fluids, hence, methanol, acetone, and water, on their thermal performance. Based on their theoretical analysis, centrifugal acceleration will increase the resistance for the vapor to penetrate through the liquid slugs to form an annular flow, which was supported by slow-motion visualization. The phase change occurs, and vapor moves to the condenser to release heat by condensing into liquid.

Sartor and Dickes [5] validated numerical results obtained from a heat transport model with experiments in large solar thermal plants at the Plataforma Solar de Almeria in Spain. They argued that the previous work done had limitations in the assessment of temperatures and computational time required for simulating large pipe networks. They proposed to model the dynamic behavior of the whole system based on a few data inputs. Some atmospheric conditions, such as local clouds, could have significant influence on the outlet temperature and other dynamic behavior of the solar field. An alternative method was used to validate a solar thermal plant considering the thermal solar gain and the inertia of the pipes in their investigation. The accuracy of the model was found to be similar to those of the one-dimensional finite volume method with a reduced simulation time.

Alexopoulos et al. [6] validate design procedure from a simulation model with an experimental study of an air finned tube CO2 gas cooler. Based on the model, the evaluation of various physical parameters such as length and diameter of tubes as well as ambient temperature was conducted. The researchers attempted to identify the most suitable design in terms of pressure losses and required heat exchange area for selected operational conditions. Hence, a simulation model of the gas cooler was developed and validated experimentally by comparing the overall heat transfer coefficient. The comparison between the model and the experimental results showed a satisfactory convergence for selected operational conditions.

Barrella et al. [7] present a feasibility study which analyzed the use of a centralized electrically driven air source heat pump for space heating. Two models were developed to obtain variables in the hourly thermal energy demand and the off-design heat pump performance. The proposed heat pump is driven by a motor with variable rotational speed to modulate the heating capacity in an efficient way. An eco-friendly refrigerant (R290 or propane) was selected for the heat pump. A back-up system was used to meet the peak demand. Renewable energy used via the heat pump showed significant reduction in CO2 which would otherwise have been produced via normal fossil fuel consumption. The researchers claimed that these results showed that the proposed technology was among the most promising measures for addressing energy demand in vulnerable households.

Kim et al. [8] propose an alternative method of swaging which is claimed to be more efficient than the traditional coating technology in the fabrication of accident-tolerant fuel cladding. In their study, it was found that the specimen exhibited a pseudo-single tube structure with higher thermal stability. They reported that the specimen had a uniform and well-bonded interface structure under optical microscopy and scanning electron microscopy images. The specimen did not show significant structural collapse, even after being stored at 1200 ◦C for one hour. The experimental results show that tube process has a high potential for development of an ATF cladding with a length of several meters with their geometries calculated according to the design.

Grabowski et al. [9] present a numeric heat transfer investigation validated with an experimental study in flow boiling of water through an asymmetrically heated, rectangular and horizontal mini-channel, with transparent side walls. The mathematical model assumed the heat transfer process in the measurement module to be steady-state with temperature independent thermal properties of solids and flowing fluid. Grabowski et al. applied laminar characteristic flow in the Reynolds numbers study. The experimental data taken were temperatures at strategic points, volume flux of flowing water, inlet pressure and pressure

*Energies* **2021**, *14*, 175

drop, current, and the voltage drop in the heater power supply. They defined two inverse heat transfer problems which were solved by the meshless Trefftz method with two sets of T-functions.

Kim et al. [10] demonstrate the use of vortex tube in an air conditioning system with the objective to get rid of the use of refrigerant gas. The success of the eco-friendly technology will avoid environmental impact due to refrigerant. The vortex tube is a temperature separation system capable of separating air at low and high temperatures with compressed air. In their experimental study, both direct and indirect heat exchange were investigated to test low-temperature air flow rate according to temperature and pressure. The direct heat exchange method was found to have low flow resistance, and ease in control of temperature and flow-rate. As a result, it is judged to be a more feasible method for use in air-conditioning system by the authors.

The papers in this special issue reveal an exciting area, namely the "Thermal Systems" that is continuing to grow. The pursuit of work in this area requires expertise in thermal and fluid dynamics, system design, and numerical analysis as well as experimental validation. We are extremely delighted to be invited as the Guest Editors of this "Special Issue". We have received great support from many colleagues and top researchers of prestigious universities and research institutions. We are heartened to see such a contribution with the aim of tackling the environmental impact or providing low-cost energy options to the humble communities. We firmly believe that, with the continuing collaboration of all researchers, we can enhance our contribution to tackling the numerous challenges faced by global society. We hope that this Special Issue helps to bring the research community into closer contact with each other. Finally, we would like to thank all our authors, reviewers, and editorial staff who have contributed to this publication. I am sure all readers of this Special Issue of Energies will find the scientific manuscripts interesting and beneficial to their research work in the years to come.

**Funding:** This research received no external funding.

**Conflicts of Interest:** The authors declare no conflict of interest.

### **References**


*Energies* **2021**, *14*, 175


© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

### *Article* **Thermo-Economic Assessment of a Gas Microturbine-Absorption Chiller Trigeneration System under Di**ff**erent Compressor Inlet Air Temperatures**

### **Guillermo Valencia Ochoa 1,\*, Carlos Acevedo Peñaloza <sup>2</sup> and Jorge Duarte Forero <sup>1</sup>**


Received: 5 November 2019; Accepted: 3 December 2019; Published: 6 December 2019

**Abstract:** This manuscript presents a thermo-economic analysis for a trigeneration system integrated by an absorption refrigeration chiller, a gas microturbine, and the heat recovery steam generation subsystem. The effect of the compressor inlet air temperature on the thermo-economic performance of the trigeneration system was studied and analyzed in detail based on a validated model. Then, we determined the critical operating conditions for which the trigeneration system presents the greatest exergy destruction, producing an increase in the costs associated with loss of exergy, relative costs, and operation and maintenance costs. The results also show that the combustion chamber of the gas microturbine is the component with the greatest exergy destruction (29.24%), followed by the generator of the absorption refrigeration chiller (26.25%). In addition, the compressor inlet air temperature increases from 305.15 K to 315.15 K, causing a decrease in the relative cost difference of the evaporator (21.63%). Likewise, the exergo-economic factor in the heat exchanger and generator presented an increase of 6.53% and 2.84%, respectively.

**Keywords:** thermo-economic assessment; exergy analysis; trigeneration system; gas microturbine; absorption chiller

### **1. Introduction**

The increase in global warming, adding to the scarcity of fossil fuels, has motivated the development of new technologies to improve the efficiency of existing processes in power plants [1,2]. Among the available options, multi-generation processes such as the trigeneration cycle have been widely used as they allow for greater efficiency, lower costs, and reduced emissions [3]. Therefore, researchers have been working to increase the potential of this type of energy generation process through heat recovery under the steam generator, organic Rankine cycles [3], and absorption chillers [4,5].

Absorption chiller cooling technology is increasingly used because it utilizes refrigerants and absorbents that do not have a negative effect on the environment. In addition, it is possible to feed this type of device with waste heat or some other renewable energy source such as solar energy [6]. Therefore, they are systems widely used in the industrial sector because of the lower energy cost production and potential gas emission reduction [7].

Several studies have developed relevant contributions to the thermo-economic analysis and optimization of absorption refrigeration systems [8], but few are related to the trigeneration system working at the different operating conditions. These studies mainly involve the application of the thermodynamic second law to conduct the evaluation and thermal analysis of the system, which is based on the exergy approach [9]. This method allows us to measure the work potential or quality of different forms of energy with respect to environmental conditions [10,11]. Therefore, the environmental condition plays a key role in the thermo-economic performance of thermal cycles.

Kaynakli and Kilic [12] analyzed the effect of an H2O-LiBr absorption refrigeration system (ARS) on operating conditions by means of the first and second laws of thermodynamics. It was observed that there is an increase in system performance with the increase in generator temperatures and a decrease in condenser and absorber temperatures. However, the effect of the integration of the ARS with the exhaust gases of a thermal prime mover was not studied, and the value of the generator temperature should be determined. In contrast, Martinez and Rivera [13] conducted an energy and exergy analysis for a dual absorption system using the H2O-LiBr as a working fluid and also concluded that higher generator and evaporator temperatures and lower absorber temperatures lead to improved system performance. Then, a change in inlet air temperature means different operation conditions on the ARS, and a substantial effect on the thermo-economic indicator of the trigeneration systems, such as the relative cost difference and exergo-economic factor.

Consequently, Kaushik and Arora [14] developed an energetic and exergetic analysis of the single and double effect of a cooling absorption system with parallel free water flow. According to the presented results, the coefficient of performance (COP) presented for the single effect ARS was ranging from 0.6 to 0.75, while in the case of the double effect the COP increased from 1 to 1.28, as a result of different operating temperatures of the heat source and evaporator.

In order to identify the exergetic improvement potential in the H2O-LiBr double effect, ARS, Gomri and Hakimi [15] conducted an energetic and exergetic analysis, calculating the exergy destruction of the system components. They concluded that the absorber and the high-pressure generator are the components that most influence the total exergy destruction of the system. On the other hand, renewable energy had been used as a heat source of refrigeration systems to increase global thermal efficiency. Hence, Rosiek [16] studied a cooling system integrated into a flat plate solar collector, and the results demonstrated that it is possible to obtain the best results from the exergetic viewpoint supplying water to the absorption cooler in a temperature range of 70–80 ◦C.

A novel configuration was proposed by Pourfayaz et al. [17], by means of an exergetic analysis to increase the overall performance of the ARS, for a fuel cell cooling system in which nanofluids were used as absorbers.

There are different trigeneration systems, which can be classified mainly according to their driving force, the amount of energy used, and the size of the plant. Each of these classifications has a series of classifications, which have certain advantages and disadvantages regarding the acquisition cost, installation, maintenance, operation ranges, necessary conditions, among others [18].

The availability of sources for electricity generation and global warming are alarming factors that lead to concern about the sustainability of energy production in the future, which brings with it the transcendental impact to design more efficient energy systems [19]. Combined Heat, Cold, and Power (CHCP) are some of the alternative technologies to address problems such as growing energy demand, rising energy costs, the security of energy supply, and large environmental impact [20]. Thus, it is presented as a solution with relevant technical potential, economic, and ecological benefits, which allow reducing the use of primary energy sources to energy generation [21]. The trigeneration system is composed of five main elements: primary engine, electric generator, waste heat recovery system, thermal activation equipment, and a control unit [22].

A promising alternative to trigeneration systems that address the energy problem is based on the use of low-capacity primary sources, also called small-scale technologies, which deliver power between 28 and 200 kW [23], such as the gas microturbine considered in this study. These systems are particularly suitable for applications in commercial buildings, hospitals, schools, local industries, office blocks, and single or multi-family residential buildings [24].

Although some research results based primarily on exergetic analysis show an increase in the COP of the ARS, it is not a complete enough analysis to design a thermal system and ignores the economic part of the system. Therefore, the exergo-economic aspect is necessary to incorporate both exergetic and economic analysis into the system. In this way, it is possible to have a better guide for the thermal study of the components [25,26]. Therefore, the optimization of the ARS performance by means of the thermo-economic assessment was applied [27].

Some trigeneration systems had been studied in industrial and commercial applications. The thermo-economic potential of a trigeneration biomass plant was studied [28], using different configurations, parameters, both economic and operational. The exergetic simulation allowed to determine a 72.8% of the energy efficiency, and the exergetic efficiency ranging from 20.8% to 21.1%, but a parametric case studied is not presented to determine the relevant parameters of the trigeneration process. Also, a complete study was conducted considering some performance energetic, economic and environmental indicators, where the performance of a steam turbine trigeneration system for large buildings based on the energy demands of the facility was calculated, and the results were compared with conventional power generation systems [29]. The results show a decrease in the primary energy saving of 12.1%, CO emission reduction of 2.6%, and CO2 emission reduction of 2.6%. However, a thermo-economic model was not proposed in this research to identify exergy destruction opportunities.

On the other hand, some thermo-economic studies using a chiller in the trigeneration system have been considered, but the use of a gas microturbine as prime mover operating in a trigeneration system is not reported in the literature. Therefore, the integration of an absorption chiller to a trigeneration system was proposed to generate the required energy, and thus assess energy costs and savings, obtaining an annual cost of \$US 384,300 per year, and a payback period of 1.8 years [30]. In addition, an economic analysis of a trigeneration system based on a LiBr chiller was developed. The results were compared with respect to other heat and cold generation systems, and the primary energy consumption decreased by 26.6% with respect to cogeneration.

In the case of the trigeneration system using gas turbines as a prime mover, Ahmadi et al. [31] presented energy and exergetic analysis in a trigeneration system with a combined gas turbine cycle. The results showed a greater exergetic destruction in the combustion chamber, in addition to the environmental impact assessment, where the thermal energy efficiency increase 75.5%, the thermal exergetic efficiency increase 47.5%, and the emission of the CO2 decrease to 158 kg/MWh.

To increase the performance of the trigeneration system, an exergo-economic optimization was conducted using an evolutionary algorithm, where the economic indicators used to optimize the systems are the total revenue requirement and the total cost of the system [32]. The optimization result of the system allows an improvement of 0.207 \$/s in the objective function studied, which is 15% lower than the value in the base case.

From this literary review, the main contribution of this paper is to present a parametric study conducted in a trigeneration system integrated by a Li-Br ARS, a gas microturbine, and a waste heat recovery, to study the effect of inlet gas compressor temperature on the energy, exergy and termo-economic indicator. In addition, the analysis includes the application of the energy, exergy balance, exergy destruction calculation, cost balances application, and the thermo-economic modeling by components, considering in detail the acquisition, maintenance, and operating costs.

### **2. Methodology**

### *2.1. Description of the System*

The physical structure of the trigeneration system is presented in Figure 1. Starting with the gas power cycle, where ambient air at atmospheric pressure (state 1) enters to the compressor, where it is compressed, and its pressure rises, from which it goes out to the preheater (state 2) where it interacts with the exhaust gases of the turbine to increase its temperature and obtain a better combustion.

**Figure 1.** Physical structure of the trigeneration system.

In the combustion chamber, the air flow (state 3) enters from the preheater at high temperature and pressure, and the methane flow (state 4) enters, which will be burned during mixing with excess air. The modeling of the combustion chamber is obtained, assuming an expansion process of the air, which corresponds to an isobaric process inside the system. The resulting combustion gases (state 5) move the turbine in which the hot gases expand and cool rapidly through an adiabatic expansion, generating the power required for the compressor and the net power of the system.

The output gases of the turbine (state 6) are directed to the preheater equipment where its temperature decreases, and then in the Heat Recovery Steam Generator (HRSG) the heat transfer process allows us to generate the steam (state 10), from the water at ambient temperature (state 9). The exhaust gases (state 8), as an energy source, enters the generator where it separates the solution resulting in an H2O-LiBr mixture with a low concentration (state 11) and the generation of refrigerant saturated steam (state 12). The mixture is modeled as a sub-cooled liquid type (state 19), which expands through a flow valve and arrives at the absorber as a low concentration H2O-LiBr mixture (state 20). On the other hand, the heat is removed inside the condenser heat exchanger from the refrigerant to the environment (state 23), going from a gaseous to a liquid phase (state 13).

In the evaporator heat exchanger, the fluid takes heat from the refrigerated space or room, and induces a phase change in the refrigerant producing a pressure difference between the evaporator and the absorber, where the refrigerant exits as saturated steam (state 15) directly to the absorber, in which there is an energy change between the external water (state 27) and lithium bromide (state 20). As a result, a loop of the lithium bromide mixture is obtained to give a saturated liquid solution (state 16). The pressure of this solution is increased and entered into the heat exchanger by the flow energy supplied by the motor of the pump (state 17) through a counter-current configuration, which increases the temperature to improve efficiency. Finally, the fluid arrives at the generator to continue with the system cycle (state 18).

### *2.2. Thermodynamic Modeling*

In the thermodynamic modeling of the trigeneration system [33], the components of the system are considered as open systems where a steady-state mass balance is applied according to Equation (1). For the case of constant flow systems, such as the generator and absorber, this balance results as shown in Equation (2).

$$
\sum \dot{m}\_{\text{out}} = \sum \dot{m}\_{\text{in}} \tag{1}
$$

$$
\sum \dot{m}\_{\text{out}} \cdot \mathbf{x}\_{\text{out}} = \sum \dot{m}\_{\text{in}} \cdot \mathbf{x}\_{\text{in}} \tag{2}
$$

where *<sup>x</sup>* is the concentration, . *mout* and . *min* are the output and input mass flows to the system in kg/s.

Also, the energy balance applied to each component of the trigeneration system based on the first law of thermodynamics is expressed in Equation (3).

$$
\sum \dot{Q} - \sum \dot{W} = \sum \dot{m}\_{\text{out}} \cdot l\_{\text{out}} - \sum \dot{m}\_{\text{in}} \cdot l\_{\text{in}} \tag{3}
$$

where *<sup>h</sup>* is the specific enthalpy in kJ/kg, . *<sup>Q</sup>* is heat flow rate in kW, and . *W* is the power rate in kW.

The performance coefficient of the ARS (*COPARS*) is expressed by Equation (4), which is defined as the ratio of the heat transfer of the evaporator ( . *QEvaporator*) in kW, and the amount of heat transfer in the generator ( . *QGener*) plus the energy rate of the pump ( . *WP*), both in kW.

$$\dot{\bf{COP}}\_{\rm{ARS}} = \frac{\dot{\dot{Q}}\_{\rm{Expuator}}}{\dot{Q}\_{\rm{Gewer}} + \dot{\dot{W}}\_{\rm{P}}} \tag{4}$$

Applying the energy balance to each of the components of the trigeneration system gives the equations shown in Table 1.

**Component Energy Balance** Compressor . *<sup>m</sup>*1·*h*<sup>1</sup> <sup>+</sup> . *Wcomp* <sup>−</sup> . *m*2·*h*<sup>2</sup> = 0 Combustion Chamber . *<sup>m</sup>*3·*h*<sup>3</sup> <sup>+</sup> . *m*4·*h*<sup>4</sup> ·*ncc* <sup>−</sup> . *m*5·*h*<sup>5</sup> = 0 Turbine . *<sup>m</sup>*5·*h*<sup>5</sup> <sup>−</sup> . *Wturb* <sup>−</sup> . *m*6·*h*<sup>6</sup> = 0 Pre-heater . *<sup>m</sup>*7·*h*<sup>7</sup> <sup>−</sup> . *<sup>m</sup>*6·*h*<sup>6</sup> <sup>=</sup> . *<sup>m</sup>*3·*h*<sup>3</sup> <sup>−</sup> . *m*2·*h*<sup>2</sup> HRSG . *<sup>m</sup>*7·*h*<sup>7</sup> <sup>−</sup> . *QHRSG* <sup>+</sup> . *m*8·*h*<sup>8</sup> = 0 Generator . *<sup>m</sup>*18·*h*<sup>18</sup> <sup>+</sup> . *Qgener* <sup>−</sup> . *<sup>m</sup>*12·*h*<sup>12</sup> <sup>−</sup> . *m*11·*h*<sup>11</sup> = 0 Condenser . *<sup>m</sup>*12·*h*<sup>12</sup> <sup>−</sup> . *<sup>m</sup>*13·*h*<sup>13</sup> <sup>+</sup> . *Qcond* = 0 Evaporator . *<sup>m</sup>*13·*h*<sup>13</sup> <sup>+</sup> . *Qevap* <sup>−</sup> . *m*15·*h*<sup>15</sup> = 0 Absorber . *<sup>m</sup>*15·*h*<sup>15</sup> <sup>+</sup> . *<sup>m</sup>*20·*h*<sup>20</sup> <sup>−</sup> . *<sup>m</sup>*16·*h*<sup>16</sup> <sup>+</sup> . *Qabs* = 0 Heat exchanger . *<sup>m</sup>*17·*h*<sup>17</sup> <sup>+</sup> . *<sup>m</sup>*11·*h*<sup>11</sup> <sup>−</sup> . *<sup>m</sup>*18·*h*<sup>18</sup> <sup>+</sup> . *m*20·*h*<sup>20</sup> = 0

**Table 1.** Energy balance equations by components of the trigeneration system.

To calculate the specific physical exergy ( . *E Ph* ) was not considered the kinetic and potential energy, resulting in the Equation (5).

$$\dot{E}^{\text{Ph}} = (h - h\_0) - T\_{\text{0}} \text{(s-so)}\tag{5}$$

where *h* is the specific enthalpy in kJ/kg, *s* is the specific entropy in kJ/kg · K of the working fluid flow, *h*<sup>0</sup> and *s*<sup>0</sup> are the state enthalpy and entropy at reference condition (*T*<sup>0</sup> = 298.15 K and *P*<sup>0</sup> = 101.325 kPa).

On the other hand, the chemical exergy for water ( . *E Ch water*) was calculated using Equation (6), while for the microturbine exhaust gases (states 6, 7, 8, and 21) was used the Equation (6) since the change of chemical exergy for lithium bromide was not considered.

$$
\dot{E}\_{\text{water}}^{\text{Ch}} = \dot{m} \frac{\text{(}\tau\_{\text{water}}\text{)}}{\text{(}\dot{M}\_{\text{water}}\text{)}} \dot{E}\_{\text{Ch, water}}^{\text{0}}\tag{6}
$$

$$\dot{E}^{ch} = \sum\_{k=1}^{n} \mathbf{x}\_{k} \dot{E}^{ch\_{k}} + \mathcal{R} \cdot T\_{0} \sum\_{k=1}^{n} \mathbf{x}\_{k} \cdot \ln \mathbf{x}\_{k} \tag{7}$$

where ( . *E* 0 *<sup>Q</sup>*, *water*) is the standard chemical exergy of the water, *xk* is the molar fraction, and *exchk* is the exergy per mol unit for the *k* gas.

The exergy balance was applied to each component of the trigeneration system according to Equation (8) [34].

$$
\sum \dot{m}\_{\text{in}} \cdot \dot{E}\_{\text{in}} - \sum \dot{m}\_{\text{out}} \cdot \dot{E}\_{\text{out}} + \dot{Q} \cdot \left(1 - \frac{T\_0}{T}\right) - \dot{W} - \dot{E}\_D = 0 \tag{8}
$$

where . *min*· . *Ein* is the inflow exergy, . *mout*· . *Eout* is the outflow exergy, and . *ED* is the destroyed exergy.

The exergetic efficiency (η*ex*) based on the second law of thermodynamics, is expressed by the Equation (9).

$$
\eta\_{c\alpha} = \frac{\dot{E}\_P}{\dot{E}\_F} \tag{9}
$$

where the amount of fuel exergy ( . *EF*) to the system, and the exergy produced ( . *EP*) per system are related to the destroyed exergy ( . *ED*), and the lost exergy ( . *EL*), as shown in Equation (10).

$$
\dot{E}\_F = \dot{E}\_P + \dot{E}\_D + \dot{E}\_L \tag{10}
$$

The Fuel and Product structure in each component of the trigeneration system was calculated, as shown in Table 2.


**Table 2.** Fuel and Product exergy equations.

#### *2.3. Thermo-Economic Analysis*

To calculate the total production cost, it is considered the capital investment costs ( . *ZCI*), operation and maintenance ( . *ZOM*), as shown in Equation (11).

$$
\dot{Z} = \dot{Z}\_{\text{CI}} + \dot{Z}\_{\text{OM}} \tag{11}
$$

The equations used to calculate the Purchase Equipment Costs (PEC) for the components of the ARS were: heat exchangers (Equation (12)), pump (Equation (13)), motor (Equation (14)), where the sub-index "0" represents the reference of the studied component [35–37].

$$PECK = \ \_{PEC0, \&} \cdot \left(\frac{A\_k}{A\_0}\right)^{0.6} \tag{12}$$

where the reference area (*A*0) is 100 m2, the reference costs (PEC0, K) considered are Evaporator (16,000 USD), Condenser (8000 USD), Absorber (16,500 USD), and Heat Exchanger (12,000 USD) [26]. Also, the PEC for the pump is calculated based on Equation (13).

$$PEC\_{pump} = \text{PEC}\_{0,pump} \cdot \left(\frac{\dot{W}\_{pump}}{\dot{W}\_{0,pump}}\right)^{m\_B} \cdot \left(\frac{1 - \eta\_{pump}}{\eta\_{pump}}\right)^{n\_{pump}} \tag{13}$$

where the pump efficiency (η*pump*) is 75%, the pump size power ratio (*mB*) is 0.26, and the reference cost (*PEC*0,*pump*) is 2100USD. In addition, the model used to estimate the PEC of the pump motor is presented in Equation (14).

$$PEC\_{mat} = \;PEC\_{0,mt} \cdot \left(\frac{\dot{W}\_{mat}}{\dot{W}\_{0,mt}}\right)^{m\_{mat}} \cdot \left(\frac{1 - \eta\_{mat}}{\eta\_{mat}}\right)^{m\_{mat}} \tag{14}$$

where the motor size power ratio (*mmot*) is 0.87, the motor reference power ( . *W*0, *mot*) is 10 kW, the motor efficiency *(*η*mot)* is 90%, the efficiency ratio of motor size (*nmot*) is 1, and the reference cost (*PEC*0,*mot*) is 500 USD [26].

The components of the gas microturbine were used some well-known models [26]. For the PEC of the compressor was used the Equation (15), combustion chamber (Equation (16)), and turbine (Equation (17)).

$$PEC\_{comp} = \left(\frac{C\_{11} \cdot \dot{m}\_{air}}{C\_{12} - n\_{comp}}\right) \cdot \left(\frac{P\_{a2}}{P\_{a1}}\right) \cdot \ln\left(\frac{P\_{a2}}{P\_{a1}}\right) \tag{15}$$

where the compressor coefficients *C*<sup>11</sup> and *C*<sup>12</sup> are 71.10 and 0.9 USD/(kg/s), respectively. In addition, the PEC of the combustion chamber was calculated according to Equation 16.

$$PEC\_{cc} = \left(\frac{C\_{21} \cdot \dot{m}\_{CH4}}{C\_{22} - \frac{P\_{gv4}}{P\_{a3}}}\right) \cdot \left[1 + e^{(C\_{23} \cdot T\_4 - C\_{24})}\right] \tag{16}$$

where *C*<sup>21</sup> is 46.08 USD/(kg/s), *C*<sup>22</sup> is 0.995, *C*<sup>23</sup> is 0.018 *K*−<sup>1</sup> and *C*<sup>24</sup> is 26.4 [26]. Also, for the case of the turbine, Equation (17) was used to calculate the PEC, which is a relevant cost of the microturbine equipment.

$$PEC\_{turb} = \left(\frac{C\_{31} \cdot \dot{m}\_{comb}}{C\_{32} - n\_{turb}}\right) \cdot \ln\left(\frac{P\_{\text{gc4}}}{P\_{\text{gc5}}}\right) \left[1 + e^{\left(C\_{33} \cdot T\_4 - C\_{34}\right)}\right] \tag{17}$$

where the model turbine coefficients are, *C*<sup>31</sup> in 479.34 USD/(kg/s), *C*<sup>32</sup> is 0.92, *C*<sup>33</sup> is 0.036 *K*−<sup>1</sup> and *C*<sup>34</sup> is 54.4 [26].

On the other hand, the values of the leveled costs (*PECL*) by components were calculated by the mean of the Equation (18), which consider the money transactions occur at the end of each year in the economic life of the trigeneration system.

$$PECL = CRF \cdot \sum\_{j=1}^{n} \frac{PEC\_{j}}{\left(1 + i\_{eff}\right)^{j}} \tag{18}$$

where *CRF* is the Capital Return Factor, *ieff* is the interest rate and *PECj* is the value of the Purchase Equipment Costs in the year *j*. Also, the *CRF* is calculated using the Equation (19).

$$\text{CRF} = \frac{i\_{eff} \left(1 + i\_{eff}\right)^n}{\left(1 + i\_{eff}\right)^n - 1} \tag{19}$$

where *n* is the lifetime of the equipment.

To obtain the value of the capital investment in term of unit cost per time ( . *Zk*), without having to calculate leveled costs, the Equation (20) can be used.

$$
\dot{Z}\_k = \frac{\text{PEC-CRF-}\varphi}{\tau} \tag{20}
$$

where τ is the total operation time in hours of the system at full load, and ϕ is the maintenance factor [38].

### *2.4. Exergy Cost Balance and Thermo-Economic Indicators*

The exergetic cost balance, as shown in Equation (21), the inefficiencies presented in the equipment are evaluated, and the intermediate and final cost of the streams of the thermal process. This analysis allowed us to estimate all exergies of stream in the trigeneration cycle, considering the total costs (acquisition costs, operating costs, and maintenance costs) [39–41].

$$\sum\_{i=1}^{n} \dot{\mathcal{C}}\_{out,i} + \dot{\mathcal{C}}\_{W,i} = \sum\_{i=1}^{n} \dot{\mathcal{C}}\_{in,i} + \dot{\mathcal{C}}\_{Q,i} + \dot{Z}\_k \tag{21}$$

where the terms . *Cin*,*i*, and . *Cout*,*<sup>i</sup>* are the exergy costs associated with the inlet and outlet flow, which are calculated using the Equation (22). The terms . *CW*,*<sup>i</sup>* and . *CQ*,*<sup>i</sup>* are the costs associated with the power and heat transfer exergy cost, which are calculated by the mean of the Equations (23) and (24), respectively.

$$
\dot{\mathbf{C}}\_i = \mathbf{c}\_i \dot{\mathbf{E}}\_i \tag{22}
$$

$$
\dot{\mathbf{C}}\_{W,i} = \mathbf{c}\_{W^\cdot} \dot{\mathbf{W}} \tag{23}
$$

$$
\dot{\mathbf{C}}\_{Q,i} = \mathbf{c}\_Q \cdot \dot{\mathbf{E}}\_Q \tag{24}
$$

where *ci*, *cW* and *cQ* are the specific costs per unit of exergy expressed in dollars per Gigajoules (USD/GJ).

The equations presented in Table 3 are obtained, applying the general cost balances to each component of the trigeneration system. The cost balance can be expressed as shown in Equation (25), as a function of the cost rates of the exergy lost ( . *CL*,*i*), product cost rate ( . *CP*,*i*) according to Equation (26), fuel cost rate ( . *CF*,*i*) attending to Equation (27), and the cost rate of exergy destruction ( . *CD*,*i*) by means of Equation (28). .

$$
\dot{\mathbf{C}}\_{P,i} = \dot{\mathbf{C}}\_{F,i} - \dot{\mathbf{C}}\_{L,i} + \dot{\mathbf{Z}}\_{i} \tag{25}
$$

$$c\_{P,i} = \frac{\dot{C}\_{P,i}}{\dot{E}\_{P,i}} \tag{26}$$

$$
\omega\_{F,i} = \frac{\dot{C}\_{F,i}}{\dot{E}\_{F,i}} \tag{27}
$$

$$
\dot{\mathcal{C}}\_{D,i} = c\_{F,i} \dot{E}\_{D,i} \tag{28}
$$


**Table 3.** Cost balance equation by components.

The relative cost difference (*r*) is calculated from the specific fuel and product cost, which shows the average relative cost increase per exergetic unit between the input power and the product, as shown in Equation (29).

$$r = \frac{c\_{F,i} - c\_{F,i}}{c\_{F,i}} \tag{29}$$

Likewise, another thermo-economic indicator studied is the exergo-economic factor (*f*), which measures the relation between the capital cost investment compared to the loss costs rate and exergy destruction, and it is calculated using Equation (30).

$$f = \frac{\dot{Z}\_i}{\dot{Z}\_i + \dot{\mathbb{C}}\_{L,i} + \dot{\mathbb{C}}\_{D,i}}\tag{30}$$

A high value of this factor reflects a decrease in investment costs by improving energy efficiency, in contrast to low rates of this factor, which suggests saving costs throughout the system for improvement in efficiency.

### **3. Results and Discussion**

This section presents the results of the energetic, exergetic, and thermo-economic analysis of the base case of the trigeneration system, and the parametric study result at different inlet air compressor temperatures. Table 4 presents the thermodynamic properties, physical and chemical exergy of the trigeneration system obtained from the energy balance.


**Table 4.** Thermodynamic properties, physical and chemical exergy of the trigeneration system.

These values allow conducting the exergy balance to calculate the exergy destruction in each component of the system, obtaining the values shown in Figure 2, which presents the percentage of the exergy destruction fraction by equipment and/or components of the system studied.

**Figure 2.** Percentage Exergy Destruction for each component of the trigeneration system.

The exergy destruction represents a loss of useful work that can be taken advantage of by the components of the trigeneration system for the improvement of the operating and thermo-economic conditions. This translates into a great inefficiency and a considerable quantity of energy that must be minimized just when it is needed to maximize the overall thermal efficiency of the process. However, it represents an opportunity to optimize and develop innovate techniques based on new proposals and alternatives of operation and manufacture to design and to rate the devices, all with the purpose of reducing the investment cost associated with the components and, therefore, to the general system.

The results show that the combustion chamber of the gas microturbine is the component with the greatest exergy destruction (29.24%), followed by the generator of the ARS (26.25%). The compressor has a contribution of 0.08% due to the low heat transfer irreversibility presented in this device as a result of the high operating temperature. The greatest amount of exergy destroyed (72%) is located in the combustion chamber, generator and absorber, which suggests a greater technological and operational effort focused on the design of these heat exchange, allowing us to reduce the temperature difference between the fluids. Likewise, from the component cost balances, the exergetic costs of the stream in the system were calculated, as shown in Table 5.


**Table 5.** Exergetic cost rates ( . *C*) in USD/s, and the costs per unit of exergy of the stream (*c*) in USD/GJ.

These costs allow us to determine the cost of raw materials, products, and destruction, besides the thermo-economic indicators shown in Table 6. The highest costs are associated with the pre-heater assembly (532.40 USD/GJ), which exceeds the cost of the products obtained by the turbine.

**Table 6.** Average costs per unit of fuel ( . *CF*,*<sup>i</sup>* ) and product ( . *CP*,*<sup>i</sup>* ), destruction cost ( . *CD*,*i*), relative cost difference (*r*), and exergo-economic factor (*f*) for each component of the trigeneration system.


### *3.1. Energy and Exergy Analysis*

This section presents the parametric study results applied in the trigeneration system, through the variation in compressor inlet temperature (*T1*) from 293.15 K to 328.15 K. The energy performance was evaluated for the microturbine subsystem, HRSG, and the evaporator of the ARS, which will be analyzed and discussed in detail below.

Figure 3a shows that the increase in temperature causes a decrease in the net power supplied by the turbine with respect to the heat absorbed by the evaporator at different air flow ratio, because of the enhancement in the heat removed in the evaporator as a consequence of the increase in the air flow temperature in the inlet compressor. For a given inlet compressor temperature of 313.15 K, the Trigeneration System delivers almost 236.2% power per unit of heat in the evaporator with an air-fuel ratio of 0.7 with respect to 0.5. This means that an increase in the air-fuel ratio causes the microturbine to deliver more power, given the higher airflow. However, the decreasing trend among the studied energy coefficient is preserved, so the inlet compressor temperature is a basic parameter that

significantly affects the energy performance of the trigeneration system. On the other hand, the steam energy supplied by the HRGS (Figure 3b) increases with the increase of the air flow inlet compressor temperature, which is due to the improvement of the combustion chamber thermal efficiency, having a combustion with an oxidant both at higher temperature, and more air mass flow by increasing the air-fuel ratio, which increases the thermal capacitance in the HRSG.

**Figure 3.** Energy performance of the trigeneration system at different compressor inlet air temperatures, (**a**) Net Power/Evaporator heat, and (**b**) Heat in the Heat Recovery Steam Generator (HRSG).

### *3.2. Exergy Destruction*

The exergy destruction analysis by component represents an opportunity for the improvement of the operative and thermo-economic of the thermal system. With this parameter, some alternatives can be developed to reduce the effects of both economic and operational losses. It is important to know the behavior of the minimum and maximum values that can be presented in a device with the purpose of obtaining the possible exergy losses in the trigeneration system and the effect of the increase of the inlet air compressor temperature in the devices of the system that is studied here. The exergy destroyed in the trigeneration system at different compressor inlet air temperatures is shown in Figure 4.

**Figure 4.** Exergy destroyed in the trigeneration system at different compressor inlet air temperatures.

For the exergy destroyed fraction in terms of components, a tendency was observed. We used an accumulation of the percentage losses to obtain a graphical model that describes the behavior of the exergy losses of the general system. With the exclusion of the combustion chamber in the study, we observed a 12% reduction in the exergy destroyed in the turbine until 300 K. Then, it remained at its lowest values while the inlet air compressor temperature increased. Another relevant result is the behavior of the exergy destruction in the generator and condenser, which both present the same range of variation of exergy destroyed (4% to 5%). In this case, below the temperature of 300 K, the exergy destroyed in the generator is greater than the values presented. After this, an increase in temperature caused a significant enhancement of the exergy efficiency of this device at higher compressor inlet air temperature. Therefore, the highest values of exergy destroyed in the devices of the trigeneration system happen when the air temperature is lower than 300 K. However, there are exceptions such as the HRSG and the evaporator, because a higher air temperature limits their functionality in the system.

### *3.3. Cost Rate of Exergy Destruction*

The costs related to the exergy destruction represent the economic losses in dollars with respect to the operating time of the system and, therefore, of each one of the components. A comparative analysis of the costs per component of the microturbine is proposed for the variation of the inlet air compressor temperature, while simultaneously analyzing the behavior of the costs related to the exergy destruction in the assembly and different components both in the ARS and microturbine. Figure 5 presents the trends of the cost rate of exergy destruction of the trigeneration system at different compressor inlet air temperatures.

**Figure 5.** Cost rate of exergy destruction of the trigeneration system at different compressor inlet air temperatures, (**a**) microturbine components, and (**b**) absorption chiller components.

The results show that for the ARS (Figure 5b), both the generator and the heat exchanger decrease their costs as the temperature increases, while in the evaporator assembly there is an unusual behavior which suggests that the variation in the compressor air temperature does not infer exergy destruction or, therefore, related costs. In addition, it is evident that both the evaporator assembly and the heat exchanger have lower costs in the microturbine than in the ARS, which can be explained by the operating conditions of the refrigerant fluids, which directly influence the efficiency of the general system.

### *3.4. Relative Cost Di*ff*erence and Exergo-Economic Factor*

The exergo-economic factor results evidenced that, in the systems, there are values of 100% for the compressor and this remains constant for changes of temperature as shown in Figure 6 because the

flow energy input in microturbine is free and the destroyed exergy cost rate of the compressor is equal to zero.

**Figure 6.** Exergoeconomic factor of the trigeneration system at different compressor inlet air temperatures, (**a**) microturbine components, and (**b**) absorption chiller components.

Therefore, in the microturbine, the compressor is the equipment in the microturbine with the lowest exergo-economic values, which implies that destruction costs and maintenance costs are relevant in the system. On the other hand, in the ARS system, both the generator and the heat exchanger have high values compared to the microturbine system, which is the same in the relative cost behavior of the evaporator assembly with a very high exergo-economic factor of around 26%. Thus, the ARS has a higher exergo-economic factor than the microturbine system, which is a consequence of high acquisition costs, which are more relevant than the costs for exergy destruction and maintenance, which suggests a reduction of non-fixed costs that can be modified.

In this case, the relative costs increase as the compressor air temperature increases, which indicates a thermodynamic limitation in the system due to the temperature limits allowed. However, a temperature increase of 10 K provides improvements in the overall performance of the microturbine components, while the other subsystem suffers increases of 6.53% and 2.84% of the exergo-economic factor in the heat exchanger and generator, respectively, as shown in Figure 6b.

Figure 7 enables the analysis of the relative cost difference of the main components of the trigeneration systems in this study. Thus, we obtained the tendency of this thermo-economic indicator in a wide range of operations, which allows us to determine the critical equipment that represents the main construction cost of the system and, subsequently, to reduce these costs through design strategies or operational changes of the general thermal system.

In the microturbine subsystem (Figure 7a), the results show a very little effect of the inlet air compressor temperature on the relative costs associated with the different components, which is a consequence of the law variation presented in the exergy destroyed by component and similar values of product and fuel cost in the range of the evaluated temperature.

The results show that, in the case of the ARS, the variation of the compressor inlet temperature causes a decrease in the relative cost of the evaporator assembly (21.63%) with increasing temperature from 305.15 K to 315.15 K, as shown in Figure 7b. However, this behavior does not occur in the entire range studied, which allows us to predict the thermo-economic indicator of the trigeneration system under different ambient temperatures.

**Figure 7.** Relative cost difference of the trigeneration system at different compressor inlet air temperatures, (**a**) Microturbine components, and (**b**) Absorption chiller components.

### **4. Conclusions**

This study has been carried out considering a trigeneration system with defined limits and considerations. The modification of these considerations and extension of the limits of this process will require an evaluation of some exergetic costs in the exergo-economic model that were not considered, in addition to the exergy required in the condenser, evaporator and absorber currents. The exergetic costs not considered, if evaluated, could compromise the thermo-economic viability of the system, which must be evaluated in detail by means of thermo-economic indicators such as the recovery period of the investment, the specific investment cost and the leveled cost of the energy since this was not defined in the scope of the present study.

The study allowed us to undertake a thermo-economics approach to evaluate energy conversion systems from the energy perspective, and in a broad way, in complement with economic considerations. It also enabled us to verify the real viability of the trigeneration system when it operates with different air temperatures at the inlet of the compressor.

The thermo-economic parametric was developed based on a gas microturbine-ARS–HRSG trigeneration system model, evaluating the different exergy destroyed cost, exergo-economic factor, and relative cost according to each component. The result also showed the opportunities for improvement in components, and the amount of useful energy available that can be recovered from the exhaust gases of the gas microturbine. However, there is an operation condition where the exergy is highly destroyed due to the exergy inefficiencies of the equipment, and a greater purchase equipment cost is presented with a high exergo-economic factor.

In addition, it is concluded that the system is technically and economically viable, which represents a considerable alternative for the implementation in the secondary energy sector, where the cooling and power generation is required, such as the operation of shopping malls, supermarkets, and hotels.

**Author Contributions:** Conceptualization: G.V.O.; Methodology: G.V.O. and C.A.P.; Software: G.V.O., J.D.F. and C.A.P.; Validation: G.V.O., J.D.F. and C.A.P.; Formal Analysis: G.V.O., J.D.F. and C.A.P.; Investigation: G.V.O.; Resources: G.V.O. and C.A.P.; Writing-Original Draft Preparation: G.V.O.; Writing-Review & Editing: J.D.F. and C.A.P.; Funding Acquisition: G.V.O.

**Funding:** This work was supported by the Universidad del Atlántico and Universidad Francisco de Paula Santander.

**Acknowledgments:** The authors are grateful to Universidad del Atlántico, Universidad Francisco de Paula Santander, and Jhonatan De la Cruz for his collaboration received in the case study simulation.

**Conflicts of Interest:** The authors declare no conflict of interest.

*Energies* **2019**, *12*, 4643

### **Abbreviations**

The following abbreviations are used in this manuscript:


### **References**


© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

### *Article* **Advance Exergo-Economic Analysis of a Waste Heat Recovery System Using ORC for a Bottoming Natural Gas Engine**

### **Guillermo Valencia Ochoa 1,\*, Jhan Piero Rojas <sup>2</sup> and Jorge Duarte Forero <sup>1</sup>**


Received: 19 November 2019; Accepted: 2 January 2020; Published: 5 January 2020

**Abstract:** This manuscript presents an advanced exergo-economic analysis of a waste heat recovery system based on the organic Rankine cycle from the exhaust gases of an internal combustion engine. Different operating conditions were established in order to find the exergy destroyed values in the components and the desegregation of them, as well as the rate of fuel exergy, product exergy, and loss exergy. The component with the highest exergy destroyed values was heat exchanger 1, which is a shell and tube equipment with the highest mean temperature difference in the thermal cycle. However, the values of the fuel cost rate (47.85 USD/GJ) and the product cost rate (197.65 USD/GJ) revealed the organic fluid pump (pump 2) as the device with the main thermo-economic opportunity of improvement, with an exergo-economic factor greater than 91%. In addition, the component with the highest investment costs was the heat exchanger 1 with a value of 2.769 USD/h, which means advanced exergo-economic analysis is a powerful method to identify the correct allocation of the irreversibility and highest cost, and the real potential for improvement is not linked to the interaction between components but to the same component being studied.

**Keywords:** advanced exergo-economic analysis; waste heat recovery system; ORC; endogenous exergy; exogenous exergy

### **1. Introduction**

In thermodynamics systems, the irreversibility in the components of the system produce exergy destruction, and the continuous improvement of the performance of energy conversion systems based on exergy analysis has been a priority among researchers in this field of study [1], from energy analysis of the thermal systems [2] and trigeneration systems [3]. However, from all the information available with these traditional analyses, there is no relevant development on the analysis carried out with this methods to make further improvements to the components of a system [4], so advanced analyses and economic analysis are proposed to facilitate the thermal and economic improvement of the system.

In recent years, new concepts of exergy, such as endogenous/exogenous and avoidable/ unavoidable exergy destruction, have been employed to obtain relevant information for the identification of irreversibilities and thermodynamic inefficiencies in the systems [5]. In any thermodynamic system, the exergy destruction in the components can be generated in two ways. The first is due to the irreversibilities of the component under study, which is called endogenous exergy destruction, while the second is due to the irreversibilities of the other components that affect the component under study; this is called exogenous exergy destruction [6]. Thus, its optimization process will depend on the technical and economic limitations of the system, so there will only be a part of the exergy destruction

and avoidable/unavoidable investment costs for each component. By uniting these two concepts, it is possible to identify improvements to the system [7].

Long et al. [8] evaluated the importance of the working fluid in the thermal performance of an organic Rankine cycle (ORC) by means of an external and internal exergetic analysis, and an optimization analysis based on a genetic algorithm with exergetic efficiency as an objective function. The results of the exergy analysis showed that the organic working fluid affects the exergetic efficiency of the cycle, with the opposite case in the internal part, where the efficiency did not present changes. The results of the optimization showed that the selection of the working fluid depends, to a greater degree, on the optimal evaporation temperature, which increases the exergetic efficiency of the cycle. Long et al. [9] performed an exergy analysis to evaluate the impacts of the evaporation pressure and ammonia fraction on the ammonia–water mixture of the system performance Kalina, obtaining that the evaporation pressure plays an important value in the internal and external exergetic efficiency. Additionally, optimal values are obtained from these in their ideal operation, as well as the ammonia fractions increasing the exergetic efficiency depending on the evaporation pressure. However, the exergetic efficiency of the cycle depends on the input temperature of the heat source, evaluating the impact of this parameter on internal and external exergetic efficiency.

Tian et al. [10] developed a techno-economic analysis of a system consisting of an ORC and an internal combustion engine operating with 235 kW diesel, in order to study the performance of 20 organic fluids. The results showed that the highest energy generated per unit of mass flow and the highest energy efficiency are obtained for refrigerant R-141b and refrigerant R-123, respectively. The study is limited to a single engine operating condition, and a traditional exergetic analysis where the real opportunities for both endogenous and exogenous component savings are not shown.

On the other hand, Zare V. [11], in order to find savings opportunities, added economic criteria to the thermal performance studies, applied to three configurations of an ORC. However, this application was limited to binary geothermal power plants, where the RORC presented better energy results, while from the economic point of view, the simple ORC was the best option because it is integrated by a smaller amount of equipment, which implies a lower acquisition cost. The results do not consider the evaluation of costs by components but of a global system. In addition, studies from the exergetic point of view have been developed in a traditional way, and thermal-economical studies for waste heat recovery systems of gas generation engines through ORC have not been widely integrated. Thus, the literature reports the results of the modeling developed by Kerme and Orfi [12], who studied the effect of the temperature of the organic fluid at the entrance of an ORC turbine on the energy and exergy performance, obtaining that the increase of the temperature increases the efficiency while total exergy destruction decreases it.

The combination of traditional and advanced exergetic analysis can provide significant thermodynamic information, such as the source and the amount of exergy destroyed by each component [13], and how much this destruction can be avoided [14], as in the case of solar energy collectors with a flat plate and a flat plate with a thin plate, resulting in the exergy destruction in the absorbent plate being greater than the rest of components, but according to the advanced exergetic analysis performed, this exergy destruction is endogenous and unavoidable, which means that the irreversibilities of this component are inherent in its operation mode [15].

Mohammadi et al. [16] studied a combination of conventional and advanced exergetic analyses in a supercritical CO2 recompression cycle to determine the potential for improving the thermal cycle performance, where the overall exergetic efficiency reached 17.13%, the system's maximum best potential was 106.85 MW, and approximately 35% of exergy destruction could be avoided by focusing on components, such as the heat exchanger, turbine, and main compressor. These investigations can be complemented with the help of the combination of exergetic analysis [17] and economic analysis to obtain thermo-economic costs based on the irreversibilities of the components [18].

In addition, comparative studies have been carried out on different configurations of waste heat recovery cycles integrated to gas engines [18] and applications, such as Petrakopoulou et al. [19], where the first application of an exergo-economic analysis in a CO2 capture power plant was evaluated, revealing that the costs associated with exergy and investment analyses are endogenous for most components, where it proposed a suggestion for improving some components, such as the reactor, expander, and compressor. The literature review shows the case of a polygeneration plant operating in a geothermal cascade system coupled to an organic Rankine cycle that produces 40 kWe, where improvement potentials were found in the ORC cycle (10.61 kW) and heat exchanger (2.28 kW), while the exergo-economic analysis revealed an electricity production cost of 7.78 \$/h and the advanced exergo-economic analysis suggests that the plant heat exchanger is the component with the greatest opportunity to reduce the exergy destruction of the heat exchanger equipment [20].

Another application was in a combined steam-organic Rankine cycle to recover waste heat from a gas turbine, where an exergo-economic analysis was performed using three different organic fluids (R124, R152a, and R34a), obtaining that the maximum exergy efficiency and the minimum rate of product costs are 57.62% and 396 \$/h, respectively. In addition, the parametric study was complemented with genetic algorithm optimization, where it was obtained that the combined cycle with R152a has the best performance from the thermodynamic and exergo-economic point of view among the fluids analyzed [21].

Advanced exergetic analyses have focused on the ORC cycle, taking into account the advantage of adapting this cycle to another thermal system for different applications, such as waste heat recovery [22], thermodynamic optimization [23], and emergy analysis [24]. Also, several works have combined these studies to obtain improvement potentials. In applications in turbocharged combustion engines, conventional exergetic analysis gives the evaporator and the expander priority improvement potential while advanced exergy analysis suggests the expander and pump as a priority, and the cycle exergy destruction can be reduced by 36.5% [25]. For applications of advanced exergo-economic analysis taking into account waste heat recovery in geothermal applications, low-temperature solar applications, and waste heat recovery from engine gases, the exergetic efficiency of the ORC improves by 20%, optimizing the system through advanced exergetic analysis and proposing the expander, evaporator, condenser, and pump as improvement potentials. Different organic fluids have been tested in the ORCs to improve their performance, obtaining that pentane, cyclohexane, iso-butene, iso-pentane, and cyclohexane have the highest avoidable endogenous cost corresponding to the heat sources evaluated. In addition, the avoidable endogenous cost is sensitive to the heat source temperature, and it is possible to reduce the heat source temperature increase from 100 to 150 ◦C by 28% [26]. Therefore, it has been identified that advanced exergetic and thermo-economic analysis is one of the alternatives to achieve technically and economically favorable operating conditions, and to achieve its application in real conditions.

In response to the inadequate management of energy resources in industrial processes, there is a need to improve the efficiency of equipment and processes, in addition to reducing the environmental impact. Thus, the energy recovery of the exhaust line of the natural gas generation engines is one of the alternatives to increase the thermal efficiency of these systems [27]. However, this issue has been approached from different approaches but not articulated with alternative generation systems, which leads to an enormous scientific impact since if it is true that different ORC configurations have been studied, these have not been studied from an advanced exergetic point of view and integrated with thermo-economic modeling in real contexts of operation of stationary high-power natural gas turbocharged engines as a means of heat recovery, in order to obtain technically and economically viable solutions that allow their commercial application [28].

Thus, the main contribution of this work was to perform an advanced thermo-economic analysis of an organic Rankine cycle for a bottoming natural gas engine, and its respective comparison with the results obtained with conventional exergetic and exergo-economic analyses. The analysis of the irreversibilities of each component is presented, and the possible improvements to the cycle are found using the concepts of endogenous/exogenous and avoidable/unavoidable exergy destruction, combinined with the exergo-economic analysis, thus finding the advance cost rate improvement opportunities for each component based on the irreversibilities of the thermal system.

### **2. Methodology**

### *2.1. Description of the Cycle*

The cycle to be analyzed can be seen in Figure 1. The natural gas generation engine operates with an air/natural gas mixture, which is compressed before it enters the cylinders to improve the engine's thermal performance. The exhaust gases are expanded by means of a turbo compressor flow S1 (708 K, 102 kPa), where energy is transferred by means of the heat exchanger 1 (HX1) to the thermal oil in stream S5, and discharged to the environment in stream S2. The thermal oil circulates through the energy supplied by the thermal oil pump (P1), and the hot fluid coming out of HX1 in stream S3 (616 K, 101.4 kPa) operates as a thermal source to evaporate and reheat through the heat exchanger 2 (HX2) the organic fluid, which is toluene in this case study. The maximum values of thee pressure and temperature of the organic Rankine cycle are presented in the turbine inlet (546 K, 675 kPa), where the organic fluid then expands into the S7 stream (475 K, 22 kPa) in the turbine (T1), generating additional energy without increasing the fuel consumption. To complete the thermal cycle, the organic fluid decreases the pressure to its lowest point, passing to the condensation stage from S7 to S8 (338 K, 675 kPa).

**Figure 1.** The organic Rankine cycle waste heat recovery system.

The 2 MW Jenbacher engine JMS 612 GS-N. L was modeled and studied, as shown in Figure 2, with its technical specifications and nominal operating conditions [29]. This engine operates with natural gas as fuel, since its high robustness allows it to better adapt to variable load regimes. This engine is widely used for self-generation purposes worldwide and is installed in a company of the plastic sector in the city of Barranquilla, Colombia without any waste heat recovery system. The engine regulates fuel consumption to operate between a minimum load of 1000 kWe and a maximum load of 1982 kWe, with an excess air number (lambda) of 1.79 and 1.97, respectively, generating unused exhaust gases in each of its 12 cylinders with a temperature ranging from 580 to 650 ◦C.


**Figure 2.** Jenbacher JMS-612 GS-N.L technical specification.

### *2.2. Energy and Exergy Analyses*

The exergy analysis is defined from the second law of thermodynamics, which means that, unlike the analysis of energy, it depends on the ambient temperature and pressure in which the process in study operates, which allows any system to be investigated in changing environmental conditions. The following assumptions were considered to develop thermodynamic modeling of the RORC:


The global equation of exergy balance, valid for any volume control system, is shown in Equation (1) [21]:

$$\dot{X}\_{\text{heat}} + \sum\_{i=-1}^{n} \left( \dot{m}\_{i} \cdot \varepsilon\_{i} \right)\_{\text{IN}} = P + \sum\_{i=-1}^{k} \left( \dot{m}\_{i} \cdot \varepsilon\_{i} \right)\_{\text{OUT}} + \dot{E}\_{\text{ex, D}} \tag{1}$$

where . *Xheat* is the exergy of heat transfer in kW, . *m* is the mass flow in kg/s, ε is the specific entropy in kJ·<sup>K</sup> kg , *<sup>P</sup>* is the power in kW, and . *Eex*, *<sup>D</sup>* is the exergy destruction [30]. The exergy by heat transfer at temperature T is defined according to Equation (2) [31]:

$$
\dot{X}\_{\text{heat}} = \sum \left( 1 - \frac{T\_0}{T} \right) \dot{Q}. \tag{2}
$$

The exergy of a fluid flow stream is defined as the energy power of the fluid flow, with the mass flow ratio of the fluid flow and the pressure and temperature of the fluid flow being necessary, as well as knowing the environmental conditions (pressure and temperature) in which the fluid flow operates [32]. Therefore, the exergetic power of the fluid flow stream was calculated according to Equation (3): .

$$
\dot{E}\_{\text{ex}\_r \text{ i}} = \dot{m}\_{\text{i}} \varepsilon\_{\text{i} \prime} \tag{3}
$$

where the specific exergy (ε*i*) was calculated according to Equation (4) [33]:

$$
\varepsilon\_i = (h\_i - h\_0) - T\_0 (s\_i - s\_0). \tag{4}
$$

Therefore, the exergy efficiency (η*ex*) for a thermal system was calculated according to Equation (5), as a function of the exergy output ( . *Eout*) and exergy input ( . *Ein*) to the system or device:

$$
\eta\_{\rm ex} = \frac{\dot{E}\_{\rm out}}{\dot{E}\_{\rm in}}.\tag{5}
$$

*Energies* **2020**, *13*, 267

Some energy and exergy performance indicators were calculated for the waste heat recovery system based on the ORC [34], cycle thermal cycle efficiency (η*I*, *<sup>c</sup>*), calculated according to Equation (6); the heat recovery efficiency (ε*hr*) as shown in Equation (7); and the overall energy conversion efficiency <sup>η</sup>*I*, *overall* , given by Equation (8) [20]:

$$
\eta\_{l,C} = \frac{\dot{W}\_{\text{net}}}{\dot{Q}\_G},
\tag{6}
$$

$$\varepsilon\_{\rm lr} = \frac{\dot{Q}\_G}{\dot{m}\_{10} \mathbb{C}\_{P10} (T\_{10} - T\_0)} \,\tag{7}$$

$$
\eta\_{\text{l.overall}} = \eta\_{\text{l.C.E.}} \varepsilon\_{\text{l.r.}} \tag{8}
$$

In addition, to measure the thermal efficiency improvement, the increase in thermal efficiency was calculated through Equation (9), as a function of the net power generated by the ORC ( . *Wnet*), and the heat supplied by the fuel mass rate ( . *mf uel*):

$$
\Delta\eta\_{thermal} = \frac{\dot{\mathcal{W}}\_{net}}{\dot{m}\_{fuel} \cdot LHV}.\tag{9}
$$

The specific fuel consumption (BSFC) was calculated by the mean of Equation (10) [20], and the absolute reduction in the specific fuel consumption as a consequence of the waste heat recovery was calculated as presented in Equation (11):

$$BSFC\_{ORC-engine} = \frac{\dot{m}\_{fuel}}{\dot{W}\_{engine} + \dot{W}\_{net}} \,' \tag{10}$$

$$
\Delta BSFC = \frac{\left| BSFC\_{ORC-engine} - BSFC\_{engine}\right|}{BSFC\_{engine}} \cdot 100.\tag{11}
$$

### *2.3. Advanced Exergetic Analysis*

In advanced exergetic analysis, the values of exergy destruction are divided into four basic parts: Endogenous, exogenous, avoidable, and unavoidable exergy destruction. Avoidable and unavoidable exergy destruction refers to the system improvement potentials. The destruction of avoidable exergy, . *ED*,*<sup>c</sup> AV* , represents the potential for improvement, during the destruction of unavoidable exergy, . *ED*,*<sup>c</sup> UN* , which represents the limitations. The avoidable part of exergy destruction is described in Equation (12):

$$
\dot{E}\_{D,\varepsilon}^{\quad AV} = \dot{E}\_{D,\varepsilon} - \dot{E}\_{D,\varepsilon}^{\quad UN} \tag{12}
$$

where the destruction of unavoidable exergy can be calculated with Equation (13):

$$
\dot{E}\_{D,\varepsilon}^{\dots lIN} = \dot{E}\_{P,\varepsilon} \left(\frac{\dot{E}\_{D,\varepsilon}}{\dot{E}\_{P,\varepsilon}}\right)^{lIN}.\tag{13}
$$

The destruction of endogenous exergy, . *ED*,*<sup>c</sup> EN* , and exogenous . *ED*,*<sup>c</sup> EX* are related to the operational relation between the components of the system. The endogenous part of the exergy destruction is associated only with the irreversibilities that occur in component *c*, where all the other components operate ideally, and component *c* operates with its real conditions. On the other hand, the exogenous part of the exergy destruction is produced by the other components. This part can be determined by subtracting the endogenous exergy destruction from the real exergy destruction of component *c*, as shown in Equation (14): .

$$
\dot{E}\_{D,\varepsilon}^{\quad EX} = \dot{E}\_{D,\varepsilon} - \dot{E}\_{D,\varepsilon}^{\quad EN} \,. \tag{14}
$$

In addition, the destruction of unavoidable endogenous exergy, . *ED*,*<sup>c</sup> UN*, *EN* , was calculated by the Equation (15), the destruction of unavoidable exogenous exergy, . *ED*,*<sup>c</sup> UN*,*EX* , through Equation (16), the destruction of avoidable endogenous exergy, . *ED*,*<sup>c</sup> AV*, *EN* , with Equation (17), and the destruction of avoidable exogenous exergy, . *ED*,*<sup>c</sup> AV*,*EX* , by mean of Equation (18) [7]:

$$
\dot{E}\_{D\mathcal{L}} \stackrel{lIN, EN}{=} \dot{E}\_{D\mathcal{L}} \stackrel{EN}{\left(\dot{E}\_{D\mathcal{L}}\right)} \stackrel{lIN}{\quad} \tag{15}
$$

$$
\dot{E}\_{D,\varepsilon}^{\text{LIN, EX}} = \dot{E}\_{D,\varepsilon}^{\text{LIN}} - \dot{E}\_{D,\varepsilon}^{\text{LIN, EN}},\tag{16}
$$

$$
\dot{E}\_{D,\varepsilon}^{\quad AV,EN} = \dot{E}\_{D,\varepsilon}^{\quad EN} - \dot{E}\_{D,\varepsilon}^{\quad UN,EN} \tag{17}
$$

$$
\dot{E}\_{D,\varepsilon}^{\quad AV,EX} = \dot{E}\_{D,\varepsilon}^{\quad AV} - \dot{E}\_{D,\varepsilon}^{\quad AV,EN}.\tag{18}
$$

### *2.4. Conventional Exergo-Economic Analysis*

Exergetic analyses are used to determine the location, type, and magnitude of thermodynamic inefficiencies in system components. On the other hand, exergo-economic analyses combine the concept of exergy and economic analyses to obtain a tool for the optimization of energy systems [4]. In addition, the economic model takes into account the components' cost, including amortization, maintenance, and fuel costs. To define a cost function that depends on interest optimization parameters, the cost of components must be expressed in terms of thermodynamic design parameters. The cost balance equations applied to component *c* of the system under study show that the sum of rates associated with all outgoing exergy flows equals the sum of cost rates of all incoming exergy flows, plus those corresponding to charges due to capital investment and operating and maintenance costs, as shown in Equation (19) [12]:

$$
\sum\_{\epsilon} \dot{\mathcal{C}}\_{\epsilon,\varepsilon} + \dot{\mathcal{C}}\_{w,\varepsilon} = \dot{\mathcal{C}}\_{q,\varepsilon} + \sum\_{i} \dot{\mathcal{C}}\_{i,\varepsilon} + \dot{Z}\_{\varepsilon,\prime} \tag{19}
$$

where the cost rates of input exergy flow ( . *Ci*) are defined in Equation (20), the cost rates of output exergy flow in Equation (21), the heat transfer cost rate in Equation (22), and the cost rate related to energy transfer by work in Equation (23) [35,36]:

$$
\dot{\mathbf{C}}\_i = \ c\_i \dot{E}\_i = \ c\_i \Big[ \dot{m}\_i \mathbf{e}\_i \Big],\tag{20}
$$

$$
\dot{\mathbf{C}}\_{\varepsilon} = \mathbf{c}\_{\varepsilon} \cdot \dot{\mathbf{E}}\_{\varepsilon} = \mathbf{c}\_{\varepsilon} \begin{bmatrix} \dot{m}\_{\varepsilon} \mathbf{c}\_{\varepsilon} \end{bmatrix} \tag{21}
$$

$$
\dot{\mathbf{C}}\_{q} = \mathbf{c}\_{q} \dot{\mathbf{E}}\_{q\nu} \tag{22}
$$

$$
\dot{\mathbf{C}}\_{\text{av}} = \mathbf{c}\_{\text{av}} \dot{\mathbf{W}}\_{\text{\prime}} \tag{23}
$$

where *ci*, *ce*, *cq*, and *cw* are the costs per unit of exergy in \$/GJ, and . *Zc* is the sum of the cost rates associated with the cost of capital investment, . *Zc CI* , and operation and maintenance costs, . *Zc OM* , as shown in Equation (24) [37]:

$$
\dot{Z}\_{\varepsilon} = \dot{Z}\_{\varepsilon}^{\text{CI}} + \dot{Z}\_{\varepsilon}^{\text{CI}} = \text{CRF} \cdot \left[ \frac{q\rho\_r}{N} \cdot 3600 \right] \text{PEC} \tag{24}
$$

where the *PECc* is the purchase equipment cost of component *c*, which is given for all components of the system; *N* is the number of annual hours that the unit operates; and ϕ*<sup>r</sup>* is the maintenance factor, which is generally approximately 1.06 [38]. The modeling and sizing of the plate heat exchanger equipment was done for each of the evaporator, condenser, and recovery zones [39]. In addition, thermoeconomic modeling and cost balances for the integrated configuration with the engine were developed by the authors [40]. Also, the CRF is the capital recovery factor, which depends on the interest rate and the estimated lifetime of the equipment, which was calculated based on Equation (25):

$$\text{CRF} = \frac{i[1 + i]^n}{[1 + i]^n - [ \\']} \tag{25}$$

where *i* is the interest rate, and *n* is the total period of operation of the system in years.

### *2.5. Advanced Exergo-Economic Análisis*

### 2.5.1. Unavoidable and Avoidable Costs

The avoidable and unavoidable cost ratios associated with the exergy destruction within each component of the system were calculated by Equations (26) and (27), respectively:

$$
\dot{\mathcal{L}}\_{D,\varepsilon} \stackrel{lIN}{=} \mathcal{c}\_{\mathcal{F},\varepsilon} \dot{\mathcal{E}}\_{D,\varepsilon} \stackrel{lIN}{\quad} \tag{26}
$$

$$
\dot{\mathcal{C}}\_{D\mathcal{E}}{}^{AV} = \mathcal{c}\_{\mathcal{F}\mathcal{E}} \dot{E}\_{D\mathcal{E}}{}^{AV}{}\_{\mathcal{I}} \tag{27}
$$

where the sum of the avoidable and unavoidable costs of exergy destruction is equal to the total cost associated with exergy destruction, as shown in Equation (28) [5]:

$$
\dot{\mathbf{C}}\_{\rm D\mathcal{L}} = \mathbf{c}\_{\rm F\mathcal{L}} \cdot \dot{\mathbf{E}}\_{\rm D\mathcal{L}} = \dot{\mathbf{C}}\_{\rm D\mathcal{L}}^{\rm lIN} + \dot{\mathbf{C}}\_{\rm D\mathcal{L}}^{\rm AV} \,. \tag{28}
$$

The unavoidable investment cost was calculated considering an extremely inefficient version of component *c* [41]. Therefore, for the calculation of the unavoidable investment cost rate in the components [42], some operational conditions were proposed, as shown in Table 1.

**Table 1.** Main assumptions for real conditions, ideal conditions, unavoidable exergy, and unavoidable investment cost.


The avoidable investment cost was calculated by subtracting the unavoidable cost rate from the total investment cost, as shown in Equation (29):

$$
\dot{Z}\_{\mathfrak{c}}^{\text{AV}} = \dot{Z}\_{\mathfrak{c}} - \dot{Z}\_{\mathfrak{c}}^{\text{LIN}}.\tag{29}
$$

2.5.2. Endogenous and Exogenous Cost Rates

The endogenous cost rates . *CD*,*<sup>c</sup> EN* and exogenous . *CD*,*<sup>c</sup> EX* are defined according to Equations (30) and (31), respectively: .

$$
\dot{\mathbf{C}}\_{\rm D\_{\mathcal{E}}}^{\rm EN} = \mathbf{c}\_{\rm F\_{\mathcal{E}}} \dot{\mathbf{E}}\_{\rm D\_{\mathcal{E}}}^{\rm EN},\tag{30}
$$

$$
\dot{\mathcal{C}}\_{D\mathcal{E}}{}^{EX} = \mathcal{c}\_{\mathcal{F}\mathcal{E}} \dot{E}\_{D\mathcal{E}} \, ^{EX} \, \text{\,} \tag{31}
$$

where the sums of the endogenous and exogenous cost rates of exergy destruction are equal to the total cost rate associated with exergy destruction, as defined in Equation (32):

$$
\dot{\mathcal{C}}\_{D\mathcal{L}} = \mathcal{c}\_{\mathcal{F}\mathcal{L}} \dot{E}\_{D\mathcal{L}} = \dot{\mathcal{C}}\_{D\mathcal{L}}{}^{EN} + \dot{\mathcal{C}}\_{D\mathcal{L}}{}^{EX} \,. \tag{32}
$$

Therefore, the exogenous investment rate was calculated with Equation (33) as follows:

$$
\dot{Z}\_{\mathfrak{c}}^{\mathit{AV}} = \dot{Z}\_{\mathfrak{c}} - \dot{Z}\_{\mathfrak{c}}^{\mathit{UV}}.\tag{33}
$$

### 2.5.3. Splitting Cost Rates

The cost rates with respect to exergy desegregation can be calculated with Equations (34)–(37) [6]:

$$
\dot{\mathbf{C}}\_{D\boldsymbol{\omega}}{}^{ENAV} = \mathbf{c}\_{\mathbb{F}\boldsymbol{\omega}} \dot{\mathbf{E}}\_{D\boldsymbol{\omega}}{}^{ENAV} \, \tag{34}
$$

$$
\dot{\mathcal{C}}\_{D,\varepsilon}^{EN,IN} = \left. c\_{F,\varepsilon} \dot{E}\_{D,\varepsilon}^{EN,IN} \right. \tag{35}
$$

$$
\dot{\mathbf{C}}\_{\mathrm{D}\mathcal{L}}^{\mathrm{EX},AV} = \mathfrak{c}\_{\mathrm{F}\mathcal{L}} \dot{\mathbf{E}}\_{\mathrm{D}\mathcal{L}}^{\mathrm{EX},AV} \tag{36}
$$

$$
\dot{\mathcal{C}}\_{D,\varepsilon}{}^{EX,lIN} = \varepsilon\_{F,\varepsilon} \dot{E}\_{D,\varepsilon}{}^{EX,lIN} \,, \tag{37}
$$

where . *CD*,*<sup>c</sup> EN*,*AV* represents the unavoidable cost rate without component *c*, associated with the operation of the same component, and the value can be reduced by optimizing the component through technological improvements. Also, the cost . *CD*,*<sup>c</sup> EX*,*AV* is the avoidable exogenous cost rate that can be reduced by optimizing other components of the cycle while the costs . *CD*,*<sup>c</sup> EN*,*UN* and . *CD*,*<sup>c</sup> EX*,*UN* are the unavoidable endogenous and unavoidable exogenous cost rates, respectively.

The endogenous cost rate of component *c* can be calculated with Equation (38), and the rate of unavoidable endogenous investment costs can be calculated using Equation (39) [7]:

$$
\dot{Z}\_c^{\quad EN} = \dot{E}\_{P,c} \left(\frac{\dot{Z}}{\dot{E}\_P}\right)\_c^{lIN} \tag{38}
$$

$$
\dot{Z}\_{\mathcal{L}}^{EN,LIN} = \dot{E}\_{P,\mathcal{L}} \begin{bmatrix} \dot{Z}\_{\mathcal{L}} \\ \dot{E}\_{P,\mathcal{L}} \end{bmatrix} \tag{39}
$$

Therefore, the equations used to calculate the desegregated investment costs are presented from Equations (40) to (42): .

$$
\dot{Z}\_{\text{c}}^{\text{EN,AV}} = \dot{Z}\_{\text{c}}^{\text{EN}} - \dot{Z}\_{\text{c}}^{\text{EN,IN}},\tag{40}
$$

$$
\dot{Z}\_{\text{c}}^{\text{EX,LIN}} = \dot{Z}\_{\text{c}}^{\text{LIN}} - \dot{Z}\_{\text{c}}^{\text{EN,LIN}},\tag{41}
$$

$$
\dot{Z}\_{\varepsilon}^{EX,AV} = \dot{Z}\_{\varepsilon}^{EX} - \dot{Z}\_{\varepsilon}^{EN,IN}. \tag{42}
$$

### **3. Results and Discussion**

In this section, the influence of the engine load on the heat recovery system energy and exergy performance integrated into the natural gas engine was studied, as an alterntive to reduce the global operational cost and increase the thermal efficiency [43]. The engine power control system adjusts internal engine variables, such as the pressure and temperature of the air-fuel mixture before entering the cylinders, and the recirculation percentage, to provide high efficiency in partial load operation of the gas engine. Some energy indicators were proposed to study the performance of the waste heat recovery system based on ORC, as shown in Figure 3, while the evaporating pressure was set to 675.8 kPa, and toluene was selected as the working fluid [33]. For safety restriction, all feasible operating points of the proposed system at different engine loads ensured that toluene evaporates completely at the outlet of the evaporator to prevent corrosion of the liquid in the expander, in addition to a gas temperature at the outlet of the evaporator (state 11) being higher than the acid dew temperature (200 ◦C) to avoid acidic corrosion of the exhaust [34].

**Figure 3.** Energy and exergy indicator of the waste heat recovery system with different organic fluids and engine load, (**a**) met power, (**b**) Absolute increase in thermal efficiency, (**c**) specific fuel consumption, (**d**) global energy conversion efficiency, (**e**) ORC thermal efficiency, and (**f**) global exergetic efficiency.

The results show that the absolute increase in thermal efficiency (Figure 3b) decreases for toluene and cyclohexane, as does the overall energy conversion efficiency (Figure 3d) with an increasing engine load, while the net power output (Figure 3a) presents its maximum value with toluene (89.4 kW—97.9%), cyclohexane (73.2 kW—97.9%), and acetone (53.2 kW—91.81%), respectively. However, in an engine operating interval, acetone presents a slight increase in thermal efficiency, and then decreases, presenting a maximum at an 82.68% engine load.

These results are due to a higher engine load, implying an increase in the exhaust gas flow according to the first and second laws of thermodynamics [44] while a greater energy loss is presented in the recuperator heat exchanger 1 (HX1) because of the evaporation pressure and thermal oil temperature have been limited. As the engine load increases, there is an increase in the fluid evaporating temperature at the evaporator. Therefore, the power increases, which is the main factor for thermal and exergetic efficiency. However, the isentropic turbine efficiency decreases slightly as a consequence of the increase in the thermal oil temperature, causing a decrease in the energy indicators at high engine loads. Likewise, the tendency to increase the power with the engine load is a consequence of both the increase in the inlet thermal oil temperature to the evaporator, which leads to an increase in the toluene mass flow, and the enthalpy difference between the outlet and the inlet of the pump and turbine, but this is more relevant in the turbine.

In addition, the results obtained from the traditional exergetic and exergo-economic analysis are shown in Table 2, where the exergy and fraction of exergy destroyed, *yD*,*c*, shows that the greatest values are present in the heat exchanger 1 (shell and tube heat exchanger) with 32.54%, the evaporator (28.32%), and the condenser with 27.97%. The component with the highest destroyed exergy value (41.95 kW) is heat exchanger 1, being one of the components with the lowest exergetic efficiency of the cycle, due to the large heat exchanger area required and the high temperature difference. The greater the investment and the cost of exergy destroyed, the greater the influence of the component in the system, therefore, the component with the greatest improvement in cost efficiency of the total plant can be defined. In the case study, the components with the greatest opportunities for improvement in this ratio are the condenser and HX 1. Therefore, these components are the most important from a thermodynamic point of view.


**Table 2.** The results of conventional analysis for all components in the waste heat recovery system.

The exergo-economic factor, *fc*, is the effective parameter that allows us to compare and evaluate the components that make up the system. A high value for this parameter indicates that for the component under study, acquisition costs predominate over operation and maintenance costs. For example, in the case of the condenser, which is the component with the lowest value of the exergo-economic factor, it can be concluded that expenses are mostly related to operating and maintenance costs compared to acquisition costs.

By means of the advanced exergetic analysis, the exergy destruction can be disaggregated for each one of the components. In this way, the real possibilities of improvement can be determined both through the operational and design point of view of the component, and the global consideration of the thermal system. From the solution of Equations (12)–(18), and the unavoidable operation conditions described in Table 1, the disaggregation of the exergy can be found in its endogenous, exogenous, avoidable, and unavoidable part, as well as avoidable and unavoidable endogenous and its avoidable and unavoidable exogenous counterpart. The determination of the avoidable part of the destroyed exergy is a significant step because it allows identification of opportunities for improvement in the component and its interaction with the rest of the components. Also, this result allows knowledge of which is the optimal way to increase the thermal efficiency of the system, besides providing valuable information about how the components operate together as a global system. Figure 4 presents a graphical version of the improvement opportunities in each component from the exergetic point.

**Figure 4.** Advanced exergy analyses for each component in the ORC cycle, (**a**) HX 1, (**b**) P1, (**c**) turbine, (**d**) P2, (**e**) evaporator, and (**f**) condenser.

The results of the advanced exergetic analysis and economic exergetic analysis are presented in Table 3, where the disaggregation of the destroyed exergy was calculated as a function of the endogenous, exogenous, avoidable, and unavoidable for each of the components under study. The results show that most of the destroyed exergy is endogenous (78.53% of the total destroyed exergy), emphasizing that the interaction between components does not have a significant effect on the overall exergetic performance of the cycle. Similarly, it is noted that the component with the greatest avoidable exergy destroyed in the system is the turbine, with a value of 11.075 kW, where 69.625% is endogenous and 30.374 is exogenous, which means that in the turbine, there is a real great opportunity for improvement. On the other hand, in the unavoidable part, the components with the greatest technological limitations are the HX 1 and evaporator, representing 96.1% of the total unavoidable exergy of the cycle.


**Table 3.** Splitting exergy destruction for each component.

The equations presented in Sections 2.5.2 and 2.5.3 were used to calculate the advance exergy destruction costs as shown in Table 4, which is based on the result of the advanced destroyed exergy.

It can be observed that the endogenous exergy destruction is higher than the exogenous cost in the components of the thermal cycle, which is a consequence of the high endogenous investment costs values for all components of the system with respect to the exogenous investment cost, as shown in Table 5. Therefore, it can be established that the interaction between components in terms of investment costs is not very relevant in the system; however, for the component under study, it is a parameter of vital importance. Also, it can be observed that the rates of unavoidable investment costs for the components studied showed an inclination in the unavoidable part.


**Table 4.** Advanced exergy destruction cost rates for all components in the waste heat recovery system.

**Table 5.** Advanced investment costs for all components in the WHR (Wast Heat Recovery).


Negative values of exogenous investment cost rates (ZEX, ZAV,EX, ZUN,EX) revealed that investment costs within these components might decrease if investment costs within the other components are increased.

In order to show a comparison of the destroyed exergy relationship between traditional and advanced exergetic analysis, Figure 5 is shown. In Figure 5A, a slight difference of the parameter under study is denoted because the one that was calculated by means of the advanced exergetic analysis only emphasizes the exergy that is destroyed by each component, that is to say, the endogenous part (without the interaction that this one has with its surroundings).

**Figure 5.** Contribution of each component to the overall exergy destruction of the cycle based on (**A**) traditional and (**B**) advanced exergy analysis.

However, the relationship of exergy destroyed by the component was also calculated by emphasizing which fraction is borne by the component itself or by the interaction of the component with its surroundings, as shown in Figure 5B. From this graph, what has been mentioned before is supported, that is, that the interaction between each of the components of the system is not significant in comparison to the exergy that destroys the component under its own operating conditions. A comparative analysis was performed by implementing a new exergo-economic factor calculated by advanced exergetic analysis, as shown in Figure 6.

As a percentage, it can be seen that the exergo-economic factor, as well as the traditional and advanced approach, presents a similarity in the components that make up the system. So, the main efforts should concentrate on designing the most efficient heat exchangers [35], with a smaller heat transfer area and less exergy destruction, without increasing the purchase equipment costs.

**Figure 6.** Comparison of the calculated exergo-economic factor based on traditional and advanced exergy analysis.

### **4. Conclusions**

In this paper, the benefit offered by developed traditional and advanced exergetic analysis in thermal systems was shown, in particular in the organic Rankine cycle systems. Exergetic analysis allows determination of the sources of irreversibility in a thermal system, and therefore indicates the starting points of an optimization procedure, and contributes to the rational use of the energetic resources. In the study carried out, it was possible to determine which equipment that resulted in greater destruction of exergy introduced in the waste heat recovery system based on the organic Rankine cycle. The equipment in which the design or operational improvements can be made was also determined, since the implementation of some recommendations is not practical for optimizing the cycle due to operational or design limitations. Therefore, traditional exergy, advanced exergy, and exergo-economic analysis were applied to gain a better understanding of the system performance. Moreover, a comprehensive comparison was conducted to further assess the system from various points of view.

The conventional exergy showed that the heat exchanger 1 had the largest exergy destruction and exergy destruction, and highest investment costs (41.95 kW, 32.54%, and 2.67 USD/h). The results of the energetic and exergetic analysis of the system showed that the exergy destroyed is a measure of the degree of process irreversibility. Thus, in the case of heat exchanger 1, the causes of the irreversibility were due to the heat transfer through a finite temperature difference higher than 100 ◦C. Similarly, the results of exergy destructions appeared to be in accordance with the exergy efficiencies. That is, a smaller exergy efficiency implies greater exergy destruction in the system components.

Also, the highest exergo-economic factor was found in the pump 2, turbine, and pump 1, with 90.78%, 89.20%, and 85.24%, respectively. These results were a consequence of the high effect of the purchased equipment cost, and the low thermodynamic efficiency in the aforementioned devices, where the probable solution could be the implementation of low-cost components, which are usually characterized by a lower energy efficiency.

Most of the exergy destruction calculated was endogenous (78.53%), emphasizing that the interaction between components does not have a significant effect on the overall exergetic performance of the cycle. The maximum unavoidable exergy was found for the heat exchanger 1, with 90.44%. This indicates that there are not too many ways to improve this component. Nevertheless, other components, such as pump 2, pump 1, and turbine 1, have the minimum unavoidable exergy destruction, with 0.028, 0.03, and 2.819 kW. In addition, the component with the highest cost rate was the condenser with 7.188 USD/h, followed by the heat exchanger 1 with 2.3 USD/h, but the highest avoidable cost rate was found for the turbine with a value of 0.701 USD/h.

On the other hand, the advanced exergo-economic analyses showed that the turbine is the component with the major purchase equipment cost in the system, with a value of 7.898 USD/h, which is 54.09% of the total equipment cost of the system. For all components studied, the endogenous investment cost was higher than the exogenous part, showing the weak relation between them. A comparison was realized between the traditional and advanced exergo-economic factor, which resulted in a similar effect in each component, but the advanced exergy approach presented a slightly higher value, implying that the advanced exergetic analysis gives greater precision in terms of results without ignoring the really great opportunities for improvement.

**Author Contributions:** Conceptualization: G.V.O.; Methodology: G.V.O. and J.D.F.; Software: G.V.O., J.P.R. and J.D.F.; Validation: G.V.O., J.P.R. and J.D.F.; Formal Analysis: G.V.O., J.P.R. and J.D.F.; Investigation: G.V.O.; Resources: G.V.O.; Writing—Original Draft Preparation: G.V.O. and J.P.R.; Writing—Review and Editing: G.V.O. and J.D.F.; Funding Acquisition: G.V.O. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Acknowledgments:** Acknowledgments to Universidad del Atlántico, Universidad Francisco de Paula Santander, and the E2 Energía Eficiente S.A E.S. P company by the support received to conduct this research.

**Conflicts of Interest:** The authors declare no conflict of interest.

### **Abbreviations**

The following abbreviations are used in this manuscript:


*Energies* **2020**, *13*, 267

### **Subscripts**


### **References**


© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

### *Article* **A Li-Ion Battery Thermal Management System Combining a Heat Pipe and Thermoelectric Cooler**

### **Chuanwei Zhang 1, Zhan Xia 1,\*, Bin Wang 2, Huaibin Gao 1, Shangrui Chen 1, Shouchao Zong <sup>1</sup> and Kunxin Luo <sup>1</sup>**


Received: 29 November 2019; Accepted: 11 February 2020; Published: 14 February 2020

**Abstract:** The temperature of electric vehicle batteries needs to be controlled through a thermal management system to ensure working performance, service life, and safety. In this paper, TAFEL-LAE895 100Ah ternary Li-ion batteries were used, and discharging experiments at different rates were conducted to study the surface temperature increasing characteristics of the battery. To dissipate heat, heat pipes with high thermal conductivity were used to accelerate dissipating heat on the surface of the battery. We found that the heat pipe was sufficient to keep the battery temperature within the desired range with a midlevel discharge rate. For further improvement, an additional thermoelectric cooler was needed for a high discharge rate. Simulations were completed with a battery management system based on a heat pipe and with a combined heat pipe and thermoelectric cooler, and the results were in line with the experimental results. The findings show that the combined system can effectively reduce the surface temperature of a battery within the full range of discharge rates expected in the battery used.

**Keywords:** thermal management; Li-ion battery; heat pipe; thermoelectric cooler

### **1. Introduction**

With increasing environmental pollution and fossil fuel depletion, electric vehicles are gradually appearing in the automotive market, as they produce no pollution and use electricity instead of fuel [1,2]. As the power source of electric vehicles, Li-ion batteries are characterized by high energy density, low self-discharge rate, no memory effect, and a small size, so they are favored by automobile manufacturers [3–5]. However, electric vehicles are at risk of spontaneous combustion, attributed to the high surface temperature of the Li-ion battery used in electric vehicles, causing thermal runaway where the temperature of the batteries rises constantly, even causing fire. As a result, the battery does not work normally [6]. During the Li-ion battery charging and discharging process, complex chemical reactions occur in its interior and heat is generated, so the temperature of the battery continuously rises [7]. If the heat cannot be dissipated in time, thermal runaway of the Li-ion battery may occur [6]. Therefore, a battery thermal management system (BTMS) to control the temperature of the battery is essential to ensure the normal working and the safety of electric vehicles.

Some studies have been conducted on BTMS, and various cooling methods for Li-ion batteries have been proposed. For an air-cooling system (ACS), where air is used as the cooling medium, the battery pack is primarily cooled by heat convection. However, the ACS is unable to meet the cooling demands of the battery [8]. A liquid cooling system (LCS) for batteries was also developed to cool or heat the battery pack. Rao et al. [9] designed a BTMS based on liquid cooling for a cylindrical Li-ion battery pack and researched

the change of in-battery contact surface to determine the length of pack block when the flowing rate of cooling liquid is at the optimal efficiency. Under a battery heat generation of 30 W, Pesaran et al. [10] tested oil and air as cooling media. Their results showed that the surface temperature of the battery in the oil cooling system is nearly 10 ◦C lower than in ACS under the same working conditions.

Versatile materials have gradually been developed. Researchers explored materials with special characteristics that can be used in BTMS and applied phase change heat transfer media to BTMS. The principle involves using the characteristics of phase change materials (PCMs) to store or release energy through a phase change without significantly changing the temperature. At high discharge rates, Sabbah et al. [11] conducted experiments using 18,650 batteries to demonstrate that the thermal management effect of PCM cooling is better than that of air. Rao and Zhang et al. [12–14] studied different models and concluded that by allowing batteries to be in direct contact with the PCM, the surface temperature of the battery can be reduced by 18 ◦C. Wu et al. [15] constructed a new type of reinforced composite phase change material (CPCM) based on copper mesh, paraffin, and expanded graphite for BTMS. The heat dissipation and temperature uniformity of the CPCM were better than that of the traditional air convection and PMC.

In addition to the above studies, a heat pipe (HP) with high thermal conductivity in the longitudinal direction and isothermal properties in heat conductive direction, and a thermoelectric cooler (TEC) characterized by refrigeration ability, rapidly transferring heat, small size, and high reliability, were both investigated [16–18]. Joshua et al. [19] applied the HP to BTMS, and found it can better control the temperature of the Li-ion battery than the traditional liquid cooling system. Deng et al. [20] combined an L-style HP with an aluminum plate to build a BTMS and showed that with ambient temperature increasing, the heat dissipation from HP increases and the increasing rate of battery temperature reduces. Wang [18] selected the optimal working current of the TEC and built a BTMS based on it. The results showed the advantage of controlling the temperature of the TEC. High-temperature environments have no significant effect on the performance of the TEC. Lu et al. [21] verified that the TEC is highly effective for heating or cooling the battery pack. Shashank et al. [22] arranged battery cells in an inverted position within a battery pack filled with a PCM, which improved the performance of the BTMS. Further research found that a TEC connected with one side of the battery pack can improve the reliability of the system.

However, at a high discharge rate, the HP cannot sufficiently dissipate battery heat. Though the heat of the battery can be dissipated quickly by the cold side of TEC, the heat produced on the hot side is hard to be dissipated under natural conditions. Therefore, this paper proposes a BTMS combining an HP with TEC. With this system, the heat produced on the hot side of the TEC can be quickly transferred to the device dissipating heat by the HP. Both numerical simulations and discharging experiments were conducted to research the surface temperature rising characteristics of the battery with different discharge rates. To reduce the energy consumption but maintain the surface temperature of the battery within the optimal temperature, we performed experiments during which the TEC began to work at different temperatures.

### **2. Theoretical Analysis**

A TAFEL-LAE895 100 Ah ternary Li-ion battery, which was produced by TAFEL company of Jiangsu province, China, was selected for analysis in this research. Its main parameters are provided in Table 1.


**Table 1.** TAFEL-LAE895 battery specifications.

The heat generated by the battery includes electrochemical reaction heat, joule heat, polarized heat, and side reaction heat [23], which can be calculated by

$$Q = Q\_\mathbf{r} + Q\_j + Q\_p + Q\_s \tag{1}$$

where *Q* is the total heat generated by the battery; *Qr* is the electrochemical reaction heat, which refers to the heat generated by the electrochemical reaction when lithium ions are intercalated and deintercalated between the anode and cathode materials in the process of charging and discharging; *Qj* represents joule heat generated by joule internal resistance; *Qp* represents polarized heat, which refers to the heat generated by the polarized internal resistance due to the polarization of the battery during the charging and discharging process; and *Qs* is side reaction heat, which is the heat generated in the process of the decomposition and reaction of the electrolyte, the thermal decomposition of the anode and cathode materials, and the decomposition of the separator.

Due to the influence of the environment and heat loss of transmission, directly measuring the heat generated by a Li-ion battery in the experiment was difficult. In this work, the calculation method proposed by Bernardi et al. [24] was adopted. Under the assumption that the heat generation inside the Li-ion battery is uniform, a theoretical formula of the heat generation rate was proposed [25,26]. The formula is simplified as

$$Q = I(E - \mathcal{U}I) - IT\frac{\partial E}{\partial T} \tag{2}$$

where *I* represents the charging and discharging current of the Li-ion battery; *E* and *U* are the open and closed circuit voltages, respectively; and *T* represents the surface temperature of the battery, the unit of which is kelvins. The first term on the right side of the equation represents the irreversible reaction heat, such as joule heat and electrochemically polarized heat, which can be represented by *I* <sup>2</sup>*R*, where *R* is the total internal resistance of the Li-ion battery. The second term represents the heat generated by the reversible reaction, and the Li-ion battery generally takes a reference value of 0.042 [25]. After the values are inserted into the formula, we have

$$Q = I^2 R - 0.042I$$

The value of *R* can be tested using a Hybrid Pulse Power Characteristic (HPPC) experiment. The value of *R* is calculated in Table 2.


**Table 2.** The value of *R* (total internal resistance of the Li-ion battery).

The heat generated by the Li-ion battery at different discharge rates can be calculated as shown in Figure 1.

Note that the battery discharge rate is measured in C, which means that battery discharges its rated capacity within the specified time (C = current/rated capacity).

As shown in Figure 1, with increasing discharge rate, the heat-generating rate rises in an accelerated manner, from 11.5 W ata1C discharge rate to 80.7 W ata3C discharge rate.

**Figure 1.** Theoretical rate of battery heat generation.

### **3. Numerical Simulations**

Analyzing the temperature of a battery is a prerequisite for designing a highly efficient battery thermal management system. Numerical simulations were used to examine battery surface temperature change at different discharge rates to verify the effect of the theoretical model and help design the experimental BTMS model.

### *3.1. Simulation Model Design*

3D design software, that is, Solidworks of Dassault Systemes S.A company, was used to build models of the battery, the BTMS based on an HP, and the combined HP and TEC system. The sizes of the models were designed as experimental devices. To achieve a more uniform temperature distribution and to accelerate heat dissipation, aluminum sheets were added to the surface of the battery, and grooved aluminum sheets were fixed to the HP. Heat sinks and fans were added. The model is shown in Figure 2.

**Figure 2.** Battery thermal management system (BTMS) model: (**a**) experimental; (**b**) simulation models.

The 3D model was imported into Fluent, the computational fluid dynamics (CFD) software. For simulation, the following assumptions were made:


(3) The complex internal structure of the HP was regarded as a solid piece with equivalent thermal conductivity [27].

In the transient simulation model, first-order upwind and laminar flow method were adopted, and then the energy equation was turned on. The environmental temperature was set to ambient room temperature, which is 298 K. The boundary conditions of the battery were set to natural convection. The convective heat transfer coefficient was set to 3.5 Wm−2K−1, adding a wind speed of 7 m s−<sup>1</sup> to the heat sink to simulate the airflow speed generated by the fan. The heat sink was composed of aluminum, so the parameters of the heat sink, aluminum sheet, and grooved aluminum sheet were the same: density of 2719 kg m−3, specific heat capacity of 871 J kg−1K−1, and thermal conductivity of 202.4 Wm<sup>−</sup>1K−1.

During the calculation, the residuals of the conservation equations, temperature of characteristics, and the characteristic plane were detected. When the residual was less than 10−<sup>6</sup> and the temperature change at the detection point was less than 0.1%, the calculation was considered to reach convergence.

### *3.2. Simulation Results and Analysis*

### 3.2.1. Simulation of Surface Temperature Rise of Battery in Natural State

For the natural working state of the battery pack, without any additional device to help dissipate heat, the model was meshed with 216,448 elements. The rate of generating heat is provided in Figure 1 and the simulation results of the battery pack are shown in Figure 3.

**Figure 3.** The surface temperature of the battery pack at different discharge rate in natural state: (**a**) temperature distribution at 1 C discharge rate; (**b**) temperature distribution at 3 C discharge rate and (**c**) surface temperature and ratio of difference.

Figure 3c shows the highest and lowest surface temperatures of the battery and the percentage of difference in surface temperature at the end of discharging. In Figure 3, with the increase in the discharge rate, the surface temperature shows an upward trend, as does the temperature difference ratio. The difference in temperature was moderate at less than 1.2%; thus, it was regarded as uniformly distributed. The temperature was 316.1 K at a 1 C discharge rate, which is within the range of optimal working temperature between 293 K to 318 K. However, at 2.5 and 3 C discharge rates, the surface temperatures all exceeded 323 K. Working for a long time at such a high temperature will have an irreversible impact on battery performance and cycle life. The heat generated by the battery will accumulate, causing thermal runaway and may lead to battery fire or explosion.

### 3.2.2. Simulation with an Added HP

With an HP, a heat sink, and a fan attached to the surface of the battery, the model was modeled with 292,605 elements with the same assumptions applied. The heat generation rate of the Li-ion battery selected in this paper was 80.7 W at the maximum 3 C discharge rate. Therefore, the HP with a power of 80 W, a length of 240 mm, a width of 11 mm, and a thickness of 2 mm was selected in this study. The wick material of the selected HP was copper, and the heat pipe was regarded as special copper with high thermal conductivity. The thermal conductivity of the HP is shown in Table 3, which was measured by an experimental method previously reported [27,28].

**Table 3.** The thermal conductivity of the heat pipe (HP).


The density and specific heat capacity of the HP were 4000 kg m−<sup>3</sup> and 400 J kg−1K−1, respectively [27]. The results of the simulation are shown in Figure 4.

**Figure 4.** The surface temperature with different discharge rates based on the heat pipe HP model: (**a**) 1.5 C; (**b**) 2 C; (**c**) 2.5 C; (**d**) 3 C.

Figure 4 shows that the surface temperatures were significantly lower than in the previous natural state at different discharge rates. The maximum temperatures were 337.1 K at a 3 C discharge rate and 315 K at 1.5 C, which are 9.1 and 7.7 K lower than those in the natural state shown in Figure 3, respectively. However, the difference in surface temperature was substantially larger than that under the natural state because the width of the HP was narrower than that of the battery, and the surface heat of the battery near the HP was more quickly dissipated. The effectiveness of the HP was limited. At a

3 C discharge rates, the surface temperature still exceeded 323 K, indicating that additional measures were needed to further lower the temperature.

3.2.3. Simulations for Combined HP and TEC

As shown in Figure 5, the TEC was mainly composed of a series of N- and P-type thermo-elements, which were connected to each other using a metal conductor. When current passes through the unit, the energy transfers through the circuit and forms the cold and hot sides of the TEC [16].

**Figure 5.** The thermoelectric cooling (TEC) structure.

When a TEC is used for refrigeration, the cold side is attached to the power source and the hot side is attached to the dissipating device. Therefore, TEC was regarded as a cuboid for simulation. The side of the cuboid with negative power simulated the cold side of the TEC to absorb heat, and the other side, with positive power, simulated the hot side to dissipate heat [29]. The power was calculated by the optimal working current, the corresponding voltage, and the coefficient of performance (COP). The optimal working current of the TEC selected was 3.96 A, and the corresponding voltage was 7.8 V. The power was 17.92 W = 3.96 A × 7.8 V × 0.58. The other part of the cuboid was regarded as particular material, whose thermal conductivity was 1.5 Wm<sup>−</sup>2K−1. With the added TEC, the model was meshed into 296,569 elements. The simulation results are shown in Figure 6.

**Figure 6.** Surface temperatures at different discharge rates based on HP and TEC: (**a**) 2 C; (**b**) 2.5 C; (**c**) 3 C.

Figure 6 shows that the surface temperatures of the battery were significantly lower than that with only the HP at different discharge rates. The maximum temperature was 325.4 K at a 3 C discharge rate and 316.4 K at a 2 C discharge rate, which are 11.7 and 9 K lower, respectively, than that in the HP alone case shown in Figure 4. This combined BTMS can reduce the surface temperature of the battery pack and basically keep the battery working within a reasonable temperature range.

### **4. Experiment**

The three discussed BTMS models analyzed by the simulations were tested to verify the practicality of the design and the simulation results.

### *4.1. Experimental Design*

To avoid being influenced by the environmental temperature during the experiment, the battery pack was placed in a thermotank. The thermotank was set to 298 K to simulate the ambient environmental temperature; the battery pack was charged and discharged by a charging tester. The first constant current and then constant voltage (CC-CV) method was adopted to charge the battery. High accuracy and low-cost Omega T-type thermocouples were used as temperature sensors and were connected to a data acquisition instrument to acquire the temperature of the battery over time. A micro temperature control module was used to operate the TEC according to the surface temperature of the battery pack. Four Li-ion batteries were connected in series into a battery pack, and its working voltage was 11.2–16.8 V. Four temperature measuring points were selected on the front and back sides of each battery, and two more measuring points were arranged on opposite sides of each battery. Data obtained from the six measuring points were then averaged to represent the temperature characteristics of this battery. The measuring points of the battery are shown in Figure 7.

**Figure 7.** Placement of temperature sensors on a Li-ion cell.

### *4.2. BTMS Experiment and Results Analysis*

### 4.2.1. Surface Temperature Rise of Battery in Natural State

Four batteries were connected in series with a distance of 2 Cm and placed in the thermotank for 1 h to stabilize the surface temperature of the battery at 298 K. Then, the experiments were conducted at different discharge rates. To ensure the same discharging capacity each time and avoid damage to the battery, a discharge capacity of 90 Ah was used. After each charging and discharging cycle, the battery was kept for 7 h in the thermotank to stabilize the surface temperature of the battery at 298 K again. The experimental results are shown in Figure 8.

As shown in Figure 8, the surface temperature of the battery continued to rise during the discharging period. The drop near the end of each curve indicates that the discharge was over and the battery had begun to cool down. The maximum surface temperature of the battery was 313.2 K at a 1 C discharge rate, and the battery worked within the optimal temperature range throughout the process. As the discharge rate increased, the discharging time shortened and the surface temperature of the battery increased. The maximum surface temperature of the battery reached 343 K at a 3 C discharge rate. A battery operating at this high temperature for a prolonged duration will undergo thermal runaway, which reduces battery capacity and shortens its service life. Therefore, this is not normally permitted and is the key reason to apply a BTMS to battery packs.

**Figure 8.** Surface temperature rise at different discharge rates in the natural state.

### 4.2.2. Surface Temperature Rise of Battery with BTMS and HP

Figure 9 shows the structure of the BTMS with the HP. One side of the aluminum sheet was in contact with the surface of the battery; the other side was in contact with the evaporation section of HP, which was fixed by a grooved aluminum sheet. The HP condensation section was connected with one side of the heat sink, and the other side of the heat sink with the fan.

**Figure 9.** The BTMS based on HP.

The environmental temperature was set to 298 K, and the fan started to work at 298 K. The experimental results are shown in Figure 10. With this system, we did not study the surface temperature rising characteristics of the battery ata1C discharge rate as the temperature remained within the optimal range without the HP.

As seen in Figure 10, the temperature for discharging at 1.5 C remained within the optimal range, and the maximum temperature was 319.7 K ata2C discharge rate and 328.2 K at a 2.5 C discharge rate, above the upper limit, but 5.6 and 4.6 K lower, respectively, than that in the natural state shown in Figure 8. Thus, it was necessary to take additional measures to further reduce temperature.

**Figure 10.** Surface temperature rise at different discharge rates in a BTMS based on an HP.

4.2.3. Surface Temperature Rise of Battery in BTMS Using Both HP and TEC

Figure 11 shows the structure with a TEC installed between the aluminum sheet and the HP. The cold side of the TEC was connected to the aluminum sheet, which can accelerate heat loss. The hot side was connected to the HP, and the heat generated by the hot side of TEC was transmitted to the heat sink through the HP, and then the fan accelerated the loss of heat.

**Figure 11.** The BTMS based on an HP and TEC.

The two kinds of limits on the working conditions of the TEC are the maximum cooling capacity and the maximum COP. The maximum cooling capacity is when the TEC works at maximum current, but then it consumes a larger amount of power and the obtained COP is also smaller. Although the working condition at maximum COP is better economically as it consumes less power, the cooling capacity is very low. Therefore, the optimal working current must be selected to increase the cooling capacity of the TEC and consume less power. The TEC used in this study was TCE1-12705, which has an optimal working current of 3.96 A according to Wang [19]. The environmental temperature was set to 298 K, and the TEC and the fan begin to work at 298 K. The experimental results are shown in Figure 12.

Figure 12 shows that at 2 C and 2.5 C discharge rates, the maximum surface temperature of the battery reached 313.5 K and 316.6 K, respectively, within the optimal working temperature range. At 3 C, the surface temperature of the battery reached 318 K in about 800 s, and 322 K at the end of discharging, which is 11.8 K lower than that of the BTMS based on an HP alone, shown in Figure 10. According to the experimental results analyzed above, this combined BTMS can effectively control the surface temperature of the battery pack.

**Figure 12.** Surface temperature rise at different discharge rates in BTMS based on HP and TEC.

#### *4.3. Optimal Working Temperature of TEC*

When the BTMS based on HP and TEC was applied to the battery, the surface temperature of the battery was lower than 318 K at the 2 C and 2.5 C discharge rates. To reduce energy consumption, the optimal working temperature of the TEC was selected, and the surface temperature of the battery was still under 318 K. The environmental temperature was 298 K, and the upper limit of the optimal working temperature of the battery was 318 K. Therefore, the TEC worked when the surface temperature of battery reached 298 K, 303 K, 308 K, 313 K, and 318 K. At the 3 C discharge rate, the surface temperature of the battery exceeded 318 K, and we did not study the surface temperature rise characteristics of the battery at the 3 C discharge rate.

As the environmental temperature was set to 298 K, the fan was turned on at 298 K, and the TEC began to work at 298 K, 303 K, 308 K, 313 K, and 318 K. The experimental results of the discharge rates of 2 C and 2.5 C are shown in Figure 13.

**Figure 13.** Surface temperature rise with the TEC beginning to work at different temperatures: (**a**) 2 C; (**b**) 2.5 C.

As shown in Figure 13, with increasing TEC working temperature, the surface temperature of the battery also rose. When the TEC began to work at 308 K, 313 K, and 318 K, the temperature curve tended to decrease at the start of TEC working. At that moment, the cold side of the TEC accelerated the loss of surface temperature of the battery, and the heat generated from the hot side was transferred to the heat sink through the HP.

At the 2 C discharge rate, when the TEC began to work at 318 K, the surface temperature of the battery remained near 318 K until discharge was over. At a 2.5 C discharge rate, when the TEC began to work at 303 K, the surface temperature of the battery remained near 318 K until discharge was over. Therefore, the optimal working temperature for TEC to begin working is 318 K at a 2 C discharge rate and 303 K at a 2.5 C discharge rate.

### **5. Discussion and Conclusions**

### *5.1. Discussion*

### 5.1.1. Comparison of Simulation and Experimental Results

Based on the different BTM models, simulations and experiments were conducted to research temperature increase characteristics of a battery at different discharge rates. The maximum battery surface temperatures at different discharge rates are shown in Figure 14.

**Figure 14.** Comparison of maximum temperatures between simulation and experimental results.

As seen in Figure 14, with increasing discharge rate, the battery surface temperature increased. For the same model under the same discharge rate, the surface temperature was the highest in the natural state and the lowest with the combined HP and TEC, and the maximum surface temperature measured in the experiment was lower than in the simulation, which is reasonable. The maximum difference in surface temperature between the experiment and simulation was within 5 K, showing a very good agreement between the two approaches, verifying the simulated results.

### 5.1.2. Comparison of Battery Surface Temperature Increase Rate

The rising rates of battery temperature directly reflect the change in temperature. The temperature rising rates of the battery are shown in Figure 15.

Under the same discharge rate, the rising rate of the battery surface temperature was the highest in the natural state and the lowest with the BTMS based on combined HP and TEC. The progressive improvement in reducing battery surface temperature demonstrated the effectiveness of the designed model.

**Figure 15.** Battery surface temperature rising rates.

### *5.2. Conclusions*

As the discharge rate increases, the surface temperature of a battery rises. We designed a battery thermal management system that combined HP and TEC and investigated its performance. Simulations and experiments were conducted to study the surface temperature rising characteristics of the battery with the system at different discharge rates. Based on TAFEL-LAE895 100 Ah ternary lithium ion batteries, the conclusions were drawn as follows.


**Author Contributions:** Conceptualization, C.Z.; methodology, H.G.; software, S.Z. and K.L.; formal analysis, Z.X. and S.C.; data curation, B.W. and Z.X.; writing—original draft preparation, Z.X.; writing—review and editing, H.G. and Z.X. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the Natural Science Foundation of China, grant number 51974229" and "The APC was funded by Xi'an University of Science and Technology.

**Acknowledgments:** This study was supported by a research grant from the Natural Science Foundation of China (Grant No. 51974229).

**Conflicts of Interest:** Declare no conflicts of interest.Authors must identify and declare any personal circumstances or interest that may be perceived as inappropriately influencing the representation or interpretation of reported research results. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

### **References**


© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

### *Article* **Heat Transport Capacity of an Axial-Rotating Single-Loop Oscillating Heat Pipe for Abrasive-Milling Tools**

### **Ning Qian 1, Yucan Fu 1,\*, Marco Marengo 2, Jiuhua Xu 1, Jiajia Chen <sup>3</sup> and Fan Jiang <sup>1</sup>**


Received: 25 March 2020; Accepted: 7 April 2020; Published: 30 April 2020

**Abstract:** In order to enhance heat transfer in the abrasive-milling processes to reduce thermal damage, the concept of employing oscillating heat pipes (OHPs) in an abrasive-milling tool is proposed. A single-loop OHP (SLOHP) is positioned on the plane parallel to the rotational axis of the tool. In this case, centrifugal accelerations do not segregate the fluid between the evaporator and condenser. The experimental investigation is conducted to study the effects of centrifugal acceleration (0–738 m/s2), heat flux (9100–31,850 W/m2) and working fluids (methanol, acetone and water) on the thermal performance. Results show that the centrifugal acceleration has a positive influence on the thermal performance of the axial-rotating SLOHP when filled with acetone or methanol. As for water, with the increase of centrifugal acceleration, the heat transfer performance first increases and then decreases. The thermal performance enhances for higher heat flux rises for all the fluids. The flow inside the axial-rotating SLOHP is analyzed by a slow-motion visualization supported by the theoretical analysis. Based on the theoretical analysis, the rotation will increase the resistance for the vapor to penetrate through the liquid slugs to form an annular flow, which is verified by the visualization.

**Keywords:** oscillating heat pipes; heat transfer; milling cooling; abrasive-milling processes

### **1. Introduction**

During abrasive-milling processes, a large amount of heat is generated in the contact zone between the abrasive-milling tool and workpiece. Consequently, when abrasive-milling hard-to-machine materials with low heat conduction (e.g., nickel-based superalloy, carbon fiber reinforced composites or ceramic matrix composites), the heat gathers in contact zone and leads to critical hot spots and serious thermal damage for both the workpiece and tool [1–4]. To solve this problem, the heat should be transferred out through the tool in time to reduce the stored heat, lowering the temperature in the contact zone, with the result decreasing the risk of thermal damage.

Heat pipes are passive heat transfer devices with excellent heat transport capacity. They are widely used in the microelectronics, manufacturing, and automobile and aerospace industries for thermal management. Lately, heat pipes have been applied to the machining process to enhance heat transfer [5–9]. The heat is absorbed by the working fluid in the evaporator of the heat pipe, which, in the case of the tool, is near the contact zone. The phase change occurs, and vapor moves to the condenser to release heat by condensing into liquid. After this, the liquid returns to the evaporator either by capillary wick or unbalanced pressure distribution. In the case of applying heat pipes in the cutting tools, such fluid motion continues to transfer out the heat, with high heat transport capacity during the machining processes. The combination of heat pipes and cutting tools was conducted by several researchers [5–7], and the heat pipe cutting tool shows a great advantage in enhancing the heat and cooling the tool in the processes. He et al. [8,9] and Chen et al. [10,11] designed and produced axial-rotating and radial-rotating heat pipe grinding wheels. A titanium alloy and nickel-based superalloy were successfully dry-grinded by heat pipe grinding wheels without grinding burnout. An in-depth numerical analysis was conducted, where results show that the wickless heat pipe can assist heat transfer under a low rotational speed. Nucleate boiling was found to be totally suppressed while increasing the centrifugal acceleration. Laminar convection heat transfer is the heat transfer mode in the evaporator when the centrifugal acceleration is more than 1000 times the gravity acceleration.

The oscillating heat pipe (OHP, also called the pulsating heat pipe) is a novel heat transfer device, which has been developed since the 1990s. It is made of a meandering capillary tube which is partially filled with working fluid [12]. In recent years, experimental and theoretical research has been carried out on the thermal performance of OHPs and their heat transfer mechanism. The heat transport capacity of OHPs is affected by geometric, physical and operational parameters, such as inner diameters, working fluid properties, filling ratio, etc. It was found that within the critical inner diameter, an OHP with a larger inner diameter has a better heat transfer capability [13]. The influence of the working fluid on the heat transport capacity of the OHP has been described by several researchers [14,15]. The thermal performance of an OHP is also determined by its flow regimes, which in turn are affected by the working fluids. In general, fluids with lower dynamic viscosity and smaller surface tension lead to better performance, due to their lower friction and the higher wettability, which lead to a higher response to the pressure variation and the presence of a liquid film on the wall (decreasing the probability of dry spot occurrence). Fluids with a lower latent heat (such as acetone) bring powerful oscillations into the OHP, while they have a higher possibility of dry-out. OHPs charged with fluids with higher latent heat (e.g., water) may have better heat transfer performance in terms of maximum heat flux, but the start-up power is higher. The influence of operational parameters on the heat transfer was investigated by some researchers, such as Stevens et al. [16]. Operational parameters, which include a filling ratio (FR), inclination angle and heat load, etc., have significant effects on the heat transfer. In most of the experiments, there exists an optimized filling ratio which is within 40% to 60%, and the OHP will have a better heat transfer ability under the inclination of 50◦ to 90◦ for the proposed four turns closed-loop OHP [13,14,17,18]. Below the maximum heat flux, increasing the heat load, the flow regime develops from a bubbly flow to a slug flow to a transient flow, and then to an annular flow. As a result, the heat transfer performance varies depending on the flow patterns [19,20].

Due to the high heat transfer capability and the simple wickless structure, OHPs have the potential for application in the machining process for the enhancement of heat transfer. The OHP cutting tool was first proposed by Wu et al. [21,22]. A Ti-6Al-4V alloy was dry-cut by an OHP cutting tool. Compared with conventional cutting tools, the cutting tool with the OHP reduced the process temperature, and the tool wear-rate became slower by 20%. As for the grinding process, Qian et al. [19] introduced the concept of the OHP grinding wheel. The heat transfer analysis was studied to the point of the proof of concept. The heat transfer prediction model was built, and the start-up performance was also analyzed [23,24].

In order to dissipate the heat generated in the contact zone during the abrasive-milling of hard-to-machine materials, the concept of combining the OHP with the abrasive-milling tool is proposed, as illustrated in Figure 1. Several single-loop OHPs are formed by the independent channels machined in the tool matrix. The generated heat in the contact zone is transferred into the evaporator through the abrasive layers, and is then brought to the condenser to transfer out.

In this case, a rotation is involved. Experimental investigations of the heat transfer and flow patterns of "flower-shaped" radial-rotating OHPs were conducted by Aboutalebi et al. [25], Ebrahimi et al. [26], Liou et al. [27] and On-ai et al. [28]. In their studies, the maximum rotational speed reached 800 revolutions per minute (rpm), and the maximum centrifugal acceleration reached around 20 times gravity acceleration. For such a range of rotational speed and centrifugal acceleration, they all found

that employing rotational speed was a way to enhance the internal flow, and the centrifugal acceleration has a positive influence on the increase of thermal performance. These oscillating heat pipes rotate around the central axis of the tool with the condenser at the centre and the evaporator at the external side. In this case, centrifugal accelerations segregate the fluid between the evaporator and condenser, pushing the fluid to the evaporator.

**Figure 1.** Illustration of an OHP abrasive-milling tool.

In our case, the OHP, whose central axis is aligned to the rotational axis of the tool, is positioned on the plane parallel to the rotational axis of the tool. In this case, centrifugal accelerations do not segregate the fluid between the evaporator and condenser. To the best of our knowledge, there is no previous paper on the influence of rotation on such positioning of OHPs, especially for a single-loop OHP.

In this paper, experiments are conducted to investigate the effects of centrifugal acceleration, heat flux and working fluids on the thermal performance of OHPs under axial-rotating conditions. The flow of liquid slugs and vapor plugs inside the axial-rotating SLOHP is also analyzed by the theoretical method. Moreover, the flow is also observed through a slow-motion camera.

### **2. Experimental Preparation and Data Processing**

### *2.1. Description of Experimental Apparatus*

The oscillating heat pipe (OHP) abrasive-milling tool (OHP tool) comprises several axial-rotating single-loop oscillating heat pipes (SLOHPs). For simplification, in this paper, one axial-rotating SLOHP is designed and made for the experiment, as shown in Figure 2. The axial-rotating OHP includes a holder, a copper SLOHP, a slip ring, thermocouples and Ni-Cr heating wires. The SLOHP is mounted on the holder 30 mm away from its rotational axis, which is the same as the radius of the OHP tool. The outer and inner diameters of a SLOHP are 5 mm and 3 mm, respectively. The evaporator, adiabatic section and condenser are 30 mm, 40 mm and 30 mm long, respectively. The SLOHP is evacuated to the pressure of 10−<sup>2</sup> Pa, and then filled with acetone (CH3COCH3), methanol (CH3OH) or deionized water (DI water) as working fluids, with a filling ratio of 55%, through the vacuum and injection tube. The axial-rotating SLOHP is bottom-heated by a Ni-Cr wire, which is connected to the KIKUSUI PZB40-10 power supply (limits of error: ±0.1 W) through the slip ring and copper brush. The axial-rotating SLOHP is wrapped with an aerogel insulation to avoid heat leak. The condenser is cooled by the 0.4 MPa cold air jet at a temperature under 5 ◦C. Temperatures of the evaporator and condenser are measured by Omega type-K thermocouples (limits of error: 0.4%). Signals are transferred through the slip ring and copper brush to the NI-USB 6366 card under the sampling rate of 10 Hz, and signals are processed by NI LabView and NI DIAdem software (National Instruments,

Shanghai, China). The environment temperature is maintained at 25±1◦ C. The setup of the experiment is illustrated in Figure 2.

**Figure 2.** Experimental preparation of axial-rotating oscillating heat pipe: (**a**) illustration of axial-rotating OHP; and (**b**) experimental setup (aerogel insulation was removed for better illustration).

In the abrasive-milling process, the rotational speed of the tool varies from several hundred to a thousand rpm, and the effective heat flux generated in the contact zone is around 10<sup>5</sup> W/m2. In this case, in order to study the thermal performance of the axial-rotating SLOHP under the same conditions with the abrasive-milling process, the rotational speeds in the experiment are 0, 300, 950, and 1500 rpm, with a relative centrifugal acceleration of 0, 30, 296, and 738 m/s2. The heat flux used in the experiment varies from 9100 W/m2 to 31,850 W/m2, with an increment of 4550 W/m2. The detailed conditions are shown in Table 1.



The thermophysical properties of working fluids (i.e., methanol, acetone and DI water) are shown in Table 2.



### *2.2. Data Processing*

Thermal resistance is used as the measurement of the thermal performance of the axial-rotating SLOHP. The thermal resistant is defined as follows

$$R = \frac{Q}{T\_{\text{evap}} - T\_{\text{cond}}},\tag{1}$$

where *Q* is the heat load, and *T*evap and *T*cond are the time-averaged temperatures in evaporator and condenser. Lower thermal resistance means a better thermal performance of the axial-rotating SLOHP. The uncertainty of the thermal resistance is evaluated based on the method in [19]. The uncertainty is less than 5%.

### **3. Flow Visualization and Modeling**

### *3.1. Flow Inside Axial-Rotating Single-Loop Oscillating Heat Pipe (SLOHP) Under Rotation*

A glass axial-rotating SLOHP is made for visualization. The flow of working fluid inside the axial-rotating SLOHP is captured by a portable camera at a 240-frames-per-second speed (eight times slower). Due to the limitation of the apparatus and in the light of safety, the rotational speed is set up to 300 rpm in this case.

Under the static state and heat flux of 22,750 W/m2, the flow pattern of the working fluid, acetone, is a veering circulation with an annular flow and slug flow. When the vapor bubbles and plugs generate, expand and move up to the condenser, the vapor easily penetrates the liquid slugs to form the annular flow (see Figure 3a(i–iii)), while the flow in the neighboring tube remains a train of liquid slugs and vapor plugs (see Figure 3a(ii,iii)). Besides, the direction of circulation changes from time to time, as shown in Figure 3a(ii,iii).

When the axial-rotation is applied, as mentioned in the theoretical analysis in Section 3.2, the resistance for the vapor to penetrate the liquid slug becomes larger, and as a result it is more difficult for the vapor to penetrate through the liquid slug to form the annular flow. The flows inside the axial-rotating SLOHP at the speeds of 60, 150 and 300 rpm are shown in Figure 3b–d. At the rotating speed of 60 rpm, the flow is a veering slug flow, which generates a train of liquid slugs and vapor plugs (see Figure 3b(ii)). Sometimes, when a train of liquid slugs and vapor plugs moves from evaporator to condenser, it will slow down and stop, and a new train of liquid slugs and vapor plugs is generated in the neighboring tube and moves upwards, causing the flow direction to change (see Figure 3b(ii,iii)). The whole process is illustrated in Figure 3b. At the rotating speeds of 150 and 300 rpm, the unidirectional circulation forms. Since the SLOHP rotates clockwise, the inertia causes the working fluid to flow in a one-way, clockwise direction inside the SLOHP. In this case, liquid slugs generate and move up to the condenser in the left tube, while in the right tube, liquid slugs and vapor plugs oscillate at the adiabatic section (see Figure 3c(i,iii) and Figure 3d(ii,iv)). Though flows at rotational speeds of 150 rpm and 300 rpm are similar, new trains of liquid slugs and vapor plugs are generated at the speed of 150 rpm (Figure 3c(i,iii)), while just one liquid slug forms and moves up at a time at the speed of 300 rpm (Figure 3d(ii,iv)).

As mentioned before, under the static state, vapor easily penetrates through the liquid slugs to form an annular flow. Under axial-rotation, especially when the speed exceeds 150 rpm, unidirectional circulation forms. It takes a relatively long time to form a long train of liquid slugs and vapor plugs at the speed of 150 rpm. Nevertheless, one short liquid slug is generated at a shorter time at the speed of 300 rpm, as shown in Figure 4.

**Figure 3.** *Cont*.

**Figure 3.** Flow inside the axial-rotating SLOHP filled with acetone under 22,750 W/m2: (**a**) at static state; (**b**) at 60 rpm rotating speed; (**c**) at 150 rpm rotating speed; and (**d**) at 300 rpm rotating speed.

**Figure 4.** Comparison of liquid slug lengths under different rotating speeds (heat flux of 22,750 W/m2).

### *3.2. Modelling of Flow Under Rotation*

For a well-working oscillating heat pipe (OHP), a train of liquid slugs and vapor plugs must form and exist in the system. At the moment that the evaporator is heated, the liquid turns into vapor, causing the expansion of the vapor volume. When the vapor velocity exceeds some critical value, the vapor plug penetrates all liquid plugs and generates an annular flow. At this moment, the mass-spring system, consisting of a train of liquid slugs and vapor plugs, vanishes. The oscillating motion stops, and the OHP reaches the highest heat transport capacity [29]. The illustration of the flow inside the OHP is shown in Figure 5.

When the vapor penetrates the liquid slug, the momentum that is generated by the vapor plugs will be used to overcome the total forces applying on the liquid plug, i.e.,

$$\int\_{0}^{\tau\_{0}} 2\pi \rho\_{v} u\_{v}^{2} r dr = (-p\_{1} + p\_{2})\pi r\_{0}^{2} + 2\pi r\_{0} \sigma (\cos a\_{r} + \cos a\_{\mathfrak{t}}) + 2\pi r\_{0} \left(\int\_{0}^{\mathcal{L}\_{\parallel}} \tau\_{\mathfrak{u}} d\mathbf{x}\right) \tag{2}$$

where ρ<sup>v</sup> is the vapor density, *p*<sup>1</sup> and *p*<sup>2</sup> are the vapor pressure, *r*<sup>0</sup> is the radius of the tube, σ is the surface tension, α<sup>a</sup> and α<sup>r</sup> are the advancing and receding contact angles and τ<sup>w</sup> is the shear stress. The left term is the momentum produced by the vapor due to the vapor volume expansion. The right of the equation is marked as the total resistant force. The first term of the right part is the pressure differences acting on the liquid slug and vapor plug, while the second term of the right part is due to the surface tensions and the last term is due to the shear stress between the wall and fluid.

**Figure 5.** Illustration of liquid slug moving in the OHP: (**a**) right before the vapor penetrating through the liquid slug; and (**b**) right after the penetration of vapor through the liquid slug.

The pressure of vapor and liquid are

$$p\_1 = p\_{0,v} + \frac{\pi}{2} \rho\_v r\_0 \omega^2 R\_\prime \text{ } \omega \ge 0\text{,}\tag{3}$$

$$p\_2 = p\_{0,l} + \frac{\pi}{2} \rho\_l r\_0 \omega^2 \mathbb{R}\_\prime \text{ or } \ge 0,\tag{4}$$

where *p*0,v and *p*0,l are the pressure of the vapor and liquid phases under a static state, *R* is the distance between the OHP and the rotational axis and ω is the angular velocity.

The shear stress on the wall can be expressed as:

$$
\pi\_w = \mu \left( -\frac{d\mu\_l}{dr} \right)\_w \tag{5}
$$

According to the solution proposed by Ma [29], Equation (1) can be solved as follows

$$\rho\_l r\_0 u\_{l, \text{max}}^2 \int\_0^1 u\_l^{\alpha 2} r^\* dr^\* = \frac{\pi}{4} \omega^2 r\_0^3 R \left(\rho\_l - \rho\_0\right) + \frac{\epsilon\_0}{2} \left(p\_{0,l} - p\_{0,p}\right) + o\left(\cos \alpha\_l + \cos \alpha\_4\right) + \left(\frac{1}{\epsilon\_0} \mu\_l \left(-\frac{4\omega^\*}{4r}\right)\_{l'=-1} \int\_0^{L\_l} u\_{l, \text{max}} dx\right) \tag{6}$$

where *u*l,max is the maximum velocity of the liquid phase and μ<sup>l</sup> is the dynamic viscosity of the liquid.

The following parameters are defined as follows:

$$A = \int\_0^1 u\_l^{\*2} r^\* dr^\*,\tag{7}$$

$$B = \left(-\frac{d\mu^\*}{dr^\*}\right)\_{r^\*=1} \tag{8}$$

$$C = \frac{r\_0}{2}(p\_{0J} - p\_{0J}) + \sigma(\cos \alpha\_V + \cos \alpha\_d) \tag{9}$$

Therefore, Equation (6) will be rewritten as:

$$
\rho\_v r\_0 u\_{l, \text{max}}^2 A = \frac{\pi}{4} \omega^2 r\_0^3 R (\rho\_l - \rho\_v) + \mathcal{C} + \frac{B}{r\_0} \mu\_l \int\_0^{L\_l} u\_{l, \text{max}} dx \tag{10}
$$

Under the static state, the driving force of vapor and the total resistant force to overcome are:

$$
\rho\_{\upsilon} r\_0 u\_{l, \text{max}}^2 A = C + \frac{B}{r\_0} \mu\_l \int\_0^{L\_l} u\_{l, \text{max}} d\mathbf{x} \tag{11}
$$

Meanwhile, under the axial-rotating condition, the driving force of vapor and the total resistant force to overcome is shown in Equation (10). It is clear that total resistant force is increased in case of an axial-rotating condition by πω2*r*<sup>3</sup> <sup>0</sup>*R*(ρ*<sup>l</sup>* − ρ*v*)/4, so the vapor plug is more difficult to penetrate the liquid slug under axial rotation, and the annular flow will turn into slug flow under axial-rotation.

### **4. Results and Discussion**

### *4.1. E*ff*ects of the Centrifugal Acceleration*

The influence of rotation is described through the rotational speed and radius. Figure 6 illustrates the effects of centrifugal acceleration on the thermal performance of the axial-rotating single-loop oscillating heat pipe (SLOHP) filled with methanol, acetone and deionized water (DI water). It is clear that the centrifugal acceleration has a positive influence on the heat transport capacity of the axial-rotating SLOHP filled with methanol and acetone. When the centrifugal acceleration increases up to 30 m/s2, the thermal resistance decreases rapidly for all the fluids. In the range of the centrifugal acceleration from 30 m/s2 to 738 m/s2, the thermal performance enhances for most cases of the axial-rotating SLOHP filled with methanol and acetone. In the range of 0 m/s<sup>2</sup> to 738 m/s2, the thermal resistance decreases by 22% and 137% for the axial-rotating SLOHP filled with methanol and acetone, respectively. According to previous studies [19,25,28,30], the inside flow patterns can be estimated through the temperature curves. In the evaporator, the longer the vapor plug is, the longer it passes through the temperature measuring point, and the slower the vapor moves, the longer it passes through the measuring point. This leads to a higher amplitude of the temperature, as the hot long vapor plug will increase the temperature of the damped measuring point gradually until a liquid slug passes the measuring point (see Section 3.1). Moreover, the temperature differences decrease when the fluid circulation is faster and the vapor plugs and liquid slugs become shorter.

In order to emphasize the influence of the centrifugal acceleration on the flow pattern, a rotational speed of 950 rpm (296 m/s2) is applied on the SLOHP filled with methanol under the heat flux of 27,300 W/m2. Figure 7a shows the evaporator temperature curve. Under static conditions (*a* = 0 m/s2), the value of the evaporator temperature is high, with the peak-to-valley value being around 20 ◦C, and the evaporator temperature fluctuating at low frequencies. On the contrary, under the centrifugal acceleration of 296 m/s2, the temperature curve is obviously smoother, with a smaller amplitude and high oscillation frequencies. The peak-to-valley amplitude is around 10 ◦C under these conditions. In terms of the methanol and acetone, the centrifugal acceleration likely changes the annular flow to slug flow, as it is qualitatively proven by the theoretical analysis in Section 3.2, and increases the circular motion velocity. As a result, for methanol and acetone as working fluids, the thermal performance of the axial-rotating SLOHP enhances with the increase of centrifugal acceleration.

**Figure 6.** Effects of centrifugal acceleration: (**a**) methanol; (**b**) acetone; and (**c**) DI water.

**Figure 7.** Evaporator temperature under rotational speed of 0 rpm and 950 rpm (filled with methanol): (**a**) temperature curve; and (**b**) and (**c**) FFT analysis of the temperature curve in areas (I) and (II).

On the contrary, as shown in Figure 6c, a different behavior is observed when DI water is the working fluid. Up to *a* = 30 m/s2, the heat transfer performance increases sharply, the thermal resistance drops by up to 74%, from 2.5 K/W to 1.43 K/W, under the heat flux of 13,650 W/m2. As the centrifugal acceleration exceeds 30 m/s2, the thermal resistance rises gradually, which means the thermal performance deteriorates. Specifically, when the centrifugal acceleration is higher than 296 m/s2, the thermal performance is lower than that of the static state at the heat flux range from 18,200 W/m2 to 31,850 W/m2. Therefore, in terms of the axial-rotating SLOHP filled with DI water, the thermal performance first enhances and then decreases along with the increase of centrifugal acceleration, which peaks at 30 m/s2. This result is consistent with previous studies [22,25]. The reasons for such different behavior of water can be found in its higher viscosity (Table 2), which increases the flow friction [28], the higher thermal conductivity, which can produce local transient dry-out [25], and the lower value of (d*p*/d*T*)sat, which brings it to a weaker expansion, i.e. a weaker flow circulation (see also Section 3.1). From the visualization, it is observed that under 300 rpm, the flow of water is a discontinuous train of vapor bubble and liquid slug. Meanwhile, when the speed exceeds 300 rpm, it is difficult to form effective circulation. An intermittent train of bubble and slug forms, and sometimes moves to the condenser. Most of the time, no obvious motion is observed—just local oscillation.

### *4.2. E*ff*ects of the Heat Flux*

When the heat flux increases from 9100 W/m2 to 31,850 W/m2, the thermal performance of the axial-rotating SLOHP enhances at all centrifugal accelerations for all working fluids (Figure 8). In detail, for the axial-rotating SLOHP filled with methanol, the thermal resistance decreases from around 2.67 K/W to around 0.64 K/W by up to 4.5 times. As for the acetone as the working fluid, the thermal resistance decreases from around 3.01 K/W to around 0.34 K/W by a maximum of 8.6 times. Similarly, in terms of DI water, the thermal resistance decreases from around 3.33 K/W to around 0.63 K/W by four times. Based on our previous research [19], as the heat flux increased, changes of flow pattern and motion modes—which develops from bubbly flow to vapor plug flow, then to annular flow, and from oscillation to circulation—enhance the OHP's thermal performance.

### *4.3. E*ff*ects of the Working Fluid*

From Table 2, acetone and methanol have a high (d*p*/d*T*)sat value, low dynamic viscosity and low surface tension, which lead to a higher driving force and the lower resistance of flow motion. On the other hand, DI water has a higher thermal conductivity, latent heat and specific heat. Due to differences in thermophysical properties, the working fluid has an important impact on thermal performance under the axial-rotating condition. The influence of the working fluids is different under low and high centrifugal accelerations, as illustrated in Figure 9. The axial-rotating SLOHP filled with acetone has the highest heat transport capacity, regardless of the heat flux and centrifugal acceleration. At the low centrifugal acceleration (e.g., 30 m/s2), the effective heat transfer coefficient of the axial-rotating SLOHP filled with DI water is higher than that of methanol in general. However, when the centrifugal acceleration is high (e.g., 738 m/s2), the thermal performance of the axial-rotating SLOHP filled with methanol is better than that of DI water. This is due to different flow patterns and motion modes, which are caused by the differences of thermophysical properties at varied centrifugal accelerations. In the light of this—for the given temperatures and heat load ranges, and even if the deeper reasons are not completely clear—as a conclusion, acetone is the recommended working fluid for this axial-rotating SLOHP, which is supposed to work in the abrasive-milling tool to enhance heat transfer in abrasive-milling processes.

**Figure 8.** Effects of heat flux: (**a**) methanol; (**b**) acetone; and (**c**) DI water.

**Figure 9.** Effects of working fluids under centrifugal acceleration of: (**a**) 30 m/s2; and (**b**) 738 m/s2.

#### **5. Conclusions**

The influences of centrifugal acceleration, heat flux and working fluid on the thermal performance of the axial-rotating SLOHP are investigated through visualization, theoretical analysis and experiments. The main conclusions are drawn as follows:

The flow inside axial-rotating SLOHP is analyzed through the theoretical and visualization method. Based on theoretical analysis, the rotation will increase the resistance for the vapor to penetrate through the liquid slugs to form an annular flow, which is verified by the visualization. Besides, under the static state and heat flux of 22,750 W/m2, the flow pattern of working fluid, acetone, is a veering circulation with the annular flow and slug flow. Meanwhile, under the axial-rotating conditions, the flow turns

into slug flow. When the rotating speed increases from 60 rpm to 300 rpm, the flow changes from veering circulation to unidirectional circulation, and the liquid slug becomes shorter with a higher generating frequency.

For the acetone and methanol as working fluids, the thermal performance is enhanced by increasing the centrifugal acceleration, which is demonstrated through the change of internal flow motions. Nevertheless, the thermal performance of the axial-rotating SLOHP filled with DI water enhances first, and then decreases when the centrifugal acceleration increases. The acetone is recommended since the axial-rotating SLOHP filled with acetone has the best thermal performance. Furthermore, the heat transport capacity of the axial-rotating SLOHP improves with the increase of heat flux for all cases.

**Author Contributions:** Conceptualization, N.Q., Y.F. and J.X.; methodology and analysis, N.Q., M.M. and F.J.; writing—original draft preparation, N.Q., J.C.; writing—review and editing, Y.F., M.M. and J.X.; supervision, Y.F., M.M. and J.X.; funding acquisition, Y.F., M.M. and J.C. All authors have read and agreed to the published version of the manuscript.

**Funding:** The authors gratefully acknowledge the financial support for this work by the National Natural Science Foundation of China (Grant Nos. 51175254 and 51905275), National Natural Science Foundation of China for Creative Research Group (Grant No. 51921003) and Natural Science Foundation of Jiangsu Province (Grant No. BK20190752). Marco Marengo has been supported by UK EPSRC through the grant EP/P013112/1 (HyHP project, https://blogs.brighton.ac.uk/hyhp/).

**Acknowledgments:** The authors would like to address a great thank to Jiusheng Xu for his assistant to build the apparatus.

**Conflicts of Interest:** The authors declare no conflict of interest.

### **Nomenclature**


*Greek*


*Subscripts*


### **References**


© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

### *Article* **Experimental Validation of Heat Transport Modelling in Large Solar Thermal Plants**

### **Kevin Sartor \* and Rémi Dickes**

Thermodynamics Laboratory, University of Liège, 4000 Liège, Belgium; rdickes@ulg.ac.be **\*** Correspondence: kevin.sartor@uliege.be

Received: 16 April 2020; Accepted: 6 May 2020; Published: 8 May 2020

**Abstract:** Solar thermal plants are often considered as a convenient and environmentally friendly way to supply heat to buildings or low temperature industrial processes. Some modelling techniques are required to assess the dynamic behaviour of solar thermal plants in order to control them correctly. This aspect is reinforced while large plants are considered. Indeed, some atmospheric conditions, such as local clouds, could have significant influence on the outlet temperature of the solar field. A common modelling approach to assess the heat transport in pipes is the one-dimensional finite volume method. However, previous work shows limitations in the assessment of the temperatures and in the computational time required for simulating large pipe networks. In this contribution, a previous alternative method developed and validated in a district heating network is used and extended to a solar thermal plant considering the thermal solar gain and the inertia of the pipes. The present contribution intends to experimentally validate this model on an existing solar plant facility available at the Plataforma Solar de Almeria in Spain.

**Keywords:** solar network; dynamic modelling; plug flow; control

### **1. Introduction**

Solar thermal plants are one of the current solutions to change the world energy sources and to contribute to a green energy transition by recovering solar energy. The use of this kind of energy system has grown significantly for several years, especially to produce power or to supply heat to industrial processes or, in some cases, to residential buildings [1]. Moreover, among these systems, the concentration thermal plants tend to be economically competitive compared to conventional energy sources (60–100 EUR/MWh) [2], increasing their future use.

However, some modelling techniques are required to assess the dynamic behaviour of solar thermal plants to control them correctly. Indeed, like most renewable energy systems, a solar plant is facing to the intermittency of the energy source, namely the Sun, and requires some modelling techniques to couple them to other energy systems or energy storage facilities [3]. This aspect is reinforced while large plants are considered such as the solar park in Cape Town in South Africa, where the area of the solar plant can reach several km<sup>2</sup> [4]. In this plant configuration, some atmospheric conditions, such as local clouds, could have significant influence on the outlet temperature of the solar field and a control strategy should be used to maintain the outlet plant temperature near the defined set point as a predictive control based on climatic observations [5] or with dedicated algorithms [6]. Another solution is to instrument the network at numerous key locations to measure the operating conditions (temperatures and mass flow rates). However, this method is often expensive because of the numerous expensive sensors used.

A common modelling approach to assess the heat transport in pipes is the one-dimensional finite volume method. However, previous works [7–9] show limitations in the assessment of the temperatures and in the computational time required for simulating large pipe networks. Indeed, the finite volume method requires a high discretisation scheme to get accurate results. However, the increasing discretisation leads to a higher computational time. This computational time could be incompatible with the control strategy.

The present contribution intends to use an existing model developed by the authors dedicated to district heating networks and to extend it to assess the dynamic behaviour of solar thermal networks. Indeed, the proposed model, based on a Lagrangian approach [10], was validated considering the working fluid as water and a reduced heat transfer between the pipes and the ambient due to the large insulation of the pipes. In this contribution, the working fluid considered is oil and the heat transfer is higher due to the thermal power recovered from the Sun. Therefore, the model is extended into solar thermal plants, considering the thermal solar gain and the inertia of the piping system. The case study used to validate the dynamic model is a parabolic trough collectors (PTC)-based solar plant facility installed at the Plataforma Solar de Almeria in Spain.

The developed model allows for a further study to investigate control strategy definitions and their optimisation. Although this paper focuses on the experimental validation of a medium-size solar thermal plant, it could be extended to the modelling of larger solar plants, while the key element behaviour, namely a pipe, is studied and validated.

### **2. Problem Statement**

A typical solar thermal plant is composed of a piping system dedicated to transport the energy recovered from the Sun to a defined process. In large solar thermal plants, the total pipe network length can reach over several kilometres. Therefore, the heat propagation through the network depends on the fluid velocity, which can induce some delay. These delays could reach some minutes or even hours, depending of the operating conditions and the pipe length. Indeed, the order of magnitude of the working fluid velocity is generally some meters per second to limit the pressures losses and the related electric pump consumption and it depends on the control strategies of the plant.

Previous works of the authors [7,11] show that a one-dimensional finite-volume method is able to model the dynamic behaviour of this kind of network. However, this modelling method requires an important spatial discretisation of the pipe to get accurate results. This spatial discretisation leads to high computational time requirements, which is not compatible with the final purpose of the modelling, i.e., the investigation and the optimisation of control strategies. For lower spatial discretisation layouts, a phenomenon named "numerical diffusion" occurs [7,8,12] and reduces the accuracy of the finite volume method by anticipating the heat waves at the pipe outlet.

Therefore, it is proposed to counter these issues by considering an alternative modelling method called the "plug flow" method, which is briefly detailed in the following section. The plug flow method accuracy is of the same order of magnitude as the one of the finite-volume method, with important spatial discretisation. Moreover, the simulation speed of the plug flow method is improved and is clearly one advantage for the definition and the optimisation of control strategies, such as model-based predictive control methods. Indeed, in large thermal networks, the weather conditions, and more specifically, the solar irradiance, can have an important influence on the performance of some parts of the plant, especially if a given outlet temperature is required for the supplying process connected to the solar thermal plant.

The plug flow modelling method has already been validated experimentally on a district heating network. However, such a district heating network involves a very low heat transfer between the pipe and the ambient due to the high insulation used on the piping system. In order to go a step further, this contribution proposes to experimentally validate the "plug flow" method on a solar thermal network facility characterised by a very high transfer. More specifically, the plug flow model is used to simulate a field of PTC, as described in Section 4.

### **3. Modelling**

Regarding the flow inside the pipe, the proposed modelling method is derived on the standard TRNSYS Type 31 component. This modelling method is based on a Lagrangian approach, i.e., the properties of each fluid particle are considered along with their direction in function of time, considering the energy balance [10]. In order to simplify the resolution of the whole system, the flow is considered as incompressible, which is valid if the fluid is a liquid and for low pressure variations [13]. Moreover, the wall friction dissipation and the energy coming from the pressure drop into the pipe are neglected as the axial diffusion. This method was previously analysed, developed and validated under the Matlab platform for different operating conditions and pipe layouts [7,14], but has now been successfully ported to Modelica language [11] and the related open-source Modelica library is available in [15].

The current model is based on a simplified combination of energy and continuity equation (Equation (1)) [11]:

$$\frac{\partial \left(\rho c\_p A T\right)}{\partial t} + \frac{\partial \left(\rho c\_p A v T\right)}{\partial x} = -\dot{q}\_c \tag{1}$$

where ρ denotes the density, *cp* is the fluid specific heat, *A* is the pipe cross section (m2), *v* is the flow velocity, *<sup>x</sup>* is the spatial coordinate, *<sup>t</sup>* is the time, *<sup>T</sup>* is the temperature and . *qe* is the heat loss per unit length (which is positive for heat loss from pipe to the ground).

Considering in the Lagrangian approach that the observer is attached to a fluid parcel, there is no notion of spatial coordinate. Finally, after rearranging the previous equation, it is possible determine the outlet pipe temperature (*Toutlet*) through the Equation (2):

$$T\_{outlet} = T\_{ground} + \left(T\_{inlet} - T\_{ground}\right) \cdot \exp\left(-\frac{\Delta T}{RC}\right) \tag{2}$$

Δ*T* is the delay time of the fluid parcel travelling the pipe; *R* is the thermal resistance between the fluid temperature and the ground temperature, which can be calculated by [16] or [17] depending on the geometric characteristics of the pipe, and *C* is the heat capacity per unit length of the water into the pipe. In this case, the fluid temperature is assumed uniform throughout the cross section of the pipe.

In the Modelica modelling, a unique volume is considered for the pipe. To assess the time delay of a fluid parcel between the inlet and the outlet of a pipe, the following dynamic equation (Equation (3)) is solved:

$$\frac{\partial z\left(y,t\right)}{\partial t} + v(t)\frac{\partial z\left(y,t\right)}{\partial y} = 0,\tag{3}$$

where *z* (*y*, *t*) is the transported quantity, *y* is the normalised spatial coordinate, *t* is the time and *v*(*t*) the normalised velocity. An approximation of the one-way wave equation was successfully introduced with the spatialDistribution() operator defined in the Modelica Language Specification [18]. This modelling method considers the heat propagation, the thermal inertia of the piping system, the pressure losses and the heat transfer with the ambient.

The hydraulic behaviour is assessed by a previously developed model denoted HydraulicDiameter of the Annex 60 Library [15]. In this case, the pressure drop (Δ*P*) is coupled to the mass flow rate (m) ˙ using a quadratic relation as in (4): .

$$
\dot{m} = k \sqrt{\Delta P} \,\tag{4}
$$

where k is a constant depending on the piping system in function of the nominal conditions (nominal pressure losses for a nominal mass flow rate). A thermal capacity is added to the core pipe model at the outlet of the pipe to account for the thermal inertia of the pipe constituting material as [19]. The interested reader could refer to [11] for further details.

Regarding the heat transfer to the pipe, the radial heat balance between all the heat collection element components (see Figure 1) is modelled according to the deterministic model proposed by Forristal [20], which accounts for the most important phenomena, i.e.,:

*Energies* **2020**, *13*, 2343


**Figure 1.** Energy balance around the heat collection element. In blue, the glass envelope; in grey, the metal pipe; and in white, the vacuum between the two. Heat transfer is highlighted with red arrows. Based on [21].

It is assumed that temperatures, heat transfer coefficients and thermodynamic properties are considered uniform around the circumference of the heat collector. Moreover, thermal losses through the support brackets are neglected and solar absorption in the tube and the glass envelope is treated as a linear phenomenon [21]. For further information about the original model, the reader can refer to [21–23].

For larger network where several pipes are interconnected, energy, mass and momentum balances could be performed on the key location to assess the energy and mass repartition inside the whole system according to [18].

### **4. Experimental Apparatus**

The reference data used for the model validation are experimental measurements gathered on the Parabolic Trough Test Loop (PTTL) facility installed at the Plataforma Solar de Almería in Spain [24]. The data are recovered from a previous work presenting a finite-volume dynamic model of parabolic trough collectors [21]. As depicted in Figure 2, the facility includes three parallel lines of PTC from different manufacturers (AlbiasaTrough, EuroTrough and UrssaTrough) but only the EuroTrough collectors (ETC) are used in this work. The ETC line is composed of 6 EuroTrough modules connected in series and 18 prototype receiver tubes from a Chinese manufacturer, for a total length of 70.8 m and a net aperture area of 409.9 m2. The system is a closed loop, with an East-West orientation and it is charged with the thermal oil Syltherm 800 [25].

The system operation is quite straightforward. A centrifugal pump drives the fluid from the point (1) through one of the three parallel PTC lines of the solar field. In the collectors (i.e., from points (2) to (3)), the heat transfer fluid is heated by absorbing solar energy reflected by the collectors onto the receiver tubes. The fluid is then cooled down by two air-cooled heat exchangers characterised by a maximum thermal capacity of 400 kWth. A 1 m<sup>3</sup> expansion vessel (filled by nitrogen) is placed in between the two air coolers to regulate the loop pressure (18 bar max). Finally, two additional electric heaters are installed at the outlet of the pump to control the oil temperature at the PTC inlet port.

The temperatures at the inlet and at the outlet of the PTC are measured with temperature transmittance sensors (denoted T in Figure 2). The direct normal irradiation (DNI) is measured with a pyrheliometer model CH1 by Kipp&Zonen [26]. A weather station installed nearby the solar field is used to measure the ambient conditions (temperature, humidity, wind speed and direction). The data acquisition system is based on National Instrument devices with a sampling time of 5 seconds. LabView software is used for data visualisation. The mass flow rate is derived from a data calibration study performed on the pump to assess it in function of the working conditions.

**Figure 2.** Process flow diagram of the Parabolic Trough Test Loop (PTTL) facilities with the relative sensors and devices position [21].

In Table 1, the working conditions of the main variables and of the external ambient parameters during the experimental campaign are reported.

**Table 1.** Range of the operation of the EuroTrough collectors (ETC) main variable and of the external ambient condition during the experimental campaign.


In order to perform the dynamic validation of the modelling method, three different scenarios are tested on the system:


### **5. Results and Discussion**

In this section, the model presented in Section 3 is validated onto the data presented in Section 4. The simulation is performed under Dymola2017 using Modelica language and the Differential Algebraic System Solver [27] with a relative tolerance of 10<sup>−</sup>4. The results of [21], obtained with the one-volume finite modelling, are also used to compare the accuracy of the current model with this other common modelling method.

A good agreement is found between the experimental data and the modelling results as it can be seen in the figures available into this section. In Figure 3a,b, only TE and MFE scenarios are considered. When the mass flow rate increases (resp. decreases), the outlet pipe temperature decreases (resp. increases) a few times after, due to the quasi-constant energy transmitted to the pipe (during this experiment the DNI is quite constant). On the other hand, a heat wave propagation induced by an increase of the inlet temperature is also correctly modelled, while the trend of the outlet pipe temperature follows the experimental data.

**Figure 3.** Experimental results performing different days versus plug flow and finite volume methods for several inlet temperature and mass flow rate conditions. (**a**) Test conditions #1; (**b**) Test conditions #2.

For the sake of clarity, in Figure 4a,b, only the SBE steps are noted on the experimental trends (denoted SBE) while TE and MFE conditions can be deduced from the figures. These SBE conditions are used to determine the influence of a cloud passing over the PTC system.

**Figure 4.** Experimental results performing different days versus plug flow and finite volume methods for several inlet temperature, mass flow rate and solar beam conditions. (**a**) Test conditions #3; (**b**) Test conditions #4.

Once again, the evolution of the predicted outlet temperature follows the experimental trends for several mass flow steps, inlet temperature steps and solar energy steps or a combination of them (Figure 4).

In all the cases studied, the accuracy is of the same order of magnitude as the finite volume method, with a discretisation of 50 cells which is the reference case of the study [21].

While the aim of this study is to consider active control strategies, the last experiment is presented in Figure 5 to assess the dynamic behaviour of the system when combinations of operating conditions steps are considered. Once again, there is a good agreement between plug flow modelling and the experimental data for wide variations of the mass flow rate, inlet pipe temperature and solar beam.

**Figure 5.** Experimental results performing one day versus plug flow and finite volume methods for a combination of inlet temperature, mass flow rate and solar beam step conditions.

As shown in previous study of district heating network modelling, the thermal inertia of the system and the operating mass flow rate have a significant influence on the outlet temperature response. Indeed, the mass flow rate influences how the fluid is propagated into the pipe, while the influence of thermal inertia can be seen when SBE conditions are set and can lead to significant extra delay (over 200 s) to get an outlet temperature with a similar value as the inlet temperature (considering the heat transfer to the ambient constant).

In this contribution, the simulation time required to model the process with plug flow modelling is generally reduced from a factor 20 to 30 compared to the finite volume method. This variation depends on the operating conditions and their trends. While this factor could seem quite reduced, it can be increased drastically in a larger network due to the higher complexity of the general problem. Indeed, in large networks with some intersections between the pipes, the hydraulic part of the problem can become complex to solve due to quadratic law relationship between pressure losses and the mass flow rates.

#### **6. Conclusions**

Solar thermal networks take place in the world energetic transition to supply clean thermal energy to process or residential needs. To optimise them and develop dedicated control strategies, it is required to assess the state of the heat transfer fluid correctly in several key points of the network, depending on the operating conditions. A practical solution could be to instrument these plants, but this solution is generally expensive and is not always easy to implement due to technical constraints. A second solution consists of modelling the dynamic behaviour of the whole system based on few input data.

In this contribution, the use of a previously developed plug flow approach for the district heating network modelling is extended through an experimental validation on a solar thermal network test bench. This validation step is required to check the assessment of the dynamic behaviour of the plant. Indeed, the main difference between district heating networks and thermal networks is the magnitude of heat transfer between the piping system and the ambient. This experimental validation considers inlet temperature and mass flow rate variations, as it could be happening in some control strategies. On the other hand, the weather influence is simulated by a modified solar energy recovered to the pipes through a modified and controlled focussing and defocussing of the PTC system.

All the results show a good agreement between experimental data and simulation results for a wide range of operating conditions of the thermal plant. Moreover, the accuracy of the model is similar to those of the one-dimensional finite volume method with a reduced simulation time by a minimal factor of 20.

The next step of this work will consist of the definition and the optimisation of dedicated control strategies of the thermal plant to ensure the best control of the temperature required by the energy process.

**Author Contributions:** Conceptualization, K.S.; methodology, K.S.; software, K.S. and R.D.; validation, K.S.; formal analysis, K.S.; investigation, K.S.; resources, K.S. and R.D.; data curation, K.S. and R.D.; writing—original draft preparation, K.S.; writing—review and editing, K.S.; visualization, K.S. and R.D.; supervision, K.S.; project administration, K.S. and R.D. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Conflicts of Interest:** The authors declare no conflict of interest.

### **Nomenclature**


### **References**


© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

### *Article* **Design Evaluation for a Finned-Tube CO2 Gas Cooler in Residential Applications**

### **Charalampos Alexopoulos 1, Osama Aljolani 2, Florian Heberle 2,\*, Tryfon C. Roumpedakis 1, Dieter Brüggemann <sup>2</sup> and Sotirios Karellas <sup>1</sup>**


Received: 22 April 2020; Accepted: 10 May 2020; Published: 12 May 2020

**Abstract:** Towards the introduction of environmentally friendlier refrigerants, CO2 cycles have gained significant attention in cooling and air conditioning systems in recent years. In this context, a design procedure for an air finned-tube CO2 gas cooler is developed. The analysis aims to evaluate the gas cooler design incorporated into a CO2 air conditioning system for residential applications. Therefore, a simulation model of the gas cooler is developed and validated experimentally by comparing its overall heat transfer coefficient. Based on the model, the evaluation of different numbers of rows, lengths, and diameters of tubes, as well as different ambient temperatures, are conducted, identifying the most suitable design in terms of pressure losses and required heat exchange area for selected operational conditions. The comparison between the model and the experimental results showed a satisfactory convergence for fan frequencies from 50 to 80 Hz. The absolute average deviations of the overall heat transfer coefficient for fan frequencies from 60 to 80 Hz were approximately 10%. With respect to the gas cooler design, a compromise between the bundle area and the refrigerant pressure drop was necessary, resulting in a 2.11 m2 bundle area and 0.23 bar refrigerant pressure drop. In addition, the analysis of the gas cooler's performance in different ambient temperatures showed that the defined heat exchanger operates properly, compared to other potential gas cooler designs.

**Keywords:** supercritical carbon dioxide; experimental testing; finned-tube gas cooler

### **1. Introduction**

The use of air conditioning systems is expanding rapidly around the world. An estimated amount of 700 million air conditioners will be operating in the world by 2030 [1]. This growing demand for air conditioning systems has enormous impacts on the environment.

Currently, a number of present regulations have been applied worldwide to control the use of harmful refrigerants [2–4]. The key implications of the use of conventional refrigerants include the depletion of the ozone layer and global warming. Based on the Montreal Protocol, a complete abolishment of chlorofluorocarbons (CFCs) was decided, due to their high ozone depletion potential (ODP) [5]. In addition, the phase out of hydrochlorofluorocarbon (HCFC) refrigerants was implemented [6]. On the other hand, the F-gas Regulation, first issued in the European Union in 2006, aimed to introduce measures for the reduction of fluorinated gases—hydrofluorocarbons (HFCs) and perfluorocarbons (PFCs)—in form of a phase-down, due to their high global warming potential (GWP) [7].

Instead, natural refrigerants are proposed as the substitute for the harmful refrigerants. Carbon dioxide (CO2) is a natural, low cost, non-flammable, non-toxic refrigerant. Subsequently, it has emerged as a credible natural refrigerant to replace HFCs and HCFCs. However, its unique critical point, high critical pressure of 73.8 bar, and low critical temperature of 30.98 ◦C, remarkably affects the performance of CO2 refrigeration systems, as well as imposes special design and control challenges. When the ambient temperature is higher than the critical temperature of CO2, the system operates at supercritical conditions, and the heat rejection process occurs at a supercritical regime. In consequence, a phase change does not take place, and the heat exchanger in which this change of state occurs is called the gas cooler.

The impact of the gas cooler on the CO2 refrigeration systems plays an important role, due to its high exergy loss. Therefore, it is considered vital to be further investigated and designed properly [8]. The finned-tube type for gas coolers is well established in the heating, ventilation, and air conditioning (HVAC) and refrigeration industries, due to its compactness and manufacturing flexibility. The design of the finned-tube heat exchangers affects considerably the overall heat transfer performance and system efficiency. Particularly, fin and tube thickness and the respective materials, spacing, and dimensions of the tubes and fins are crucial parameters of the design [9]. Fundamental studies about the heat transfer characteristics during the heat rejection process in tubes have been performed theoretically and experimentally by many researchers since Lorentzen and Pettersen [10] proposed the transcritical CO2 cycle for mobile air conditioning systems.

Pitla et al. [11] conducted an investigation about heat transfer phenomena and pressure losses of CO2 at supercritical conditions into a tube. They found that the majority of the deviations between the numerical and experimental values are within ± 20%, and a new heat transfer correlation was presented. Son and Park [12] carried out an experiment in order to investigate the gas cooling process of CO2 in terms of heat transfer coefficient and pressure drop of the refrigerant. They described the variations of local heat transfer coefficient in the cooling process in the direction of the flow and proposed a more accurate heat transfer correlation. Zhang et al. [13] evaluated the performance of a printed circuit heat exchanger for cooling CO2 with water. The analysis concluded that rapid variations in the thermodynamic properties of supercritical CO2 increase entropy generation and therefore, to optimize the second law efficiency of the investigated heat exchanger, higher CO2 mass flow rates should be used. Jadhav et al. [14] evaluated, using simulations, CO2 gas coolers for air conditioning applications. For their investigation, a counter crossflow plain fin and staggered tube configurations were considered. According to the study, transverse tube spacing, gas cooler width, and air volumetric flow were the most influential parameters in the heat transfer mechanisms of the gas cooler.

Liu et al. [15] investigated experimentally the supercritical CO2 characteristics in horizontal tubes with inner diameter of 4, 6, and 10.7 mm in terms of heat transfer phenomena and pressure losses. The authors concluded that the tube diameter significantly affects the heat transfer performance, and they proposed a new heat transfer correlation for the large diameter. Jiang et al. [16] investigated the convection heat transfer of CO2 at supercritical pressures in a vertical small tube with inner diameter of 2.0 mm, experimentally and numerically. They studied the effects of various operational parameters and buoyancy on convection heat transfer in a small diameter. They concluded that when the CO2 bulk temperatures are in the near-critical region, the local heat transfer coefficients vary significantly along the tube. Chai et al. [17] investigated, using simulations, the performance of finned-tube supercritical CO2 gas coolers, combining a distributed modeling approach with the ε-NTU method. The results indicated that the performance of the gas coolers was enhanced by higher mass flow rates and lower tube diameters at the expense of higher pressure drops.

Other researchers have also investigated the performance of the air-cooled CO2 gas coolers. Cheng et al. [18] presented an analysis of heat transfer and pressure drop experimental data and correlations for supercritical CO2 cooling in macro- and micro-channels. Ge and Cropper [19] presented a detailed mathematical model for air-cooled finned-tube CO2 gas coolers. They used a distributed method in order to obtain more accurate refrigerant thermophysical properties and local heat transfer

coefficients during cooling processes. The model was compared with published test results. The comparison showed that the approach temperature and the heat capacity are simultaneously improved with the increase of heat exchanger circuit numbers. Marcinichen et al. [20] conducted simulations to optimize the working fluid charge of the gas cooler. The optimal design of the study reduced CO2 charge by 14%, compared to a reference design. Moreover, the analysis revealed the importance of the oil concentrations in the CO2 pressure drop, which is up to 2.65 times higher for oil concentrations of up to 3%. Zilio et al. [21] experimentally evaluated two different gas coolers, one with continuous, and one with separated fins, and on two different circuit arrangements for a transcritical CO2 cycle. Using a coil with fins, a heat flux improvement of up to 5.6% was identified, which corresponded to a coefficient of performance (COP) increase of up to 6.6% for a conventional CO2 refrigeration cycle.

Gupta and Dasgupta [22] applied a similar modelling method to the one from Ge and Cropper, [19] in order to evaluate the performance of the heat exchanger being affected by the airflow velocity. Here, a higher gas cooler performance is achieved at a higher air flow velocity as it decreases the refrigerant's approach temperature, and thus the heating capacity of the gas cooler is increased. Santosa et al. [23] built two CO2 finned-tube gas coolers with different structural designs and controls, connected with a test rig of a CO2 booster refrigeration system. They carried out experiments at different operating conditions while they developed models of the finned-tube CO2 gas cooler. The analysis was conducted based on the distributed and lumped methods. They concluded that the heat exchanger design can affect the performance of both the component and the integrated system.

Although the heat transfer and pressure drop characteristics of the supercritical CO2 in tubes have been investigated extensively using experimental and theoretical methods, research on the air-cooled finned-tube CO2 gas coolers is still limited.

In this paper, mathematical calculations of the finned-tube CO2 gas cooler are conducted, in order to establish a reliable design procedure. With focus on the heat transfer characteristics of the air side, the developed model is validated with an experimental setup using water as working fluid. Investigations of the effects of fan frequency, water inlet temperature, and water mass flow on the overall heat transfer coefficient are conducted, while deviations between the model and the test results are extracted according to the fan frequency. In addition, potential heat transfer correlations for the airand refrigerant-side heat transfer coefficients have been studied. Finally, the model was applied to identify a reliable and efficient finned-tube CO2 gas cooler design, as well as to evaluate its performance in different off-design conditions under varying ambient temperatures.

### **2. Materials and Methods**

This study is part of a larger project of CO2 air conditioning systems for residential applications and focuses on the gas cooler. A scheme of the considered CO2 air conditioning system is depicted in Figure 1. Particularly, an efficient and reliable air finned-tube gas cooler was designed based on the boundary conditions, which are given in Table 1.


**Table 1.** On-design specifications of the gas cooler.

**Figure 1.** Scheme of the considered CO2 air conditioning system.

In order to design the CO2 gas cooler, a script in MATLAB R2019a [24] was developed. The simulation model calculated the overall heat transfer coefficient *U* of the gas cooler, based on the mass flows, inlet and outlet temperatures, and pressures of the medium. Subsequently, the required exchange area *AR* was determined. Both parameters were based on the following equations:

$$\mathcal{U} = \left[ \frac{1}{h\_{air}} + \frac{A\_o \times \ln \frac{d\_o}{d\_i}}{2 \times \pi \times L \times k} + \frac{A\_o}{A\_i} \times \frac{1}{h\_{refri}} \right]^{-1} \tag{1}$$

where *Ao* and *Ai* represent the outer and inner surface of the tube, respectively.

$$A\_R = \frac{Q}{U \times \Delta T\_{LMTD}}\tag{2}$$

The *U*-value consists of three parts: air convection, refrigerant convection, and the conduction. Compared to the other two heat transfer contributions, the conduction plays a minor role. The heat transfer coefficients on the air- and refrigerant-side are crucial for the overall heat transfer coefficient, and the validation for them are considered necessary. In order to validate the calculations, the overall heat transfer coefficient of a defined air-cooled heat exchanger was investigated experimentally.

The experimental part was based on a test rig using water and employed a specific design of an air-cooled heat exchanger (HEX). The HEX type was a finned tube with a fan air cooling system. The *U*-value was investigated experimentally for the entire HEX for different conditions. The second part of the validation consisted of modelling the heat exchanger to simulate the finned tube of the HEX.

### *2.1. Theoretical Model*

For the model calculations, a script was created in MATLAB R2019a [24], modelling the defined gas cooler. The model overall heat transfer coefficient of the HEX was calculated from Equation (1).

### 2.1.1. Air-Side Heat Transfer

In-line arrangement and circular finned tubes were assumed, in order to model the HEX. Based on the assumption of crossflow type, the air- and refrigerant-side heat transfer coefficients can be calculated. The proposed correlation from VDI-Heat Atlas [25] was used in order to calculate the Nusselt number, using the following equation:

$$Nu = C \times Re\_{d\_o} ^{0.6} \left(\frac{A\_o}{A\_{t0}}\right)^{-0.15} Pr^{\frac{1}{3}} \tag{3}$$

with *C* = 0.22 for in-line arrangement. *Ao*/*Ato* is the ratio of the finned surface to the surface of the base tube, and for circular fins was calculated from the following equation:

$$\frac{A\_o}{A\_{to}} = 1 + 2 \times \frac{H\_f \times \left(H\_f + d\_o + t\_f\right)}{s \times d\_o} \tag{4}$$

The Reynolds number was calculated by the equation:

$$Re\_{d\_o} = \frac{\rho\_{air} \times w\_s \times d\_o}{\mu} \,\tag{5}$$

where *ws* is the velocity in the smallest cross-section and was calculated from the following equation:

$$w\_{\ $} = w\_{\inf} \times \frac{A\_{\inf}}{A\_{\$ }} \tag{6}$$

The air-side heat transfer coefficient was calculated from its definition, as the following equation shows.

$$h\_{\rm air} = \frac{Nu \times k\_{\rm air}}{d\_0} \tag{7}$$

However, the air-side heat transfer coefficient was affected by the fins. The fins should be taken into consideration, thus the following equation was used [25]:

$$h\_{air,f} = h\_{air} \times \left[1 - \left(1 - \eta\_f\right) \times \frac{A\_f}{A\_o}\right] \tag{8}$$

The fin efficiency is defined as the ratio of the heat removed by the fin to the heat removed by the fin at wall temperature. The efficiency of the fin was calculated from the following equation:

$$
\eta\_{fin} = \frac{\tanh X}{X},
\tag{9}
$$

with [25]:

$$X = \varphi \times \frac{d\_o}{2} \times \sqrt{\frac{2 \times h\_{air}}{k \times t\_f}} \tag{10}$$

and

$$\varphi = \left(\frac{d\_f}{d\_o} - 1\right) \left[1 + 0.35 \ln\left(\frac{d\_f}{d\_o}\right)\right] \tag{11}$$

for circular fins [25]. So, the equation of the overall heat transfer coefficient used the updated air-side heat transfer coefficient as follows:

$$\Delta U = \left[ \frac{1}{h\_{\text{air},f}} + \frac{A\_o \times \ln \frac{d\_o}{d\_i}}{2 \times \pi \times L \times k} + \frac{A\_o}{A\_i} \times \frac{1}{h\_{refri}} \right]^{-1} \tag{12}$$

### 2.1.2. Refrigerant-Side Heat Transfer

On the refrigerant-side, the Gnielinski correlation [26] was used to calculate the Nusselt number:

$$Nu = \frac{\frac{f}{8}(Re\_b - 1000)Pr}{12.7\sqrt{\frac{f}{8}}(Pr^{\frac{2}{3}} - 1) + 1.07} \tag{13}$$

which is valid in the range of 2300 < *Reb* < 106.

In addition, the friction factor *f* was calculated by:

$$f = \left[ 0.79 \ln(Re\_b) - 1.64 \right]^{-2} \tag{14}$$

Here, the Reynolds number was defined as

$$\text{Re}\_b = \frac{G\_{refri} \times d\_i}{\mu} \tag{15}$$

and *G* was defined as the mass velocity, and was calculated from the following equation:

$$\mathbf{G}\_{refri} = \frac{m\_{refri}/N\_t}{\pi \times d\_l^2/4} \tag{16}$$

Finally, heat transfer coefficient at the refrigerant side was calculated from:

$$h\_{nefri} = \frac{Nu \times k\_{nefri}}{d\_o} \tag{17}$$

Further investigation of potential heat transfer correlations was conducted. More specifically, comparisons between the correlations of Gnielinski [26] and Dittus–Boelter [27], and the correlations proposed by VDI-Heat Atlas [25] and by Schmidt [28] were made for the refrigerant- and air-side, respectively. The equations below illustrate the Dittus–Boelter's [27] and Schmidt's [28] correlations, respectively:

$$Nu = 0.023 \times Re^{4/5} Pr^{\mu} \tag{18}$$

where *n* = 0.3 for the fluid being cooled.

$$Nu = C \times Re^{0.625} Pr^{1/3} \left(\frac{A\_o}{A\_{to}}\right)^{-0.375} \tag{19}$$

where *C* = 0.3 for in-line arrangement.

### *2.2. Experimental Set Up*

As it is referred, the *U*-value consists of three parts: the air convection, refrigerant convection, and the conduction. In the present case, the air-side heat transfer represents the main thermal resistance. Thus, the experimental set up aimed to identify a suitable heat transfer correlation with focus on the air side. In order to validate the calculations, especially the air-side heat transfer, a well-known working medium should be used for the experiments on the refrigerant-side. Here, water with well-known thermophysical properties and reliable heat transfer correlations at single-phase regime was selected as a working medium.

The test rig consisted of the heater, the gas cooler, the measurement equipment, and controls. To enable the information to be read and recorded, the instrumentations were connected to a data logging system. The test rig is shown in Figure 2 below.

**Figure 2.** Schematic of the gas cooler test rig.

The air-cooled HEX used for the experiments was the CU-713CX2 from Panasonic. The heater was from the Single® company, model STW 150/1-18-45-KS7. The K-type thermocouples and pressure transducers used were from OMEGA company with uncertainties of ±0.2 ◦C and ±1.5%, respectively, while the mass flow valve and meter with uncertainties of ±0.5% were manufactured by Bürkert.

The experiments were carried out for different water mass flows, fan frequencies, and inlet water temperatures, as the following Table 2 shows. Figure 3a,b illustrate the air mass flow rate and fan power as function of the fan frequency.


**Table 2.** Range and interval of the experimental variables.

**Figure 3.** (**a**) Air mass flow rate depending on fan frequency; (**b**) fan power as a function of fan frequency.

Using Equations (20)–(24), the experimental overall heat transfer coefficient was calculated by:

$$Q = m\_{\text{w}} \times \overline{c}\_{p\_{\text{w}}} \times \Delta T\_{\text{w}} \tag{20}$$

$$T\_{\rm air,o} = \frac{Q}{m\_{\rm air} \times \overline{c}\_{\rm p\_{air}}} + T\_{\rm air,i} \tag{21}$$

$$
\Delta T\_{LMTD} = \frac{\Delta T\_2 - \Delta T\_1}{\left(\frac{\Delta T\_2}{\Delta T\_1}\right)}\tag{22}
$$

$$
\Delta T\_2 = T\_{w,i} - T\_{air\rho} \text{ and } \Delta T\_2 = T\_{w\rho} - T\_{air,i} \tag{23}
$$

$$\mathcal{U} = \frac{\mathcal{Q}}{A\_{\mathcal{R}} \times \Delta T\_{LMTD}} \tag{24}$$

### *2.3. Model Application*

The validated model was finally applied to define an efficient and reliable gas cooler design, as well as to evaluate its performance for different off-design conditions. Within scope, the on-design analysis investigated different designs of heat exchangers, such as the number of rows and length of the tube. The target of the off-design analysis was to evaluate the operation of the gas cooler in different ambient temperatures.

### 2.3.1. On-Design

The developed model was utilized to design an air-finned CO2 gas cooler. Therefore, three and four numbers of rows and finned tubes with inner diameters of 6.85, 16, and 22.5 mm were considered. In order to identify an efficient gas cooler design, a compromise between the bundle area and the refrigerant pressure drop of the gas cooler was aimed for. The model was provided with the inlet, outlet temperatures and pressures, medium mass flows, and the duty of the heat exchanger as input variables. In addition, the design properties of the tube were specified, so the total area of the tube was calculated. The air and refrigerant heat transfer coefficients, which were initially unknown, were assumed, and then the overall heat transfer coefficient was calculated from Equation (1). The required exchanged area was calculated from Equation (2), and so the required number of tubes was obtained. In the iteration process, the air and refrigerant heat transfer coefficients were updated using Equations (3)–(17). The iteration was continued until the relative tolerance of the overall heat transfer coefficient for two continuous iterations was equal to 0.001.

The pressure drop was calculated by the following equation:

$$
\Delta P = f \times \frac{{G\_{refn}}^2}{2 \times \rho\_{refn}} \times \frac{L}{d\_i} \tag{25}
$$

$$f = \left[1.82\ln(Re\_b - 1.64)\right]^{-2} \tag{26}$$

which Filonenko [27] applies for 10<sup>4</sup> <sup>≤</sup> *Reb* <sup>≤</sup> <sup>5</sup> <sup>×</sup> 106.

In case of *Reb* <sup>≤</sup> <sup>10</sup>4, Blasius [27] correlations was applied:

$$f = \frac{0.316}{\text{Re}\_b \text{\AA}} \tag{27}$$

The thermophysical properties of air and refrigerant like density, viscosity, specific heat capacity, and thermal conductivity were obtained from REFPROP version 10 database. The model used the mean temperature and pressure of the mediums in order to calculate the heat transfer coefficients.

### 2.3.2. Off-Design

The calculations were used to evaluate the overall performance of the gas cooler in different conditions. The defined gas cooler was investigated in different ambient temperatures. The boundary conditions and operational parameters for the off-design analysis were defined according to Table 3.


**Table 3.** Off-design boundary conditions.

Based on the off-design data, the overall heat transfer coefficient was calculated using Equations (3)–(17). In addition, investigation of different potential heat exchanger designs was conducted by comparing their overall heat transfer coefficient and the refrigerant pressure drop. The potential heat exchanger's designs were based on the on-design analysis, and were chosen in terms of pressure losses and bundle area. The selected air-cooled HEXs were designed with identical finned-tube properties.

### **3. Results and Discussion**

### *3.1. Validation of the Model*

An experimental campaign was carried out in order to validate the mathematical calculations for the performance of an air-finned CO2 gas cooler. The *U*-value from experimental results was compared with the *U*-value obtained from the model. The following results were obtained using the proposed correlation from VDI-Heat Atlas [25] and Gnielinski's [26] correlation for the air-, and refrigerant-side, respectively. Figure 4a illustrates that most of the deviations were lower than 10%. Particularly, the average absolute deviation for 45 ◦C inlet water temperature and 50 Hz was 5%, while the maximum absolute deviation was 12%. The average absolute deviation for 50 Hz for different inlet water temperature was 11%. The calculations seem to approach the experimental results from 50 to 80 Hz. Figure 4b depicts that the absolute deviations for 80 Hz were even lower than 10%. Particularly, the absolute average deviations for 45 ◦C of 80 Hz was 5%, while the maximum absolute deviation was 9%. The absolute average deviation for 80 Hz of fan frequency for different water inlet temperatures was 6%.

**Figure 4.** Deviations between the model and experimental overall heat transfer coefficient for 45 ◦C water inlet temperature at (**a**) 50 Hz; and (**b**) 80 Hz fan frequency.

### *3.2. Sensitivity Analysis*

The most satisfying results were obtained by using the correlations proposed by VDI-Heat Atlas [25] and by Gnielinski [26] for the air- and the refrigerant-side heat transfer, accordingly. Fundamental investigations regarding the dependence of the heat transfer characteristics (*U*-value) on different operational parameters and applied correlations were conducted. Particularly, the impact of the fan frequency, the water mass flow, and the water inlet temperature on the performance of the air-cooled HEX were investigated. The increase of both the water mass flow and the fan frequency enhanced the overall performance of the heat exchanger, as Figure 5a illustrates. The increase of the water mass flow rate had a strong impact on the overall heat transfer coefficient by 20% and 30% of 50 and 80 Hz, accordingly. The increase of the fan frequency showed an average increase of the overall heat transfer coefficient of 8%, 6%, and 6% in the ranges of 50–60, 60–70 and 70–80 Hz, accordingly. The increase of the water inlet temperature affected in a similar way the overall heat transfer coefficient of the heat exchanger Figure 5b. An average increase of 35% of the overall heat transfer coefficient from 35 to 75 ◦C of water inlet temperature for fan frequency of 60 Hz was revealed.

The investigation of different heat transfer correlations is considered necessary to identify the most suitable air-side heat transfer correlation. The comparison between the refrigerant-side heat transfer correlations shows that the deviations are approximately 1% for all the different fan frequencies, while the correlation by Gnielinski [26] calculates the heat transfer as lower than by that of Dittus–Boelter in Figure 6a. The comparison between the air-side heat transfer correlations shows that the deviations are highly affected by the fan frequency. Particularly, absolute average deviations for 50 and 80 Hz are 8% and 4%, accordingly, while the maximum absolute deviations are 10% and 5%. It is revealed that with the increase of the fan frequency, the deviations between the correlations become lower, as shown in Figure 6b.

**Figure 5.** Overall heat transfer coefficient of (**a**) varying mass flow rate and fan frequencies; and (**b**) varying water inlet temperatures and water mass flow rates.

**Figure 6.** Comparison for 45 ◦C water inlet temperature for (**a**) refrigerant-side heat transfer correlations; and (**b**) air-side heat transfer correlations.

### *3.3. Design Procedure*

The discussed model was applied. Thus, an on-design analysis was investigated in order to identify an efficient and reliable gas cooler in terms of pressure losses and required exchange area. The mathematical calculations were applied to on-design analysis using the correlations by Gnielinski [26] and those proposed by VDI-Heat Atlas [25] for the refrigerant- and air-side heat transfer characteristics, accordingly. The investigation of different number of rows (NR) and finned-tubes showed that the heat exchanger with four rows and the smallest diameter has a smaller bundle area (Figure 7).

**Figure 7.** Bundle area vs. tube length for 3 and 4 number of rows (**a**) outside diameter of 7.35 mm; and (**b**) outside diameter of 18.5 mm.

Taking the pressure drop into consideration, the results show that the tube with the outer diameter of 7.35 mm cannot be used, as the occurring pressure drop causes serious operation problems to the system, due to pressure drops exceeding 0.3 bar (Figure 8a). Instead, the heat exchanger designed with the tube with outside diameter of 18.5 mm can be used for lengths to 20 m (see Figure 8b). Figures 9 and 10 show the Reynolds number for air- and refrigerant-side, respectively.

**Figure 8.** Refrigerant pressure drop vs. length of the tube for 4 rows (**a**) outside diameter of 7.35 mm; and (**b**) outside diameter of 18.5 mm.

**Figure 9.** Reynolds number air side (**a**) outside diameter of 7.35 mm; and (**b**) outside diameter of 18.5 mm.

**Figure 10.** Reynolds number refrigerant side (**a**) outside diameter of 7.35 mm; and (**b**) outside diameter of 18.5 mm.

Based on the presented results, the final design of the air-finned CO2 gas cooler is defined for a finned-tube (type-*U*), according to the specifications listed in Table 4.

### *3.4. Evaluation at O*ff*-Design Conditions*

The off-design analysis of the defined gas cooler was carried out for different ambient temperatures from 20 to 35 ◦C, and a comparison between potential gas coolers was made. The potential heat exchangers were chosen in terms of pressure losses and bundle area, while they had identical finned tube properties, except for the length of the tube. Particularly, the gas cooler with bundle area of 2.39 m<sup>2</sup> was characterized by its high pressure losses of 0.3 bar, while the gas cooler with bundle area of 2.81 m2 showed relatively low pressure losses of 0.06 bar. The gas cooler with bundle area of 2.43 m<sup>2</sup> was chosen, due to its combination of comparatively low pressure losses of 0.14 bar, as well as its low bundle area.


**Table 4.** Gas cooler technical specifications.

Figure 11a below shows that the defined gas cooler operates better than the other air-cooled HEX in different ambient temperatures. The defined gas cooler's overall heat transfer coefficient was remarkably higher than the heat exchangers with the bundle area of 2.43 and 2.81 m2. On the other hand, the *U*-value of the heat exchanger of 2.39 m<sup>2</sup> bundle area is near to the defined heat exchanger, (*Ab* = 2.11 m2) but the refrigerant pressure drop is significantly higher (Figure 11b).

**Figure 11.** Comparison between potential heat exchangers in different ambient temperatures of (**a**) the overall heat transfer coefficient; and (**b**) the refrigerant pressure drop.

### **4. Conclusions**

A design procedure for a CO2 finned-tube gas cooler was developed and validated experimentally. The experimental focus was laid on the validation of the air-side heat transfer. Absolute deviations between the model and the experiments were extracted and prove the reliability of the selected heat transfer correlations. The developed simulation model was used to design an efficient gas cooler and

evaluate its performance in different ambient temperatures. Based on these results, the following conclusions are justified:

The deviations between the calculations and the experiments are highly affected by the fan frequency, the water mass flow, and the water inlet temperature. It was found that the absolute average deviations for 60–80 Hz are less than 10%. The increase in the fan frequency, the water mass flow, and the water inlet temperature caused an improvement of the overall heat transfer coefficient of the heat exchanger. The comparison between potential heat transfer correlations showed that the combination of the correlations proposed by VDI-Heat Atlas [25] and by Gnielinski [26] for the airand refrigerant-side heat transfer coefficient, respectively, approached the experimental results better.

The heat exchanger with four rows has a smaller bundle area than with three rows, while the investigation of different tubes showed that the heat exchanger designed with the smallest diameter has the smallest bundle area and the highest refrigerant pressure drop. As a compromise between the bundle area and the refrigerant pressure drop, a gas cooler of 2.11 m<sup>2</sup> and refrigerant pressure drop of 0.23 bar was defined. Comparison between potential gas coolers was made, and the results show that the defined gas cooler operates more efficiently for different ambient conditions, compared to other potential heat exchangers.

Due to the pressure limitations of the equipment, the experimental validation of the model was conducted with water as working fluid. This procedure enables a reliable validation of the applied heat transfer correlation at the air side. However, the use of CO2 as working fluid in the tubes would improve validation approach. Therefore, experiments with CO2 are suggested for further work, next to the off-design analysis for the entire air conditioning system including the developed gas cooler model.

**Author Contributions:** All authors contributed to this work by collaboration. Methodology, software, validation, investigation and writing—original draft preparation: C.A.; off-design evaluation: O.A.; conceptual design of the study and review and editing: F.H.; methodology, software and writing, and review and editing: T.C.R.; supervision: S.K. and D.B. All authors have read and agreed to the published version of the manuscript.

**Funding:** This publication is strongly related to a research project (Project-ID TEW01CO2P-73707) at the Center of Energy Technology (ZET) of the University of Bayreuth, which is financed by the Bavarian State Ministry of Environment and Consumer Protection. The main research work of this publication was conducted during a stay of Charalampos Alexopoulos in Bayreuth, which was funded by the ERASMUS+ Programme. Additionally, the fourth author would like to state that this scientific paper was supported by the Onassis Foundation; Scholarship ID: G ZO 025-1/2018-2019.

**Conflicts of Interest:** The authors declare no conflicts of interest.

#### **Nomenclature**



### **References**


© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

*Article*
