Energy Storage Analysis of a Mixed R161/MOF-5 Nanoparticle Nanofluid Based on Molecular Simulations

The thermal properties of refrigerants can be modified by adding porous nanoparticles into them. Here, molecular simulations, including molecular dynamics and grand canonical Monte Carlo, were employed to study the thermal energy storage properties of an R161/MOF-5 nanofluid. The results show that the thermodynamic energy change of MOF-5 nanoparticles is linear to the temperature. The adsorption heat calculated by grand canonical Monte Carlo is close to that calculated by the Clausius–Clapeyron equation. Additionally, a negative enhancement of the thermal energy storage capacity of the R161/MOF-5 nanofluid is found near the phase transition area.


Introduction
With the rapid development of modern society, serious problems, such as energy shortage and environmental pollution, have become increasingly prominent [1,2]. Energy conservation and emissions reduction are effective ways to alleviate the energy and environmental problems. The thermodynamic cycle is the main approach for energy conversion. Additionally, the working fluid is the energy carrier of the thermodynamic cycle. Therefore, using various means to improve the thermophysical properties of working fluid can improve the efficiency of energy conversion.
Dispersing nanoparticles into traditional working fluids, e.g., water, alcohol, oil, and refrigerants, to produce uniform and stable nanofluids can effectively improve the thermal conductivity of working fluids [3][4][5][6]. Nanofluids have great prospects for application in the fields of energy, chemical engineering, automobiles, building and construction, etc. The enhancement of the thermal conductivity of nanofluids is mainly due to the effective medium theory [3]. Also, the layering phenomenon, Brownian motion, clustering, ballistic phonon motion, thermal boundary resistance, and mass difference scattering will impact on the thermal conductivity of nanofluids [7][8][9][10].
Besides this, a nanofluid with porous materials can be used for thermal energy storage [11]. Carbon nanotubes, zeolites, and Metal-Organic Frameworks (MOFs) are common nanoporous materials. These materials have a very large specific surface area compared to conventional materials. Thus, extra energy can be stored and output by the fluid desorbing/adsorbing in the nanoporous materials. Chen et al. [12] performed comprehensive research on the energy storage of carbon nanotube nanofluids under the action of heat, force, and electric coupling. McGrail et al. [13] proposed to use metal-organic heat carrier nanofluids (MOHCs, a refrigerant to which has been added MOF nanoparticles to form the nanofluid) to improve the efficiency of the Organic Rankine cycle (ORC, using refrigerants as the working fluid for the Rankine cycle). They reported [14] noticeable adsorption and reversible desorption of refrigerant molecules in MOFs. Therefore, MOHCs have great potential in low-grade energy utilization and refrigeration cycles.
Theoretically, the heat adsorption energy of MOHCs (h MOHCs ) is mainly composed of three parts: (i) the enthalpy change of the working fluid (∆h Fluid ); (ii) the thermodynamic energy change of MOF particles (( C p dT) MOFs ); and (iii) the desorption heat of the fluid in MOFs (∆h desorption ), i.e., [11,13], where x is the mass fraction of MOFs in the MOHC. Additionally, Equation (1) can be written as It can be inferred from Equation (2) that MOHCs are able to store more thermal energy than pure fluid when the sum of the thermodynamic energy change of MOF particles and the desorption heat of the fluid in MOFs is larger than the enthalpy change of the working fluid.
MOFs are a class of organic-inorganic hybrid materials constructed from organic ligands and inorganic metal units. They own a large specific surface area, large porosity, a diverse structure, and thermal stability, which have provided MOFs with a wide range of applications in the fields of energy, chemicals, materials, and so on [15]. In the numerous MOF structures, MOF-5 [16] is an MOF that consists of Zn 4 O clusters connected by 1,4-benzodicarboxylate linkers as shown in Figure 1. It has a specific surface area of 2000 m 2 /g and a good thermostability whose structure remains stable below 573 K [17]. Fluoroethane, also known as R161, is one of the most used organic refrigerants with zero ozone depression potential, low global warming potential, and good thermal properties. Additionally, the adsorption properties of R161 in MOF-5 have been rarely reported [14,18]. Thus, the thermal energy storage property of an MOHC based on R161 and MOF-5 is investigated in the present paper.  [13] proposed to use metal-organic heat carrier nanofluids (MOHCs, a refrigerant to which has been added MOF nanoparticles to form the nanofluid) to improve the efficiency of the Organic Rankine cycle (ORC, using refrigerants as the working fluid for the Rankine cycle). They reported [14] noticeable adsorption and reversible desorption of refrigerant molecules in MOFs. Therefore, MOHCs have great potential in low-grade energy utilization and refrigeration cycles. Theoretically, the heat adsorption energy of MOHCs (hMOHCs) is mainly composed of three parts: (i) the enthalpy change of the working fluid (∆hFluid); (ii) the thermodynamic energy change of MOF particles ((∫CpdT)MOFs); and (iii) the desorption heat of the fluid in MOFs (∆hdesorption), i.e., [11,13], where x is the mass fraction of MOFs in the MOHC. Additionally, Equation (1) can be written as It can be inferred from Equation (2) that MOHCs are able to store more thermal energy than pure fluid when the sum of the thermodynamic energy change of MOF particles and the desorption heat of the fluid in MOFs is larger than the enthalpy change of the working fluid.
MOFs are a class of organic-inorganic hybrid materials constructed from organic ligands and inorganic metal units. They own a large specific surface area, large porosity, a diverse structure, and thermal stability, which have provided MOFs with a wide range of applications in the fields of energy, chemicals, materials, and so on [15]. In the numerous MOF structures, MOF-5 [16] is an MOF that consists of Zn4O clusters connected by 1,4-benzodicarboxylate linkers as shown in Figure 1. It has a specific surface area of 2000 m 2 /g and a good thermostability whose structure remains stable below 573 K [17]. Fluoroethane, also known as R161, is one of the most used organic refrigerants with zero ozone depression potential, low global warming potential, and good thermal properties. Additionally, the adsorption properties of R161 in MOF-5 have been rarely reported [14,18]. Thus, the thermal energy storage property of an MOHC based on R161 and MOF-5 is investigated in the present paper. Besides this, due to the length scale of the nanopore structure in MOFs being small, it is difficult to investigate the thermal energy storage mechanism of the organic working fluid in MOFs by using conventional experimental and theoretical methods. Molecular simulation [19] is an ideal tool to reveal the micro-mechanisms of materials. Therefore, molecular simulation is employed to investigate the thermal energy storage of R161/MOF-5-based MOHCs. The work is also expected to provide useful insights on adsorption refrigeration and heat pumps [20]. Besides this, due to the length scale of the nanopore structure in MOFs being small, it is difficult to investigate the thermal energy storage mechanism of the organic working fluid in MOFs by using conventional experimental and theoretical methods. Molecular simulation [19] is an ideal tool to reveal the micro-mechanisms of materials. Therefore, molecular simulation is employed to investigate the thermal energy storage of R161/MOF-5-based MOHCs. The work is also expected to provide useful insights on adsorption refrigeration and heat pumps [20].

Method and Computational Details
According to the thermodynamic principles discussed in Equation (2), the thermodynamic parameters of a pure working fluid are classic, and can be obtained by an experimental or a theoretical method. Since MOFs are a novel material, the thermal properties of MOFs and the adsorption properties of the refrigerant in an MOF need further research. Therefore, here, the thermodynamic parameters of R161 were obtained from the National Institute of Standards and Technology (NIST [21]). The thermodynamic energy change of MOF-5 can be obtained by molecular dynamics (MD) simulation [22]. Additionally, the desorption heat of R161 in MOF-5 can be easily calculated from grand canonical Monte Carlo (GCMC) simulations [23].

Simulation Models
The calculation model of MOF-5 comprising 2 × 2 × 2 unit cells (X: 51.788 Å, Y: 51.788 Å, Z: 51.788 Å) is shown in Figure 1 with 1536 carbon atoms, 832 oxygen atoms, 768 hydrogen atoms, and 256 Zn atoms. The structure of R161 (CH3CH2F) is shown in Figure 2. The MD simulations and GCMC simulations were performed using Materials Studio [24]. In the simulations, the COMPASS force field [25] was used to describe the intra and inter molecular interactions. The Ewald method [26] was used to deal with the long-range Coulombic interactions. Periodic boundary conditions were applied in the X, Y, and Z directions in all of our simulations.

Method and Computational Details
According to the thermodynamic principles discussed in Equation (2), the thermodynamic parameters of a pure working fluid are classic, and can be obtained by an experimental or a theoretical method. Since MOFs are a novel material, the thermal properties of MOFs and the adsorption properties of the refrigerant in an MOF need further research. Therefore, here, the thermodynamic parameters of R161 were obtained from the National Institute of Standards and Technology (NIST [21]). The thermodynamic energy change of MOF-5 can be obtained by molecular dynamics (MD) simulation [22]. Additionally, the desorption heat of R161 in MOF-5 can be easily calculated from grand canonical Monte Carlo (GCMC) simulations [23].

Simulation Models
The calculation model of MOF-5 comprising 2 × 2 × 2 unit cells (X: 51.788 Å , Y: 51.788 Å , Z: 51.788 Å ) is shown in Figure 1 with 1536 carbon atoms, 832 oxygen atoms, 768 hydrogen atoms, and 256 Zn atoms. The structure of R161 (CH3CH2F) is shown in Figure 2. The MD simulations and GCMC simulations were performed using Materials Studio [24]. In the simulations, the COMPASS force field [25] was used to describe the intra and inter molecular interactions. The Ewald method [26] was used to deal with the long-range Coulombic interactions. Periodic boundary conditions were applied in the X, Y, and Z directions in all of our simulations.

MD Simulation Details
The MD simulations were performed in the Forcite module of Materials Studio. The timestep was set as 1 fs and the system was equilibrated for 100 ps to calculate the thermodynamic energies of MOF-5 nanoparticles at 280, 300, 320, 340, 360, 380, 400, and 420 K. The simulations were computed in the canonical ensemble (NVT) ensemble with the Berendsen method [27] to control the temperature.

GCMC Simulation Details
The GCMC simulations were performed in the Sorption module of Materials Studio. The adsorption isotherms (280, 300, 320, 340, 360, 380, 400, and 420 K) of the adsorption of R161 in MOF-5 nanoparticles were calculated from 0 to 6000 kPa. The fugacity was calculated by the Peng-Robinson equation. For each point of the adsorption isotherms, the equilibration time was 10,000 steps with another 100,000 steps for statistics.

Thermodynamic Energy of MOF-5
The thermodynamic energies of MOF-5 at different temperatures are plotted in Figure 3. It shows that the thermodynamic energy of MOF-5 increases linearly with the temperature rise. The increment of internal energy per unit of temperature is Cp of MOF. Here, the

MD Simulation Details
The MD simulations were performed in the Forcite module of Materials Studio. The timestep was set as 1 fs and the system was equilibrated for 100 ps to calculate the thermodynamic energies of MOF-5 nanoparticles at 280, 300, 320, 340, 360, 380, 400, and 420 K. The simulations were computed in the canonical ensemble (NVT) ensemble with the Berendsen method [27] to control the temperature.

GCMC Simulation Details
The GCMC simulations were performed in the Sorption module of Materials Studio. The dsorption isotherms (280, 300, 320, 340, 360, 380, 400, and 420 K) of the adsorption of R161 in MOF-5 nanoparticles were calculated from 0 to 6000 kPa. The fugacity was calculated by the Peng-Robinson equation. For each point of the adsorption isotherms, the equilibration time was 10,000 steps with another 100,000 steps for statistics.

Thermodynamic Energy of MOF-5
The thermodynamic energies of MOF-5 at different temperatures are plotted in Figure 3. It shows that the thermodynamic energy of MOF-5 increases linearly with the temperature rise. The increment of internal energy per unit of temperature is C p of MOF. Here, the ρc p = 1.24(J/cm 3 ·K), which is close to the value of a typical MOF, 1.44 J/cm 3 ·K [13,28].

Adsorption Isotherms and Enthalpy of Desorption
The adsorption isotherms of R161 adsorbed in MOF-5 are shown in Figure 4. The adsorption isotherms are also fitted by the Langmuir equation, details of which are shown in Appendix A. The adsorption snapshot in Figure 4 shows that R161 mainly distributes around the O and Zn atoms of MOF-5 because these atoms are more attractive than the benzene ligand. Also, the adsorption isotherms decrease as the temperature rises. This is because the raised temperature results in the increased motion of molecules. Thus, the kinetic energy of the fluid molecules can overcome the adsorption energy of the MOF structure. This is the basis for thermal energy storage in MOHCs.

Adsorption Isotherms and Enthalpy of Desorption
The adsorption isotherms of R161 adsorbed in MOF-5 are shown in Figure 4. The adsorption isotherms are also fitted by the Langmuir equation, details of which are shown in Appendix A. The adsorption snapshot in Figure 4 shows that R161 mainly distributes around the O and Zn atoms of MOF-5 because these atoms are more attractive than the benzene ligand. Also, the adsorption isotherms decrease as the temperature rises. This is because the raised temperature results in the increased motion of molecules. Thus, the kinetic energy of the fluid molecules can overcome the adsorption energy of the MOF structure. This is the basis for thermal energy storage in MOHCs.

Adsorption Isotherms and Enthalpy of Desorption
The adsorption isotherms of R161 adsorbed in MOF-5 are shown in Figure 4. The adsorption isotherms are also fitted by the Langmuir equation, details of which are shown in Appendix A. The adsorption snapshot in Figure 4 shows that R161 mainly distributes around the O and Zn atoms of MOF-5 because these atoms are more attractive than the benzene ligand. Also, the adsorption isotherms decrease as the temperature rises. This is because the raised temperature results in the increased motion of molecules. Thus, the kinetic energy of the fluid molecules can overcome the adsorption energy of the MOF structure. This is the basis for thermal energy storage in MOHCs.  The enthalpy of desorption can be calculated by the Clausius-Clapeyron (C-C) equation based on the adsorption isotherm, It can be rewritten as: Thus, the desorption heat q st can be obtained by calculating the slope relationship between lnP and 1/T, which is shown in Figure 5. Additionally, the results of the desorption heat calculated by the C-C equation and GCMC simulations are listed in Table 1. The desorption heat calculated by GCMC simulations is close to that calculated by the C-C equation. With consideration of the calculation error of the C-C equation, the desorption heat calculated by GCMC will be used to calculate the thermal energy storage in the following section. It can be rewritten as: Thus, the desorption heat qst can be obtained by calculating the slope relationship between lnP and 1/T, which is shown in Figure 5. Additionally, the results of the desorption heat calculated by the C-C equation and GCMC simulations are listed in Table 1. The desorption heat calculated by GCMC simulations is close to that calculated by the C-C equation. With consideration of the calculation error of the C-C equation, the desorption heat calculated by GCMC will be used to calculate the thermal energy storage in the following section.

Thermal Energy Storage
Then, the thermal energy storage capacity of MOHCs can be computed according to Equation (1) or Equation (2). As aforementioned, the enthalpy difference of MOHCs at a different temperature is the thermal energy storage capacity. Thus, the thermal energy storage capacity of MOHCs with a different mass ratio of MOF-5 is presented as Figure 6. The temperature of the cold source was set as 280 K. Here, the enthalpy difference of pure R161 was considered to be the reference working fluid. Also, the thermodynamic energy change of MOF particles and the desorption heat of the fluid in MOFs are calculated and plotted in Figure 7. Note that the R161 experiences a phase transition when the temperature difference is about 100 K. The enthalpy difference of pure R161 experiences a sharp increase near this temperature range. The thermodynamic energy change of MOF particles and the desorption heat of the fluid in MOFs linearly increase as the temperature difference rises. It can be concluded that the thermal energy storage property of MOHCs is enhanced when the R161 is in a liquid state. Additionally, the more MOF-5 nanoparticles are added, the more thermal energy MOHCs can store. However, the MOHCs performed a negative enhancement of their thermal energy

Thermal Energy Storage
Then, the thermal energy storage capacity of MOHCs can be computed according to Equation (1) or Equation (2). As aforementioned, the enthalpy difference of MOHCs at a different temperature is the thermal energy storage capacity. Thus, the thermal energy storage capacity of MOHCs with a different mass ratio of MOF-5 is presented as Figure 6. The temperature of the cold source was set as 280 K. Here, the enthalpy difference of pure R161 was considered to be the reference working fluid. Also, the thermodynamic energy change of MOF particles and the desorption heat of the fluid in MOFs are calculated and plotted in Figure 7. Note that the R161 experiences a phase transition when the temperature difference is about 100 K. The enthalpy difference of pure R161 experiences a sharp increase near this temperature range. The thermodynamic energy change of MOF particles and the desorption heat of the fluid in MOFs linearly increase as the temperature difference rises. It can be concluded that the thermal energy storage property of MOHCs is enhanced when the R161 is in a liquid state. Additionally, the more MOF-5 nanoparticles are added, the more thermal energy MOHCs can store. However, the MOHCs performed a negative enhancement of their thermal energy storage properties near the phase transition point of R161. This is because the enthalpy difference of the phase transition of R161 is large, which is more than the sum of the thermodynamic energy change of MOF-5 and the desorption heat of R161 in MOF-5 at this temperature range. Finally, the thermal energy storage property of MOHCs is enhanced when the temperature difference is over 140 K. Obviously, the sum of the thermodynamic energy change of MOF-5 and the desorption heat of R161 in MOF-5 is larger than the enthalpy difference of R161 when the temperature difference is over 140 K. This indicates that the MOHCs cannot be used for enhancing the thermal energy storage of a working fluid in certain conditions. Materials 2018, 11, x FOR PEER REVIEW 6 of 8 storage properties near the phase transition point of R161. This is because the enthalpy difference of the phase transition of R161 is large, which is more than the sum of the thermodynamic energy change of MOF-5 and the desorption heat of R161 in MOF-5 at this temperature range. Finally, the thermal energy storage property of MOHCs is enhanced when the temperature difference is over 140 K. Obviously, the sum of the thermodynamic energy change of MOF-5 and the desorption heat of R161 in MOF-5 is larger than the enthalpy difference of R161 when the temperature difference is over 140 K. This indicates that the MOHCs cannot be used for enhancing the thermal energy storage of a working fluid in certain conditions.

Conclusions
In this paper, molecular simulations (MD and GCMC methods), as well as thermodynamic calculation, are employed to investigate thermal energy storage by the adsorption of R161 in MOF-5. The R161/MOF-5 mixed MOHCs noticeably enhanced the thermal energy storage property of R161 in the liquid state. The thermal energy storage capacity of the studied MOHCs is less than that of storage properties near the phase transition point of R161. This is because the enthalpy difference of the phase transition of R161 is large, which is more than the sum of the thermodynamic energy change of MOF-5 and the desorption heat of R161 in MOF-5 at this temperature range. Finally, the thermal energy storage property of MOHCs is enhanced when the temperature difference is over 140 K. Obviously, the sum of the thermodynamic energy change of MOF-5 and the desorption heat of R161 in MOF-5 is larger than the enthalpy difference of R161 when the temperature difference is over 140 K. This indicates that the MOHCs cannot be used for enhancing the thermal energy storage of a working fluid in certain conditions.

Conclusions
In this paper, molecular simulations (MD and GCMC methods), as well as thermodynamic calculation, are employed to investigate thermal energy storage by the adsorption of R161 in MOF-5. The R161/MOF-5 mixed MOHCs noticeably enhanced the thermal energy storage property of R161 in the liquid state. The thermal energy storage capacity of the studied MOHCs is less than that of

Conclusions
In this paper, molecular simulations (MD and GCMC methods), as well as thermodynamic calculation, are employed to investigate thermal energy storage by the adsorption of R161 in MOF-5. The R161/MOF-5 mixed MOHCs noticeably enhanced the thermal energy storage property of R161 in the liquid state. The thermal energy storage capacity of the studied MOHCs is less than that of Conflicts of Interest: The authors declare no conflict of interest.

Nomenclature ∆h MOHCs
The enthalpy of MOHCs (kJ/kg) ∆h Fluid The enthalpy of pure organic fluid (kJ/kg) ∆h desorption The enthalpy of desorption (kJ/mol) T Temperature (K) C p The heat capacity of MOFs (kcal/(mol·K)) q st The adsorption heat (kcal/(mol·K)) x The mass fraction of MOF in MOHCs ρ Density (g/cm 3 )