Molecular Simulations of Adsorption and Energy Storage of R1234yf, R1234ze(z), R134a, R32, and their Mixtures in M-MOF-74 (M = Mg, Ni) Nanoparticles

The refrigerant circulation heat can be enhanced through the mutual transformation between thermal energy and surface energy during the adsorption and separation process of fluid molecules in porous materials. In this paper, the adsorption and energy storage of R1234ze(z), R1234yf, R32 and R134a, as well as their mixed refrigerants in Mg-MOF-74 and Ni-MOF-74 nanoparticles were investigated by means of molecular dynamics simulations and grand canonical Monte Carlo simulations. The results suggested that, in the case of pure refrigerant adsorption, the adsorption quantities of R32 and R134a in MOFs were higher than those of R1234yf and R1234ze(z). However, in the case of saturation adsorption, the desorption heat of R32 was lower than that of R1234yf and R1234ze(z). The addition of MOF-74 nanoparticles (NPs) could enhance the energy storage capacity of the pure refrigerant; besides, R1234yf and R1234ze(z) nanofluids had superior enhancement effect to that of R32 nanofluid. In mixed refrigerant adsorption, the adsorption quantities of R1234ze(z) and R1234yf were lower than those of R32 and R134a; with the increase in temperature, the adsorption of R1234ze(z) and R1234yf showed a gradually increasing trend, while that of R32 was gradually decreased.

www.nature.com/scientificreports www.nature.com/scientificreports/ nanofluid shows broad application prospect in the field of energy utilization; nonetheless, study on the adsorption and energy storage of porous medium nanofluid should be further carried out in the face of various refrigerants and porous materials.
There are multiple options regarding the refrigerants of ORC based on the different heat source environments. With the increasing enhancement in people's environmental awareness, the fourth generation refrigerants with low global warming potential (GWP) and zero ozone depletion potential (ODP) will gradually replace the currently used third generation refrigerants. At present, R1234yf, R1234ze(z) and R32, as well as their mixed refrigerants are one of the current research hotspots 11,12 .
Among the numerous currently available porous materials, MOFs have been extensively recognized as the most promising materials used for adsorption and energy storage 13,14 . MOFs are the porous materials with periodic network structure formed by metal ions or clusters with the organic ligands by means of coordinate self-assembly, which are associated with the advantages of high specific surface area, large porosity, low density and adjustable structure. They have exhibited huge application potentials in the fields of gas storage, separation, heterogeneous catalysis and drug sustained release 15 . Among various MOFs structures, MOF-74 is one of the MOFs with the highest unsaturated metal bit density, and the type of metal ion can be adjusted to modify the adsorption property of MOF-74.
MOFs have nanoscale pore structure, so it is difficult to examine the adsorption property of refrigerant in MOFs through conventional experiment and theoretical method. With the rapid development of computer technology, molecular simulation technique has been extensively applied in scientific research [16][17][18][19][20][21][22][23] . Plenty of literature suggests that, molecular simulation has become the third research means apart from experiment and theory 12,[24][25][26][27][28][29] . Therefore, this study had adopted molecular simulation methods to investigate the adsorption and energy storage properties of R1234yf, R1234ze(z), R134a and R32, as well as their mixtures in M-MOF-74 (M = Mg, Ni).

Materials and Methods
energy storage thermophysical model. The addition of MOFs NPs into the organic refrigerant can obtain the metal organic heat carriers (MOHCs). Theoretically, the energy (∆h MOHCs ) during the endothermic process of MOHCs is mainly constituted by three parts 10,30 , including a. enthalpy of phase change of organic refrigerant (∆h Fluid ); b. thermodynamic energy change ((∫CpdT) MOFs ) of MOFs particles, and c. desorption heat of fluid refrigerant in MOFs (∆h desorption ), namely, Alternatively, it can be expressed as: where x is the mass fraction of MOHC in MOF. It can be discovered from formula (2) that, when the sum of (∫CpdT) MOFs and ∆h desorption is greater than ∆h Fluid , MOHCs can enhance the energy storage property of the refrigerant. At present, research on the thermophysical property of pure organic refrigerant is relatively mature, which can be obtained through experimental and theoretical methods. In this paper, the thermophysical data of R1234y, R1234ze(z), R134a and R32 were retrieved by the National Institute of Standards and Technology (NIST) 31 . However, MOFs have complex and changeable structures and components, and their adsorption performance and thermophysical property required to be further implemented. In this paper, (∫CpdT) MOFs could be calculated through the thermodynamic energy change 32 of MOF-74 structure with the increase in temperature obtained through simulation of molecular dynamics (MD). Meanwhile, the desorption heat 33 of R1234yf, R1234ze(z), R134a and R32, and their mixed refrigerants in MOF-74 would be calculated according to the grand canonical Monte Carlo (GCMC) method. During the GCMC simulations, the chemical potential, volume and temperature of the system are set to be constant value. computational model. The MOFs calculation model in this paper was constituted by the 1 × 1 × 4 M-MOF-74 unit cells. The simulation box was Z: 27.0804 Å in size, and the bottom was an equiside parallelogram with one-side angle of 120 degrees and the side length of 25.7856 Å. M-MOF-74 was constituted by the self-assembly of divalent metal ions with the ligand 2,5-dihydroxyl terephthalic acid, which had formed the 3D alveolate spatial network structure with 2D hexogonal channel, and the metal atoms were formed in the manner of octahedral coordination (occupied by 5 oxygen atoms and 1 water molecule). There were 648 atoms in the model, including 288 C, 72 H, 216 O, and 72 Mg/Ni, as presented in Fig. 1. The molecular configuration of R1234yf, R1234ze(z), R134a and R32 is shown in Fig. 2.
The MD and GCMC simulations in this paper were realized on the Materials Studio molecular simulation software 34,35 , which offers a rich set of features for materials modeling. In the simulation, MD and GCMC had adopted the COMPASS force field 36,37 to describe the inter-atom interactions, the Ewald method was employed to manage the long-range electrostatic interaction between particles, and the periodic boundary conditions were applied in the X, Y and Z directions of the simulation box [38][39][40] . The cut off radius is 10 Å in the simulations.

MD Simulation details.
The Forcite module in Materials Studio was selected for MD simulation, and the thermodynamic energy changes of MOF-74 particles under different temperatures were calculated in the NVT system. The time step was set at 1 fs. The systems were run 1,000 ps to reach the equilibrate state first. Then, system www.nature.com/scientificreports www.nature.com/scientificreports/ run another 1,000 ps to analysis the data with information of atoms stored every 1,000 fs to analyze the results. The calculation temperatures were selected at 293 K, 313 K, 333 K, 353 K, 373 K and 393 K. The Berendsen thermostat heat bath method was employed to control the system temperature 41 . All the MD simulation work repeats three times. The uncertainties of the MD simulations is the statistical error of the data with information of atoms stored every 1,000 fs after the balance.
GCMC Simulation Details. The Sorption module in Materials Studio was selected for GCMC simulation, so as to calculate the isothermal adsorption processes of R1234yf, R1234ze(z), R134a, R32, and their mixed refrigerants in MOF-74 at different temperatures (293K, 313K, 333K, 353K, 373K and 393K). Among them, the pure refrigerant simulation pressure range was 1-5,000 kPa, while that of mixed refrigerants (R1234yf/R32, R1234ze(z)/R32, R1234yf/R134a, R1234ze(z)/R134a) was 1-16,000 kPa. The fugacity was calculated using the Peng-Robinson equation. Each state point was balanced after 1,000,000 cycles and averaged with 2,000,000 cycles. The data with information of atoms stored every 2,000 cycles to analyze the results. And the uncertainties of the

Results and Discussions
In order to check the reliability of the simulation work, the adsorption of R134a in M-MOF-74 was simulated at 298K and compared with the experimental work of Zheng et al. 42,43 , as shown in Fig. 3. It would better reveal the reliability of the simulation work in this investigation to some extent. pure refrigerant adsorption. Herein, Fig. 4 presents the adsorption isotherms of R1234yf, R1234ze(z), R134a and R32 pure refrigerants in M-MOF-74 at 6 calculation temperatures. It could be discovered that, the adsorption quantities of R134a and R32 in M-MOF-74 were slightly higher than those of R1234yf and R1234ze(z) in corresponding M-MOF-74. This was mainly related to the molecular structure size. As shown in Fig. 2, the molecular sizes of R134a and R32 were smaller than those of R1234yf and R1234ze(z), thus more R134a and R32 were adsorbed into the MOF-74 pore. Additionally, the adsorption quantities of R1234yf, R1234ze(z), R134a and R32 in Mg-MOF-74 were greater than those in Ni-MOF-74, which was because that the ion radius of Mg 2+ was smaller than that of Ni 2+ , so that Mg-MOF-74 had stronger attracting effect 44,45 .
In order to calculate the energy difference (∆h MOHCs ), the desorption heat and enthalpy values of all the refrigerants are calculated. Figure 5 has displayed the enthalpy values of R1234yf, R1234ze(z), R134a and R32 under different pressures retrieved through NIST, as well as the desorption heat of R1234yf, R1234ze(z), R134a and R32 in M-MOF-74 calculated by GCMC. At saturation adsorption, the desorption heat of R32 in MOF-74 was smaller than that of R1234yf, R1234ze(z) and R134a in MOF-74. The desorption heat of R1234yf, R1234ze(z), R134a and R32 in Mg-MOF-74 was greater than that in Ni-MOF-74, which was because that, Mg-MOF-74 had greater adsorption capacity than Ni-MOF-74. However, the enthalpies of phase change of R32 within the researched temperature and pressure ranges were greater than those of R1234yf, R1234ze(z) and R134a.
Adsorption and energy storage of pure refrigerant. It had been proved in previous work 46 that, the MD method could obtain reliable thermodynamic energy changes. In this paper, MD calculation could obtain the thermodynamic energy values of M-MOF-74 particles under different temperatures, as shown in Table 1. The thermodynamic energy value of MOF-74 particles was increased with the increase in temperature, which conformed to the general variation rule of solid material under heating. The Cp represented the slope of curve of thermodynamic energy with temperature, and the calculated values in this paper were consistent with those reported in literature 44 .
As mentioned in section 2.1 above, the energy storage of MOHCs at the time of heat adsorption could be calculated by formula (1) or formula (2). Among them, the enthalpy of phase change of organic refrigerant (∆h Fluid ) and the desorption heat of fluid refrigerant in MOFs (∆h desorption ) were calculated from Fig. 5(a,b), respectively. And the thermodynamic energy change ((∫CpdT) MOFs ) of MOFs particles was calculated from Table 1. Additionally, the energy difference (∆h MOHCs ) can be calculated by the three parts of ∆h Fluid , ∆h desorption and (∫CpdT) MOFs . After that, the energy storage enhancement (∆E enhancement ) of MOHCs was calculated by, where ∆E enhancement represents the increase in energy storage. ∆h MOHCs,x=0,293K is the energy difference (∆h MOHCs ) of the discussed refrigerant with the mass fraction of MOHCs x = 0 at 293 K. Figure 6 demonstrates the relationship of energy storage enhancement of pure refrigerants R1234yf, R1234ze(z), R134a and R32 in M-MOF-74 particles at different mass fractions and temperatures during the heat www.nature.com/scientificreports www.nature.com/scientificreports/ adsorption process, which calculated at the reference pressures (2 MPa, 3 Mpa and 4 MPa). It could be discovered in Fig. 6 that, the addition of M-MOF-74 NPs in the pure refrigerants R1234yf, R1234ze(z), R134a and R32 could enhance the energy storage properties of refrigerants; besides, the energy storage effect was enhanced with the increase in the mass fraction of M-MOF-74 NPs. R1234yf shared similar molecular structure with R1234ze(z), and the M-MOF-74 NPs showed almost consistent effect on enhancing the heat storage of R1234yf and R1234ze(z), which was superior to that on R32, but was lower than that on R134a. Moreover, in the MOHCs property curve constituted by R32 and Ni-MOF-74, negative enhancement effect occurred at the temperature difference of 20 K, as presented in Fig. 6(h). This was mainly because that R32 had a too large enthalpy of phase change ∆h Fluid , which had surpassed the sum of the variation value of Ni-MOF-74 thermodynamic energy (∫C p dT) MOFs with temperature and the desorption heat of R32 ∆h desorption in Ni-MOF-74.    Fig. 7 shows the competitive adsorption behaviors of mixed refrigerants R1234yf/R32, R1234yf/R134a, R1234ze(z)/R32 and R1234ze(z)/R134a in M-MOF-74 under different temperatures and pressures. The adsorption quantities of R32 and R134a were far higher than those of R1234yf and R1234ze(z); besides, the adsorption quantities of R32 and R134a were reduced with the increase in temperature, while those of R1234yf and R1234ze(z) were increased with the increase in temperature. Compared with R32 and R134a, R1234yf and R1234ze(z) could reach the saturation adsorption status in a faster rate, which was because that the small molecule structures of R32 and R134a could effectively utilize the pore structure in MOF-74 to increase the adsorption quantity. With the increase in temperature, the molecular thermal motion was www.nature.com/scientificreports www.nature.com/scientificreports/ enhanced, some R32 could be desorbed from MOF-74, the empty points within MOFs were increased, and a small amount of R1234yf and R1234ze(z) would enter into the MOFs for adsorption. At the same time, in the mixed refrigerant adsorption, the adsorption quantity of each refrigerant was lower than that of pure refrigerant, which was because that the combined adsorption of mixed refrigerants would occupy the pore space in MOFs. On the other hand, the adsorption quantities of R1234yf and R1234ze(z) in the R1234yf/R134a and R1234ze(z)/R134a mixtures under the same temperature and pressure were greater than those of R1234yf/R32 and R1234ze(z)/ www.nature.com/scientificreports www.nature.com/scientificreports/ R32. In the meantime, the adsorption quantities of R1234yf, R1234ze(z), R134a and R32 in Mg-MOF-74 in the mixture were greater than those in Ni-MOF-74.
Then the adsorption selectivity was calculated to assess the adsorption capacity of a sorbent in mixture competitive adsorption simulations. The adsorption selectivity of sorbent A over sorbent B in M-MOF-74 was defined as,  where x i and y i are the mole fraction of sorbent i in the adsorbed phase and bulk phase, respectively. The adsorption selectivity of R32/R134a over R1234yf/R1234ze(z) in M-MOF-74 was calculated and the results were shown in Fig. 8. All selectivity values were greater than one, which indicated that the adsorption capacity of R32/R134a was superior to R1234yf/R1234ze(z) in M-MOF-74. Additionally, this superiority would decrease with increasing temperature. The adsorption selectivity of R32 over R1234ze(z) was the maximum value among the binary adsorptions in both MOF-74. Thereafter, the adsorption selectivity of R32 over R1234ze(z) in Mg-MOF-74 was higher than that in Ni-MOF-74, which indicated that Mg-MOF-74 had better adsorption selectivity than that of Ni-MOF-74. Figure 9 shows the adsorption enthalpies of mixed refrigerants in M-MOF-74 under different temperatures calculated through GCMC. In Mg-MOF-74, the desorption enthalpy of R1234ze(z)/R32 was higher than that of R1234yf/R32, while that of R1234ze(z)/R134a was lower than that of R1234yf/R134a. The desorption enthalpies of the mixed refrigerants in Mg-MOF-74 were higher than those in Ni-MOF-74.

conclusions
The efficiency of thermodynamic cycles can be improved by using the optimized working fluids. Using MD and GCMC simulation methods, this paper has studied the adsorption and energy storage performance of R1234yf, R1234ze(z), R32 and their mixed refrigerants in Mg-MOF-74 and Ni-MOF-74. The adsorption quantities of R32 and R134a in MOF-74 are higher than those of R1234yf and R1234ze(z), while at saturation adsorption state, the desorption heat of R32 in MOF-74 is lower than that of R1234yf and R1234ze(z). The adsorption quantities and desorption enthalpies of refrigerants in Mg-MOF-74 are higher than those in Ni-MOF-74. Noteworthily, the addition of MOF-74 NPs in the pure refrigerant can enhance its energy storage capacity; besides, the addition of MOF-74 NPs in R1234yf and R1234ze(z) to form MOHCs can result in superior heat storage enhancement effect to that of R32 MOHCs. Furthermore, in mixed refrigerant competitive adsorption, the adsorption quantities of R1234ze(z) and R1234yf are lower than those of R32 and R134a; with the increase in temperature, the adsorption quantities of R1234ze(z) and R1234yf are gradually increased, while those of R32 and R134a are gradually decreased. And Mg-MOF-74 has better adsorption selectivity than that of Ni-MOF-74. In general, this study is of great significance for theoretical analysis of adsorption characteristics of refrigerant molecules in porous media. And it can provide reliable basis for selecting appropriate refrigerant working mediums and MOFs to form optimal MOHCs in engineering.