Journal of Colloid and Interface Science 342 (2010) 445–454
Contents lists available at ScienceDirect
Journal of Colloid and Interface Science www.elsevier.com/locate/jcis
Atomistic calculation of adsorption in activated carbon with pore-size distribution Hadas Abir *, Moshe Sheintuch Department of Chemical Engineering, Technion, Israel Institute of Technology, Haifa 32000, Israel
a r t i c l e
i n f o
Article history: Received 19 July 2009 Accepted 14 October 2009 Available online 20 October 2009 Keywords: Adsorption isotherm Phenol para-Bromophenol m-Cresol Activated carbon Molecular mechanics
a b s t r a c t We use molecular mechanics universal force field parameters to calculate single- and multicomponent adsorption of phenol, para-bromophenol, and m-cresol from the gas or liquid phase on activated carbon (AC). The carbon pores are modeled by shallow carbon nanotubes of various pore diameters. The effect of pore length is studied for the case of phenol. This calculation yields the Gibbs free energy change (DG) of adsorption which is used to predict the adsorption isotherm following reaction rate theory. The singlepore adsorption can be integrated to predict the overall adsorption if pore-size distribution is known. We show that the Freundlich adsorption isotherm may result from pore-size distribution coupled with a linear decline of DG with pore diameter. When the pore-size distribution of the carbon under study is not known, the adsorption of one species can be used to predict that of another on the same AC. We extend this methodology for the case of a bisolute adsorption system. Reasonable prediction of measured adsorption isotherms is demonstrated for single solute adsorption when both pore-size distribution and adsorbate–adsorbate interaction are taken into accounted. The suggested methodology is highly sensitive to the numerous parameters it requires. Ó 2009 Elsevier Inc. All rights reserved.
1. Introduction Activated carbon (AC) is a widely used adsorbent for separation of pollutants from the gas and liquid phases, but despite its wide use there is little atomic level understanding of adsorption of large organic molecules on AC. Knowing the adsorption isotherm is essential for many technological applications of activated carbon, and its prediction from basic principles may be used to replace experimentation and to optimize the adsorbent. Proper predictions of the adsorption properties of the system requires modeling of the adsorption isotherm for single-component as well as for competitive adsorption. The activated carbon performance is influenced both by its surface chemistry and by the functional groups of the surface as well as by its porosity and pore distribution. The AC surface structure is highly complex and depends on the raw material used to produce it, the method of production, and pretreatment. Commercial active carbons commonly have a surface area of 800–1500 m2/g, reaching as high as 2500 m2/g [1]. This surface area is supplied mainly by micropores (with diameter <2 nm), while mesopores (2–50 nm) and macropores (>50 nm) contribute mainly to transport within the material. Several surface functional groups were identified including carboxyl, carbonyl, phenol, quinone, lactone, and lactol [2,3]. Nanoporous carbon is primarily made up of tightly curved individual carbon layers, including pores of typical size of 1 nm. * Corresponding author. E-mail address:
[email protected] (H. Abir). 0021-9797/$ - see front matter Ó 2009 Elsevier Inc. All rights reserved. doi:10.1016/j.jcis.2009.10.032
Completely closed carbon particles are also present. The carbon may have a fullerene-related structure in which pentagons and heptagons are distributed randomly throughout a hexagonal network [4,5]. Most previous atomistic simulations of adsorption on AC described the binding of small molecules from the gas phase on a single-sized slit-shaped pore model using the Grand Canonical Monte Carlo (GCMC) approach with spherical molecules described by the Lennard-Jones (LJ) potential. Yet several numerical studies of adsorption energies on carbon nanotubes (CNT) using molecular mechanics (MM) simulations were also reported [6–9]. In the present paper we use MM to characterize the adsorption on pores modeled as carbon nanotubes (NT) of various diameters and lengths. We use this information to predict the adsorption isotherm. We described these concepts and justify our approach below. We represent the complex AC structure by a simple model composed only of micropores, represented by single wall carbon nanotubes (SWNT) of various diameters. The SWNT are graphene sheets wrapped into cylinders with 1–2 nm in diameter. There are three typical configurations for the SWNT: the ‘‘armchair,” ‘‘zigzag,” and ‘‘chiral”. Carbon nanotubes are molecular tubes which now attract considerable attention. Though a slit-shape model has been commonly used to represent an AC structure, its structural stability is usually not demonstrated computationally; it usually underestimates the molecule–AC interaction and various corrections have been suggested to correct it [10]. Several studies indicated that the adsorption on AC can be represented using a nanotube model [11–14]. Moreover,
446
H. Abir, M. Sheintuch / Journal of Colloid and Interface Science 342 (2010) 445–454
several studies on AC structure found that AC is composed of curved carbon layers, containing pentagonal rings which explain the porosity, due to curvature, and hardness while parallel grapheme layer are absent [4,5]. Based on these results Harris et al. claim ‘‘that the idea that microporous carbons have fullerene-like structure has important implications for the modeling of adsorption”. In these cases a NT model may be a better representation for the real AC structure [5]. In addition, the adsorption in very narrow NT allows interaction from all directions and implies higher adsorption energy and allows adsorption even in very low concentrations. While AC may be more complex than the simple NT geometry, we show that the narrow pores used here are essential for predicting high loadings observed in adsorption from water, and for predicting high surface area. This simple geometry represent an AC with a very high surface area of about 2500 m2/g. Active carbon has a lower surface area either since the pores should be modeled by double-walled NT or since they incorporate meso and macropores as well. The Freundlich adsorption isotherm has been widely applied in describing experimental adsorption of organic pollutants in aqueous solutions, including that of phenol and phenolic compounds [15–18]. The Freundlich isotherm has been originally accounted for by Langmuir adsorption while an exponential distribution of adsorption energies is assumed. Our atomistic calculations show how DG varies with pore diameter. This approach can be used to integrate the adsorption over all pore sizes, when their distribution is known. The Freundlich isotherm can approximately describe the calculated isotherm when a certain pore-size distribution is accounted for, while the behavior in a single pore diameter follows the Langmuir adsorption isotherm. This description is exact only when the distribution accounts for all pore sizes, while in practice the smaller pores are excluded. While adsorbate–adsorbate repulsion can account for this behavior [15,19,20] this explanation does not apply here since for phenol, the model adsorbate used here, the adsorbate–adsorbate interaction is attractive. Previous work by this group used atomistic calculations to study the single-component adsorption of phenol on shallow single-size SWNT with various diameters and functional groups [21]. This work showed that the adsorption of phenol on most functional groups is thermodynamically unstable at room temperature and most of the adsorption is contributed by micropores and graphite planes [21]. The present work extends this study to active carbon with poresize distribution as well as to the case of unknown pore-size distribution. We also extend the study to long SWNT, so that the adsorption is 3-D in nature and to carbon with multicomponent adsorption. This work suggests a computational chemistry methodology to predict adsorption properties on well-defined structures. The first part of this work deals with single-component adsorption and is incorporates the following steps: 1. Characterizing the structure of the adsorbent and calculating the maximal capacity of adsorption of a certain component. 2. Calculation of the change in the thermodynamic properties (enthalpy, entropy, and Gibbs free energy) of adsorption, and when necessary the energetic barrier to adsorption. We show that in general DG may vary with fractional occupancy. 3. Using this information to calculate the adsorption equilibrium of a NT in contact with a gas phase of known (partial) pressure of adsorbent, using known kinetic theory for calculation of adsorption and desorption probabilities. This will lead to several problems in 3-D adsorption as discussed and addressed below. When DG varies with fractional occupancy we should incorporate this dependence to the calculation or choose to show the effect of the maximal and the minimal DG.
4. Evaluating the adsorption isotherm by assuming a simple Langmuir isotherm for each pore size and then averaging for a certain pore-size distribution The gas–solid adsorption equilibrium is translated to liquid–solid equilibrium, using gas–liquid solubility data (i.e., Henry constant). The estimation of adsorption of two or more competing species requires more data than those of a single component: current practice for such isotherms is to apply the multicomponent Langmuir isotherm, which applies to a single-sized pore with no adsorbate–adsorbate binary interactions, or the multicomponent Freundlich (see derivation by Sheindorf et al. [15,16]). Methods that attempt to predict the multicomponent adsorption from singlecomponent data often fail to predict the isotherm in cases where adsorbate–adsorbate interaction occurs [15,16,18]. The Freundlich multicomponent isotherm requires single-component data as well as a competition coefficient which had to be determined experimentally. In the second part of this work we apply the molecular mechanics methodology to predict the multispecies competitive adsorption on AC with a certain pore-size distribution. The structure of this work is the following: We study first the single-component gas-phase adsorption of phenol, for which we describe the energy dependence on pore size and length in detail. We then study m-cresol and para-bromophenol (PBP) adsorption, before studying coadsorption of phenol with each of the other molecules. The work was motivated by our interest in adsorption on AC [15,16] and regenerative procedures using catalytic agents that were deposited on the AC [22,23].
2. Computational methods The adsorption of the model adsorbents, phenol, m-cresol, and para-bromophenol, on activated carbon is governed by van der Waals and Coulomb interactions. To represent van der Waals interactions, we applied molecular mechanics (MM) modeling using the universal force field (UFF) parameters of Rappe and co-workers [24–27]. This approach was proven to accurately represent both short- and long-range interactions in the absence of chemical interactions. While adsorption is usually studied using GCMC simulations with simplified Lennard–Jones potentials [28], there are several numerical studies of adsorption energies on CNT using MM simulations [6–9]. The MM approach yields a better description of the interaction of large (phenol) molecules with other molecules and with the walls, at the expense of the system size. GCMC can handle larger samples, but use relatively simpler potential description (LJ). All calculations were performed using the Gaussian 03 package from Gaussian Inc. and Arguslab 4.0.1 by Planaria software. The geometries were fully optimized in all calculations in order to reach the lowest energy configuration. It should be emphasized here that throughout this work it was assumed that there is infinite time for equilibration and that the minimum energy configuration of the system was reached. The minimum energy geometries were determined to be true minima by the absence of imaginary frequencies in the calculated vibrational spectrum. Thermodynamic parameters of the adsorption were calculated, using the Gaussian 03 package, as the difference between the corresponding values of the system with optimal position of a molecule inside a nanotube and those of the corresponding separated components. For example, for calculating DG for adsorption of a single phenol molecule on a NT the Gibbs free energy values were calculated separately for phenol, for the ‘‘empty” NT, and for the corresponding system of ‘‘phenol inside the NT” and then the DG for the system was calculated as DG = G(phenol+NT) G(phenol) G(NT). Similarly for
H. Abir, M. Sheintuch / Journal of Colloid and Interface Science 342 (2010) 445–454
the addition of a molecule to NT with n other molecules we calculated DG = G(n+1phenol+NT) G(phenol) G(nphenol+NT). In adsorption of multiple molecules, the thermodynamic parameters for the adsorption of the nth molecule were calculated with respect to the ‘‘adsorbent + (n 1) adsorbed molecules” and the separate adsorbate molecule. All calculations were conducted at a temperature of 298.15 K. The nanopores were modeled by shallow (7 Å) zigzag carbon nanotubes with various diameters (8–18 Å, or 12–24 C atoms in perimeter) and lengths (7–76 Å). For the derivation of the Freundlich isotherm we study shallow NTs, over a wide domain of radii. In order to avoid an influence of boundary conditions the outer C atoms were saturated with hydrogen in all the models.
3. Phenol adsorption 3.1. Effect of pore diameter and pore occupancy Single-component phenol adsorption is studied on NTs with diameters that vary between 8 and 18 Å [denoted as (12, 0)– (24, 0)] and of 7–76 Å in length (Fig. 1). The 7 Å long pore is composed of 3 rings in length, which corresponds to one layer of adsorption. Such pores can adsorb from a single molecule per layer in the (12, 0) NT and up to six molecules per layer in the (24, 0) NT. For the (15, 0), (16, 0), and (17, 0) NT, which contribute the most adsorption, we also employ phenol adsorption on pore lengths of 38 Å (i.e., with 17 rings in length, which corresponds to five layers of adsorption) and the (16, 0) NT is later studied up to 78 Å pore (corresponds to 10 layers of adsorption). We find that all the molecules are oriented parallel to the pore wall to maximize the interaction energy. The thermodynamic parameter computed varies with the pore diameter and pore occupancy. Analysis of the adsorption thermodynamic parameters indicates attraction between adsorbed molecules; i.e., the second molecule will find its optimal position next to the first and so on, filling up to six molecules in layer (depending on the pore diameter) before the next layer, which is also adjacent, starts to fill up (see Fig. 1). This course of ‘‘layer” adsorption follows from the results, showing that the optimal position of a new molecule is to be arranged in ‘‘layers”. This leads to the wavy DG variation with occupation. The data in Fig. 2 report the change of adsorption en-
447
ergy when molecules are added to the pore at the energetically optimal position. The arrangement of the adsorbed molecules in layers causes a repeating sinusoidal variation in the Gibbs free energy with pore occupancy. Comparison of the adsorption on (15, 0)–(17, 0) NT, that can adsorb up to 3 molecules in a layer, shows that this behavior is strongest for the (16, 0) nanotube (Fig. 2a and b), as adsorption on the wider or narrower pore becomes weaker. The decline in Gibbs free energy for an adsorbed molecule is strongest when it can create a ‘‘complete layer”. For the (15, 0) model the most stable arrangement is two molecules in a layer; when a third n molecule is adsorbed, it will arrange in the same layer (forming a three molecule layer), but with some increase in DG (Fig. 2b). These trends are evident in the plot of DG for the first, second, and third molecule vs pore diameter (Fig. 2c), showing that these effects are strong at the smaller pores, which in turn contribute most of the adsorption. 3.2. Effect of pore length The effect of pore length on gas-phase adsorption of phenol was studied on three models of zigzag (16, 0) nanotube of lengths of 7, 38, and 76 Å which adsorb up to 3, 15, and 30 molecules, respectively. A sinusoidal behavior can be observed—especially in the longest investigated pore. There is some decrease in DG with increasing occupancy of the pore, the extent of which depends on the pore length. This dependence should be further investigated and was not further addressed in the scope of this work. It is especially interesting since the decline of DG with occupancy leads to hysteresis of the isotherm. We have not found experimental support for this behavior yet (Fig. 3). 3.3. Adsorption isotherms in a nanotube To translate the thermodynamic DG value to an adsorption isotherm we use the Langmuir isotherm, which is the simplest adsorption isotherm on a surface
hA ¼
K A PA : 1 þ K A PA
ð1Þ
Here PA is the partial pressure in the AC environment and KA(DG) is the adsorption equilibrium coefficient. KA can be predicted by equating the adsorption and desorption rates, as derived from reaction rate theory and kinetic theory of gases [29].
a2 S0 ðTÞ expðDEdes =RTÞPA pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi m 2pmkT kT a2 kn0 S0 ðTÞ expðDEdes =RTÞ; hm
K A PA ¼
ð2Þ
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi where a2 is adsorption site area, k ¼ h= 2pmkT is the thermal de0 Broglie wavelength, and n ¼ PA =RT is the gas-phase density. S0 is the sticking coefficient which expresses entropy effects, S0 ¼ expðDS=RÞ. The adsorption energy is given by DG ¼ DH T DS, for nonactivated adsorption DEdes ¼ DH and Eq. (2) can be written as
KA ¼
Fig. 1. Visualization of adsorption: three phenol molecules in a (16, 0) nanotube model of length 38 Å. The three molecules are arranged in one layer due to attraction.
a2 k expðDG=RTÞ K A0 eDG=RT : RT
ð3Þ
We argue now why this procedure, which was derived for surfaces, should also apply for a 3-D pore: a pore can be viewed as a stack of layers. For the most exposed layer we can equate the adsorption and desorption rates and in that case Eq. (1) applies for each layer, where h gets the discrete values corresponding to addition of every molecule (h = 0, 1/3, 1/2, 2/3, 1, for our example of (15, 0)–(17, 0) 7 Å long pore). Since the filling order completes the layer before opening a new layer we record for each molecule
448
H. Abir, M. Sheintuch / Journal of Colloid and Interface Science 342 (2010) 445–454
Fig. 2. Calculated DG kcal/mol vs number of adsorbed phenol molecules (expressed as fraction of occupancy), for various NT diameters, in pores of 38 Å (a) and 7 Å (b) in length, and calculated DG (kcal/mol) vs pore diameter (Å) for the adsorption of the first, second, and third molecule in a layer for 7 Å length NT (c). (The lines were added to guide the eye.)
0 -5 0
0.2
0.4
0.6
0.8
1
-10 -15 -20 7
-25 -30
38 76
-35 Fig. 3. Calculated DG (kcal/mol) vs fraction of pore occupancy for various pore lengths for (16, 0) nanotubes (pore diameter is about 13 Å). The lines were added to guide the eye.
the corresponding PA at which it adsorbs with hA recording the pore partial occupancy at this stage. This will yield a multivalued diagram as in Fig. 4c. When the pressure is increased during adsorption the molecules will adsorb as a group (usually a layer), so an arrangement with minimal energy will be provided, since DG is decreasing within the group (as in Fig. 2). Thus, we can denote the negative slope branches as the unstable branches of this curve. With increasing partial pressure we move through the points that correspond to the minima in the DG plot. When decreasing PA we follow a different path, suggesting that hysteresis may emerge. Thus, since DG varies with pore diameter, length, and occupancy, and since we wish to integrate over various pore sizes, we want to characterize the adsorption in a simpler way. Fig. 4a and b presents several Langmuir-type adsorption isotherm of shallow 7 Å pores that were calculated based on the above suggested methodology. Note that large pores contribute very little to the overall adsorption. However, since DG varies with occupancy
Fig. 4. (a and b) Several calculated Langmuir-type adsorption isotherms of shallow 7 Å pores and various diameters based on the calculated DG (Dp) of the first molecule (a) or of DG of the arrangement with the minimal adsorption energy (b). Also shown are the maximal adsorption weighted fractions for each pore expressed as (g phenol/g C). (c) Calculated phenol adsorption isotherm—fraction of pore occupancy vs partial pressure for a nanotube of (16, 0) and length 76 Å. (The lines were added to guide the eye.)
H. Abir, M. Sheintuch / Journal of Colloid and Interface Science 342 (2010) 445–454
we check the adsorption energy for the first molecule (Fig. 4a) as well as for the arrangement with minimal energy (Fig. 4b). 3.4. Extension to systems of pore-size distribution We show now how the isotherm can be extended to a system with pore-size distribution. Recall that a Freundlich-type adsorption isotherm can be derived using an exponential distribution of energies, which under certain conditions implies exponential distribution of pore sizes. Suppose that the fraction of pores, N, having adsorption energy of DG follows:
NðDGÞ ¼ a expðnDG=RTÞ:
ð4Þ
Taking the above expression as a continuous approach, with
a expðnDG=RTÞ as pore energy distribution, and if in every pore diameter the occupancy follows the Langmuir isotherm, then according to Eq. (4) we can say that the fraction of adsorption sites (occupancy) having energy between DG=RT and ðDG þ dDGÞ=RT is given by
K A0 expðDG=RTÞ PA a expðnDG=RTÞdðDG=RTÞ: 1 þ K A0 expðDG=RTÞ PA
ð5Þ
The total occupancy (hT ) is the integral of Eq. (5), from DG*, the adsorption energy of the narrowest pore that can host a molecule, to DG = 1 [15,29]. Then,
R1 hT ¼
K A expðDG=RTÞP A 0 DG =RT 1þK A expðDG=RTÞPA expðnDG=RTÞdðDG=RTÞ R 10 : expðnDG=RTÞdðDG=RTÞ DG =RT
a
a
ð6Þ
The denominator is simply the total number of pores, which is integrated here for the same domain of pores that are effective for adsorption. Note also that in determining the average occupancy we should have weighted the various pores according to the number of molecules that occupy its cross section. In fact, as we later translate this value into adsorption weight capacity (q, g/g C), the adsorption should be weighted according to the maximal capacity (qm, denoted in Fig. 4), which varies somewhat with the pore diameter (Dp). For n < 0, which implies that there are more pores corresponding to weaker adsorption (larger DG), and when integrated over all DG (i.e., DG* = 1), this is the classical derivation of the Freundlich equation [15,29]:
hT ¼
aK nA0 n
PnA AP nA :
ð7Þ
Note that n reflects the energy distribution and should be specific to the surface as well as to the adsorbate molecule. For n > 0, which implies more pores of stronger adsorption, Eq. (6) does not converge and was solved numerically for specific values of n and DG*. A large n implies a narrow distribution and the system can approximately be described by a Langmuir isotherm, while n = 0 implies a flat distribution. To analyze the isotherm for n > 0, which describes many actual solids (see below), let us express u K A0 expðDG=RTÞ PA ; the integral limits then are u ¼ K A0 expðDG =RTÞ P A , which corresponds to DG*, and u* = 0 which corresponds to DG = 1. Then,
R0 h¼
u
u 1þu
R0
ðK
un n A PA Þ 0
un u ðK A PA Þn 0
¼ n ðu Þn
Z 0
a du u
a du u
u
R u ¼
un du: 1þu
0
R0
u
449
adsorption becomes stronger, up to DG*. The integrand, f = un/ (1 + u), passes through a maximum with increasing u, at umax = n/ (1 n) (or approximately umax = n when n 1). Thus for u* 1 the isotherm is linear, h = nu*/(n + 1), while for large u*, u* 1, h approaches 1. For the intermediate domain the integral was solved numerically using the following procedure. For every discrete pressure value in the range of study, a corresponding value of u* was calculated by using K A0 evaluated according to reaction rate theory (Eq. (3)) and DG* value evaluated from our model. Next the integral in Eq. (8) was solved for every point and the results were plotted as h vs PA. The results may also be presented as q (adsorption capacity) vs C (concentration) in order to compare them to experimental results, as explained in the next section (Fig. 5d). Eq. (8) can also be written in the following differential form that will be used below:
dh h 1 : ¼ n þ u 1 þ u du
ð9Þ
To apply this approach we employ data of phenol adsorption on AC in a study which provides data of the pore-size distribution. Wu et al. [17] measured an almost exponential pore-size distribution in a Fir wood AC, activated with KOH: NðDp Þ ¼ 0:355 expð0:073Dp Þ, where the pore diameter, ðDp Þ, is expressed in Å (Fig. 5a). To predict their measured adsorption results, the effect of pore diameter on gas-phase adsorption of phenol was studied on shallow 7 Å long models of zigzag nanotubes of 8–18 Å in diameter (nanotubes (12, 0)–(24, 0)). The pore diameter was limited to this range, since no phenol adsorption is possible below 8 Å and the effect of pores greater than 18 Å was considered to be negligible due to the weak adsorption (Fig. 5c) and the small contribution to the pore-size distribution. The pore length was limited to 7 Å. Since DG varies with pore diameter, length, and occupancy we had to limit our study. As Fig. 3 shows the differences in DG of the first molecule for the 7, 38, and 76 Å long pores are quite small. The behavior of DG with pore occupancy is quite similar for the first layer and for occupancy <0.6 (short of one data point). This justifies our choice of using the 7 Å pore as the main structure for further studies. The pore diameter is corrected to take into account the diameter of the carbon atom (c = 0.7 Å). (This implies taking into account the inner diameter of the pore.) Since the adsorption energy varies with occupancy we study two limiting cases where we employ either the energy of placing a single molecule in every pore or the energy of placing several molecules in a layer using the arrangement that will bring about the minimal DG (Fig. 5b). Using the curves in Fig. 5a and b, the energy distribution was evaluated for these two cases. To translate the gas–solid equilibrium to liquid–solid equilibrium we use Henry’s law for phenol with [30,31] H ¼ 2:37 106 atm m3 =mol. There is a significant uncertainty in the Henry constant value and reported measured values vary between 2 106 and 1:8 107 atm m3 =mol [30,31]. The chosen values will affect the calculated adsorption isotherm. The total weight adsorption was calculated as
qT ¼ h
X
qmi NðDpi Þ;
ð10Þ
i
un du 1þu
un1 du ð8Þ
Again, for 1 < n < 0 and bi-infinite DG domain, we will find the Freundlich isotherm. Here, we are interested in 1 > n > 0, which corresponds to a situation with a larger fraction of pores as the
where h was calculated according to Eq. (8), by summing over the pore diameters explained above, and qmi is the maximal weight adsorption capacity of pore of size Dpi. Solving Eq. (8) over a range of pressures leads to a curve that can be described well by the Freundlich-type adsorption isotherm, but the exponent is sensitive to the concentration range chosen. Using phenol concentration (C) in the range of 0–30 g phenol/m3 (corresponding to pressure range of 0.01–0.07 Pa, or u* in the range
450
H. Abir, M. Sheintuch / Journal of Colloid and Interface Science 342 (2010) 445–454
a
b
0.25
25
0.2
20
0.15
15
0.1
10
0.05
5
0 0
20
0 5
40
1
10
15
20
c 0.2
1 0.15 0.1 0.05 0 0
20
40
60
Fig. 5. Measured pore-size distribution (fraction of pores of each diameter vs inner pore diameter) after Wu et al. [17] (a), calculated DG (kcal/mol) vs inner pore diameter (Dp) for the first molecule adsorption and for the arrangements of molecules with minimal DG on a shallow (7 Å) pore (b), and calculated weight adsorption isotherm according to Eq. (8) using the two extreme DG values, along with their comparison with experimental data of Wu et al. [17] and Sheindorf et al. [16] (c).
0.5–4) and DG data of the first molecule leads to a dependence that can be described well by the isotherm q ¼ 4:3 103 C 0:583 . The variation of DG with pore occupancy cannot be ignored in an attempt to predict experimental measurements. When using DG of multiple molecules arrangement this procedure yields q ¼ 5:4 103 C 0:664 . The latter shows much higher adsorption, which is still smaller than values measured by Wu et al. [17] (q = 0.049 C0.295 in the range 0–200 g/m3), who provided the pore-size distribution used here, and by Sheindorf et al. [16] (q = 0.037 C0.371). The gap between experimental and calculated values can be partially bridged by using the lowest literature H value reported, as was done here, (H = 2.37 106 atm m3/mol) as opposed to the H = 6.3 107 [30,31]: note that adsorption on the smaller pores (12, 0) is quite close to the measured one (Fig. 5, H = 2.37 106 atm m3/mol) but this strong adsorption is diluted by the much weaker adsorption on the larger pores. Assuming a double-wall CNT, as opposed to the single-wall assumption used here, will lead to q values that are half of those above, but the adsorption will be somewhat stronger. This discussion suggests that more accurate values of Henry constants are required and that the DG values determined here are underestimates, even when we use values corresponding to high occupancy. One of the main messages of this work is the need to account for adsorbate–adsorbate interaction. Thus, the use of DG of the first molecule, while being simpler and reduces the required computation, is not representative. Due to the layer structure of packing and the wavy DG pattern with occupancy we can use the DG value of the most stable arrangement as the other limiting value. Evidently, the commonly
used slit-pore model will yield lower values. Assuming that DG ascends linearly with occupancy will yield a hysteresis-shaped adsorption curve [32].
4. Predicting adsorption isotherm with unknown PSD The evaluation of the adsorption isotherm using the approach outlined requires knowing the AC pore-size distribution. This is an experimental value which is often not reported. To cope with this problem we suggest now a methodology where we predict a single solute adsorption isotherm of a new component using the known adsorption isotherm of another component (learning component). To apply this method we will employ by the following steps 1. Using the measured isotherm of the component used for learning (q vs C) we may evaluate an isotherm for the pressure vs the fraction of occupancy (h) by translating.
P A ¼ C H;
h¼
q ; qm
where qm is the maximal adsorption capacity (maximal number of molecules that can adsorb MWadsorbent/MWCNT adsorbate). 2. Calculating the values of u* vs h and the evaluated value of DG* from our model. 3. Calculating dh/du* (by evaluating Dh/Du* for every point) and using Eq. (9) to find the value of n for the first component applying successive substitution.
451
H. Abir, M. Sheintuch / Journal of Colloid and Interface Science 342 (2010) 445–454
0.3
4. Using the estimated value of ‘n’ and calculated values of adsorption energy for the learning component we can estimate the pore-size distribution of the active carbon: using Eq. (4) we can evaluate NðDGÞ=a1 (for the learning component) and plot NðDGÞ=a1 vs DG2 =RT (for the new component), this will yield the curve NðDGÞ=a1 ¼ ða2 =a1 Þ expðn2 DG2 =RTÞ—allowing the evaluation of the parameter ‘n’ for the new component. 5. The adsorption isotherm for the new component may be evaluated now by the method outlined above (Eq. (8)).
0.25 0.2 0.15 0.1
To verify this approach we employ it for estimating m-cresol adsorption with data learned from phenol adsorption, and vice versa. Single-component adsorption of these two components was experimentally studied by Leitao and Serrao on Nuchar WA-activated carbon [18]. We later study the competitive adsorption of these two components. To predict the isotherms we study the adsorption of m-cresol on a shallow 7 Å long model of zigzag nanotubes of 8–18 Å in diameter (nanotubes of (12, 0)–(24, 0), Fig. 6a). The pore diameter was limited to the range of 8–18 Å as was done for phenol on the previous section. Realizing that the energy varies with occupancy we study DG for the case of the most stable arrangement, but limit it to the range (15, 0)–(17, 0) (Fig. 6b). The calculated adsorption capacity of m-cresol (based on phenol data) is somewhat higher than the experimental one (Fig. 7), while for phenol the calculated adsorption capacity is lower than the experimental one. The calculations are extremely sensitive to the evaluation of DG* as well as for uncertainty in the values of Henry constant and the pressure range taken. To translate the gas–solid equilibrium to liquid–solid equilibrium we used Henry’s law [30]
0.05 0 0.00
60.00
5. Competitive adsorption We outline here a formalistic calculation approach for a twocomponent adsorption that follows the single-component steps: a simple Langmuir isotherm is assumed for each pore size, then it is averaged for a certain pore-size distribution, and finally adsorbate–adsorbate interactions are incorporated. We study single-component and competitive gas-phase adsorption of phenol and m-cresol on (15, 0), (16, 0), and (17, 0) NT (11– 13 Å diameter) of pore of 7 Å in length. All models can adsorb up to three molecules. The molecules are oriented parallel to the pore wall to maximize the interaction energy. Analysis of the adsorption thermodynamic parameters indicates attraction between adsorbed molecules for a single-component adsorption as well as between the adsorbed molecules in the multicomponent adsorption; i.e., adsorption of the second molecule will find its optimal position next to the first causing a reduction in the system DG, filling three molecules in layer. For (15, 0) NT adding the third molecule causes an increase in the system DG due to steric interference, as was observed for the phenol single-component system (Fig. 2b). The adsorption of m-cresol is stronger than that of phenol and, therefore, on competitive adsorption the replacement of phenol
Hm-cresol ¼ 1:29 106 atm m3 =mol: Other inaccuracies, such as the evaluation of the PSD, also contribute to the inaccuracy of the evaluation of the isotherm. This evaluation is sensitive to the pressure range used. Within the pressure range used here, 0.01–0.13 Pa (corresponding to concentrations of 4–50 g/m3), the value estimated from phenol isotherm is n(phenol) = 0.08, which implies for m-cresol that n(m-cresol) = 0.06 while learning n from m-cresol yields n(m-cresol) = 0.055, which implies n(phenol) = 0.075. Using a wider pressure range of 0.001–0.1 Pa will lead to a higher value of ‘n’ (n(phenol) 0.133,n(m-cresol) 0.1). Another method that may be used in order to predict the adsorption isotherm of a new component is to base our calculations on known PSD of another AC, where the adsorption isotherm of the first component is known on both active carbon.
25 23 21 19 17 15 13 11 9 7 5
40.00
Fig. 7. Calculated weight adsorption isotherm for m-cresol, with PSD determined from experimental isotherm of phenol (nm-cresol = 0.059), and calculated weight adsorption isotherm for phenol with PSD determined from experimental isotherm of m-cresol (nphenol = 0.075) and measured weight adsorption isotherm for phenol and m-cresol according to Leitao and Serrao [18]. The calculated data are based on adsorption energy of the first molecule.
Hphenol ¼ 2:37 106 atm m3 =mol;
a
20.00
b 20
7
12
17
22
18 16 14 12 10 8 6 4
5
10
15
20
Fig. 6. Calculated DG (kcal/mol) vs inner pore diameter (Dp) for the first molecule adsorption in single solute adsorption of phenol, m-cresol, and para-bromophenol (PBP), (a) or for the most stable position in single solute adsorption of phenol or m-cresol (b), both in shallow (7 Å) pores. (The lines were added to guide the eye; the broken line denotes extrapolation.)
452
H. Abir, M. Sheintuch / Journal of Colloid and Interface Science 342 (2010) 445–454
with m-cresol is thermodynamically favorable and phenol may be rejected and replaced by m-cresol. 5.1. Multicomponent adsorption isotherms in a nanotube Consider a system where several components compete on the same adsorption sites. We will assume that each component in a single pore follows the Langmuir adsorption isotherm. We will also assume a linear dependence of DG on pore diameter for the singlecomponent adsorption, which is reasonable in the range of importance (Fig. 6a).
DGi ¼ DGi;0 ci Dp or DGi DGi;0 ¼ ci Dp ¼ Dg i
ð11Þ
Moreover, for the system studied here we assume c1 = c2 and that the PSD is exponential (Eq. (4)), so that
NðDgÞ ¼ a exp
nDGi0 ncDp nDg exp ¼ a0i exp RT RT RT
ð12Þ
and that in every pore diameter the occupancy follows the competitive Langmuir isotherm
hi ¼
K i Pi ; k P 1 þ K j Pj
ð13Þ
j¼1
where Ki is the coefficient of the single-component adsorption, Pi is the partial pressure, and Dp is the pore diameter. From the kinetic theory
a2i k expðDGi =RTÞ a2i k expðDGi0 =RTÞ expðDg i =RTÞ ¼ RT RT 0 ¼ K i expðDg i =RTÞ;
Ki ¼
ð14Þ
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi where a2i is adsorption site area, and ki ¼ h= 2pmi kT is the thermal de-Broglie wavelength of species i. The multicomponent multipore isotherm can be found in analogy to the monocomponent isotherm from solving the equation:
R1
Dg i =RT
hTi ¼
K 0 expðDg=RTÞP
i a expðnDg i =RTÞdðDg=RTÞ Pik 0 1þ K expðDg=RTÞPj j¼1 j R1 : 0 Dg =RT ai expðnDg=RTÞdðDg=RTÞ
ð15Þ
i
Dg i is the adsorption energy difference of the narrowest pore that can host a molecule i, this holds for all components if c is constant. For the case of two-component adsorption, assuming that c1 ¼ c2 then n1 ¼ n2 , this isotherm can be analyzed in a similar approach to the single-component isotherm. For n > 0 let us express v K 01 P1 þ K 02 P2 expðDg=RTÞ; the integral limits then are 0 0 v ¼ K 1 P1 þ K 2 P2 expðDg =RTÞ, which corresponds to Dp , and v* = 0, which corresponds to DG = 1. Then
R0 hT1 ¼
K 01 P1
v ½K 01 P1 þK 02 P2
R0
v
a1
v a 1
1þv
n n K 01 P1 þK 02 P 2
v
½
v
½
n n K 01 P 1 þK 02 P2
K 0 P1 nðv Þn ¼ 01 K 1 P1 þ K 02 P2
Z v 0
K1 K 01
vn 1þv
n
K1 K 01
n
dðv Þ v
dðv Þ v
dðv Þ:
Z v 0
vn 1þv
ð17Þ
Therefore, translating measured q1/q2 into h1h2, using known maximal capacity (qmi ), and plotting it vs C1/C2 should yield a straight line that should correspond to the competition coefficient between the two components as we test below. 5.2. Competitive adsorption of phenol and m-cresol We apply this approach to the bisolute system of phenol and m-cresol. Leitao and Serrao [18] studied this adsorption, as well as the individual component adsorption, by placing varying amounts of Nuchar WA AC into batch systems, using a constant initial concentration of 200 g/m3 phenol and 200 g/m3 m-cresol. After reaching equilibrium the concentrations were measured and weight adsorption isotherms were measured. To represent this form of experiment in our calculation we represented the data as having equilibrium in a solution containing the reported concentration of both components. Since the PSD is not reported in their paper, the method described in Section 4, was used in order to evaluate the PSD and the average value of the parameter ‘n’ from the single-component isotherm was used (nphenol = 0.08, nm-cresol = 0.059). The competitive adsorption was calculated using Eq. (16), over a range 0–200 g/m3 of concentration of each component (corresponding to pressures of 0–0.05 Pa phenol and 0–0.2 m-cresol), using the single solute adsorption energy as a function of pore diameter (Fig. 6). While the calculated DG dependence on pore diameter is similar (but not identical) we assume here that the slopes in Fig. 6 are identical so that c1 c2 . Two approaches were employed using for every pore size the DG of either the first molecule of each component or that of the most stable arrangement. While the evaluation of m-cresol adsorption, in the first-molecule approach, gave a rather good agreement with experimental results the evaluation of phenol adsorption is two orders of magnitude lower than the experimental results (Fig. 8a and c). The pore diameter was limited to the range of 8–18 Å, for both components. The approach based on the most stable position yields a much better description but still differences with experimental results are evident (Fig. 8a and b). Note that the data for m-cresol in this case was based on three points only and was linearly extrapolated to the entire pressure range (Fig. 6b). The differences between the calculated and the experimental results are due to the high sensitivity of competitive adsorption. To demonstrate this sensitivity let us consider differences in DG0 of 1.8 kcal/mol (DDG0 = DG2,0 DG1,0) and an inaccuracy of 0.1 kcal/mol in the evaluation of DDG; this will translate to an inaccuracy of about 38 mol/kcal in the evaluation of K 0i . 5.3. Competitive adsorption of phenol–PBP
ð16Þ
This reduces, of course, to the same integral for the single solute case (P2 = 0). Similarly, for the case of K 01 P1 K 02 P2, m reduces to v ffi K 01 PA1 expðDg=RTÞ and
hT1 ¼ nðv Þn
K 01 P 1 nðv Þn R v v n dðv Þ hT1 ½K 01 P1 þK 02 P2 0 1þv K 0 P1 ¼ 10 : ¼ K 0 P nðv Þn R v n v dðv Þ hT2 2 2 K 2 P2 ½K 01 P1 þK 02 P2 0 1þv
dðv Þ;
which is, again, the single-component case. For a two-component system we can evaluate a competition coefficient:
Before studying the competitive adsorption we estimated the adsorption isotherm of PBP, using the data for phenol adsorption; both components were studied by Sheindorf et al. [16] on granular activated carbon (Calgon Filtrasorb 400). Sheindorf et al. measured for phenol the isotherm q = 0.037 C0.371 which is in rather good agreement with the phenol adsorption isotherm measured by Wu et al. on (Fir wood) AC: q = 0.049 C0.295. We will assume that the same pore-size distribution that was measured by Wu et al. applies also for the Calgon Filtrasorb 400 AC. To predict the isotherms the adsorption of PBP was studied on the same AC model as was described in the previous section: a shallow 7 Å long model of zig-
453
H. Abir, M. Sheintuch / Journal of Colloid and Interface Science 342 (2010) 445–454
a
0.12
b
0.1 0.08 0.06 0.04 0.02 0 0
50
100
150
c1
200
0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0
0
100
200
0
1 1 1 2 1 3 1 4 1 5 0
50
100
150
200
Fig. 8. Measured bisolute weight adsorption isotherm of phenol and m-cresol according to Leitao and Serrao [18] (a), calculated bisolute isotherm according to Eq. (16), based on the adsorption energy of the most stable position (b) and based on the adsorption energy of the first molecule (c) (in b and c the parameter ‘n’ was calculated as average for the two components from the single-component isotherm of phenol; nphenol = 0.07).
zag nanotubes of 8–18 Å in diameter (nanotubes of (12, 0)–(24, 0)). The adsorption was limited for the case of a single molecule in every pore in order to simplify the calculations; the molecules are oriented parallel to the pore wall to maximize the interaction energy (Fig. 6a). Using PBP concentration in the range of 0–160 g PBP/m3 and DG data of the first molecule and the procedure described in 3.4 leads to a dependence that can be described well by the Freundlich-type adsorption isotherm q ¼ 0:267 C 0:117 g PBP/grC. This is a somewhat higher than the values measured by Sheindorf et al. [16] q ¼ 0:13 C 0:225 which was determined for a much narrower range of concentrations (Fig. 9). We examine now phenol adsorption in the presence of various constant concentrations of PBP. This system was studied by Sheindorf et al. [16]. Again we assume that slopes in Fig. 6 are close enough so that c1 c2 . The value of c was taken as an average value between that of phenol and PBP. The pore diameter was limited to the range of 8–18 Å, for both components. We will assume here that the PSD that was measured by Wu et al. can be applied in this case. For the single solute case we showed that using the PSD provided by Wu et al. [17] gave a rather good agreement with the experimental results of Sheindorf et al. [16], using adsorption data of the first molecule in every pore. Eq. (17) was solved over a range of phenol pressures (0.001– 0.02 Pa) for four different PBP concentrations (0, 0.1, 1, and 10 g/ m3 corresponding to 0, 9 106, 9 105, and 9 104 Pa). The results are shown in Fig. 10. The competition coefficient indeed varies linearly with Cphenol/ CBPB (Fig. 10) but the slope differs by several orders of magnitude from the experimental results of Sheindorf et al., and predicts a much weaker phenol adsorption. Much of the difference may be attributed to the use of the first-molecule energies and ignoring
0.5 0.45 0.4 0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 0
.
50
100
150
200
Fig. 9. Calculated weight adsorption isotherm for PBP according to Eq. (8) using the pore-size distribution reported by Wu et al. [17] for Fir wood AC, along with a comparison to experimental data for PBP adsorption of Sheindorf et al. [16].
adsorbate–adsorbate interaction and the use of parameters approximation. 6. Conclusions We use a molecular mechanics approach with universal force field parameters to calculate single- and multicomponent adsorption of phenol, para-bromophenol, and m-cresol from the gas or liquid phase on activated carbon. The carbon pores are modeled by carbon nanotubes of various pore diameters and different lengths. This calculation yields the Gibbs free energy change on
454
H. Abir, M. Sheintuch / Journal of Colloid and Interface Science 342 (2010) 445–454
b
0.25
q(Phenol) /q (PBP)
a
0.2 0.15 0.1 0.05 0 -35
-15
1. 02 1. 03 1. 04 1. 05 1. 06 1. 07 1. 08 1 02
1 00
1 02
1 04
Fig. 10. Evaluated energy distribution for single-component phenol and PBP assuming the PSD is like that measured by Wu et al. [17] (a), and calculated competition coefficient, defined as the ratio between qT1/qT2 and C1/C2 evaluated according to Eq. (17) of phenol/PBP, for different, constant PBP concentrations (b).
adsorption which is used to predict the adsorption isotherm following reaction rate theory. We provide a heuristic atomistic calculation explanation to the phenomena of Freundlich isotherms that are due to pore-size distribution and specifically due to the coupling of DG dependence on pore diameter and exponential distribution of pore sizes. That results in weaker adsorption as it progresses. We provide a methodology for the prediction of the adsorption isotherm in several cases:
Predicting a single solute adsorption isotherm based on a known (experimental) pore-size distribution of the active carbon.
Predicting a single solute adsorption isotherm based on a known adsorption isotherm of another component.
Predicting a multicomponent adsorption isotherm in the case of known pore-size distribution of the AC or when the pore-size distribution is evaluated from the single-component adsorption isotherm. The prediction of the adsorption isotherm by this method gives a rather good agreement with experimental results for the case of single-component data, when accounting for adsorbate–adsorbate interaction, especially when the pore-size distribution is known. As more uncertainty is added by the evaluation of the pore-size distribution the prediction becomes less accurate. This method for prediction of the adsorption isotherm is extremely sensitive to inaccuracy in the evaluation of the required parameters especially for the evaluation of DG0 and its dependence of adsorption due to exponential dependence on this parameter. The results are also sensitive to the range of pressure taken and to accuracy in the provided values of the pore-size distribution and Henry constant. Acknowledgment Work supported by the Grand Institute of Water Research.
References [1] R.C. Bansal, J.B. Donnet, F. Stoeckli, Active Carbon, Dekker, New York, 1988. [2] H.P. Boehm, Carbon 40 (2002) 145–149. [3] J.L. Figueiredo, M.F.R. Pereira, M.M.A. Freitas, J.J.M. Orfao, Carbon 37 (1999) 1379–1389. [4] P.J.F. Harris, A. Burian, S. Duber, Philos. Mag. Lett. 80 (2000) 381–386. [5] P.J.F. Harris, L. Zheng, S. Kazu, J. Phys. Condens. Matter 20 (2008) 362201. [6] A. Ricca, C.W. Bauschlicher, Chem. Phys. 324 (2006) 455–458. [7] J. Zhao, A. Buldum, J. Han, J.P. Lu, Nanotechnology 13 (2002) 195. [8] G.A. Sznejer, I. Efremenko, M. Sheintuch, AIChE J. 50 (2004) 596–610. [9] I. Efremenko, M. Sheintuch, Langmuir 21 (2005) 6282–6288. [10] T.J. Bandosz, M.J. Biggs, K.E. Gubbins, Y. Hattori, T. Iiyama, K. Kaneko, J. Pikunic, K.T. Thomson, Chem. Phys. Carbon (2003) 41–228. [11] B.H. Kim, G.H. Kum, Y.G. Seo, Korean J. Chem. Eng. (2003) 104–109. [12] E. Pantatosaki, D. Psomadopoulos, Th. Steriotis, A.K. Stubos, A. Papaioannou, G.K. Papadopoulos, Colloids Surf. 241 (1–3) (2004) 127–135. [13] E. Pantatosaki, A. Papaioannou, A.K. Stubos, G.K. Papadopoulos, Appl. Surf. Sci. 253 (13) (2007) 5606–5609. [14] E. Pantatosaki, G.K. Papadopoulos, J. Chem. Phys. 127 (16) (2007) 164723. [15] Ch. Sheindorf, M. Rebhun, M. Sheintuch, J. Colloid Interface Sci. 79 (1981) 136. [16] Ch. Sheindorf, M. Rebhun, M. Sheintuch, Water Res. 16 (1982) 357. [17] F.C. Wu, R.L. Tseng, R.S. Juang, J. Colloid Interface Sci. 283 (2005) 49–56. [18] A. Leitao, R. Serrao, Adsorption 10 (2005) 167–179. [19] R.J. Sips, Chem. Phys. 16 (1948) 490. [20] C.H. Yang, J. Colloid Interface Sci. 208 (1998) 379–387. [21] I. Efremenko, M. Sheintuch, Langmuir 22 (8) (2006) 3614–3621. [22] Yu.I. Matatov-Meytal, M. Sheintuch, Ind. Chem. Eng. Res. 36 (1997) 4374. [23] Yu.I. Matatov-Meytal, M. Sheintuch, Ind. Chem. Eng. Res. 39 (2000) 18. [24] A.K. Rappe, W.A.J. Goddard, Phys. Chem. 95 (1991) 3358–3363. [25] A.K. Rappe, C.J. Casewit, K.S. Colwell, W.A. Goddard, J. Am. Chem. Soc. 114 (1992) 10024–10035. [26] A.R. Leach, Molecular Modeling Principles and Applications, Pearson, HarlowEngland, 1996. [27] A.K. Rappe, C.J. Casewit, Molecular Mechanics across Chemistry, University Science Books, Sausalito, CA, 1997. [28] A. Heyden, T. Duren, F.J. Keil, Chem. Eng. Sci. 57 (2002) 2439–2448. [29] I. Chorkendorff, J.W. Niemantsverdriet, Concepts of Modern Catalysis and Kinetics, Wiley-VCH, New York, 2003. [30] V. Feigenbrugel, S. Calvé, P. Mirabel, F. Louis, Atmos. Environ. 38 (2004) 5577– 5588. [31] R.H. Perry, D.W. Green, Perry’s Chemical Engineer’s Handbook, McGrew-Hill, New York, 1997. [32] H. Abir, Research Thesis, Technion – Israel Institute of Technology, 2009.