High-Pressure Adsorption Capacity and Structure of CO2 in Carbon

We present new simulation results for the packing of single-center and three-center models of carbon dioxide at high pressure in carbon slit pores. Th...
0 downloads 0 Views 516KB Size
9612

Langmuir 2004, 20, 9612-9620

High-Pressure Adsorption Capacity and Structure of CO2 in Carbon Slit Pores: Theory and Simulation S. K. Bhatia,* K. Tran, T. X. Nguyen, and D. Nicholson Department of Chemical Engineering, The University of Queensland, Brisbane, Queensland 4072, Australia Received June 9, 2004. In Final Form: August 19, 2004

We present new simulation results for the packing of single-center and three-center models of carbon dioxide at high pressure in carbon slit pores. The former shows a series of packing transitions that are well described by our density functional theory model developed earlier. In contrast, these transitions are absent for the three-center model. Analysis of the simulation results shows that alternations of flat-lying molecules and rotated molecules can occur as the pore width is increased. The presence or absence of quadrupoles has negligible effect on these high-density structures.

1. Introduction The separation of carbon dioxide from gaseous streams is an important industrial problem, the requirements for which are projected to increase manifold due to tightening controls on greenhouse gas emissions. In particular, the large volumes of CO2 emitted from thermal power plants pose a serious challenge, demanding also measures for sequestration of the captured carbon dioxide. While various alternatives for the latter case are being pursued,1 among the most promising appears to be the high-pressure geological storage of carbon dioxide in deep coal seams, with simultaneous displacement of the inherent methane,2-4 or in abandoned coal mines. Effective process development for such sequestration requires an understanding of the factors affecting carbon dioxide adsorption capacity of coals and carbons, including its dependence on structure. Such understanding will also aid in developing and optimizing processes for the separation of carbon dioxide from flue gases and other industrial streams that are based on adsorption and storage in carbon. Besides its obvious importance as a physical performance limit for any process, the capacity is a key parameter in many classical adsorption models, such as the Langmuir or Dubinin models and their many variants,5 or others based on multisite or vacancy solution theory6-9 that have a kinetic, a thermodynamic, or a statistical thermodynamic starting point. In such models, the capacity is generally considered to be dependent only on temperature and independent of pore size and determined by fitting macroscopic adsorption data, though carbons invariably possess a distribution of pore sizes. This seems reasonable for an adsorbed fluid in large pores, such as * Author to whom correspondence should be addressed. Email: [email protected]. (1) Voormeij, D. A.; Simandl, G. J. Geosci. Can. 2004, 31, 11. (2) Busch, A.; Kroos, B. M.; Gensterblum, Y.; van Bergen, F.; Pagnier, H. J. M. J. Geochem. Explor. 2003, 78-79, 671. (3) Karacan, C. O ¨ ; Mitchell, G. D. Intl. J. Coal Geol. 2003, 53, 201. (4) Ramanathan, C.; Bencsik, M. Magn. Reson. Imag. 2001, 19, 555. (5) Yang, R. T. Gas Separation by Adsorption Processes; Imperial College Press: London, 1997. (6) Nitta, T.; Shigetomi, T.; Kuro-Oka, M.; Takashi, K. J. Chem. Eng. Jpn. 1984, 17, 39. (7) Bering, B. P.; Serpinskii, V. V.; Yakubov, T. S. Izv. Akad. Nauk SSSR Ser. Khim. 1977, 4, 727. (8) Bhatia, S. K.; Ding, L. P. AIChE J. 2001, 47, 2136. (9) Ding, L. P.; Bhatia, S. K. AIChE J. 2003, 49, 883.

macropores or mesopores, where capacity is independent of the pore volume. However, this is not the case for micropores whose pore widths are only a few adsorbate molecular diameters. Here, the close-packed configuration of the adsorbed phase can vary with pore width, passing sequentially through hexagonal and square geometries, in a manner similar to that observed experimentally when colloid suspensions are confined between parallel glass plates.10-13 In earlier work, Schmidt and Lo¨wen14 confirmed the dependence of close-packing density on plate separation distance, by means of theoretical studies and computer simulations of frozen hard-sphere fluids between two parallel plates. Several studies also show pore-sizedependent density at high pressure due to packing transitions in all geometries.15-17 The transitions, occurring continuously with increasing pore width, allow rearrangement of adsorbed molecules in pores in order to attain a minimum free energy configuration.15 Although much theoretical work has been done on high-pressure layering transitions in confined geometry, the studies have invariably involved only hard-sphere systems, such as colloidal particle systems, where the capacity can be calculated purely from the geometry of the particles.14,18 Little theoretical study of high-pressure layering transitions of soft-sphere systems has been attempted due to their complexity. This problem has, however, been successfully examined in our laboratory using density functional theory (DFT),19 with methane adsorption capacities in good agreement with simulation. Because of the special importance of CO2 in the context of sequestration and the greenhouse issue, a study of the (10) Pieranski, P.; Strzelecki, L. Phys. Rev. Lett. 1983, 50, 331. (11) Murray, C. A. Phys. Rev. B 1990, 42, 688. (12) Murray, C. A. In Bond-Orientational Order in Condensed Matter Systems; Strandburg, K. J., Ed.; Springer: New York, 1992. (13) Weis, J. A.; Oxtoby, D. W.; Grier, D. G.; Murray, C. A. J. Chem. Phys. 1995, 103, 1180. (14) Schmidt, M.; Lowen, H. Phys. Rev. E 1997, 55, 7228. (15) Davies, G. M.; Seaton, N. A. Carbon 1998, 36, 1473. (16) Nicholson, D.; Stubos, T. In Recent Advances in Gas Separation by Microporous Ceramic Membranes; Kanellopoulos, N. K., Ed.; Elsevier: New York, 2000. (17) Cracknell, R. F.; Nicholson, D.; Tennison, S. R.; Bromhead, J. Adsorption 1996, 2, 193. (18) Pansu, B.; Pieranski, P.; Pieranski, P. J. Phys. (Paris) 1984, 45, 331. (19) Nguyen, T. X.; Bhatia, S. K.; Nicholson, D. J. Chem. Phys. 2002, 117, 10827.

10.1021/la048571i CCC: $27.50 © 2004 American Chemical Society Published on Web 09/30/2004

CO2 High-Pressure Adsorption Capacity and Structure

carbon dioxide adsorption capacity in carbons has been conducted in our laboratory using our recent model,19 and the results are reported in the present article. Since CO2 is a linear molecule, the Lennard-Jones (LJ) model assumption and the associated packing transitions may not necessarily describe the correct high-pressure behavior and may lead to inaccurate predictions for the capacity. Consequently, it is important to consider also a more representative potential model that better describes the molecular shape. Here, we consider also the accepted three-center potential model with associated charges20 and find qualitative as well as quantitative differences in behavior in comparison to the single-site LJ model while validating our theory in the latter case with simulation. Earlier studies with the three-center model for carbon dioxide have investigated isotherms and transport in slit pores,21-24 with only passing attention to the microstructure of the adsorbed phase. Further, the difference in packing phenomena compared with those using the singlesite LJ model has not been seriously explored. Here, we address both of these concerns. In the next section, we review the DFT theory for packing transitions of singlesite molecules adsorbed in slit pores. After a brief review of our simulation method, we describe new results for the adsorption capacity of single-center and three-center carbon dioxide in Section 4. 2. Packing Theory for LJ Fluids at High-Density Limit A. General Approach. DFT has proved to be a powerful tool for the investigation of adsorption equilibrium of single-site molecules in inhomogeneous systems25-29 and, in our laboratory, has also been used in the investigation of dense packing.19 The approach minimizes the grand potential

Ω[F] ) F[F] +



F(r)[vex(r) - µ] dr

(1)

Here, F(r) is the density profile in the presence of external potential, vex. Further, µ is the chemical potential of the adsorbate, and F[F] is the intrinsic Helmholtz free energy, decomposed in terms of an ideal gas part, Fid[F], and an excess part, Fex[F]

F[F] ) Fid[F] + Fex[F]

(2)

Given suitable models for these free energies, the general procedure is to obtain the optimal profile, F(r), that minimizes the grand potential, Ω[F], for the chosen system. B. Application to Close-Packed Limits in 2D and 3D. In the close-packed limit, the adsorbed phase is considered to be a set of identical unit cells, which are detailed in a later section. The local density can be then be treated as the density in the unit cell. Consequently, (20) Murthy, C. S.; O’Shea, S. F.; McDonald, I. R. Mol. Phys. 1983, 50, 531. (21) Samios, S.; Stubos, A. K.; Kanellopoulos, N. K.; Cracknell, R. F.; Papadopoulos, G. K.; Nicholson, D. Langmuir 1997, 13, 2795. (22) Nicholson, D. Langmuir 1999, 15, 2508. (23) Nicholson, D.; Gubbins, K. E. J. Chem. Phys. 1996, 104, 8126. (24) Papadopoulos, G. K. J. Chem. Phys. 2001, 114, 8139. (25) Tarazona, P. Phys. Rev. 1985, 31, 2672. (26) Evans, R. In Fundamental of Inhomogeneous Fluids; Henderson, D., Ed.; Marcel Dekker: New York, 1992. (27) Kierlik, E.; Rosingberg, M. L. Phys. Rev. A 1991, 44, 5025. (28) Denton, A. R.; Ashcroft, N. W. Phys. Rev. A 1991, 44, 8242. (29) Bhatia, S. K. Langmuir 1998, 14, 6231.

Langmuir, Vol. 20, No. 22, 2004 9613

eq 1 can be rewritten in terms of the unit cell potential as follows

Ω ) Fid + Fex + N(vex - µ) ) Fid + Fex + Nvex - Nµ (3) where N is the number of adsorbed molecules in a unit cell. For the two-dimensional unit cell, the ideal gas free energy is given by

Fid ) NkT[ln(N*λ2) - 1]

(4)

Here, λ is the de Broglie wavelength, λ ) xh2/2πmkT, where h is Planck’s constant, m is the mass of the adsorbate particle, k is Boltzmann's constant, and T is the temperature. The (2D) density of the unit cell, N*, is given by

N* )

N S

(5)

where N is the number of molecules in a unit cell and S is the surface area occupied on an adsorbent plane. The chemical potential may be expressed in terms of the bulk activity, a, as

µ ) kT ln(λ3) + kT ln(a)

(6)

In eq 3, the excess part of the Helmholtz free energy may be obtained from the canonical partition function, Q, through

Fex ) -kT lnQ

(7)

Q ) e-βE

(8)

where

since the degeneracy is unity in the close-packed limit. Here, E is the total interaction potential energy of the unit cell and β ) 1/kT. The excess Helmholtz free energy is therefore given by

Fex ) E

(9)

To proceed further, it is necessary to specify the potential models and unit-cell geometry that will be employed for minimization of the grand free energy, Ω. C. Potential Models. The adsorbate-adsorbate interactions are described by a (LJ) 12-6 pair potential

Φff(r) ) 4ff

[( ) ( ) ] σff r

12

-

σff r

6

(10)

The parameters ff and σff for the potential were determined by fitting the bulk isotherm for CO2 as described later. For the solid-fluid interaction, we use the 10-4-3 Steele potential,30 applied to an infinitely wide slit pore having infinitely thick walls

Φsf(z) ) A

[(

) ( )

2 σsf 5 z

10

-

σsf z

4

-

(

σ4sf

)]

3∆(0.61∆ + z)3

(11)

Here, A ) 2πFsstσ2st∆, z is the distance between an adsorbate molecule and the solid surface. ∆ is interplanar distance in the carbon, and Fs is carbon site density. (30) Steele, W. A. Surf. Sci. 1973, 36, 317.

9614

Langmuir, Vol. 20, No. 22, 2004

Bhatia et al.

Figure 1. Evolution of monolayer to form a second layer. (a) The rows containing empty circles move down and those having filled balls move up when width increases. The rectangle represents the selected unit cell projected in different directions. (b) The hexagonal configuration is divided into two distinct types of rows. (c) Configuration in a layer, continuously evolving to the square geometry.

The parameters sf and σsf are determined using the Lorentz-Berthelot mixing rules. The carbon parameters used in the study have the values30 σss ) 0.34 nm, ss/k ) 28.0 K, ∆ ) 0.335 nm, and Fs ) 114 nm-3. The external potential in a pore of width H (center-to-center distance of carbon atoms on opposing planes) is then given by

vex ) Φsf(z) + Φsf(H - z)

(12)

D. Unit-Cell Energy. Using the above fluid-fluid and fluid-solid potential models, the total interaction energy, E, of a unit cell can now be evaluated for any chosen transition. We consider close-packed LJ particles, forming a two-dimensional structure, comprising unit cells. To determine the total interaction potential for the unit cell, the following assumptions are applied (i) Pairwise additivity. (ii) The adsorbed phase consists of an integral number of identical unit cells. (iii) Structural transitions of the adsorbed phase correspond to a series of following crystalline structures as the pore width increases, following

14 f 20 f 24 f 30 f 34 ... The symbol 4 represents hexagonal packing, and the symbol 0 represents square packing. The numbers preceding the symbols represent the number of layers in the adsorbed phase. This sequence of crystalline structure has been studied for hard-sphere systems,14,31 as well as for frozen phases of an LJ fluid in slit-shaped pores.32-33 In particular, two main transition types may be identified. (1) The transitions n0 f n4, or rhombic transitions. In this case, only the structure within each layer changes. (2) The transitions n0 f (n + 1)4, or buckling transitions. In this case, a new layer is formed together with structural variation of the adsorbed phase. On the basis of the above assumptions, the unit cells corresponding to each transition are presented in turn below. 1. Evolution of the Monolayer and Formation of a Second Layer. a. Buckling transition (14 f 20). Figure 1 depicts the transition occurring with increase in pore size when the number of adsorbed layers increases from one to two. At the smaller pore size when only one layer is formed, the particles are organized in hexagonal close-packed geometry. As the pore size increases, a (31) Johnson, E. Science 2002, 296, 477. (32) Vishnyakov, A.; Neimark, A. V. J. Chem. Phys. 2003, 118, 7585. (33) Ayappa, K. G.; Ghatak, C. J. Chem. Phys. 2002, 117, 5373.

Figure 2. Continuous deformation of a unit cell. The square configuration is continuously deformed to transform into hexagonal configuration as pore width enlarges. Empty circles represent adsorbate molecules in the first layer and filled ones represent the adsorbate molecules in the second layer.

continuous transition occurs to the square geometry attained when exactly two layers can be formed. From Figure 1, each adsorbate molecule has the same number of nearest neighbors. Thus,the interaction potential, E, of the unit cell is given as follows

E ) N[Φff(d) + Φff(d1) + Φff(d2) + 2Φff(xd12 + d22)] (13) where d, d1, and d2 are the unit cell parameters identified in Figure 1c and N is the number of particles in the unit cell. Here, N ) 2. The relationship between the unit cell dimensions, d, d1, and d2, and the pore width, H, and density can be derived from the geometry of the unit cell.19 b. Rhombic Transition (24 f 20). Figure 2 describes the transition after the buckling transition, which ends up with exactly two layers, as the pore width widens. In this transition, a continuous deformation of the square configuration into a hexagonal configuration occurs without the formation of a new layer. For the unit cell shown in the Figure 2, the total interaction, E, of the cell for the rhombic transition is given as follows

E ) N[2Φff(d1) + Φff(d42) + Φff(d31) + 1.5Φff(d) + 0.5Φii(d53)] (14) where N is the number of molecules in the unit cell. In this case, N ) 2. The variables d42, d31, and d53 are distances between vertexes of the unit cell. 2. Evolution of the First Two Layers and Formation of a Third Layer. a. Buckling Transition (24 f 30). Figure 3 demonstrates the buckling transition

CO2 High-Pressure Adsorption Capacity and Structure

Langmuir, Vol. 20, No. 22, 2004 9615

Figure 3. (a) Adsorbed molecules from the first two layers contribute to the formation of a third layer. At the pore width where the third layer starts to form, these adsorbed molecules do not lie on the same plane. (b) For sufficiently large pore widths, the adsorbed molecules lie on the same plane. The selected unit is body-centered cubic (BCC) instead of face-centered cubic (FCC).

occurring with increase in pore width when the number of adsorbed layers increases from two layers to three layers. This transition can be modeled in stages. In the early stage of the transition at the smaller pore size shown in Figure 3a, adsorbates belonging to the third layer do not lie on the same plane since there is a stronger attraction between the nearest pore wall and the adsorbate and between adsorbates within each layer than between the more-distant pore wall and the adsorbates and between adsorbates in different layers. In this stage, second-nearest neighbors, or higher order, rather than only nearest neighbors, should be taken into account for an accurate solution, as the molecules are relatively close to each other. However, this complexity may be avoided if the unit cell parameter d is assigned to be 21/6σff at this stage. With this simplification, it was found that the predicted results matched the simulated ones reasonably well, as will be discussed subsequently. In the final stage of the transition, at the larger pore width shown in Figure 3b, adsorbates belonging to the third layer locate at the pore center. Furthermore, the adsorbates in this stage are relatively far apart. Consequently, it is necessary to consider only nearest neighbors to obtain a reasonable match to grand canonical Monte Carlo (GCMC) simulation. From Figure 3, the total interaction energy is given by

E ) 8Φff(d) + 2Φff(d1) + 2Φff(d2) + 2Φff(d1x2) + 2Φff(d2x2) (15) where d, d1, and d2 are as shown in Figure 3. b. Rhombic Transition (30f 3∆). In this transition, depicted in Figure 3b, the first and second layers are relatively far apart and the third layer is assumed to be located at the pore center. Thus, if only nearest neighbors have been taken into account, the body-centered cube can be split into two identical cubes,similar to that shown in Figure 2. These cubes are symmetrical through the

centerline of the pore. The total interaction energy, E, of the unit cell is given

E ) N[2Φff(d1) + Φff(d42) + Φff(d31) + 1.5Φff(d) + 0.5Φff(d53)] (16) E. Free Energy Minimization. Equations 3, 4, and 6 provide the grand potential of the unit cell

[ (1λ N*a ) - 1]

Ω ) E + Nvex - NkT ln

(17)

where the external potential, vex, follows eq 12. The unit cell potential energy, E, corresponding to any structure, is obtained as discussed in Section D. Through these, the grand potential can be expressed in terms of the characteristic parameters related to the geometry of unit cell and its position to relative to pore surfaces. The equilibrium structures will correspond to a minimum in free energy with respect to the characteristic parameters of the unit cell whose energetically optimized values are then employed to determine the equilibrium density of the adsorbed phase corresponding to each transition, as detailed in our earlier work.19 3. Simulation GCMC simulations have been conducted to verify the predictions of the above theory for LJ fluids and to study the capacity variation with pore size for a three-center model of CO2. These simulations mimick the µ,V,T ensemble used the established Metropolis sampling scheme34 for moving (including rotation for the three-center model), creating, or deleting molecules. Throughout a simulation, the numbers of attempted deletions and creations were kept equal to maintain microscopic reversibility. To ensure reliable values for the capacity, a range of high fugacities between 1000 and 6000 bar was investigated. It was found that larger fugacities were needed to attain saturation at (34) Metropolis, N.; Rosenbluth, A. W.; Rosenbluth, M. N.; Teller, A. N.; Teller, E. J. Chem. Phys. 1953, 21, 1087.

9616

Langmuir, Vol. 20, No. 22, 2004

Bhatia et al.

Figure 4. Fit of LJ model isotherm to the three-center model isotherm for bulk carbon dioxide at 308 K. higher temperature for a given pore size or for a larger pore size at a given temperature. At high pressure, the acceptance rate for creating or destroying particles can be low, and 10 creation or deletion attempts to one move attempt were therefore performed. In the simulation, the LJ potential is truncated at 2.5 nm. For the three-center model, the potential parameters used were the same as those established by Murthy et al.,20 with the values σoo ) 0.3026 nm, σcc ) 0.2824 nm, oo/k ) 75.2 K, and cc/k ) 26.3 K, and point charges on the atom centers of qo/e ) -0.332 and qc/e ) 0.664. The O-O distance, doo, is 0.2324 nm and the C-O distance is 0.1162 nm. The Lorentz-Berthelot rules were used to estimate parameters for the atom-atom cross interactions. The parameters of the single LJ site CO2 were determined by matching the bulk isotherm for this fluid with that of the three-center fluid between 0 and 2500 bar at 308 K. This yielded σf ) 0.37 nm and f/k ) 220 K. Figure 4 depicts the bulk isotherms obtained using the three-center and single-center models with the above parameters. The interaction potential between an adsorbate LJ center and each wall in the graphitic slit-pore was represented by the Steele30 10-4-3 potential. The size of the simulation box was set to 6 nm, and periodic boundary conditions were employed. The dependence of adsorbed density on the box size was insignificant under these conditions, particularly since the aim is to compare spherical and linear CO2 models. The number of configurations varied from 1 × 107 to 6 × 107. GCMC simulations were run at 308, 318, and 328 K to investigate the effect of temperature on adsorbate density at the high-pressure limit over a wide range of pore size.

Figure 5. Variation of adsorbed high-pressure capacity with pore size at (a) 308, (b) 318, and (c) 328 K.

4. Results and Discussion Lennard-Jones Fluid. Initially, computations were done for the one-center CO2 using both simulation and theory as discussed above. As in the earlier work19 with methane, good agreement was found between theory and simulation. Figure 5a-c depict the variations of capacity with slit width for the single-center fluid at three different temperatures, with the open circles representing simulation and solid lines representing the theoretical predictions. As in the case of methane,19 the packing structure and transitions considered in the theory provide a satisfactory quantitative explanation of the capacity variation with slit width. Thus, the branch AB represents the monolayer with LJ molecules in two-dimensional hexagonal configuration, with the pores in this range having a single potential minimum. In the range of pore width covered by this branch, the capacity reduces with increase in pore size because of reducing packing density, as the interparticle separation is nearly constant at ro ≈ 21/6σff. The presence of hexagonal packing in the monolayer is also confirmed by the pair distribution function shown in Figure 6, which yields peaks at ro, x3ro, 2ro, and x7ro,

Figure 6. Pair distribution function of LJ carbon dioxide obtained by GCMC simulation illustrating the different structure in hexagonal and square configuration.

corresponding to the locations of increasingly distant neighbors in this configuration. On branch BC, the structure adjusts gradually from the hexagonal packing to a square packing through a

CO2 High-Pressure Adsorption Capacity and Structure

Langmuir, Vol. 20, No. 22, 2004 9617

Figure 7. Density profiles at three different pore sizes depicting the formation of the third layer at the center.

continuous sequence of intermediate configurations, illustrated in Figure 1. In this case, there are two solutions for the value of the layer position corresponding to the inception of a layer near each wall. The gradual establishment of the second layer leads to an increase in packing density, evident on the branch BC. The transition on the branch BC is three-dimensional, as particles separate from the monolayer, and corresponds to a buckling transition. At C, the two layers are fully established in square configuration, which is confirmed by the pair distribution in Figure 6, showing peaks at approximately ro, x2ro, 2ro, and x5ro. On further increase in pore size on the branch CD, no new layer is formed, but the two layers undergo a gradual transition from square at C to a perfectly hexagonal configuration at D, as illustrated in Figure 2, corresponding to a rhombic transition. In this range of pore widths, the packing density and, therefore, capacity drops, as no new layer is formed while the interparticle spacing increases from 1.164σsf to 1.202σsf. On further increase in pore width beyond D, a third layer separates out with variation in the structure of each layer (i.e., a buckling transition) from hexagonal to square, as depicted in Figure 3. The unit cell structure is modeled as body-centered cubic (BCC), and the agreement between theory and simulation along branch DE supports this treatment. The third layer is, however, initially not at the pore center, as depicted by the solid line showing two shallow peaks on either side of the center for a pore width of 1.19 nm in Figure 7. With an increase in pore width, these two peaks coalesce into an increasingly narrow peak at the center, where the third layer finally becomes established, as seen in Figure 7 for pore widths of 1.225 and 1.3 nm. At the end of the transition, at E, each of the three layers is in perfect square configuration, with a BCC unit cell representing the 3D packing structure for the particles in the three layers. The formation of this third layer leads to a gradual increase in density on branch DE. On branch EF, again a rhombic transition occurs, similar to branch CD, with the square configuration in each layer being continuously deformed toward the hexagonal one. Three-Center Fluid: Capacity. Simulations of the capacity were also performed using the three-site model, with a charge at each site, as discussed earlier. In this case, no clear packing transitions are found and the capacity increases continuously with pore size except for very low pore sizes below H ) 2σsf, for which a decrease is noted with an increase in pore width. Here, σsf is the one-center LJ parameter, used here for scaling purposes only. The open triangles in Figure 5a-c depict the results

Figure 8. Profiles of (a) density and (b) 〈cos2(θ)〉 for pores of size 0.68 and 0.75 nm dominated by monolayer adsorption. In the latter case, the initial stages of the formation of two layers is apparent from the shallow peaks, one on either side of the center.

obtained, showing the above trend, in sharp contrast to the oscillatory behavior for the LJ fluid resulting from the packing transitions discussed above. It is clear that the one-center LJ model does not correctly represent the true packing behavior and leads to inaccurate values of the capacity at any pore size. Nevertheless, since the LJ model overpredicts in some ranges of pore size, while underpredicting in others, it is likely that in the presence of a wide pore size distribution (PSD), as exists in most carbons, the error in the overall capacity following

Cmt )

∫o∞Cmh(H)f(H) dH

(18)

may not be large. Here Cmh(H) is the capacity for a pore of size H, f(H) is the pore volume distribution, and Cmt the overall capacity. For a narrow PSD, however, unacceptably large error may be expected. It may be noted that the above error in the capacity of the LJ approximation is generally in the range of 5-8% and is not serious, except in the region ABC (i.e., pore size between 0.64 and 0.95 nm), where a transition occurs from one to two layers. In this region, the error can be as high as 30%. The inaccuracy of the LJ approximation in this region has been previously noted by Vishnyakov et al.,35 who also note the adequacy at other pore sizes for the low and moderate densities considered by them. Simulations with the three-center fluid without charges at the sites were also conducted and yielded similar trends as given by the fluid with charges, with only a small (35) Vishnyakov. A.; Ravikovitch, P. I.; Neimark, A. V. Langmuir 1999, 25, 8736.

9618

Langmuir, Vol. 20, No. 22, 2004

Bhatia et al.

Figure 9. Snapshots of adsorbate microstructure at different pore sizes. For clarity in viewing the (dense) packing, the molecules are displayed as lines, with the red portion representing oxygen and the central gray part representing carbon.

reduction in capacity. Thus, the difference in behavior compared to the LJ fluid is clearly related to the different molecular shape represented by the three-center fluid and not to electrostatic effects. To explain this difference, it is necessary to examine in some detail the microstructure of the three-center fluid at various pore sizes, and this is considered in the next section. Three-Center Fluid: Structure. While the variation in capacity with pore size, as discussed above, does not reveal any distinct phase transitions or structural transformations, studies of the adsorbate structure do reveal a pattern of change with respect to pore size, associated with layer formation for the three-center fluid. For very small pore size below H ) 0.71 nm (i.e., H* ) 2) a monolayer exists, with CO2 molecules aligned parallel to the walls. The solid line in Figure 8a depicts only a single central density peak for a small pore width of 0.68 nm at 308 K, consistent with monolayer formation, with a mean value of cos2(θ) of only about 0.06 (cf. Figure 8b). This corresponds to an effective orientation angle, θ, of 76°, where θ is the angle between the molecule axis and a normal through the pore walls, consistent with the nearly parallel alignment. Here, the density distribution corresponds to that of the molecule center of mass (i.e., carbon atom center) location. Figure 9 depicts snapshots of the adsorbate structure at different pore sizes at 308 K, showing the nearly parallel alignment of the molecules at a pore width of 0.68 nm. With an increase in pore width, there is a tendency for the molecules to rotate from the flat position, thereby permitting additional molecules to adsorb. For pores below a size of about 0.71 nm, the density still decreases with an increase in pore width, as seen in Figure 5, due to the persistence of a monolayer. However, at larger pore widths, an additional layer begins to form near the walls, as seen by the incipient shallow peaks about 0.075 nm on each side of the center for the pore width of 0.75 nm in Figure 8a. The formation of the additional layers is also evident in Figure 9, showing the appearance of some flat-lying molecules near the pore walls with the molecules at the

Figure 10. Profiles of (a) density and (b) 〈cos2(θ)〉 for pores of size 0.85, 0.923, and 1.065 nm illustrating the broadening of the central density peak and the separation of new layers of rotated molecules on each side of the center, as well as subsequent formation of a central layer of flat-lying molecules.

centerline showing a tendency to rotate. This is confirmed by the increase in mean value of cos2(θ) at the center to about 0.31 for H ) 0.75 nm, indicating an effective orientation angle of 56° from the normal. Such rotation of linear molecules adsorbed on surfaces has recently been reported by Zhao et al.36 on the basis of simulation of the adsorption of propane on graphite. These authors have determined orientation-angle distributions, showing peaks near 90° and 0° from the normal, with enhancement of the peak corresponding to the normal orientation at high density. This rotation provides an extra degree of freedom that serves to stabilize the adsorbate, preventing packingstructure transitions similar to those for the LJ fluid. The rotation of the molecules can be justified from a purely energetic viewpoint alone. As a first approximation, consider only fluid-solid interactions and assume two possible orientationssflat and normal. It is estimated that the interaction energy of a CO2 molecule lying flat on a 10-4-3 surface is -13.9 kJ/mol, while it is -8.6 kJ/mol for the molecule orientated along the normal. The addition of a CO2 molecule requires some flat-lying molecules to rotate and release space. Assuming a CO2 molecule to be cylindrical, with diameter σoo ()0.3026 nm) and length σoo + 2dco () 0.535 nm), the area occupied by the flat molecule is estimated as 0.1619 nm2, while for the vertical orientation it is 0.0916 nm2 (σ2oo). It is readily estimated that addition of a normally oriented molecule must effectively result in 1.3 flat-lying molecules to reorient along the normal. In terms of energy, the interaction energy of 1.3 mol of flat-lying molecules is -18.1 kJ, while for 2.3 mol of normally oriented molecules, it is -19.8 kJ, (36) Zhao, X.; Kwon, S.; Vidic, R. D.; Gorquet, E.; Johnson, J. K. J. Chem. Phys. 2002, 117, 7719.

CO2 High-Pressure Adsorption Capacity and Structure

Langmuir, Vol. 20, No. 22, 2004 9619

Figure 11. Profiles of (a) density and (b) 〈cos2(θ)〉 for pores of size 1.207, 1.491, and 1.5975 nm illustrating the broadening of the central density peak and the separation of new layers of rotated molecules on each side of the center, as well as subsequent formation of a central layer of flatter-lying molecules.

permitting the spontaneous addition and reorientation. Nevertheless, it may be noted that this rudimentary analysis overlooks intermolecular interactions in the adsorbate, which will affect the result quantitatively, although the qualitative trends are probably unaffected. Further, entropic considerations arising from the rotational degrees of freedom will also favor the rotation as the pore width increases while a monolayer is maintained. As discussed above, for the pore size of 0.75 nm, incipient formation of layers near the pore walls, with relatively flat-lying molecules, is evident. These layers become established as the pore size increases, as seen in Figure 10a for H ) 0.85 nm, which depicts prominent density peaks near the walls in addition to the central peak. Figure 10b shows that the molecules at the center are increasingly rotated toward the normal, with an effective deviation of 35° from the normal (corresponding to 〈cos2(θ)〉 ) 0.67) for H ) 0.85 nm. These features are corroborated in Figure 9 by the snapshot for this pore size. With further increase in pore size to 0.923 nm, the central density peak becomes broad and shallow, with incipient formation of a shoulder with a weak peak about 0.06 nm on each side of the center. Thus, the central layer begins to separate into two shoulders of rotated molecules, corroborated by the broadening of the central peak for the rotated molecules (〈cos2(θ)〉 ) 0.67) at H ) 0.923 nm in Figure 10b. These shoulders become established at larger pore size, as seen by the density peaks for H ) 1.065 nm in Figure 10a. The shoulders occur at about 0.14 nm from the center at H ) 1.065 nm and are comprised largely of rotated molecules, as seen from the peak value of 〈cos2(θ)〉 of about 0.53 in Figure 10b, with the molecules in the wall layers being nearly flat-lying with 〈cos2(θ)〉 ≈ 0. Their appearance is due to the dislocation of the center of mass on rotation of

Figure 12. Pair distributions at various pore sizes and for bulk fluid, using a three-center model for carbon dioxide.

the molecule, which is confirmed by their separation of about 0.075 nm from the peak for the flat surface layer, that is close to the projected root-mean-square length of the C-O bond of about 0.085 nm for 〈cos2(θ)〉 ) 0.53. Further incipient formation of a central layer with nearly flat-lying molecules is evident from the shallow peak beginning to appear at the center for H ) 1.065 nm in Figure 10a and the relatively small value of 〈cos2(θ)〉 of about 0.1 at the center in Figure 10b. The above features are corroborated by the snapshots for H ) 0.923 and 1.065 nm in Figure 9. On further increase in pore width, the above pattern is repeated, with the central relatively flat layer rotating and yielding a broader density peak, before again separating into two distinct layers of rotated molecules, evident from the results in Figure 11a and b for H ) 1.207, 1.491,

9620

Langmuir, Vol. 20, No. 22, 2004

Bhatia et al.

on the structure, although they enhance the density quantitatively, as discussed previously. At low densities, simulations of carbon dioxide adsorption on graphite surfaces do show the formation of structures mediated by electrostatic interactions.37 However, in the present case for the small intermolecular spacing arising at highdensity packing, the repulsive part of the intermolecular potential overshadows the electrostatic effects and no distinct structure is found. This is evident in Figure 13, which depicts the relatively flat adsorbed monolayer in a pore of width 0.68 nm, showing some clustering of CO2 molecules in nearly parallel alignment but no apparent ordering of the molecules in T-like configurations mediated by electrostatic effects that are found37 at low densities. 5. Conclusions

Figure 13. Central adsorbed monolayer for a pore of width 0.68 nm. For clarity in viewing the (dense) packing, the molecules are displayed as sticks, with the red portion representing oxygen and the central gray part representing carbon.

and 1.591 nm, and the snapshots in Figure 9 corresponding to these sizes as well as for H ) 1.3 nm. All of the above results reveal a clear and distinct pattern of formation of a central relatively flat layer followed by rotation and subsequent separation of two distinct rotated layers, with increase in pore size. In this way, additional layers are created as the pore widens. Near the pore surface, the peak for the layer of rotated molecules occurs as a shoulder to the main peak of flat-lying molecules. Contrary to the case of the LJ fluid, no packing transformations occur. This is confirmed by the pair distribution functions (PDFs) for the various pore widths investigated above, depicted in Figure 12a-c. These figures reveal a consistent pattern in the PDF with similar peak positions at all pore sizes, confirming the absence of transitions. The short-range structure matches that for the bulk fluid, as seen in Figure 12c, while the long range part of the PDF does not because of fluid molecule absence due to the effect of the pore walls. Clearly, the LJ approximation for a linear molecule such as CO2 leads to a distorted picture of the adsorbate structure and cannot be relied upon for meaningful interpretations and predictions, at least at high densities. However, this conclusion may not necessarily hold at low and moderate densities, where intermolecular repulsions are not dominant, as seen in the work of Vishnyakov et al.35 For the high-density system investigated here, Coulombic interactions do not have any qualitative influence

The capacity and adsorbate microstructure of CO2 in carbon slit pores is investigated here, using a LennardJones approximation as well as a more realistic threecenter model incorporating local charges. For the LJ fluid, a recent theory from our laboratory explains the adsorbate structure adequately and predicts capacities in agreement with simulation. It is found that the theory considering successive continuous transitions between hexagonal and square packing successfully predicts the oscillatory variations of high-pressure capacity with pore size for the LJ fluid. However, simulations with the three-center potential model for CO2 yield significantly different results from the LJ model for the filled pores with a capacity that increases monotonically with pore size beyond the monolayer width. At monolayer widths (pores smaller than 0.71 nm), reduction in capacity with increase in pore size occurs. Investigation of the adsorbate microstructure has shown that for the three-center fluid-packing transitions are absent, and instead, rotation of flat-lying molecules occurs at the pore center to accommodate additional adsorbed molecules as pore size increases, starting from monolayer pore widths. At sufficiently large pore size, however, the central layer splits into two layers, followed by the appearance of a new layer of flat-lying molecules at the center at larger pore size. With continued increase in size, this process repeats itself. The results suggest that the single-center LJ approximation is inadequate to account for the structure of linear molecules such as CO2 in narrow pores. Thus, simulation of competitive adsorption of CO2 and methane for high-pressure CO2 sequestration applications that employ the single-center LJ model, widely applied for CO2, may yield misleading results. Further studies along these lines will be reported in due course. Acknowledgment. This research has been supported by a grant from the Australian Research Council under the Linkage scheme. LA048571I (37) Johnson, J. K. Personal Communication, May 2004.