Next Article in Journal
Star-Shaped Fe3-xO4-Au Core-Shell Nanoparticles: From Synthesis to SERS Application
Next Article in Special Issue
Two-Phase Equilibrium Conditions in Nanopores
Previous Article in Journal
Metal-Based Nanoparticles as Antimicrobial Agents: An Overview
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Gibbs Ensemble Monte Carlo Simulation of Fluids in Confinement: Relation between the Differential and Integral Pressures

1
Engineering Thermodynamics, Process & Energy Department, Faculty of Mechanical, Maritime and Materials Engineering, Delft University of Technology, Leeghwaterstraat 39, 2628CB Delft, The Netherlands
2
PoreLab, Department of Chemistry, Norwegian University of Science and Technology, 7031 Trondheim, Norway
*
Author to whom correspondence should be addressed.
Nanomaterials 2020, 10(2), 293; https://doi.org/10.3390/nano10020293
Submission received: 12 January 2020 / Revised: 3 February 2020 / Accepted: 5 February 2020 / Published: 9 February 2020
(This article belongs to the Special Issue Nanoscale Thermodynamics)

Abstract

:
The accurate description of the behavior of fluids in nanoporous materials is of great importance for numerous industrial applications. Recently, a new approach was reported to calculate the pressure of nanoconfined fluids. In this approach, two different pressures are defined to take into account the smallness of the system: the so-called differential and the integral pressures. Here, the effect of several factors contributing to the confinement of fluids in nanopores are investigated using the definitions of the differential and integral pressures. Monte Carlo (MC) simulations are performed in a variation of the Gibbs ensemble to study the effect of the pore geometry, fluid-wall interactions, and differential pressure of the bulk fluid phase. It is shown that the differential and integral pressure are different for small pores and become equal as the pore size increases. The ratio of the driving forces for mass transport in the bulk and in the confined fluid is also studied. It is found that, for small pore sizes (i.e., < 5 σ fluid ), the ratio of the two driving forces considerably deviates from 1.

Graphical Abstract

1. Introduction

The widespread application of nanoporous materials in several fields, such as chromatography, membrane separation, catalysis, etc., has lead to a growing interest in the accurate description of the thermodynamic behavior of fluids confined in nanopores [1,2,3,4,5,6]. The pressure of a nanoconfined fluid is one of the most important thermodynamic properties which is needed for an accurate description of the flow rate, diffusion coefficient, and the swelling of the nanoporous material [7,8,9,10]. Various approaches for calculating the pressure of a fluid in a nanopore have been proposed [1,11,12,13]. The main difficulty of the pressure calculation arises from the ambiguous definition of the pressure tensor inside porous materials due to the presence of curved surfaces and confinement effects [7,14,15]. Traditional thermodynamic laws and concepts, such as Gibbs surface dynamics, Kelvin equation, etc., may not be applicable at the nano-scale [16]. In the past decade, several methods were reported using different simulation techniques i.e., classical density functional theory [17,18], equation of state modeling [19], etc., to model the behavior of fluids in confinement.
Recently, Galteland et al. [11] reported a new approach for the calculation of pressure in nanoporous materials using Hill’s thermodynamics for small systems [20]. In this approach, two different pressures are needed to account for confinement effects in nanoporous materials: the differential pressure P, and the integral pressure p ^ . In an ensemble of small systems, the differential pressure times the volume change of the small systems is equal to the work exerted on the surroundings by the volume change. The differential pressure corresponds to the macroscopic pressure and does not depend on the size of the system. The addition of a small system to an ensemble of small systems exerts work on the surroundings which is equal to the integral pressure times the volume of the added small system. The differential and integral pressures are different for small systems and become equal in the thermodynamic limit. For a system with volume V, the two pressures are related by [11,20]:
P ( V ) = ( p ^ ( V ) V ) V T , μ = p ^ ( V ) + V ( p ^ ( V ) ) V T , μ .
As shown in Equation (1), the two pressures are different only when the integral pressure p ^ depends on the volume of the system. Galteland et al. [11] performed equilibrium and non-equilibrium molecular dynamics simulations of Lennard–Jones (LJ) fluids in a face-centered lattice of spherical grains representing a porous medium. Using Hill’s thermodynamics of small systems [20] and the additive property of the grand potential, Galteland et al. [11] defines the compressional energy ( p ^ V ) of the representative elementary volume (REV) in the simulations as follows:
p ^ V = p ^ f V f + p ^ r V r γ ^ fr ω fr ,
where p ^ is the integral pressure of volume V, p ^ f and V f are the integral pressure and volume of the fluid, p ^ r and V r are the integral pressure and volume of the grain particles, and the γ ^ fr and ω fr are the integral surface tension and surface area between the fluid and grain particles. Based on the obtained results, it was concluded that the definition of two pressures is needed to calculate the pressure of the fluid in nanoporous medium.
In this study, the relation between the differential and integral pressure is investigated by performing Monte Carlo (MC) simulations of LJ fluids. The simulations are carried out in a modified Gibbs ensemble, using two simulation boxes in equilibrium with each other. One box represents the bulk fluid, while the other simulation box represents the nanoconfined system, including walls that interact with the fluid particles, as shown in Figure 1. The effect of confinement on the integral pressure is investigated by considering different fluid-wall interaction strengths and pore geometries, namely, a cylinder and a slit pore. To investigate the relation between the differential and integral pressure, the difference of the two pressures, P p ^ , and the ratio of driving forces for mass transport, d p ^ d P , are computed. In Section 2, the devised ensemble and the equations used to calculate the pressure and energy of the system are presented. In Section 3, a rigorous derivation of the used expression to compute the ratio of driving forces is shown. In Section 4, the results for the different differential pressures and pore geometries are shown. In Section 5, our conclusions are summarized.

2. Simulation Details

All MC simulations are carried out using an in-house simulation code. The MC simulations consist of two simulation domains, Simulation Box 1 and Simulation box 2 (see Figure 1). Throughout the manuscript the terms Simulation Box 1 and 2 are used to refer to the two domains, however, simulation domain 2 does not correspond to an actual box. The total number of particles in the system, N T , is fixed and particles can be exchanged between the simulation boxes. Box 1 is cubic and has periodic boundary conditions imposed in all directions. Box 1 represents a bulk fluid. The differential pressure, P, and temperature, T, in Simulation Box 1 are imposed, while the volume, V 1 , and the number of particles, N 1 , can fluctuate. Simulation Box 2 is a cylinder or a slit pore with a fixed volume V 2 . The size of the cylinder and the slit pore are defined by the radius, R, and by the distance between the two parallel planes, 2 R , respectively (Figure 1). In Box 2, periodic boundary conditions are applied only in the axial direction for the cylinder and in the x and y directions for the slit pore (see Figure 1). Box 2 represents the confined fluid which has an integral pressure p ^ . In Simulation Box 2, the volume, V 2 , and temperature T are imposed, while the number of particles N 2 can fluctuate by exchanging particles with Box 1. The instantaneous integral pressure fluctuates, and by definition [20] its ensemble average, p ^ , will be equal to P only for macroscopic systems ( R ). The ensemble used is a variation of the N P T -Gibbs ensemble [21]. The main difference between the ensemble used in this work and the conventional N P T -Gibbs ensemble is that in our simulations the volume of Box 2 is fixed [21]. Essentially, Box 2 corresponds to the grand canonical ensemble with the reservoir explicitly modeled in Box 1. This computational setup was also used in other studies, i.e., A. Z. Panagiotopoulos et al. [21], P. Bai et al. [22], etc.
In the MC simulations, three types of trial moves are used: translation, volume change, and particle exchange. The translation and particle exchange trial moves are used in both simulation boxes, while volume change trial moves are only performed in Simulation Box 1. The acceptance rules of the trial moves can be found elsewhere [23]. The particle exchange trial move, ensures that Box 1 and Box 2 are in chemical equilibrium, i.e., the chemical potentials of the two boxes are equal ( μ 1 = μ 2 ). The chemical potentials of the boxes are defined as the sum of the ideal and excess chemical potentials of the fluid in the respective box ( μ 1 = μ 1 id + μ 1 ex , μ 2 = μ 2 id + μ 2 ex ). The ideal gas chemical potentials ( μ 1 id , μ 2 id ) are calculated based on the density and temperature of the fluid. The excess chemical potential ( μ 1 ex , μ 2 ex ) can be calculated using different methods, i.e., Widom’s test particle insertion method [23], Contionous Fractional Component Monte Carlo method [24,25], Bennet acceptance ratio method [26], etc. Although the Bennet acceptance ratio method is computationally more efficient than the Widom’s test particle insertion method, the Widom’s test particle method is sufficient for the systems considered in this study. Separate simulations are carried out using different cylinder radii R, while imposing P = 0.2 (in reduced units) in Box 1. In Figure 2, the chemical potential of Box 1 and Box 2 (Figure 2a), as well as the density of Box 1 and Box 2 (Figure 2b), are shown as a function of the cylinder radius at P = 0.2 . As can be seen from Figure 2a, the chemical potentials of Box 1 and Box 2 are equal within the uncertainties of the simulations.
In all simulations, the total potential energy, U, is calculated using the 12-6 LJ interaction potential:
U = U fluid fluid + U fluid wall ,
where U fluid fluid is the potential energy due to interaction between the fluid particles, and U fluid wall represents the potential energy contribution from the interactions between the fluid particles and the wall of Box 2. U fluid fluid in both simulation boxes is calculated according to:
U fluid fluid = i < j 4 ϵ fluid σ fluid r i j 12 σ fluid r i j 6 U shift r i j < 2.5 σ fluid 0 , otherwise ,
where r i j is the distance of particle i and j, U shift makes the interaction potentially continuous at the cut-off distance ( 2.5 σ fluid ), and ϵ fluid , σ fluid are the LJ parameters. In the past, several studies were reported using different types of interaction potentials to model the fluid-solid interactions in confined spaces [22,27,28]. Since the aim of this study is to show the difference between the differential and integral pressures and not to simulate some specific adsorption system, only two types of interaction potentials are considered for the interaction of fluid particles with the wall in Simulation Box 2. In the first case, the wall has only repulsive interactions with the fluid particles. The potential energy contribution of this type is calculated based on the Weeks–Chandler–Andersen potential [29]:
U fluid wall = i = 0 N 2 4 ϵ fw σ fw r w i 12 σ fw r w i 6 + 1 r w i < 2 1 / 6 σ fw 0 , otherwise ,
where N 2 is the number of particles in Box 2, r w i is the closest distance of particle i from the walls, and ϵ fw and σ fw are the LJ parameters for the interaction between the wall and the fluid particles. In the second case, the attractive interactions between the wall and fluid particles are taken into account using the traditional form of the 12-6 LJ interaction potential:
U fluid wall = i = 0 N 2 4 ϵ fw σ fw r w i 12 σ fw r w i 6 U shift fw r w i < 2.5 σ fw 0 , otherwise ,
where U shift fw makes the interaction potential continuous at the cut-off distance, and ϵ fw and σ fw are the LJ parameters.
The expression for calculating the integral pressure in Box 2 is as follows [30]:
p ^ = N 2 V 2 k B T + p ^ fluid fluid + p ^ fluid wall ,
where p ^ fluid fluid represents the contribution of the fluid-fluid interaction to the integral pressure, k B is the Boltzmann constant, and p ^ fluid wall represents the contribution of the fluid-wall interaction to the integral pressure. The p ^ fluid fluid , and p ^ fluid wall terms represent the virial contributions of the integral pressure. The terms are calculated based on the virial theorem [7], i.e., using the derivative of the potential energy function with respect to r. The pressure term p ^ fluid fluid is calculated as follows [7]:
p ^ fluid fluid = 48 3 V 2 ϵ fluid i < j σ fluid r i j 12 0.5 σ fluid r i j 6 .
The third term in Equation (7) represents the pressure contribution related to the interactions of the LJ particles with the wall of Box 2. The term p ^ fluid wall for the repulsive wall potential is derived based on Equation (5) [7]. The following expression is obtained:
p ^ fluid wall = i = 1 N 2 24 3 V 2 ϵ fw σ fw r w i 12 0.5 σ fw r w i 6 r w i < 2 1 / 6 σ fw 0 , otherwise .
The following expression is used to calculate the p ^ fluid wall when the wall has also attractive interactions with the fluid [7]:
p ^ fluid wall = i = 1 N 2 24 3 V 2 ϵ fw σ fw r w i 12 0.5 σ fw r w i 6 r w i < 2.5 σ fw 0 , otherwise .
By comparing Equations (8)–(10), it can observed that the multiplication factor 48 in Equation (8) is replaced by the factor 24 in Equations (9) and (10). This difference means that only 50 % of the fluid-wall interactions are taken into account at the calculation of the p ^ fluid wall term [7]. In the MC simulations, Equations (7)–(10) yield instantaneous values from which ensemble averages are computed, i.e., p ^ .
In all simulations, the LJ parameters of the fluid particles are ϵ fluid = 1 , σ fluid = 1 , and the cut-off radius is r cut = 2.5 σ fluid for the fluid-fluid interactions. Regardless of the type of interaction potential used to calculate the interaction of the fluid particles and the wall in Simulation Box 2, the LJ parameter, σ fw = 1 is used. In the case of the purely repulsive wall potential, the LJ parameter ϵ fw = 1 is used. In case of the attractive wall potential, five different values for the ϵ fw LJ parameters are considered: ϵ fw = 0.3 , 0.5 , 0.7 , 1.0 , 1.5 . In this study, all of the reported parameters are in dimensionless units. The LJ interactions are truncated and shifted (i.e., no tail corrections are applied). To avoid phase transitions, the temperature is fixed at T = 2 [31,32].

3. Theory

To investigate the difference between the differential and integral pressures, the ratio d p ^ d P is calculated. The term d p ^ d P is the ratio of the pressure gradient for mass transport in the bulk phase, d P d L , and in the confined space d p ^ d L (see Figure 3). Essentially, the ratio of driving forces, d p ^ d P , equals the ratio of transport coefficients when either P or p ^ is used as driving force for mass transport in the corresponding transport equation. d p ^ d P is referred to as the ratio of driving forces throughout this work. One possible approach to compute the ratio of driving forces is to perform simulations at different imposed differential pressures and calculate the difference in the differential and integral pressures. To avoid the necessity of performing several simulations, in this study the ratio of driving forces is calculated based on the fluctuation theory. Using this approach, the ratio of driving forces can be obtained by performing a single simulation.
To obtain an expression for d p ^ d P , the partition function of the system is needed [23]:
Q = C N 1 = 0 N T V 2 N T N 1 N 1 ! ( N T N 1 ) ! 0 d V 1 V 1 N 1 e β P V 1 d r N T e β U ,
where C is a constant, β = 1 k B T , N 1 is the number of particles in Box 1, N T is total number of particles in the simulation, V 1 is the volume of Box 1, V 2 is the volume of Box 2, and U is the potential energy. The ensemble average of a thermodynamic property X can be obtained using:
X = N 1 = 0 N T V 2 ( N T N 1 ) N 1 ! ( N T N 1 ) ! 0 d V 1 V 1 N 1 e β P V 1 d r N T e β U X N 1 = 0 N T V 2 ( N T N 1 ) N 1 ! ( N T N 1 ) ! 0 d V 1 V 1 N 1 e β P V 1 d r N T e β U .
Therefore, to obtain the expression for d p ^ d P , the following relation is used:
d p ^ d P = d d P N 1 = 0 N T V 2 ( N T N 1 ) N 1 ! ( N T N 1 ) ! 0 d V 1 V 1 N 1 e β P V 1 d r N T e β U p ^ N 1 = 0 N T V 2 ( N T N 1 ) N 1 ! ( N T N 1 ) ! 0 d V 1 V 1 N 1 e β P V 1 d r N T e β U .
By switching the order of the integration and differentiation, the following expression is obtained:
d p ^ d P = N 1 = 0 N T 0 d V 1 d r N T V 2 ( N T N 1 ) V 1 N 1 N 1 ! ( N T N 1 ) ! ( β V 1 ) e β P V 1 e β U p ^ N 1 = 0 N T 0 d V 1 d r N T V 2 ( N T N 1 ) V 1 N 1 N 1 ! ( N T N 1 ) ! e β P V 1 e β U 2 × N 1 = 0 N T 0 d V 1 d r N T V 2 ( N T N 1 ) V 1 N 1 N 1 ! ( N T N 1 ) ! e β P V 1 e β U N 1 = 0 N T 0 d V 1 d r V 2 ( N T N 1 ) V 1 N 1 N 1 ! ( N T N 1 ) ! ( β V 1 ) e β P V 1 e β U N 1 = 0 N T 0 d V 1 d r N T V 2 ( N T N 1 ) V 1 N 1 N 1 ! ( N T N 1 ) ! e β P V 1 e β U 2 × N 1 = 0 N T 0 d V 1 d r N T V 2 ( N T N 1 ) V 1 N 1 N 1 ! ( N T N 1 ) ! e β P V 1 e β U p ^ = β V 1 p ^ ( β V 1 p ^ ) = β ( V 1 p ^ V 1 p ^ ) = V 1 p ^ V 1 p ^ k b T .
The final expression in Equation (14) is essentially the cross correlation between V 1 and p ^ . In this study, this expression is used to calculate the ratio of driving forces.

4. Results and Discussion

In this work, the difference between the differential and integral pressure, P p ^ , and the ratio of driving forces d p ^ d P , are investigated. Two different pore geometries, a cylinder and a slit pore, are studied with varying fluid-wall interaction potentials. The effect of confinement is investigated for gas ( ρ 0.1 ) and liquid ( ρ = 0.58 , 0.8 ) phases, corresponding to P = 0.2, 2.0, 6.0, respectively.

4.1. Difference between the Differential and Integral Pressure

In Figure 4, the difference between the differential, P, and the ensemble average of the integral pressure, p ^ , is shown as a function of the inverse radius, R 1 , of Box 2 for the cylindrical and slit pore cases for the two types of wall potentials. As can be seen in Figure 4, as R 1 decreases, the difference between the differential and integral pressure decreases in all cases. The data are fitted to A R 1 +B, where A and B are constants. The coefficient of determination of the fitted lines is above 0.99 showing that the relation between the R 1 and P p ^ is indeed linear. For large radii ( R > 30 σ ), where the fluid in Box 2 behaves like in the bulk, P p ^ approaches 0, which is also indicated by the fitted lines. In Figure 4a,c,e, the difference between the two pressures are shown for the cylindrical pore. In Figure 4b,d,f, the difference in the pressures are shown for the slit pore. By comparing the magnitude of P p ^ for the cylindrical and slit pore cases, it can be seen that the pressure difference is larger in the cylindrical pore than in the slit pore. The larger value of P p ^ for the cylindrical pores can be attributed to the stronger confinement effects compared to the slit pore. The stronger confinement effects are also indicated by the steeper slopes (constant A in the fitted lines) of the cylindrical pore compared to the slit pore. It can also be observed that by increasing the interaction strength between the wall of Box 2 and the fluid, ϵ wf , the difference between the differential and integral pressure at the same pore size decreases. By comparing the calculated values of P p ^ for P = 0.2 and P = 6.0, it can be seen that the effect of the interaction strength between the wall in Box 2 and the fluid, ϵ wf , is considerably larger at the lower differential pressure. For example, in the case of the slit pore, the ratio of slopes of the fitted line for ϵ wf = 1.0 and ϵ wf = 1.5 at P = 0.2 is 0.88, as shown in Figure 4b, and at P = 6.0 is 0.96, as shown in Figure 4f. This can be caused by the different number of particles inside the pore at the two differential pressures. In the case of P = 0.2 , by increasing the interaction strength, more particles can enter the pore which results in larger integral pressures, while at P = 6.0 , the pore is practically saturated for all ϵ wf ; therefore, the interaction strength has a lower effect on the integral pressure.
Based on the work of Galteland et al. [11], the slope of the fitted line can also be related to the effective surface tension between the fluid particles and the wall in Box 2, i.e., P p ^ γ ^ effective fr R . The linear relation between P p ^ and R 1 shows that the effective surface tension does not depend on the curvature of the wall. In Figure 4, it can be observed that the effective surface tension decreases as the fluid-wall interactions become more attractive. In Figure 4a–f, it can be seen that the ratio of the effective surface tensions (slopes of the fitted lines) for the same ϵ wf with the slit and cylindrical pore is nearly constant. For example, at ϵ wf = 0.5 and P = 2.0 , the ratio of the slope of the fitted lines for the slit and cylindrical pore is 0.54 , as shown in Figure 4c,d, and at ϵ wf = 0.3 and P = 0.2 , the ratio of the two slopes is also 0.54 , as shown in Figure 4a,b. The ratio of the effective surface tensions between the slit and cylindrical pores are in the range of ∼0.52–0.56 and considered constant since it is within the uncertainties of simulations. The constant ratio of effective surface tension between the slit and cylindrical pore also indicates the larger confinement effects in the cylindrical pore.

4.2. Ratio of Driving Forces

In Figure 5, the ratio of driving forces for mass transport, d p ^ d P , is shown as a function of the inverse radius of Box 2, R 1 , at three different pressures, P = 0.2 , 2.0 , 6.0 , for both the slit and cylindrical pores. From Figure 5, it can be observed that d p ^ d P is considerably smaller than 1 for small pore sizes. This means that a change in the differential pressure, P, of Box 1 results in a smaller change in the integral pressure, p ^ , of Box 2. This difference underlines the effect of the confinement on the pressure of the fluid in nanopores and shows the difference in driving forces in the bulk and confined fluid. As R increases, d p ^ d P approaches 1, i.e., the fluid in the pore behaves more like a bulk fluid. As can be seen in Figure 5, d p ^ d P is larger for the slit pore than for the cylindrical pore at the same conditions. This means that in case of the slit pore a change in the differential pressure, P, results in a larger difference in the integral pressure, p ^ , than for the cylindrical pore. The larger change in the integral pressure indicates that the confinement effects are weaker in the slit pore. In Figure 5a,b, it can be observed that by increasing the interaction strength between the wall of the pore and the fluid particles, the ratio of the driving forces decreases. At higher differential pressures the decrease in the ratio of driving forces due to the increasing fluid-wall interaction strength becomes less pronounced, as shown in Figure 5c–f. The smaller influence of the interaction strength, ϵ fw , on the ratio of driving forces may be caused by the increasing contribution of the fluid-fluid interactions to the integral pressure due to the larger number of fluid particles in the pore at higher differential pressure. In Figure 5a–f, the ratio of driving force is shown for the slit and cylindrical pores at the same differential pressure. It can be seen that the ratio of the slopes of the fitted lines for the same value of ϵ wf with the slit and cylindrical pore is constant within the uncertainties of the simulations. The ratio of the slopes is in the range of ∼0.52–0.58, which indicates that the confinement effects in the cylindrical pore are almost twice as strong as in the slit pore.

5. Conclusions

In this study, the new approach reported by Galteland et al. [11] is used to investigate the effects of confinement on a fluid in a nanopore by performing MC simulations. The simulations are carried out in a variation of the Gibbs ensemble with two simulation boxes in chemical equilibrium. One of the simulation boxes represents the bulk fluid with differential pressure P, and the other a slit or cylindrical pore with repulsive or attractive wall interaction potential. In case of the attractive wall potentials, several scenarios are considered for the strength of the interaction between the wall and the fluid particles. The effect of confinement is investigated for three differential pressures, P = 0.2, 2.0, 6.0, corresponding to gas ( ρ 0.1 ) and liquid phases ( ρ = 0.58 , 0.8 ). It is concluded that the difference between the differential and integral pressure P p ^ , for all studied cases, approaches 0 when R→∞. It is shown that the increase in the interaction strength between the wall and the fluid particles has smaller effect on the difference in the pressures, P p ^ , as the differential pressure increases. Based on the work of Galteland et al. [11], the difference of the differential and integral pressure is related to the effective surface tension between the fluid particles and wall of the pore. It is shown that the effective surface tension does not depend on the curvature of the wall. It is found that by considering a bulk fluid in equilibrium with a cylindrical or slit nanopore, the ratio of driving forces for mass transport in the bulk phase is larger than in the nanopore ( d p ^ d P < < 1 ) for small pore sizes. As R increases, d p ^ d P approaches 1, i.e., the fluid in the pore behaves more like a bulk fluid. This clearly shows that the approximation that p ^ P does not hold on the nanoscale.

Author Contributions

T.J.H.V. and M.E. wrote the software, M.E. and O.G. carried out the simulations. S.K., D.B., and O.G. contributed to the conceptual design and the development of the work. The project was supervised by O.A.M., S.K., D.B., and T.J.H.V. All authors contributed to the writing of the manuscript and the scientific discussions. All authors have read and agree to the published version of the manuscript.

Funding

This research was funded by the Center of Excellence funding scheme of the Research Council of Norway, project no 262644, PoreLab. This work was also sponsored by NWO Exacte Wetenschappen (Physical Sciences) for the use of supercomputer facilities, with financial support from the Nederlandse Organisatie voor Wetenschappelijk Onderzoek (Netherlands Organisation for Scientific Research, NWO).

Acknowledgments

T.J.H.V. acknowledges NWO-CW (Chemical Sciences) for a VICI grant.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Gubbins, K.E.; Long, Y.; Śliwinska Bartkowiak, M. Thermodynamics of confined nano-phases. J. Chem. Thermodyn. 2014, 74, 169–183. [Google Scholar] [CrossRef]
  2. Furukawa, H.; Cordova, K.; O’Keeffe, M.; Yaghi, O. The chemistry and applications of metal-organic frameworks. Science 2013, 341, 1230444. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  3. Gu, Z.Y.; Yang, C.X.; Chang, N.; Yan, X.P. Metal-organic frameworks for analytical chemistry: From sample collection to chromatographic separation. Acc. Chem. Res. 2012, 45, 734–745. [Google Scholar] [CrossRef]
  4. Yilmaz, B.; Müller, U. Catalytic applications of zeolites in chemical industry. Top. Catal. 2009, 52, 888–895. [Google Scholar] [CrossRef] [Green Version]
  5. Glaser, R.; Weitkamp, J. The application of zeolites in catalysis. Springer Ser. Chem. Phys. 2004, 75, 159–212. [Google Scholar]
  6. Erdős, M.; de Lange, M.F.; Kapteijn, F.; Moultos, O.A.; Vlugt, T.J.H. In Silico Screening of Metal–Organic Frameworks for Adsorption-Driven Heat Pumps and Chillers. ACS Appl. Mater. Interfaces 2018, 10, 27074–27087. [Google Scholar] [CrossRef] [Green Version]
  7. Irving, J.H.; Kirkwood, J.G. The Statistical Mechanical Theory of Transport Processes. IV. The Equations of Hydrodynamics. J. Chem. Phys. 1950, 18, 817–829. [Google Scholar] [CrossRef]
  8. Férey, G.; Serre, C. Large breathing effects in three-dimensional porous hybrid matter: Facts, analyses, rules and consequences. Chem. Soc. Rev. 2009, 38, 1380–1399. [Google Scholar] [CrossRef]
  9. Song, H.; Yu, M.; Zhu, W.; Wu, P.; Lou, Y.; Wang, Y.; Killough, J. Numerical investigation of gas flow rate in shale gas reservoirs with nanoporous media. Int. J. Heat Mass Transf. 2015, 80, 626–635. [Google Scholar] [CrossRef]
  10. Huber, P. Soft matter in hard confinement: Phase transition thermodynamics, structure, texture, diffusion and flow in nanoporous media. J. Phys. Condens. Matter 2015, 27, 103102. [Google Scholar] [CrossRef] [Green Version]
  11. Galteland, O.; Bedeaux, D.; Hafskjold, B.; Kjelstrup, S. Pressures Inside a Nano-Porous Medium. The Case of a Single Phase Fluid. Front. Phys. 2019, 7, 60. [Google Scholar] [CrossRef] [Green Version]
  12. Todd, B.D.; Evans, D.J.; Daivis, P.J. Pressure tensor for inhomogeneous fluids. Phys. Rev. E 1995, 52, 1627–1638. [Google Scholar] [CrossRef]
  13. Ikeshoji, T.; Hafskjold, B.; Furuholt, H. Molecular-level Calculation Scheme for Pressure in Inhomogeneous Systems of Flat and Spherical Layers. Mol. Simul. 2003, 29, 101–109. [Google Scholar] [CrossRef]
  14. Walton, J.; Tildesley, D.; Rowlinson, J.; Henderson, J. The pressure tensor at the planar surface of a liquid. Mol. Phys. 1983, 48, 1357–1368. [Google Scholar] [CrossRef]
  15. Blokhuis, E.M.; Bedeaux, D. Pressure tensor of a spherical interface. J. Chem. Phys. 1992, 97, 3576–3586. [Google Scholar] [CrossRef] [Green Version]
  16. Wang, G.M.; Sevick, E.M.; Mittag, E.; Searles, D.J.; Evans, D.J. Experimental Demonstration of Violations of the Second Law of Thermodynamics for Small Systems and Short Time Scales. Phys. Rev. Lett. 2002, 89, 050601. [Google Scholar] [CrossRef] [Green Version]
  17. Nilson, R.H.; Griffiths, S.K. Influence of atomistic physics on electro-osmotic flow: An analysis based on density functional theory. J. Chem. Phys. 2006, 125, 164510. [Google Scholar] [CrossRef]
  18. Lee, J.W.; Nilson, R.H.; Templeton, J.A.; Griffiths, S.K.; Kung, A.; Wong, B.M. Comparison of Molecular Dynamics with Classical Density Functional and Poisson–Boltzmann Theories of the Electric Double Layer in Nanochannels. J. Chem. Theory Comput. 2012, 8, 2012–2022. [Google Scholar] [CrossRef]
  19. Gjennestad, M.A.; Wilhelmsen, Ø. Thermodynamic stability of droplets, bubbles and thick films in open and closed pores. Fluid Phase Equilibria 2020, 505, 112351. [Google Scholar] [CrossRef]
  20. Hill, T.L. Thermodynamics of Small Systems, 1st ed.; Dover: New York, NY, USA, 1964. [Google Scholar]
  21. Panagiotopoulos, A.Z. Direct determination of phase coexistence properties of fluids by Monte Carlo simulation in a new ensemble. Mol. Phys. 1987, 61, 813–826. [Google Scholar] [CrossRef]
  22. Bai, P.; Siepmann, J.I. Selective adsorption from dilute solutions: Gibbs ensemble Monte Carlo simulations. Fluid Phase Equilibria 2013, 351, 1–6. [Google Scholar] [CrossRef]
  23. Frenkel, D.; Smit, B. Understanding Molecular Simulation, 2nd ed.; Academic Press: London, UK, 2001. [Google Scholar]
  24. Shi, W.; Maginn, E.J. Continuous Fractional Component Monte Carlo: An Adaptive Biasing Method for Open System Atomistic Simulations. J. Chem. Theory Comput. 2007, 3, 1451–1463. [Google Scholar] [CrossRef] [PubMed]
  25. Poursaeidesfahani, A.; Torres-Knoop, A.; Dubbeldam, D.; Vlugt, T.J.H. Direct Free Energy Calculation in the Continuous Fractional Component Gibbs Ensemble. J. Chem. Theory Comput. 2016, 12, 1481–1490. [Google Scholar] [CrossRef] [PubMed]
  26. Bennett, C.H. Efficient estimation of free energy differences from Monte Carlo data. J. Comput. Phys. 1976, 22, 245–268. [Google Scholar] [CrossRef]
  27. Stecki, J. Steele (10-4-3) Potential due to a Solid Wall. Langmuir 1997, 13, 597–598. [Google Scholar] [CrossRef]
  28. Jiménez-Serratos, G.; Cárdenas, H.; Müller, E.A. Extension of the effective solid-fluid Steele potential for Mie force fields. Mol. Phys. 2019, 117, 3840–3851. [Google Scholar] [CrossRef]
  29. Weeks, J.D.; Chandler, D.; Andersen, H.C. Role of Repulsive Forces in Determining the Equilibrium Structure of Simple Liquids. J. Chem. Phys. 1971, 54, 5237–5247. [Google Scholar] [CrossRef]
  30. Allen, M.P.; Tildesley, D.J. Computer Simulation of Liquids, 2nd ed.; Oxford University Press: Oxford, UK, 2017. [Google Scholar]
  31. Smit, B. Phase diagrams of Lennard–Jones fluids. J. Chem. Phys. 1992, 96, 8639–8640. [Google Scholar] [CrossRef] [Green Version]
  32. Potoff, J.J.; Panagiotopoulos, A.Z. Critical point and phase behavior of the pure fluid and a Lennard–Jones mixture. J. Chem. Phys. 1998, 109, 10914–10920. [Google Scholar] [CrossRef]
Figure 1. Schematic representation of the Monte Carlo (MC) simulation scheme. Simulation Box 1 represents the bulk fluid with differential pressure P. Simulation Box 2 contains the confined fluid with (average) integral pressure p ^ . (a,b) The two investigated systems are shown where the bulk fluid is in equilibrium with the nanoconfined fluid in a cylinder and in a slit pore, respectively. Due to the particle exchange, the chemical potential of the two boxes are equal, but in general P p ^ .
Figure 1. Schematic representation of the Monte Carlo (MC) simulation scheme. Simulation Box 1 represents the bulk fluid with differential pressure P. Simulation Box 2 contains the confined fluid with (average) integral pressure p ^ . (a,b) The two investigated systems are shown where the bulk fluid is in equilibrium with the nanoconfined fluid in a cylinder and in a slit pore, respectively. Due to the particle exchange, the chemical potential of the two boxes are equal, but in general P p ^ .
Nanomaterials 10 00293 g001
Figure 2. The chemical potential (a) and density (b) of the bulk, Box 1, and confined fluids, Box 2, as a function of the cylinder radius, R, at P = 0.2 (in Simulation Box 1). The blue and red colors represent Box 1 and Box 2, respectively. The temperature is fixed at T = 2 . All values are presented in reduced units. The error bars are smaller than the symbol sizes.
Figure 2. The chemical potential (a) and density (b) of the bulk, Box 1, and confined fluids, Box 2, as a function of the cylinder radius, R, at P = 0.2 (in Simulation Box 1). The blue and red colors represent Box 1 and Box 2, respectively. The temperature is fixed at T = 2 . All values are presented in reduced units. The error bars are smaller than the symbol sizes.
Nanomaterials 10 00293 g002
Figure 3. Schematic representation of a bulk fluid in equilibrium with a nanoconfined fluid in a pore. The concept of the ratio of driving forces for mass transport, d p ^ d P , can be introduced based on the definition of ratio of driving forces in the two systems, d p ^ d L and d P d L .
Figure 3. Schematic representation of a bulk fluid in equilibrium with a nanoconfined fluid in a pore. The concept of the ratio of driving forces for mass transport, d p ^ d P , can be introduced based on the definition of ratio of driving forces in the two systems, d p ^ d L and d P d L .
Nanomaterials 10 00293 g003
Figure 4. The difference between the differential, P, and the ensemble average of integral pressure, p ^ , is shown as a function of the inverse radius R 1 of Box 2 at P = 0.2, 2.0, and 6.0 for cylindrical and slit pores with fluid-wall interactions. (a,c,e) The pressure difference is shown for cylindrical pores at differential pressure P = 0.2 , 2.0, and 6.0, respectively. (b,d,f) The pressure difference is shown for slit pores at differential pressure P = 0.2 , 2.0, and 6.0, respectively. The simulation results are shown with symbols, while the lines are fits to the data points. The equation used for the fitting is A R 1 + B , where A and B are constants. The colors denote the different level of attraction between the wall of Box 2 and the fluid, repulsive wall potential (black), ϵ wf = 0.3 (blue), 0.5 (red), 0.7 (gray), 1.0 (orange), and 1.5 (cyan). The results for the cylindrical pore are shown with closed circles and for the slit pores with open rectangles. The temperature of both boxes is set to T = 2 . The average densities of Box 1 are ρ 0.10, 0.58, and 0.8 at P = 0.2, 2.0, and 6.0, respectively.
Figure 4. The difference between the differential, P, and the ensemble average of integral pressure, p ^ , is shown as a function of the inverse radius R 1 of Box 2 at P = 0.2, 2.0, and 6.0 for cylindrical and slit pores with fluid-wall interactions. (a,c,e) The pressure difference is shown for cylindrical pores at differential pressure P = 0.2 , 2.0, and 6.0, respectively. (b,d,f) The pressure difference is shown for slit pores at differential pressure P = 0.2 , 2.0, and 6.0, respectively. The simulation results are shown with symbols, while the lines are fits to the data points. The equation used for the fitting is A R 1 + B , where A and B are constants. The colors denote the different level of attraction between the wall of Box 2 and the fluid, repulsive wall potential (black), ϵ wf = 0.3 (blue), 0.5 (red), 0.7 (gray), 1.0 (orange), and 1.5 (cyan). The results for the cylindrical pore are shown with closed circles and for the slit pores with open rectangles. The temperature of both boxes is set to T = 2 . The average densities of Box 1 are ρ 0.10, 0.58, and 0.8 at P = 0.2, 2.0, and 6.0, respectively.
Nanomaterials 10 00293 g004
Figure 5. The ratio of driving forces, d p ^ d P , is shown as a function of the inverse radius, R 1 , of Box 2 at P = 0.2 , 2.0 , and 6.0 for the cylindrical and slit pore cases with repulsive and attractive wall potentials. (a,c,e) d p ^ d P is shown for cylindrical pores at differential pressure P = 0.2 , 2.0, and 6.0, respectively. (b,d,f) d p ^ d P is shown for slit pores at differential pressures P = 0.2 , 2.0, and 6.0, respectively. The simulation results are shown with symbols, while the lines are fits to the data points. The equation used for the fitting is A R 1 + B , where A and B are constants. The results for the cylindrical pore are shown with closed circles and for the slit pores with open rectangles. The colors denote the different wall potential used in Box 2, repulsive wall potential (black), ϵ wf = 0.3 (blue), ϵ wf = 0.5 (red), ϵ wf = 0.7 (gray), ϵ wf = 1.0 (orange), and ϵ wf = 1.5 (cyan). The temperature of both boxes is set to T = 2 . The average densities of Box 1 are ρ 0.10, 0.58, and 0.8 at P = 0.2, 2.0, and 6.0, respectively.
Figure 5. The ratio of driving forces, d p ^ d P , is shown as a function of the inverse radius, R 1 , of Box 2 at P = 0.2 , 2.0 , and 6.0 for the cylindrical and slit pore cases with repulsive and attractive wall potentials. (a,c,e) d p ^ d P is shown for cylindrical pores at differential pressure P = 0.2 , 2.0, and 6.0, respectively. (b,d,f) d p ^ d P is shown for slit pores at differential pressures P = 0.2 , 2.0, and 6.0, respectively. The simulation results are shown with symbols, while the lines are fits to the data points. The equation used for the fitting is A R 1 + B , where A and B are constants. The results for the cylindrical pore are shown with closed circles and for the slit pores with open rectangles. The colors denote the different wall potential used in Box 2, repulsive wall potential (black), ϵ wf = 0.3 (blue), ϵ wf = 0.5 (red), ϵ wf = 0.7 (gray), ϵ wf = 1.0 (orange), and ϵ wf = 1.5 (cyan). The temperature of both boxes is set to T = 2 . The average densities of Box 1 are ρ 0.10, 0.58, and 0.8 at P = 0.2, 2.0, and 6.0, respectively.
Nanomaterials 10 00293 g005

Share and Cite

MDPI and ACS Style

Erdős, M.; Galteland, O.; Bedeaux, D.; Kjelstrup, S.; Moultos, O.A.; Vlugt, T.J.H. Gibbs Ensemble Monte Carlo Simulation of Fluids in Confinement: Relation between the Differential and Integral Pressures. Nanomaterials 2020, 10, 293. https://doi.org/10.3390/nano10020293

AMA Style

Erdős M, Galteland O, Bedeaux D, Kjelstrup S, Moultos OA, Vlugt TJH. Gibbs Ensemble Monte Carlo Simulation of Fluids in Confinement: Relation between the Differential and Integral Pressures. Nanomaterials. 2020; 10(2):293. https://doi.org/10.3390/nano10020293

Chicago/Turabian Style

Erdős, Máté, Olav Galteland, Dick Bedeaux, Signe Kjelstrup, Othonas A. Moultos, and Thijs J. H. Vlugt. 2020. "Gibbs Ensemble Monte Carlo Simulation of Fluids in Confinement: Relation between the Differential and Integral Pressures" Nanomaterials 10, no. 2: 293. https://doi.org/10.3390/nano10020293

APA Style

Erdős, M., Galteland, O., Bedeaux, D., Kjelstrup, S., Moultos, O. A., & Vlugt, T. J. H. (2020). Gibbs Ensemble Monte Carlo Simulation of Fluids in Confinement: Relation between the Differential and Integral Pressures. Nanomaterials, 10(2), 293. https://doi.org/10.3390/nano10020293

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop