2540
J. Phys. Chem. B 2007, 111, 2540-2545
Constrained Equilibrium as a Tool for Characterization of Deformable Porous Media E. V. Vakarin,*,† Yurko Duda,‡ and J. P. Badiali† UMR 7575 LECA ENSCP-UPMC, 11 rue P. et M. Curie, 75231 Cedex 05, Paris, France, and Programa de Ingenierı´a Molecular, Instituto Mexicano del Petro´ leo, 07730 D. F., Me´ xico ReceiVed: October 2, 2006; In Final Form: January 19, 2007
A new method for characterizing the deformable porous materials with noncritical adsorption probes is proposed. The mechanism is based on driving the adsorbate through a sequence of constrained equilibrium states with the insertion isotherms forming a pseudocritical point or a van der Waals-type loop. In the framework of a perturbation theory and Monte Carlo simulations we have found a link between the loop parameters and the host morphology. This allows one to characterize porous matrices through analyzing a shift of the pseudocritical point and a shape of the pseudospinodals.
I. Introduction One of the methods for characterization of porous materials1,2 is a retrieving of pore size distributions through adsorption probes. Capillary condensation or freezing3 in macro- or mesopores permits them to be detected rather reliably because the condensation pressure and temperature are very sensitive to the pore size. A quite similar trend is reported3 for the freezing temperature. Nevertheless, it becomes more and more questionable4 whether the true equilibrium behavior is observed in experiments. One of the major problems is a well developed metastability leading to very long equilibration times. Therefore, the condensation is usually detected as an onset of hysteresis. Nevertheless, it has recently been argued5 that the hysteretic behavior is not necessarily associated with an underlying equilibrium phase transition. In particular, the recent analysis6 of the internal fluid dynamics reveals a quite complex relaxation process. Namely, the states outside the hysteresis region are characterized by an essentially diffusive relaxation, while the slow dynamics inside the hysteresis loop is found to be dominated by activated rearrangements of the guest density. In the case of micropores, the situation is even more difficult because such pores are too small to support any sharp (bulklike) phase transition. This makes one search for alternative characterization mechanisms which are not based on the critical behavior. In the case of deformable matrices,7,8 such as aerogels,9,10 amorphous intercalation compounds,11 or hydrated clays,12,13 the problem is additionally complicated by variations of the host morphology with the guest density or pressure. Effects of the host volume contraction and swelling have recently been analyzed within Monte Carlo14 and density functional15 techniques. It has been found that the matrix flexibility has a remarkable effect on the isotherms, shifting the condensation point from its rigid-host position. Nevertheless, the matrix characterization problem (at least in its traditional formulation) does not arise in these approaches because the host behavior is uniquely determined by a given molecular model of the host-guest coupling. Therefore, the matrix morphology results from the conditions of global (matrix + fluid) equilibrium. † ‡
UMR 7575 LECA ENSCP-UPMC. Instituto Mexicano del Petro´leo.
In reality the host structure is significantly affected by a preparation procedure (e.g., quenching, sol-gel technique, etc.). This results in a distribution of some relevant parameters (local porosities, site energies, etc.). Depending on a coupling to the fluid the distribution may be modified. Thus the guest response (e.g., adsorption isotherms) is a superposition of at least two contributions, the adsorbate thermodynamics under the host confinement and the adsorbate-induced host deformation. Quite often these two modes have different characteristic time scales, such that the isotherm can be represented as an average over the pore size fluctuations. In this paper we propose a method for a pore size determination through adsorption of noncritical fluids. The main idea is to find an explicit relation between the density-dependent matrix morphology and some peculiarity in the fluid thermodynamic response. The mechanism is based on our finding16 of negative compressibility states in confined dilatating systems. These states can appear as a consequence of a driving paths either a sequence of controlled lateral dilatations16 or a sequential adsorption of prescribed portions.17 This imposes an additional constraint18 on the fluid equilibrium. As a result, at certain conditions one can observe a sequence of states which form a van der Waals-type loop in the isotherm. Then the loop parameters can be taken as the required specificities for a detection of the confining geometry. II. Model The matrix is considered as a collection of slitlike pores with different widths hi distributed according to some probability density. For simplicity the pores are assumed to be uncorrelated. This allows us to focus on a single pore of fluctuating width h and surface area S. The width is known only statistically, that is, from a probability distribution f (h). The slit itself is a part of the hosting system, whose coupling to the guest species changes the slit geometry (see below). The fluid is modeled as N hard spheres of diameter σ ) 1, injected into the pore. The hard core potential allows us to mimic the excluded volume effects and ensures that the fluid is off-critical (we do not consider freezing-type transitions). For any h the total Hamiltonian is
H ) Hff + Hfw
10.1021/jp066460h CCC: $37.00 © 2007 American Chemical Society Published on Web 02/21/2007
(1)
Constrained Equilibrium as a Tool
J. Phys. Chem. B, Vol. 111, No. 10, 2007 2541
where Hff is the hard sphere Hamiltonian, Hfw is the slit potential N
Hfw ) A
∑ i)1
[
1
1
+
(h - zi)
zki
k
]
/ 2 e z i e h - 1/ 2
1
;
(2)
We do not take into account a short-ranged fluid-wall attraction, responsible for the surface adsorption or layering effects. The inverse-power shape for Hfw is chosen as a generic form of a soft repulsion, with k controlling the softness. In some sense this mimics the hydrophobic or nonwetting effects that, as we will see later, are essential. On the other hand, this simple model takes into account that the walls are impenetrable and ignores possible localized adsorption effects, which could complicate the overall picture. For technical purposes we are working with k ) 3. Moreover, our results are qualitatively insensitive to a particular choice of k. In principle, one may consider two adsorption mechanisms. One is an equilibrium between the pore and the bulk fluids. In that case N would change according to the chemical potential µ ) µb. In this paper we deal with another mechanism. Namely, an adsorption in doses,17 when N is injected in prescribed portions. In that case we have a sequence of the fluid equilibrium states corresponding to different N. The sequence is determined by an injection procedure, which imposes an additional constraint. It should be emphasized that in choosing the model we tried to retain only a minimal number of relevant ingredients: the absence of liquid-gas transition in the fluid bulk, excluded volume effects, and a density-dependent heterogeneous confinement. III. Perturbation Theory A. Conditional Insertion Isotherm. At any pore width h the free energy can be represented as
βF(h) ) βF0(h) - ln 〈e-βHfw〉0
(3)
where F0(h) is the free energy of a reference system, and 〈...〉0 is the average over the reference state. Taking a spatially confined hard-sphere system (the one with A ) 0) as a reference, we consider a first-order perturbation16 for the conditional free energy
βF(h) ) βF0(h) + βSF(h)
∫ ∏ dziHfw
(4)
i
This approximation is expected to be qualitatively acceptable in the case of rather weak wall-fluid interactions or high temperatures (A/(kT) should be small). In what follows the quality of this approach will be tested against the computer simulation results. Here β ) 1/(kT) and F(h) is the pore density in the “slab” approximation
F(h) )
N S(h - 1)
(5)
ignoring a nonmonotonic behavior of the density profile with increasing pore density. This is reasonable for wide pores and low-fluid densities. The reference part is estimated in the excluded volume approximation, while the perturbation contribution is NΨ(h), such that the total conditional free energy is
βF(h) ) -N ln
[
]
1 - bF(h) ΛT3F(h)
- N + NβΨ(h)
where b ) 2πσ3/3 is the excluded volume factor, ΛT is the thermal de Broglie length and
h Ψ(h) ) 16A (1 - 2h)2
(7)
The conditional insertion isotherm can be calculated in the standard way: βµ(h) )
( ) ∂βF(h) ∂N
) ln
β,S,h
[
ΛT3F(h)
1 - bF(h)
]
+
bF(h) + βΨ(h) (8) 1 - bF(h)
which relates the fluid density F(h) and the chemical potential µ(h). B. Matrix Dilatation. As is discussed above, the pore dilatation can be taken into account, assuming that the width h is density-dependent. Real insertion materials are usually rather complex (multicomponent and heterogeneous). For that reason one usually deals with a distribution of pore sizes or with an average size. In this context we assume that h is known only statistically and the probability distribution f(h|F) is conditional16 to the guest density F, which should be self-consistently found from
F)
∫ dh f (h|F) F(h)
(9)
This relation reflects an injection sequence F(h) and the matrix dilatation, given by f (h|F). In this way the average density F is not uniquely determined by the particles insertion but is also conditional to the matrix response. Namely this constraint makes the problem nontrivial, inducing a coupling between the insertion and the dilatation modes. Taking our result for µ(h), we can focus on the insertion isotherm, averaged over the width fluctuations:
µ(F) )
∫ dh f (h|F) µ(h)
(10)
The matrix characterization problem consists in retrieving the distribution f (h|F) from a model µ(h) and data µ(F), coupled by eq 10. In general, this is an inverse problem that may have no or many solutions. To chose an appropriate, physically relevant distribution, additional conditions are usually imposed. This makes one resort to linear or nonlinear19 regularization procedures. Quite often different procedures lead to different results and a way of avoiding this ambiguity is still an open problem, which is beyond the scope of the present paper. In this work we follow a much simpler strategy, which is well-accepted in the force-field parametrization20 and pore-size determination21 problems. Namely, we assume a functional form of the distribution with a set of parameters controlling its shape. Then the parameters are determined from a fitting of the calculated isotherm µ(F) to that obtained in the simulation. In what follows, we are trying to find the fitting conditions at which our theoretical estimations of these structural parameters approach their values taken in the simulation. Even without resorting to a concrete form for f (h|F), it is clear that the matrix reaction can be manifested as a change in the distribution width and/or the mean value. One of the simplest forms reflecting at least one of these features is a δ-like distribution, ignoring a nonzero width:
(6) f(h|F) ) δ[h - h(F)]
(11)
2542 J. Phys. Chem. B, Vol. 111, No. 10, 2007
Vakarin et al.
where δ(x) is the Dirac δ-function and the mean pore width h(F) is density-dependent. For concreteness we consider the swelling behavior.11 Therefore, the mean pore width h(F) increases with F. For illustration purposes we take the S-shaped form, considered previously7,16
h(F) ) h0 (1 + tanh[∆(F - F0)] + tanh[∆F0])
(12)
This form resembles a non-Vegard behavior, typical for layered intercalation compounds7,8 or hydrated clays.12,13 The dilatation is weak at low densities (F , F0), the most intensive response is at F ≈ F0, and then the pore reaches a saturation, corresponding to its mechanical stability limit. Here ∆ is the matrix response constant or dilatation rate, controlling the slope near F ≈ F0, and h0 is the mean pore width in the absence of insertion. This should be considered as a generic form that samples a nonlinear increase of the average pore width, resolved at a given density F. The model mimics certain heterogeneity aspects which, as we believe, are essential. On the other hand, this choice allows us to map the problem in question onto the negative compressibility problem16 in systems with fluctuating geometry. Our results are qualitatively insensitive to the choice (eq 12). In particular, a power-law dependence of h(F) has been found8 to give quite similar results for the isotherms. From eq 9 the average density is found to be
F)
1 N S h(F) - 1
(13)
By changing the surface density N/S (by increments in N), we vary the average pore density F. This allows us to eliminate N/S in the favor of F in all thermodynamic functions. Combining eqs 8, 9, and 10 we obtain
βµ(F) ) ln
[ ]
ΛT3F h(F) bF + 16βA + 1 - bF 1 - bF (1 - 2h(F))2
(14)
In the case of wide pores and a weak dilatation we expand βµ in terms of 1/h0 and ∆, obtaining a generic van der Waals form:
βµ(F) ) ln
[ ]
ΛT3F 4βA 12βA bF + + F∆ 1 - bF 1 - bF h0 h0
(15)
Therefore, as in the previous study,16 we have an interplay of several effects: the packing (first two terms), the fluidmatrix interaction (density-independent), and the dilatation (linear in density). It is seen that a coupling of the fluid-slit repulsion (A), spatial confinement (h0), and the dilatation (∆) acts as an effective infinite-range fluid-fluid attraction. The latter induces an inflection point (or even a loop) in the isotherm µ(F). Introducing a dimensionless temperature T* ) 1/(4A*), and solving
∂2(βµ) ∂(βµ) ) 0; )0 ∂F ∂F2
(16)
we find the pseudocritical parameters at which the inflection point appears.
Fc )
1 4 ∆ ; T*c) 3b 27 bh0
Note that our result (eq 15) is formally identical to an isotherm of a Lennard-Jones fluid adsorbed in a slit pore,22 with 12βA∆/ h0 playing a role of the attraction coefficient. Following this analogy, it is tempting to associate the inflection with a critical point and to cut the loop according to the stability requirement, which is equivalent to ∂µ(F)/∂F > 0. However, this way of thinking is wrong. In order to demonstrate this let us return to eq 10, calculating the inverse compressibility:
χ(F)-1 )
∂µ(F) ) ∂F
∂µ(h)
∫ dh
-1
+
∂f(h|F) µ(h) (18) ∂F
Here the first term is always positive. It involves the conditional inverse compressibility ∂µ(h)/∂F(h) that should be positive for each h because of the stability requirement. The second term includes the matrix response function R(h|F) ) ∂f (h|F)/∂F associated with the dilatation. Since -∞ < µ(h) < ∞ and the sign of R(h|F) is not fixed, the second term could be negative. Therefore, nothing prevents the average compressibility χ(F) from being negative, at least in some density domain. We thus see that the negative compressibility states do not violate the thermodynamic stability.18 The loop is closely related to the dilatation (nonzero R(h|F)). Our conclusions are quite general. They do not rely upon a specific fluid dynamics or matrix distribution. As has already been mentioned, the loop is not a signature of a phase transition. It is a sequence of equilibrium states, constrained by a coupling between the insertion and the dilatation modes. The physical mechanism is quite similar to that previously reported.16 We deal with a competition of two effects: the interparticle and wall-particle repulsions tend to increase the insertion energetic cost µ with increasing F, while the pore dilatation diminishes it. As a result, at certain parameters the isotherm becomes nonmonotonic. For a purely attractive interaction tail (A < 0), the inflection point does not appear. Therefore, the wall-particle repulsion is an essential ingredient. Nevertheless, the fluid-matrix interactions usually involve attractive parts (like in the frequently used 9-3 potential). Our preliminary estimations reveal that the presence of a softly attractive part does not discard our conclusions, but just restricts the density range where the loop feature can appear. The key factors are the relative magnitude and softness of the repulsion and attraction. This suggests that an appropriate combination of the repulsive and attractive parts in the fluidmatrix interaction is a criterion for choosing the adsorption probe, suitable for a given range of the pore width. As it should be, the inflection (loop) feature disappears in the bulk limit h0 f ∞ or in the case of rigid matrices ∆ f 0. It should be noted that the loop feature considered here is physically different from the hysteresis loops seen in experiments or simulations dealing with wetting fluids. In that case the fluid dynamics is found6 to be dominated by a complex process, involving diffusive relaxations and activated density rearrangements. In our case the loop results from a competition of the fluid-matrix repulsion and its relaxation by the matrix swelling. To improve the accuracy of the theoretical results, the reference part is replaced by the Carnahan-Starling form
(17)
Plugging these back into eq 14 we get the third parameter µc ) µc(b, h0, ∆).
[ ]
∂F ∫ dh f (h|F) ∂F(h) ∂F(h)
βµ ) ln (FΛT3) +
h(F) 8η - 9η2 + 3η3 + 16βA ; 3 (1 - η) (1 - 2h(F))2 η ) πF/6 (19)
Constrained Equilibrium as a Tool
J. Phys. Chem. B, Vol. 111, No. 10, 2007 2543
This allows us to avoid the unphysical behavior at high densities and escape from the flaws of the excluded volume approximation (eq 14). In particular, all the pseudocritical parameters (Fc, T*c, µc) appear to be related to the parameters of the pore width distribution (h0, F0, ∆). Therefore, driving the system close to the pseudocritical point allows an estimation of the density dependent pore size distribution. IV. Simulation In order to verify the existence of the loop predicted by the perturbation theory as well as to test the accuracy of the theoretical estimation, the simulation of the model in question has been carried out applying the Widom test particle insertion method.23 As in our previous work16 we applied the conventional canonical NVT MC simulations of a confined hard-sphere fluid. The simulation cell was parallelepiped in shape, with parallel walls at surface separation h, and constant surface area S ) Lx × Ly ) 12 × 12. The periodic boundary conditions were applied to the X and Y directions of the simulation box; the box length in the Z direction is fixed by the pore width. For a given adsorbed fluid density F, the pore width h and the number of particles N have been chosen according to eqs 5 and 12. The number of fluid particles was varying between 100 and 1600 for low and high densities, respectively. During this NVT MC simulation, a coordinate rN+1 has been generated at frequent intervals, randomly and uniformly over the simulation cell. For this value of rN+1, we have computed 〈exp(-βH(rN+1))〉N averaged over all generated trial positions. In such a manner the chemical potential of the system with density F has been calculated,
[
βµ ) - ln
F 〈exp(-βH(rN +1))〉N
]
(20)
where H is the system Hamiltonian. Each simulation runs 4 × 105 MC cycles, with the first half for the system to reach equilibrium and the second half for evaluating the ensemble averages. In the next section the chemical potential, calculated in the similation, will be compared with the theoretical predictions for the same confined system. V. Results Adsorption isotherms are plotted in Figure 1. It is seen that the inflection point around F ≈ 0.3 transforms into the loop with increasing repulsion A (or decreasing temperature T*). As expected,16 the perturbation theory underestimates the magnitude of A (overestimates T*) at which these effects appear. The reason is that our pseudocritical parameters are determined in a direct analogy to the critical parameters in the usual van der Wallstype equation of state. The latter, as any other approximate theory neglecting the fluctuations (such as mean field), is well known to overestimate the critical temperature. The loop gets sharper with increasing dilatation rate ∆ and becomes more pronounced with decreasing pore width h0. It is remarkable that the fitting on the pseudocritical (inflection) point (part (a)) or on the loop (part (b)) allows us to find the set of the distribution parameters (h0, ∆, F0) in agreement with the simulation. Attempts to fit under other criteria (e.g., the least-square deviations or matching the low-density region) give incorrect values of these parameters. In agreement with our theoretical prediction (eq 17), the simulation results confirm that the pseudocritical point (and the loop) is reachable at any appropriate choice of the model parameters. This is in contrast to what was reported16 for the
Figure 1. Insertion isotherms for a weakly (a) and strongly (b) repulsive slit. The other parameters are h0 ) 10, ∆ ) 15, F0 ) 0.3 (both theory and simulation).
negative compressibility states in the same model with lateral stretching (increasing S while N is fixed). There the loop has been found to appear only if the slit reaction ∆ reaches some threshold value ∆* which involves a combination of h0, F0, and A. In other words, the variation of the average density F due to the transversal dilatation should dominate its variation, induced by the changes in the surface density N/S owing to the lateral stretch. This difference suggests that the driving path (namely, changing the fluid density by a lateral stretch or by a particle injection) is really important for this sort of systems. On the other hand, in the context of the present study, the fact that the pseudo-critical parameters do not involve any additional condition, is quite favorable for making a link between these parameters and the matrix morphology. In Figure 2 the pseudocritical parameters are analyzed. It is seen that Fc grows with the dilatation rate ∆ almost linearly exhibiting a characteristic kink. The latter signals a crossover from a weak to a strong dilatation regime. The kink shifts to lower ∆ and becomes less pronounced with decreasing h0 (that is why it is almost undetectable for h0 ) 5). At low ∆ the pseudocritical density is practically independent of h0, while at high ∆ it approaches F0. In agreement with our simplified estimation (eq 17), T*c(∆) can be reasonably approximated by a straight line (the lower inset). Nevertheless, it exhibits a weak inflection point at the same ∆ (around ∆ ) 4) as Fc does. The pseudocritical temperature increases with decreasing h0. This is also coherent with the estimation (eq 17) predicting a decrease in the slope inversely proportional to h0. In contrast, the insertion energy µc decreases and becomes insensitive to h0 with increasing ∆ (the upper inset). This agrees well with the physical intuition that the pore dilatation diminishes an overlap between the repulsive wall fields, making the insertion less energetically consuming.
2544 J. Phys. Chem. B, Vol. 111, No. 10, 2007
Figure 2. Variation of the pseudocritical parameters with the dilatation rate ∆; F0 ) 0.3; h0 ) 5-squares;, h0 ) 10-circles.
Figure 3. Pseudospinodal curves at different dilatation rates ∆. The other parameters are h0 ) 7 and F0 ) 0.3.
In practice, an exact capture of the inflection point could be too demanding. Driving at rather low temperatures (T < T*c), such that the system develops the loop, should be much easier. In that case, a shape of the pseudospinodals could also be used as an indicator of the matrix morphology. These curves are plotted in Figure 3. It seen that the pseudospinodals become broader with decreasing dilatation rate and shift to higher F with increasing ∆. Comparing this with the previous figure, we can see again that Fc converges to F0 with increasing dilatation rate ∆. If necessary, other quantitative characteristics, such as the curve asymmetry with respect to Fc or µ values at the spinodal points, can also be easily extracted. These results illustrate a correlation between the pseudocritical point (or the loop) parameters and the density-dependent pore geometry. We did not try to obtain these parameters from the simulation because the procedure is time consuming, while new physical insights, compared to the analytical estimations, are difficult to expect. The only promising point is to consider narrow (1 e h0 e 2) pores, for which the first-order perturbative scheme could fail even qualitatively. This point is left for a future analysis. VI. Discussion and Conclusion A method for characterizing the deformable porous materials by adsorption of noncritical fluids is proposed. It is based on a driving the host-guest system through a sequence of constrained equilibrium states, forming a loop in the insertion isotherms. The constraint results from a coupling of the insertion and the
Vakarin et al. matrix dilatation modes, inducing a competition between the repulsive interactions and their relaxation by the host dilatation. As a result the insertion isotherms develop pseudocritical points or loops, whose parameters are shown to be correlated with the parameters of the pore-size distribution. In particular, we have found that, although the pseudocritical temperature is overestimated, the magnitudes of the geometrical parameters (h0, F0, ∆) which reproduce well the inflection point or the loop, perfectly agree with their values in the simulation. This suggests that the peculiar points (with ∂µ/∂F ) 0 or ∂F/∂µ ) 0) are extremely sensitive to the matrix morphology, while the overall isotherm shape is less important (unless one can achieve an exact isotherm coincidence). This is consistent with the recent finding20 that a matching of the inflection points (with ∂F/∂µ ) 0) correctly reproduces the remaining part of the experimental isotherm. From the theoretical point of view19 such a criterion (the inflection matching) corresponds to a maximal amount of information one can get from a fitting of a model to a given experimental or simulation result. In our case we deal with the pseudocritical point (∂µ/∂F ) 0). Fitting on this peculiarity corresponds to a singular information rate19 that allows us to update the geometrical parameters from their initial estimation in a limited range of F (around the singularity). Namely this type of singularity is used to characterize porous media through a shift of the fluid condensation1 or freezing3 point. Although we have taken a simplified δ-like distribution, our conclusions should hold in general because the characteristic loop feature is directly related to the dilatation mode. The latter translates into a shift of the mean pore size, that has been taken into account. Nevertheless, the role of a nonzero distribution width remains to be investigated. The only restriction is that the number of the detected loop parameters should be (at least) equal to the number of parameters determining the pore size distribution. This restriction is not too strong because in practice the distribution is usually characterized by only few parameters (e.g., leading mean pore sizes and dispersions around them). On the other hand, the proposed method is mainly intended for a complementary detection of the porosity range (e.g., micropores) where other techniques fail. For illustration and technical purposes we were working with rather large pores 10 e h(F) e 30. This, however, does not affect our main concept. On the other hand, the matrix has been modeled as a collection of uncorrelated slit pores. In reality the pores are correlated (in size, position, etc.), especially in high-porosity materials, like aerogels. In this respect we have to note that the averaging over a correlated matrix distribution would result in an additional correlation (through the matrix) among the fluid particles. In the simplest two-body case this correlation translates into an effective mean-force interaction for the guest species. This additional pair interaction could modify the pseudocritical parameters or the loop shape, but it should not alter our main conclusions because the proposed mechanism is essentially based on a coupling between the insertion and the dilatation modes. Nevertheless, a correct retrieving of a many-body (correlated) matrix distribution from an isotherm could be quite difficult. A knowledge of the fluid correlation functions might be required. This is, however, a general difficulty in the inverse problems, like the one considered here. Since the main features of our model are the fluid-matrix repulsion and the matrix flexibility, the approach is expected to be the most suitable for deformable materials with repulsive or weakly attractive pores. Hydrated clay minerals12,13 could serve as a prototypic example. These are negatively charged
Constrained Equilibrium as a Tool layered compounds kept together by cations. The water-clay interaction, involving Coulombic and dispersion-repulsion forces, is shown12,13 to reproduce quite well the experimental data on the pressure dependent average interlayer spacing. The latter closely resembles our h(F) dependence (eq 12). In agreement with our results, the computer simulation13 reveals (see Figures 4-6 in ref 13) that the loops (or inflections) in the isotherms appear at the insertion levels, where the interlayer spacing makes smooth steps (at F ≈ F0 in our language). This is explained13 as a consequence of the fact that a dose of water increases the interlaminar distance, promoting the entering of more water molecules. This mechanism is in agreement with our results, revealing a coupling of the insertion and matrix swelling. Note that the recent simulation24 of adsorption into completely rigid clay matrices has not found any loop traces. This is also consistent with our finding, suggesting that the loop can be considered as an indicator of the pore swelling. Concerning the experimental results, it should be noted that in the great majority of adsorption studies the isotherms are obtained from the conditions of the guest fluid equilibrium with a bulk reservoir. The only exception, to our knowledge, is reported in ref 17 where a dosen fluid injection is reported. On the other hand, in the intercalation experiments11 the techniques (such as chronopotentiometry), in which the insertion level (guest density) is controlled and the (electro-) chemical potential is measured, are quite common. We have recently demonstrated19 that for disordered matrices the isotherms, obtained by these two methods, are not the same. Their difference is directly related to a degree of disorder. The simulation results13 also show that the clay structure depends on how the equilibrium state is reached. Therefore, testing a fluid-matrix system along various insertion paths could help in resolving the matrix morphology. This resembles quite closely the adsorptiondesorption hysteresis,6 whose width serves1 as an estimator of the matrix structure. In this light we may expect that the present study might also stimulate new investigations in that domain. It should be emphasized that, since we deal with a constrained equilibrium, the loop does not violate the thermodynamic stability.18 Therefore, this feature can be readily found under suitable experimental conditions, similar to those realized in the low-pressure adsorption into microporous matrices.17 In that study quite similar loop features have been detected17 for nitrogen and argon adsorption in microporous carbons and synthetic zeolites. As in our case, the fluids have been inserted
J. Phys. Chem. B, Vol. 111, No. 10, 2007 2545 in sequential doses. Nevertheless, in contrast to our assumption, there is no evidence for an insertion-induced dilatation of these slit-shaped materials. Therefore, we cannot make a direct comparison to these experimental data, because of a difference in microscopic details. In particular, the loops have been interpreted as a signature of a monolayer-induced micropore filling.17 Nevertheless, there is a formal analogy. Namely, the authors of ref 17 claim that the S-shaped feature (the loop) signals a switching from the ultramicropore to micropore filling with increasing pressure. This is quite similar to the growth of the mean pore size with increasing density in our case. References and Notes (1) Gelb, L. D.; Gubbins, K. E.; Radhakrishnan, R.; SliwinskaBartowiak, M. Rep. Prog. Phys. 1999, 62, 1573. (2) Jaroniec, M.; Brauer, P. Surf. Sci. Rep. 1986, 6, 65. (3) Alba-Simionesco, C.; Coasne, B.; Dosseh, G.; Gubbins, K. E.; Radhakrishnan, R.; Sliwinska-Bartowiak, M. J. Phys.: Condens. Matter 2006, 18, R15. (4) Detcheverry, F.; Kierlik, E.; Rosinberg, M. L.; Tarjus, G. Phys. ReV. E 2005, 72, 051506. (5) Kierlik, E.; Monson, P. A.; Rosinberg, M. L.; Sarkisov, L.; Tarjus, G. Phys. ReV. Lett. 2001, 87, 055701. (6) Valiulin, R.; Naumov, S.; Galvosas, P.; Karger, J.; Woo, H. J.; Porcheron, F.; Monson, P. A. Nature 2006, 443, 965. (7) Vakarin, E. V.; Badiali, J. P.; Levi, M. D.; Aurbach, D. Phys. ReV. B 2001, 63, 014304. Vakarin, E. V.; Badiali, J. P. J. Phys. Chem. B 2002, 106, 7721. (8) Vakarin, E. V.; Badiali, J. P. Solid State Ionics 2004, 171, 261. (9) Escobedo, F. A.; de Pablo, J. J. J. Chem. Phys. 1997, 106, 793. (10) Herman, T.; Day, J.; Beamish, J. Phys. ReV. B 2006, 73, 094127. (11) Garcia-Belmonte, G.; Garcia-Canadas, J.; Bisquert, J. J. Phys. Chem. B 2006, 110, 4514. (12) Odriozola, G.; Aguilar, J. F. J. Chem. Phys. 2005, 123, 174708. (13) Odriozola, G.; Guevara-Rodriguez, F. de J. Langmuir 2004, 20, 2010. (14) Shen, J.; Monson, P. A. Mol. Phys. 2002, 100, 2031. (15) Ustinov, E. A.; Do, D. D. Carbon 2006, 44, 2652. (16) Vakarin, E. V.; Duda, Y.; Badiali, J. P. J. Chem. Phys. 2006, 124, 144515. (17) Amarasekera, G.; Scarlett, M. J.; Mainwaring, D. E. J. Phys. Chem. 1996, 100, 7580. (18) Everett, D. Colloids Surf., A 1998, 141, 279. (19) Vakarin, E. V.; Badiali, J. P. Phys. ReV. E 2006, 74, 036120. Vakarin, E. V.; Badiali, J. P. AIP Conf. Proc. 2006, 872, 411. (20) Dubbeldam, D.; S. Calero, S.; Vlugt, T. J. H.; Krishna, R.; Maesen, T. L. M.; Beerdsen, E.; Smit, B. Phys. ReV. Lett. 2004, 93, 088302. (21) Scaife, S.; Kluson, P.; Quirke, N. J. Phys. Chem. B 2000, 104, 313. (22) Schoen, M.; Diestler, D. J. J. Chem. Phys. 1998, 109, 5596. (23) Widom, B. J. Stat. Phys. 1978, 19, 563. (24) Hensen, E. J. M.; Tambach, T. J.; Bliek, A.; Smit, B. J. Chem. Phys. 2001, 115, 3322.
