Dynamic Mean Field Theory for Lattice Gas Models of Fluid Mixtures

Oct 8, 2013 - We present the extension of dynamic mean field theory (DMFT) for fluids in porous materials (Monson, P. A. J. Chem. Phys. 2008, 128, 084...
2 downloads 0 Views 4MB Size
Article pubs.acs.org/Langmuir

Dynamic Mean Field Theory for Lattice Gas Models of Fluid Mixtures Confined in Mesoporous Materials J. R. Edison†,‡ and P. A. Monson*,† †

Department of Chemical Engineering, University of Massachusetts, 159 Goessmann Laboratory, Amherst, Massachusetts 01003, United States S Supporting Information *

ABSTRACT: We present the extension of dynamic mean field theory (DMFT) for fluids in porous materials (Monson, P. A. J. Chem. Phys. 2008, 128, 084701) to the case of mixtures. The theory can be used to describe the relaxation processes in the approach to equilibrium or metastable equilibrium states for fluids in pores after a change in the bulk pressure or composition. It is especially useful for studying systems where there are capillary condensation or evaporation transitions. Nucleation processes associated with these transitions are emergent features of the theory and can be visualized via the time dependence of the density distribution and composition distribution in the system. For mixtures an important component of the dynamics is relaxation of the composition distribution in the system, especially in the neighborhood of vapor−liquid interfaces. We consider two different types of mixtures, modeling hydrocarbon adsorption in carbon-like slit pores. We first present results on bulk phase equilibria of the mixtures and then the equilibrium (stable/metastable) behavior of these mixtures in a finite slit pore and an inkbottle pore. We then use DMFT to describe the evolution of the density and composition in the pore in the approach to equilibrium after changing the state of the bulk fluid via composition or pressure changes.



INTRODUCTION Recently, we developed a dynamic mean field theory (DMFT), or mean field kinetic theory, of relaxation processes for fluids confined in porous materials.1−4 This theory is based on the lattice gas model for confined fluids5 and can be viewed as a mean field approximation to the dynamics averaged over an ensemble of trajectories obtained from dynamic Monte Carlo simulations of the system.6−8 It has its origins9 in theories of dynamics of Ising10 and binary alloy11 models as well as theories of diffusion.12−14 The primary application has been to the case of a fluid in a mesoporous material in contact with a bulk gas where DMFT allows the study of the evolution of the density distribution in the system in response to a step change in the bulk chemical potential.1 This is the counterpart of a dynamic uptake measurement in adsorption experiments. The theory described the dynamics in situations where changes in the chemical potential lead to capillary phase transitions and, notably, the nucleation mechanisms are emergent features of the dynamics.2 Another important feature of DMFT is that in the long time limit it reaches a state which is a solution of the mean field classical density functional equations for the density distribution.5,15 It is closely related to dynamic density functional and Cahn−Hilliard type theories of the dynamics of inhomogeneous systems. In this paper we describe the formulation and application of DMFT to the case of fluid mixtures. For mixtures the range of phenomena is greatly expanded by the competitive interactions of the components © 2013 American Chemical Society

with the pore walls and the fact that equilibration of the composition distribution is an important component of the relaxation dynamics. Our work can be viewed as building on previous studies of mixture adsorption using classical density functional theory (DFT).16−30 Several issues have been explored in this work. One focus has been on understanding the solution thermodynamics for adsorbed mixtures and, in particular, departures from the ideal adsorbed solution theory.17,21 Another emphasis has been on the effect of confinement upon the composition relative to the bulk phase (selective adsorption) due to confinement.18,19,25,27,29,30 There has also been a significant effort in understanding how confinement influences vapor− liquid18,23,24 and liquid−liquid20,22,28 equilibrium. As an illustration of the overall approach we present applications to two model mixtures with parameters appropriate to hydrocarbon mixtures in slit pores. One is close to an ideal mixture qualitatively representing ethane and methane, while the other is more nonideal with a wider difference in the volatility of the two components. We calculated the bulk vapor−liquid phase diagrams for the models and studied both the static and the dynamic behavior of the mixtures confined in two pore geometries, a slit pore and an inkbottle pore. We look Received: August 8, 2013 Revised: October 4, 2013 Published: October 8, 2013 13808

dx.doi.org/10.1021/la4030537 | Langmuir 2013, 29, 13808−13820

Langmuir

Article

where z is the coordination number of the lattice, M is the number of sites in the lattice, vs is the physical volume associated with a lattice site, and ρ(α) b is the density of species α in the bulk fluid mixture. Our interest lies primarily in studying states of the confined fluid at bulk conditions in the approach to the dew point from the undersaturated bulk vapor mixture. This makes it essential to know the phase diagram of the bulk fluid mixture, and we compute this using eqs 5 and 6 and the conditions for bulk vapor−liquid equilibrium at fixed temperature = μ(α) (Pl = Pv ; μ(α) l v ∀ α). Dynamic Behavior. Our presentation of the theory closely follows the work of Gouyet et al.9 and previous work on dynamics of pure fluids.1 We begin by defining the ensemble average density ρ(α) i (t) of species α at site i and at time t as

at the behavior of the system as the bulk vapor state is changed via isothermal, isobaric changes in composition or by isothermal changes in bulk pressure at fixed bulk composition. We pay particular attention to the evolution of the species density distributions in the system and especially the vapor− liquid interfaces encountered in the dynamics inside the pores. For the inkbottle pore we focus attention to the cavitation process on desorption.2,31−34 In the next section we present the mean field theory for binary mixture lattice gases and extension of DMFT to binary mixtures. Following that we present our studies for the two different binary fluid mixtures, including the bulk phase diagram, the static behavior of the confined mixtures, and the relaxation dynamics predicted by DMFT. Finally, we present a summary of our results and conclusions.



ρi(α)(t ) =

∑ ∑ ∑ ∑ εαδni(α)ni(+δ)a + ∑ ∑ ni(α)ϕi(α) i

a

α

δ

(1)

α

i

is the occupancy of site i (0 or 1) by species α and a is a where vector pointing to the nearest neighbors of site i. The second term in the Hamiltonian describes the external field created by the solid surface, ϕ(α) i , acting on the fluid species α (δ also denotes a species) at site i. The grand free energy of the system within the mean field approximation can be written as n(α) i

Ω = kT ∑ [∑ ρi(α)ln ρi(α) + (1 − α

i

1 − 2

∑∑∑∑ i

a

α

∂ρi(α) ∂t

+

∑∑

δ

δ

(9)

is the Metropolis transition probability for species α to jump from site i to site j and can be expressed as wij(α) = wo(α) exp(− E ij(α)/kT )

⎞ ⎛ ⎜ ∂Ω ⎟ = 0 ∀ α, i ⎜ (α) ⎟ ⎝ ∂ρi ⎠{μ},{ϕ}, T

(10)

where (3)

where {μ} and {ϕ} denote the sets of chemical potentials and surface fields for the components. This leads to

E ij(α)

∑ ρi(δ))] − ∑ ∑ εαδρi(+δ)a + ϕi(α) − μ(α) δ

a

a

For a given pore geometry, the above set of coupled nonlinear equations is solved at fixed chemical potentials of the individual species and at a fixed temperature to yield the species density distributions of the system at equilibrium (or metastable equilibrium). We used simple iteration to solve these equations starting from either a uniform low density for low pressure or a uniform high density at high pressure. The initial composition is taken as the bulk composition. Subsequent calculations were initiated from the density distributions obtained from previously studied neighboring states. For the bulk (uniform) limit we have the following equations for the pressure and chemical potentials

α

z 2

δ

δ

=

wij(α)ρi(α)(1



∑ δ

⎢ ρj(δ))⎢1

(13) where

δ

∑ ρb(δ))] − z ∑ εαδρb(δ)

⎛ μ (α) − μ (α) ⎞⎤ j i ⎟⎥ − exp⎜ ⎜ ⎟⎥ kT ⎝ ⎠⎦ ⎣



Jij(α)(t )

μi(α) = kT ln ρi(α) − kT ln(1 −

(5)

μ(α) = kT[ln ρb(α) − ln(1 −

(12)

δ

These are the extension to mixtures of the expression for pure fluids.35 w(α) is the ratio of the hopping frequency of species α to the 0 coordination number of the lattice. For simplicity we assume w0 to be the same for all species, but that assumption can be relaxed. Putting together eqs 8−12 results in the following expression for flux of species α

∑ ∑ εαδρb(α)ρb(δ) α

(11)

E i(α) = − ∑ ∑ εαδρi(+δ)a + ϕi(α)

(4)

∑ ρb(α)) −

⎧ Ej(α) < E i(α) ⎪0 =⎨ ⎪ Ej(α) − E i(α) Ej(α) > E i(α) ⎩

and

δ

= 0 ∀ α, i

Ω = − kT ln(1 − M

∑ ρj(δ)) − wji(α)ρj(α)(1 − ∑ ρi(δ)) δ

The necessary condition for a minimum in Ω is

Pvs = −

(8)

a

w(α) ij

(2)

kT[ln ρi(α) − ln(1 −

= − ∑ Ji(,αi +) a

Jij(α)(t ) = wij(α)ρi(α)(1 −

− μ(α))

α

i

(7)

is the net flux of species α between site i and j and can be where written within the mean field approximation as

α

ρi(α)(ϕi(α)

{n}

J(α) i,j

∑ ρi(α))ln(1 − ∑ ρi(α))] α

εαδρi(α)ρi(+δ)a

∑ ni(α)P({n}, t )

where P({n},t) is the probability of observing a given set of occupation numbers {n} that describe the state of the system. The ensemble envisaged here is an ensemble of independent trajectories of the lattice model dynamics, and P({n},t) is averaged over that ensemble. We choose Kawasaki dynamics or vacancy hopping dynamics to describe the microscopic dynamics of the system. The starting point of the dynamic theory is the master equation that describes Kawasaki dynamics of a binary mixture on a lattice gas. The evolution of the local density of any species is obtained by applying a mean field approximation to the master equation. This results in a conservation equation for each species of the form

THEORY

1 2

= t

Static Behavior. Our model is a single-occupancy lattice gas with nearest neighbor interactions in the presence of an external field. Any site i on the lattice can be occupied by a molecule of species α or remain vacant, and the Hamiltonian for this model is given by5,15

H=−

n i(α)

∑ ρi(δ)) − ∑ ∑ εαδρi(+δ)a δ

+

(6) 13809

ϕi(α)

a

δ

(14) dx.doi.org/10.1021/la4030537 | Langmuir 2013, 29, 13808−13820

Langmuir

Article

This latter quantity can be identified with the chemical potential in the limit where the flux vanishes, the system is at equilibrium (stable or metastable), and the chemical potential is everywhere uniform. Implementation of DMFT is similar to that described in our earlier papers.1,4 For a given pore geometry we begin by solving the static mean field equations at a certain initial state given by chemical potentials of the species. This gives the initial density distribution in the system. Next, we fix the density and chemical potentials of the individual species at the boundaries of the pore geometry to a value associated with the final state of the system. The system is then evolved by numerical solution of eq 8. In our initial studies of this approach we used Euler’s method, but we also explored using Runge− Kutta methods. We found Euler’s method to be of acceptable accuracy for time steps less than about w0Δt = 0.2. For instance, this leads to density profiles in the long time limit which match those from solving the status MFT equations. In a recent paper36 we compared the predictions of DMFT with results from dynamic Monte Carlo simulations for capillary condensation dynamics of single components in slit pores. The effect of fluctuations, which are missing in DMFT, leads to a distribution of nucleation times and locations rather than the sharp values from DMFT. Nevertheless, there is good qualitative agreement between the two methods, and we anticipate this to carry over to the mixture case also.

Figure 1. Bulk phase diagram (P−x−y) of mixture I at T* = 0.8 from MFT.

components A and B are at a subcritical state with the P−x or the bubble point curve almost linearly connecting the vapor pressures of pure fluid A and pure fluid B. As mentioned earlier, we are interested in states of the confined fluid which lie below the dew point curve in the phase diagram. Often a porous material is characterized by measuring the fluid uptake versus relative pressure at a fixed temperature, also called the adsorption isotherm. With binary mixtures we have an extra variable, yA, to characterize the bulk composition. Initially we can think of two different protocols to investigate fluid mixtures in pores: (i) fix the bulk composition and study fluid uptake versus relative pressure at a fixed temperature or (ii) fix the bulk pressure and study fluid uptake versus bulk composition at a fixed temperature. The former is the counterpart of an adsorption isotherm measured for pure fluids in pores, while the latter can be viewed as an isothermal stepwise displacement process. In both cases we study states of the confined fluid for which the bulk fluid mixture is in the vapor phase. Slit Pore. In Figure 2 we show the isotherms of total density ρA + ρB in a slit pore of width H = 6 at different bulk compositions varying from yA = 0.1 to yA = 0.5. The paths in the bulk phase diagram along which the isotherms are computed are shown in Figure 2a. All isotherms have a step in the low-pressure region associated with formation of a monolayer as well as a capillary condensation transition with an associated hysteresis loop closer to the dew point pressure. As discussed earlier,5 for an open slit or cylindrical-type pore the step on the desorption branch of the isotherm coincides with the equilibrium vapor−liquid transition for the system. All metastability is on the adsorption branch. The hysteresis arises from the nucleation barrier to forming a liquid bridge during the capillary condensation process. On desorption the vapor and liquid are in contact at the pore ends and there is no nucleation barrier to evaporation of the liquid. The width of the hysteresis loop shrinks as we reduce the bulk mole fraction of the less volatile component A, indicative of the approach to a critical point for the confined mixture. At a given T* the pure fluid of species B (more volatile component) is at a higher reduced temperature (T/Tc) than the pure fluid of species A (less volatile component). Hence, moving along the composition axis from yA = 1.0 to yA = 0.0 has an effect on the width of hysteresis loops similar to the approach to a capillary critical point by increasing the temperature for pure fluids. Figure 3 shows an illustrative example where we fix the bulk pressure at P* = Pvs/εAA = 0.036 and compute the isotherm by varying the bulk composition. The path in the bulk phase diagram along which the isotherm is computed is shown in Figure 3a. The entire range of composition considered above



SYSTEMS STUDIED We consider some illustrative applications to binary mixtures. The set of interaction parameters for the bulk fluid are εAA, εBB, and εAB, with the interaction between unlike species given by the geometric mean εAB = (εAAεBB)1/2. We present studies for two mixtures. The first was chosen so that the bulk vapor− liquid equilibria (VLE) qualitatively resembles those in a fairly ideal mixture such as ethane and methane. The second mixture was chosen to illustrate the effects of having more nonideality with a greater difference in volatility between the components. Each mixture is confined in a pore with completely wetting pore walls with a nearest neighbor interaction with the walls given by εAS = 3.0εAA. We fix εBS based on the relationship εAS/ εBS = (εAA/εBB)1/2. The set of interaction parameters used is shown in Table 1. The interaction parameters were chosen to Table 1. Interaction Parameters of the Two Different Binary Mixtures Studied system

εBB/εAA

εAS/εAA

εBS/εAA

mixture I mixture II

0.6 0.25

3.0 3.0

2.32 1.5

be a qualitative model for alkanes in a carbon-like slit pore. All quantities with dimensions of energy are made dimensionless with the interaction energy εAA. We consider two idealized pore geometries: the slit pore1 and the inkbottle pore.2 The slit pore that we consider has a height of H = 6 lattice units and a length L = 40 lattice units. The inkbottle pore that we consider has a bottle part with height H = 12, length L = 20, and it is connected to necks on either side, whose dimensions are H = 4,L = 10. We call this system the 12−4 inkbottle pore.33 Both pore systems that we consider are placed in contact with bulk fluid as shown in our earlier work.1



RESULTS AND DISCUSSION Mixture I: Static Behavior. Figure 1 shows the bulk VLE P−x−y phase diagram of mixture I at T* = kT/εAA = 0.8 from MFT. The bulk phase behavior resembles an ideal mixture where the VLE could be described by Raoult’s law. Both pure 13810

dx.doi.org/10.1021/la4030537 | Langmuir 2013, 29, 13808−13820

Langmuir

Article

Figure 2. (a) Paths in the bulk phase diagram for mixture I along which adsorption/desorption isotherms are computed. xA refers to liquid and yA to vapor. (b) Isotherms of total density (ρA + ρB) versus relative pressure for mixture I at fixed values of bulk composition at T* = 0.8 for the paths shown in a from MFT for the slit pore with L = 40, H = 6.

Figure 3. (a) Paths in the bulk- phase diagram for mixture I along which adsorption/desorption isotherms are computed. (b) Isotherms of total density (ρA + ρB) versus composition for mixture I at fixed bulk pressure of P = 0.036 for the paths shown in a from MFT for a L = 40, H = 6 slit pore.

Figure 4. (a) Path in the bulk phase diagram for mixture I along which the adsorption/desorption isotherm is computed. (b) Isotherms of total density (ρA + ρB) versus relative pressure in a 12−4 inkbottle pore for mixture I at fixed values of bulk composition at a temperature of T* = 0.8 for the paths shown in a.

the bulk fluid is in the vapor phase. A capillary transition occurs at the composition for which P* = 0.036 is the stability limit (within MFT) of the confined vapor phase. We also observe hysteresis associated with this capillary transition. Inkbottle Pore. The inkbottle pore is a simple model of connectivity effects in a pore network where the middle section, called the bottle, can access the bulk fluid only via narrow constrictions placed on either side called necks. This geometry has been used to illustrate the different mechanisms of desorption in pore networks, namely, pore blocking and cavitation.5,31,32 Pore blocking is a scenario in which the desorption of the condensed phase in the bottle is delayed until the necks empty. Cavitation is a phenomenon in which the bottle can desorb while the necks remain filled, and this happens when the fluid in the bottle reaches its stability limit. For a given adsorbate the factors that determine the mechanism of desorption are the sizes of the necks and bottle and the temperature. With increasing temperature the stability limits of

the liquid in the bottle and the equilibrium pressure (desorption pressure) of the necks change. At a given temperature along the desorption branch of the isotherm if the stability limit of the liquid in the bottle is reached ahead of the equilibrium pressure of the necks then cavitation is observed. It has been shown32,33 that for a pure fluid on increasing temperature the desorption mechanism changes from pore blocking to cavitation. As described earlier for mixture I, approaching a capillary critical point in a slit pore by varying the bulk composition has similar effects as varying the temperature of a pure fluid. Hence, if we study adsorption isotherms of mixture I at various bulk compositions we expect to see a transition in the desorption mechanism from pore blocking to cavitation. This is indeed the case for mixture I, and we show in Figure 4 adsorption isotherms at different bulk compositions for the 12−4 inkbottle pore. For bulk compositions yA = 0.5 and 0.3 the isotherms show a single-step desorption, where the 13811

dx.doi.org/10.1021/la4030537 | Langmuir 2013, 29, 13808−13820

Langmuir

Article

Figure 5. (a) Paths in the bulk phase diagram for mixture I along which adsorption/desorption isotherms are computed. (b) Isotherms of total density (ρA + ρB) versus composition for mixture I in a 12−4 inkbottle pore at pressures of (i) P* = 0.02 (black), (ii) P* = 0.036 (blue), and (iii) P* = 0.05 (red) for the paths shown in a.

Figure 6. (a) Evolution of the total density (black), density of species A (blue), and density of species B (red) in a L = 40, H = 6 slit pore caused by a step change in bulk pressure from P* = 0.0001 to P* = 0.036 at a temperature of T* = 0.8 with yA = 0.5. (Inset) Section of the uptake curve between w0t = 45 000 and w0t = 80 000. (b) Evolution of the density (blue full line) and mole fraction (blue dashed line) of species A in the pore.

the surface of the pore by species A (see Figure 7c and 7d). This effect is seen in binary adsorption dynamics experiments and is referred to as “roll-over”.37 Also interesting is the long time behavior. In the inset of Figure 6a we plot the final stages of the uptake between w0t = 45 000 and w0t = 80 000. The total density (black line) saturates at about w0t = 45 000. However, the component densities continue to vary significantly, indicating that the last phase of the uptake process involves compositional equilibration, where relative amounts of the species change until the chemical potentials of the species in the bulk match with the pore. We also show the evolution of the density (full line) and the mole fraction of A (dashed line) in Figure 6b. The mole fraction of A shows an initial increase until the monolayer forms and then reaches a plateau. It remains almost flat until the undulate state occurs and then continues to increase with formation of the liquid bridge until equilibration. Visualizations of the density distributions of species A (top) and species B (bottom) at different stages of evolution are shown in Figure 7. Visualizations are made with RASMOL,38,39 and the density of a certain species in a lattice site is represented by 10 shades of gray. Visualizations show formation of a monolayer filled mostly with species A during the initial stages (see Figure 7d). Later undulates appear near the ends of the system (see Figure 7e), followed by formation of two liquid bridges close to both ends of the pore (see Figure 7f). Liquid bridges widen to fill the pore completely. Visually the pore filling behavior of this mixture appears similar to the filling by a pure fluid.1 Pore Filling at Constant Pressure. We can also study the dynamics of pore filling at constant pressure by making a step change to the bulk composition. We fix the bulk pressure at P* = 0.036 and temperature at T* = 0.8 and initialize the system at a initial bulk composition of yA = 0.01. The density vs bulk

necks and bottle empty in one step, characteristic of pore blocking. The isotherm for composition yA = 0.1 shows a twostep desorption or a cavitation transition in which the bottle empties while the necks remain filled. In Figure 5 we show density vs composition at fixed bulk pressure for the 12−4 inkbottle pore. At the lowest pressure considered (P* = 0.02) we see a two-step adsorption and a single-step desorption characteristic of pore blocking. At the highest pressures considered P* = 0.05 we see a two-step desorption transition which is characteristic of cavitation. On desorption as we make isothermal decrements of the independent variable (yA or P) we will observe cavitation if the condensed phase in the bottle reaches its stability limit ahead of the phase equilibrium pressure of the neck. Mixture I: Dynamic Behavior in a Slit Pore. Pore Filling at Fixed Bulk Composition. We now present the dynamic behavior associated with pore filling in a slit pore of height H = 6 and length L = 40. The bulk composition of the mixture is fixed at yA = 0.5, and the temperature is fixed at T* = 0.80 (for bulk phase diagram see Figure 1). The slit pore is initialized at a dilute state (P* = 0.0001 or P/P0 = 0.00269970), and then a step change is made to the bulk pressure to a value of P* = 0.036. We can see from the isotherm shown in Figure 2b that at P* = 0.036 (P/P0 = 0.97189336) the fluid in the pore is in a condensed state. We integrate the DMFT equations to follow the relaxation of the system to this final equilibrium state. In Figure 6a we plot the evolution of the total density and density of the individual species. The first cusp in the total density curve is associated with formation of undulates in the pore, and the second cusp is associated with formation of liquid bridges. We note two interesting observations in the uptake behavior. At short times during formation of the monolayer on the surface of the pore we can observe a peak in the density of component B in the pore. Component B is then displaced from 13812

dx.doi.org/10.1021/la4030537 | Langmuir 2013, 29, 13808−13820

Langmuir

Article

Figure 7. Visualizations of the density distribution of species A (top) and species B (bottom) in a L = 40, H = 6 slit pore caused by a step change in bulk pressure from P* = 0.0001 to P* = 0.036 at a temperature of T* = 0.8 with yA = 0.5 at w0t = (a) 0, (b) 800, (c) 1600, (d) 14 000, (e) 28 000, (f) 32 000, (g) 44 000, and (h) 60 000. We use 10 shades of gray to represent the density at a given lattice site, and the density of species A is visualized at the top and B at the bottom of each panel.

Figure 8. (a) Evolution of the total density (black), density of species A (blue), and density of species B (red) in the pore at a fixed pressure of P* = 0.036 caused by a step change in bulk composition from yA = 0.01 to yA = 0.515 in a L = 40, H = 6 slit pore at a temperature of T* = 0.8 . (b) Evolution of the density (blue full) and composition (blue dashed) of species A in the pore.

composition for the slit pore at P* = 0.036 is shown in Figure 3a. We then make a step change to the bulk composition to yA = 0.515. This changes the chemical potentials of species A and B in the bulk and drives the mixture to a new equilibrium state. The density/composition evolution curves are shown in Figure 8 with visualizations of the pore filling shown in Figure 9. The initial state of the pore seen in Figure 9a has a species B rich monolayer adsorbed at the pore walls. Initial stages of the

uptake involve displacement of B from the monolayer, and this is shown in Figure 9b. Looking at the density evolution curve in Figure 8a, we can see that by w0t = 8000 much of the displacement has occurred. The composition evolution curve shown in Figure 8b suggests that the composition of the pore reaches much closer to its final value once the displacement is complete. The mechanism of filling after this point resembles the previous case. Two liquid bridges are observed close to the 13813

dx.doi.org/10.1021/la4030537 | Langmuir 2013, 29, 13808−13820

Langmuir

Article

Figure 9. Visualizations of the density distribution of species A (top) and species B (bottom) in a L = 40, H = 6 slit pore at pressure P* = 0.036 with initial bulk composition yA = 0.01 and final bulk composition yA = 0.515 at a temperature of T* = 0.8 at w0t = (a) 0, (b) 2000, (c) 8000, (d) 20 000, (e) 28 000, and (f) 44 000.

Figure 10. (a) Evolution of the total density (black), density of species A (blue), and density of species B (red) in a L = 40, H = 6 slit pore caused by a step change in bulk pressure from P* = 0.036 to P* = 0.0001 at a temperature of T* = 0.8 with yA = 0.5. (b) Evolution of the density (blue full) and composition (blue dashed) of species A in the pore.

ends of the pore, and they proceed to fill the pore, as shown in Figure 9c, 9d, and 9e. Also in Figure 9e we can observe an enrichment of species B at the interface between the liquid bridges and the trapped vapor bubble as well as at the pore ends. We shall revisit this point later in the paper. Pore Emptying at Fixed Bulk Composition. Next, we present the dynamics of pore emptying where we initialize the system at a condensed state P = 0.036 and make a step change in the bulk pressure to P = 0.0001 at fixed bulk composition of yA = 0.5. The density/composition evolution curves and visualizations of this process are shown in Figures 10 and 11. As shown in the visualizations the emptying process occurs via the meniscus receding into the pore. The meniscus recedes until it reaches the middle of the pore, where it snaps, leaving a monolayer-like state, and finally the emptying of the monolayers happens. The snapping of the liquid bridge corresponds to the cusp in the density evolution curve. Notice that at this instant the composition evolution curve in Figure 10b is very close to the upper limit of yA = 1.0, indicating that almost all of the less volatile component (B), present

predominantly in the condensed region away from the surface, has emptied. Finally, the adsorbed monolayer rich in species A evaporates into the bulk fluid. The density distributions of the two species are quite different with B tending to concentrate in regions away from the pore walls as the system evolves (see Figure 11b and 11c). Pore Emptying at Fixed Bulk Pressure. We also studied the dynamics of pore emptying at a fixed bulk pressure of P = 0.036. The slit pore is initially at equilibrium with a bulk vapor of composition yA = 0.515. The bulk composition is then changed to yA = 0.01. In other words, we begin with a pore filled with a condensed liquid-like phase and follow the evolution to a final state where a component B rich monolayer covers the surface of the slit pore. The density/composition evolution curves and visualizations of this process are shown in Figures 12 and 13. From the visualizations we can see that emptying proceeds via the liquid meniscus receding into the pore. As the meniscus recedes we see that it gets depleted of A and gets continually enriched with component B near the interface. This continues until about w0t ≈ 7500, following 13814

dx.doi.org/10.1021/la4030537 | Langmuir 2013, 29, 13808−13820

Langmuir

Article

Figure 11. Visualizations of the density distribution of species A (top) and species B (bottom) in a L = 40, H = 6 slit pore caused by a step change in bulk pressure from P* = 0.036 to P* = 0.0001 at a temperature of T* = 0.8 with yA = 0.5 at w0t = (a) 0, (b) 1200, (c) 2800, (d) 8000, (e) 20 000, and (f) 36 000.

Figure 12. (a) Evolution of the total density (black), density of species A (blue), and density of species B (red) in a L = 40, H = 6 slit pore at a fixed pressure of P* = 0.036 for a step change in bulk composition from yA = 0.515 to yA = 0.01 at a temperature of T* = 0.8 . (b) Evolution of the density (blue full line) and composition (blue dashed line) of species A in the pore.

molecular simulations and classical density functional theories, notably by Rowlinson and co-workers,41 Gubbins and coworkers,42 Evans and co-workers,43−45 and Davis and coworkers.46 Later, Oxtoby and co-workers studied nucleation of bubbles in liquid mixtures47,48 and noted that surface enrichment will be observed in the density profiles of the critical nuclei. Mixture I also exhibits surface enrichment under certain conditions as shown in the visualizations of Figure 9. However, the effect is much more pronounced in mixture II due to the vast difference in interaction strengths of the two species. In Figure 15a we plot the vapor−liquid density profile of mixture II at a bulk liquid composition of xA = 0.82 and T* = 1.0. We can clearly see the enrichment of the more volatile component at the interface. This lowers the surface tension of the liquid mixture. The surface tension of the vapor−liquid interface γvl can be computed from the excess grand potential as follows.

which the density of both components continues to decrease until the liquid bridge breaks at about w0t ≈ 11 000. Once the bridge snaps the liquid monolayer in the wall continually gets richer with component B until equilibrium is reached. Unlike the previous case (Figure 10b) the mole fraction of species A in the pore decreases monotonically as the system approaches equilibrium. Mixture II: Static Behavior. In mixture II the pure fluid interaction strengths of species A and B are fixed such that component B is substantially more volatile (see Table 1). Naturally the critical temperatures of the two pure fluids are very different and so are the saturation pressures of the two pure fluids at any fixed subcritical temperature. Figure 14a shows the P−x−y diagram of mixture II at T* = 1.0. At this temperature component B is above its critical temperature, and therefore, the P−x−y diagram is a loop detached from the P* axis. At the vapor−liquid interface of such binary mixtures an enrichment of component B is observed. This was first studied by Plesner et al.40 and has been studied extensively by

γvl = Ω vl + PvsM 13815

(15)

dx.doi.org/10.1021/la4030537 | Langmuir 2013, 29, 13808−13820

Langmuir

Article

Figure 13. Visualizations of the density distribution of species A (top) and species B (bottom) in a L = 40, H = 6 slit pore at pressure P* = 0.036 with initial bulk composition yA = 0.515 and final bulk composition yA = 0.01 at a temperature of T* = 0.8 at w0t = (a) 0, (b) 2000, (c) 5200, (d) 10 000, (e) 12 000, and (f) 32 000.

Figure 14. (a) Bulk phase diagram diagram of mixture II at T* = 1.0. (b) Isotherm of mixture II in a 12−4 inkbottle at a temperature T* = 1.0 and at a bulk composition of yA = 0.5 for the path shown in a.

Figure 15. (a) Density profiles for the vapor−liquid interface (mixture II) at a liquid composition of xA = 0.82 at T* = 1.0. (b) Surface tension of the vapor−liquid interface versus the mole fraction of species A in the liquid phase for mixtures I and II at a temperature T* = 1.0. For both mixtures at this temperature the vapor−liquid equilibrium terminates at a critical point before xA reaches zero.

Here Ωvl is the grand free energy of the system with a vapor− liquid interface. In Figure 15b we plot the surface tensions of mixtures I and II versus liquid phase composition at a temperature T* = 1.0. We can see that for mixture II the surface tension drops sharply as the composition of the more volatile component increases in the liquid phase.

Our studies of this mixture focus on the inkbottle pore geometry. The adsorption/desorption isotherm of this system at a bulk composition of yA = 0.5 is shown in Figure 14a. The adsorption branch has two steps in addition to the low-pressure monolayer transition step, and they correspond to filling of the necks and bottle. The desorption branch of the isotherm also 13816

dx.doi.org/10.1021/la4030537 | Langmuir 2013, 29, 13808−13820

Langmuir

Article

Figure 16. Visualizations of the density distribution of species A (top) and species B (bottom) showing the dynamics of pore filling in a 12−4 for a step change in bulk pressure from P* = 0.001 to P* = 0.14 in a 12−4 inkbottle pore at a temperature of T* = 1.0 for yA = 0.5 at w0t = (a) 0, (b) 5000, (c) 60 000, (d) 12 000, (e) 135 000, and (f) 160 000.

then quenched to a state (P* = 0.095, P/P0 = 0.636) where the bottle is empty but the necks remain filled. Visualizations of the states encountered in the dynamics are shown in Figure 17. We can see from Figure 17a that initially the fluid in the bottle is less dense than that in the neck part. As time progresses the fluid becomes more expanded due to loss of molecules by mass transfer (see Figure 17b), and at w0t = 9200 (see Figure 17d) a bubble forms in the middle of the pore. The bubble expands until the bottle is empty. Note that as soon as the bubble forms component B is seen to concentrate at the vapor−liquid interface. Oxtoby and co-workers48 used density functional theory to study bubble nucleation in liquid mixtures at conditions where both species were subcritical in their pure form. They speculated the possibility of cavitation being a twostep process with a liquid−liquid phase separation preceding bubble formation. We see some tendency toward this in Figure 17c.

shows two steps with the bottle desorbing at a higher pressure than the necks via cavitation. In the next section we discuss the dynamics of pore filling and cavitation in this system. Mixture II: Dynamic Behavior. We first present the dynamic behavior associated with pore filling in a 12−4 inkbottle. In the bulk we have an equimolar mixture II at T* = 1.0 (for the bulk phase diagram see Figure 14a). The pore is initialized at a state (P* = 0.001, P/P0 = 0.0067) and then quenched to a state where the bottle is filled with a liquid-like state (P* = 0.14, P/P0 = 0.938; see Figure 14b). Visualizations of the states encountered in the dynamics are shown in Figure 16. We can see from Figure 16a that initially the pore is empty. First, a monolayer is formed in the necks and the bottle (Figure 16b) and necks fill up. This traps a bubble in the bottle. The bottle fills up by squeezing the bubble out. As the bubble shrinks we observe surface enrichment of species B at the vapor−liquid interface of the bubble (Figure 16c and 16d). Finally, the pore is filled with a liquid-like state (Figure 16e), and then compositional equilibration takes place. Finally, we present the dynamic behavior associated with cavitation in a 12−4 inkbottle. In the bulk we have an equimolar mixture at T* = 1.0 (see Figure 14a). The pore is initialized at a condensed state (P* = 0.1, P/P0 = 0.670) and



SUMMARY AND CONCLUSIONS In this work we extended the DMFT for lattice gas models of fluid confined in porous materials to the case of mixtures. The theory describes the evolution of the density and composition distributions inside the pores in response to changes in the bulk 13817

dx.doi.org/10.1021/la4030537 | Langmuir 2013, 29, 13808−13820

Langmuir

Article

Figure 17. Visualizations of the density distribution of species A (top) and species B (bottom) showing the dynamics of cavitation after a step change in bulk pressure from P* = 0.1 to P* = 0.095 in a 12−4 inkbottle pore at a temperature of T* = 1.0 for yA = 0.5 at w0t = (a) 0, (b) 7000, (c) 9000, (d) 9200, (e) 13 000, and (f) 26 000.

phenomenon whereby at short times the more volatile component was able to concentrate at the pore walls due to its higher mobility but at later times was displaced by the less volatile component which interacts more strongly with the pore walls. This is similar to the roll-over seen in binary adsorption in fixed bed adsorbers.37 We studied the dynamics of adsorption and desorption for mixture II for an inkbottle geometry. The larger difference in volatility between the components in this mixture led to a more pronounced concentration of the more volatile component at the vapor− liquid interfaces in the system. Our results also provide evidence for localized liquid−liquid phase separation as a precursor to bubble nucleation in the cavitation processes in these systems.48 These initial results suggest that DMFT can provide useful information about the nanoscale processes in the relaxation dynamics of mixtures confined in porous materials. This will be of considerable interest, for example, in understanding the dynamics of adsorption separation systems. Future applications will include the study of model pore networks and applications to membrane separations. At this stage DMFT is a qualitative tool for mechanistic understanding of the relaxation processes in these systems. In future work we hope to make more

state, e.g., by changes in pressure or composition. To illustrate the results emerging from the theory we considered two model mixtures. The first is close to an ideal mixture, similar to ethane and methane, and the other is a more nonideal system. For the ideal system we presented the equilibrium adsorption isotherms from MFT where we varied the pressure at fixed bulk composition or where we varied the bulk composition at fixed pressure for both slit pore and inkbottle pore geometries. We have seen that the width of the hysteresis at fixed composition depends upon the bulk composition, narrowing as the mole fraction of the more volatile component in the bulk increases. This is a consequence of the decrease of the confined fluid critical temperature with increasing mole fraction of the more volatile component. For the inkbottle pore we see a change in the desorption mechanism from pore blocking to cavitation as the mole fraction of the more volatile component (in the bulk) increases. Using DMFT we studied the dynamics of adsorption and desorption in a slit pore for mixture I for isothermal bulk changes of state involving pressure changes at fixed bulk composition or composition changes at fixed pressure. Nucleation by liquid bridge formation is observed, in common with the behavior for single components. We also observed a 13818

dx.doi.org/10.1021/la4030537 | Langmuir 2013, 29, 13808−13820

Langmuir

Article

(14) Matuszak, D.; Aranovich, G. L.; Donohue, M. D. Thermodynamic driving force for molecular diffusion-lattice density functional theory predictions. J. Non-Equilib. Thermodyn. 2006, 31, 355−384. (15) Kierlik, E.; Monson, P. A.; Rosinberg, M. L.; Sarkisov, L.; Tarjus, G. Capillary condensation in disordered porous materials: Hysteresis versus equilibrium behavior. Phys. Rev. Lett. 2001, 8705, 055701. (16) Hadjiagapiou, I.; Evans, R. Adsorption from a binary fluid mixture. Mol. Phys. 1985, 54, 383−406. (17) Kierlik, E.; Rosinberg, M.; Finn, J.; Monson, P. Binary vapour mixtures adsorbed on a graphite surface: A comparison of mean field density functional theory with results from Monte Carlo simulations. Mol. Phys. 1992, 75, 1435−1454. (18) Tan, Z.; Gubbins, K. E. Selective adsorption of simple mixtures in slit pores: a model of methane-ethane mixtures in carbon. J. Phys. Chem. 1992, 96, 845−854. (19) Sowers, S. L.; Gubbins, K. E. Optimizing Removal of Trace Components from Nitrogen/X Mixtures Using Adsorption: Theory and Simulation. Langmuir 1995, 11, 4758−4764. (20) Kierlik, E.; Fan, Y.; Monson, P. A.; Rosinberg, M. L. Liquidliquid equilibrium in a slit pore: Monte Carlo simulation and mean field density functional theory. J. Chem. Phys. 1995, 102, 3712−3719. (21) Sweatman, M.; Quirke, N. Predicting the adsorption of gas mixtures: Adsorbed solution theory versus classical density functional theory. Langmuir 2002, 18, 10443−10454. (22) Woywod, D.; Schoen, M. Phase behavior of confined symmetric binary mixtures. Phys. Rev. E 2003, 67, 026122. (23) Bucior, K. Capillary condensation of a model binary mixture in slit-like pores. Colloids Surf., A 2003, 219, 113−124. (24) Bucior, K. Capillary condensation of a model binary mixture in slit-like pores with differently adsorbing walls. Colloids Surf., A 2004, 243, 105−115. (25) Rother, G.; Woywod, D.; Schoen, M.; Findenegg, G. H. Confinement effect on the adsorption from a binary liquid system near liquid/liquid phase separation. J. Chem. Phys. 2004, 120, 11864− 11873. (26) Sweatman, M.; Quirke, N. Gas adsorption in active carbons and the slit-pore model 1: Pure gas adsorption. J. Phys. Chem. B 2005, 109, 10381−10388. (27) Sweatman, M.; Quirke, N. Gas adsorption in active carbons and the slit-pore model 2: Mixture adsorption prediction with DFT and IAST. J. Phys. Chem. B 2005, 109, 10389−10394. (28) Woywod, D.; Schemmel, S.; Rother, G.; Findenegg, G. H.; Schoen, M. Phase behavior and local structure of a binary mixture in pores: Mean-field lattice model calculations for analyzing neutron scattering data. J. Chem. Phys. 2005, 122, 124510. (29) Pini, R.; Ottiger, S.; Storti, G.; Mazzotti, M. Prediction of competitive adsorption on coal by a lattice DFT model. Adsorption 2010, 16, 37−46. (30) Ritter, J. A.; Pan, H.; Balbuena, P. B. Adsorption of Binary Gas Mixtures in Heterogeneous Carbon Predicted by Density Functional Theory: On the Formation of Adsorption Azeotropes. Langmuir 2010, 26, 13968−13975. (31) Sarkisov, L.; Monson, P. A. Modeling of adsorption and desorption in pores of simple geometry using molecular dynamics. Langmuir 2001, 17, 7600−7604. (32) Ravikovitch, P.; Neimark, A. Experimental confirmation of different mechanisms of evaporation from ink-bottle type pores: Equilibrium, pore blocking, and cavitation. Langmuir 2002, 18, 9830− 9837. (33) Libby, B.; Monson, P. A. Adsorption/desorption hysteresis in inkbottle pores: A density functional theory and Monte Carlo simulation study. Langmuir 2004, 20, 4289−4294. (34) Rasmussen, C. J.; Vishnyakov, A.; Thommes, M.; Smarsly, B. M.; Kleitz, F.; Neimark, A. V. Cavitation in Metastable Liquid Nitrogen Confined to Nanoscale Pores. Langmuir 2010, 26, 10147− 10157. (35) Matuszak, D.; Aranovich, G. L.; Donohue, M. D. Lattice density functional theory of molecular diffusion. J. Chem. Phys. 2004, 121, 426−435.

quantitative comparisons with experiments on dynamic uptake for mixtures in porous materials.



ASSOCIATED CONTENT

S Supporting Information *

Tables of data from all graphs in this paper and gnuplot scripts for plotting them. This material is available free of charge via the Internet at http://pubs.acs.org.



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Present Address ‡

Soft Condensed Matter Group, Debye Institute for NanoMaterials Science, Utrecht University, Princetonplein 1, NL3584CC Utrecht, Netherlands. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This material is based upon work supported by the National Science Foundation under grant nos. CBET-0853068 and CBET-1158790.



REFERENCES

(1) Monson, P. A. Mean field kinetic theory for a lattice gas model of fluids confined in porous materials. J. Chem. Phys. 2008, 128, 084701. (2) Monson, P. A. Dynamic mean field theory for fluids confined in porous materials: Application to an ink bottle pore geometry.Characterization of Porous Solids VIII; Royal Society of Chemistry, Cambridge: 2009; pp 103−110. (3) Edison, J. R.; Monson, P. A. Dynamic mean field theory of condensation and evaporation in model pore networks with variations in pore size. Microporous Mesoporous Mater. 2011, 154, 7−15. (4) Edison, J. R.; Monson, P. A. Dynamic mean field theory of condensation and evaporation processes for fluids in porous materials: Application to partial drying and drying. Faraday Discuss. 2010, 146, 167−184. (5) Monson, P. A. Understanding adsorption/desorption hysteresis for fluids in mesoporous materials using simple molecular models and classical density functional theory. Microporous Mesoporous Mater. 2012, 160, 47−66. (6) Leung, K.; Luzar, A. Dynamics of capillary evaporation. II. Free energy barriers. J. Chem. Phys. 2000, 113, 5845−5852. (7) Luzar, A.; Leung, K. Dynamics of capillary evaporation. I. Effect of morphology of hydrophobic surfaces. J. Chem. Phys. 2000, 113, 5836−5844. (8) Woo, H. J.; Monson, P. A. Phase behavior and dynamics of fluids in mesoporous glasses. Phys. Rev. E 2003, 67, 041207. (9) Gouyet, J. F.; Plapp, M.; Dieterich, W.; Maass, P. Description of far-from-equilibrium processes by mean-field lattice gas models. Adv. Phys. 2003, 52, 523−638. (10) Penrose, O. A Mean-Field Equation of Motion for the Dynamic Ising-Model. J. Stat. Phys. 1991, 63, 975−986. (11) Martin, G. Atomic Mobility in Cahn’s Diffusion-Model. Phys. Rev. B 1990, 41, 2279−2283. (12) Matuszak, D.; Aranovich, G. L.; Donohue, M. D. Modeling fluid diffusion using the lattice density functional theory approach: counterdiffusion in an external field. Phys. Chem. Chem. Phys. 2006, 8, 1663−1674. (13) Matuszak, D.; Aranovich, G. L.; Donohue, M. D. Singlecomponent permeation maximum with respect to temperature: A lattice density functional theory study. Ind. Eng. Chem. Res. 2006, 45, 5501−5511. 13819

dx.doi.org/10.1021/la4030537 | Langmuir 2013, 29, 13808−13820

Langmuir

Article

(36) Edison, J. R.; Monson, P. A. Dynamics of capillary condensation in lattice gas models of confined fluids: A comparison of dynamic mean field theory with dynamic Monte Carlo simulations. J. Chem. Phys. 2013, 138, 234709. (37) Shen, C.; Worek, W. Cosorption characteristics of solid adsorbents. Int. J. Heat Mass Transfer 1994, 37, 2123−2129. (38) Sayle, R.; Milner-White, E. J. RasMol: Biomolecular graphics for all. Trends Biochem. Sci. 1995, 20, 374. (39) Bernstein, H. J. Recent changes to RasMol, recombining the variants. Trends Biochem. Sci. 2000, 25, 453−455. (40) Plesner, I. W.; Platz, O.; Christiansen, S. E. StatisticalMechanical Calculation of Surface Properties of Simple Liquids and Liquid Mixtures. II. Mixtures. J. Chem. Phys. 1968, 48, 5364−5370. (41) Chapela, G. A.; Saville, G.; Thompson, S. M.; Rowlinson, J. S. Computer simulation of a gas-liquid surface. Part 1. J. Chem. Soc., Faraday Trans. 2 1977, 73, 1133−1144. (42) Lee, D. J.; Telo da Gama, M. M.; Gubbins, K. E. Adsorption and surface tension reduction at the vapor-liquid interface. J. Phys. Chem. 1985, 89, 1514−1519. (43) Telo da Gama, M. M.; Evans, R. Theory of the Liquid-Vapor Interface of a Binary Mixure of Lennard-Jones Fluids. Mol. Phys. 1980, 41, 1091−1112. (44) Telo da Gama, M. M.; Evans, R. The structure and surface tension of the liquid-vapour interface near the upper critical end point of a binary mixture of Lennard-Jones fluids. Mol. Phys. 1983, 48, 229− 250. (45) Telo da Gama, M. M.; Evans, R. The structure and surface tension of the liquid-vapour interface near the upper critical end point of a binary mixture of Lennard-Jones fluids. Mol. Phys. 1983, 48, 251− 266. (46) Carey, B. S.; Scriven, L. E.; Davis, H. T. Semiempirical theory of surface tension of binary systems. AIChE J. 1980, 26, 705−711. (47) Zeng, X. C.; Oxtoby, D. W. Binary homogeneous nucleation theory for the gas−liquid transition: A nonclassical approach. J. Chem. Phys. 1991, 95, 5940−5947. (48) Talanquer, V.; Oxtoby, D. W. Nucleation of bubbles in binary fluids. J. Chem. Phys. 1995, 102, 2156−2164.

13820

dx.doi.org/10.1021/la4030537 | Langmuir 2013, 29, 13808−13820