Binary and Ternary Adsorption Kinetics of Gases in Carbon Molecular

5722

Langmuir 2003, 19, 5722-5734

Binary and Ternary Adsorption Kinetics of Gases in Carbon Molecular Sieves Huang Qinglin, S. Farooq,* and I. A. Karimi Department of Chemical & Environmental Engineering, National University of Singapore, 4 Engineering Drive 4, Singapore 117576 Received December 30, 2002. In Final Form: March 31, 2003

In two recent studies from this laboratory, a dual resistance model with chemical potential gradient as the driving force for transport and empirical correlations to account for the strong concentration dependence of the thermodynamically corrected transport parameters was proposed and verified for the uptake of oxygen, nitrogen, and their mixtures in the micropores of Bergbau Forschung and Takeda carbon molecular sieve (CMS) samples. In this paper, the same model has been further verified with integral uptakes of carbon dioxide, methane, their binary mixtures, as well as their ternary mixtures with nitrogen in the two CMS samples. Inadequacy of the Langmuir model and IAS-Langmuir combination for these gases has been discussed in detail. It has been demonstrated that, for these gases, the multisite Langmuir model is much more effective in fitting the unary equilibrium data and predicting the mixture isotherms in CMS.

Introduction 1

In a recent study from this laboratory, adsorption and diffusion of oxygen, nitrogen, carbon dioxide, and methane were measured in a Takeda carbon molecular sieve (CMS) sample (designated as Takeda I) over a wide range of pressures and temperatures. Similar measurements were also made for oxygen and nitrogen in a Bergbau-Forschung GmbH (BF) CMS sample and a second Takeda CMS sample (designated as Takeda II). A dual resistance model was shown to be the desirable unified approach that fitted the experimental results in the entire range covered in that study. As expected, the two transport parameters of the dual resistance model were observed to be functions of the adsorbed phase concentration. However, surprisingly, the functions were significantly stronger than those predicted from the use of a chemical potential gradient as the driving force for diffusion with constant intrinsic mobility. In other words, the thermodynamically corrected transport parameters were also strong functions of adsorbate concentration in the adsorbent. To account for the concentration dependence of the thermodynamically corrected transport parameters, an empirical but simple and effective procedure was proposed. The limiting (i.e., thermodynamically corrected) transport parameters, Dco and kbo, have generally been found to be independent of fractional coverage, θ, in zeolite adsorbents where the pore size is uniform. In CMS adsorbents, since micropore sizes are distributed and the pore connectivities are not fully understood, it appeared logical to assign each pore size its own characteristic Dco and kbo values. This argument combined with the fact that pores are filled in the order of increasing size requires Dco and kbo to be increasing functions of θ. The following forms were validated for unary diffusion * Corresponding author. Fax: 65-6779-1936. E-mail: chesf@ nus.edu.sg. (1) Huang Q. L.; Sundaram, S. M.; Farooq, S. Revisiting transport of gases in the micropores of carbon molecular sieves. Langmuir 2003, 19, 393.

( *(1 + β

Dco ) Dco* 1 + βp kbo ) kbo

b

θ 1-θ

) θ 1 - θ)

(1) (2)

As discussed earlier,1 the above forms satisfy the expected limiting behavior at θ f 0. βp and βb were obtained by fitting experimental Dc/Dco* vs θ and kb/kbo* vs θ data, respectively. Dc and kb are related to Dco and kbo by the Darken equation (Dc/Dco ) d ln c/d ln q) and its equivalent for barrier coefficient (kb/kbo ) d ln c/d ln q). The above hypothesis was experimentally verified with single component integral uptake data for oxygen and nitrogen in BF and Takeda I CMS. In a subsequent communication,2 further verification of eqs 1 and 2 were provided with more unary integral uptake and column dynamics results of oxygen and nitrogen in BF and Takeda I CMS. More importantly, the following multicomponent extensions of eqs 1 and 2 were proposed and validated with binary integral uptake and column dynamics experiments of oxygen-nitrogen mixtures in both the adsorbents

( (

(Dco)i ) (Dco*)i 1 +

(kbo)i ) (kbo*)i 1 +

θi

n

∑ i)1

βpi

n

1-

n

βbi ∑ i)1

1-

) )

(3)

θj ∑ j)1

θi n

(4)

θj ∑ j)1

where θi ) qi/qsi and i ) 1, 2, ..., n. The above multicomponent extensions are based on the assumption that the contributions of components in a multicomponent system are linearly additive. It is important to note that in both the studies, Dc/Dco* vs θ and (2) Huang Q. L.; Farooq, S.; Karimi, I. A. Prediction of pure and binary gas diffusion in carbon molecular sieves at high pressure. Submitted to AIChE J.

10.1021/la0270791 CCC: $25.00 © 2003 American Chemical Society Published on Web 06/11/2003

Kinetics of Gases in Molecular Sieves

Langmuir, Vol. 19, No. 14, 2003 5723

kb/kbo* vs θ were obtained from unary differential uptake measurements and no additional fitting parameter was involved when the proposed empirical model was applied to predict unary and binary integral uptake results. In continuation of our previous studies summarized above, the following advances are reported in this communication: (1) Effectiveness of the multisite Langmuir isotherm for representing unary equilibrium data and predicting mixture gas adsorption in carbon molecular sieves. (2) Additional verification of eqs 1 and 2 with unary integral uptake data of methane and carbon dioxide in both BF and Takeda I CMS. (3) Further verification of eqs 3 and 4 with binary uptakes of methane-carbon dioxide mixtures and ternary uptakes of methane-nitrogen-carbon dioxide mixtures in both BF and Takeda I CMS. Our choice of these mixtures was motivated by the fact that separating mixtures of methane with nitrogen and/or carbon dioxide is encountered in connection with natural gas upgrading. In the following sections, we will use Takeda CMS to mean Takeda I CMS. Experimental Section Pure component equilibrium isotherms for carbon dioxide and methane were measured on BF CMS sample up to about 8 bar pressure at several temperatures ranging from -10 to 70 °C. The temperature range was chosen to cover typical process operation for adsorption separation of these gases and to ensure reliable uptake measurements. The latter necessitated measurement at subzero temperature to slow the uptake rate for faster components like carbon dioxide and at fairly elevated temperature for slower diffusing methane to avoid unduly long duration of each run. In addition to covering typical process conditions, the pressure range also ensured significant coverage of adsorbate in the adsorbed phase. It should be highlighted that high-pressure diffusion data in CMS is not common in the literature. Differential uptakes of the two gases in BF CMS sample were also measured at several levels of adsorbent loading within the pressure and temperature ranges covered for equilibrium isotherms. In the present study, a constant volume apparatus was used to measure both pure component equilibrium isotherms and differential uptakes. A detailed account of this apparatus, experimental procedures, and data processing was provided in a previous publication,1 and hence repetition is avoided here. A differential adsorption bed (DAB) was used to measure single component, binary, and ternary integral uptakes in BF and Takeda CMS samples at 30 °C. Binary and ternary equilibrium isotherms for these gases were also measured in the DAB apparatus up to 10 bar pressure at 30 °C. A 10 ft × 1/8 in. column packed with Hayesep D 100/120 (Alltech Associates, Inc., Part no. 28301PC) was used for chromatographic analysis of nitrogen, methane, carbon dioxide, and their mixtures. Details of DAB apparatus and procedures for these measurements were essentially similar to those discussed elsewhere.2

Theoretical Section Equilibrium. The Langmuir model is used quite commonly to represent adsorption isotherm of a pure adsorbate. For a component i, it has the following form:

qi )

qsibici 1 + b i ci

(5)

where -∆Ui/RgT

bi ) bioe

(6)

The above model was used in a previous communication1 to fit individually the isotherms of oxygen, nitrogen, carbon

Figure 1. Measured equilibrium isotherms of carbon dioxide (solid symbols) and methane (open symbols) on BF CMS and fits of the Langmuir model. Solid and dotted lines are fits of the Langmuir isotherm model for the solid and open symbols, respectively.

dioxide, and methane on the Takeda CMS sample and oxygen and nitrogen isotherms on the BF CMS sample. Similar Langmuir fits for the isotherms of carbon dioxide and methane on BF CMS are given in Figure 1. The extracted Langmuir isotherm parameters and the corresponding regression residuals for all four gases are summarized in Table 1. It is quite clear that the Langmuir model is a poor fit for the carbon dioxide isotherm data compared to the other three gases. Since isotherm curvature can affect the concentration dependence of micropore transport parameters significantly, an accurate model for the isotherm data is important. Accurate representation of the individual isotherm is also important for reliable mixture equilibrium prediction. Since the Langmuir model is relatively unsatisfactory for carbon dioxide and methane, use of the extended Langmuir isotherm for modeling their binary equilibrium did not make sense. Moreover, recall that the extended Langmuir isotherm is thermodynamically consistent, only when the constituents of the mixture have the same qs.3 Since the individual qsi values for methane were ∼39% and ∼41% lower than those for carbon dioxide on BF CMS and the Takeda CMS, respectively, forcing qsi to be the same for both was clearly unadvisable. Therefore, it is clear that the extended Langmuir model is unsuitable for the binary equilibrium of carbon dioxide and methane in CMS. An alternative is to use other models such as ideal adsorbed solution (IAS) theory4 and multisite Langmuir model,5 which are thermodynamically consistent, even when the saturation capacities of the constituents are different. Multisite Langmuir Model. The multisite Langmuir model5 is a simple extension of the Langmuir model for single component and multicomponent adsorption equilibrium on microporous, homogeneous adsorbents. It assumes that an adsorbent has a fixed number (qs) of adsorption sites. An adsorbate molecule i, depending on its size and orientation in the adsorbed phase, occupies a certain number of these adsorption sites (ai). Therefore, qs ) qsiai for all i for thermodynamic consistency.3 (3) Rao M. B.; Sircar S. Thermodynamic consistency for binary gas adsorption equilibria. Langmuir 1999, 15, 7258. (4) Myers, A. L.; Prausnitz, J. M. Thermodynamics of mixed-gas adsorption. AIChE J. 1965, 11, 121. (5) Nitta, T.; Shigetomi, T.; Kruo-oka, M.; Katayama, T. An adsorption isotherm of multisite occupancy model for homogeneous surface. J. Chem. Eng. Jpn. 1984, 17, 39.

5724

Langmuir, Vol. 19, No. 14, 2003

Qinglin et al.

Table 1. Equilibrium Isotherm Parameters for Adsorption of Gases on BF CMS and Takeda CMS Langmuir isotherm adsorbent

bo × (cm3/mmol)

(-∆U) (kcal/mol)

qsi (mmol/cm3)

residual

bo × (cm3/mmol)

(-∆U) (kcal/mol)

a

qsi (mmol/cm3)

residual

O2a (92) N2a (82) CO2 (63) CH4 (41) O2 (76) N2 (70) CO2 (81) CH4 (53)

3.86 2.80 0.0067 0.51 3.48 4.02 0.22 0.73

3.90 4.26 8.71 6.22 3.98 4.01 6.73 5.97

3.80 2.87 5.30 3.21 4.34 3.61 6.30 3.74

0.112 0.058 3.36 0.204 0.064 0.048 6.03 0.298

1.13 0.804 0.078 0.337 1.75 0.740 0.556 0.502

4.26 4.61 7.68 6.31 4.08 4.72 6.51 5.98

2.89 3.36 2.99 3.47 3.04 3.51 3.04 3.55

7.20 6.20 6.96 6.00 8.55 7.40 8.55 7.32

0.030 0.045 0.626 0.043 0.032 0.050 0.599 0.082

BF CMS

Takeda CMSa

a

multisite Langmuir isotherm

adsorbate (data points)

103

103

Data taken from Huang et al.1

For multicomponent equilibrium, the multisite Langmuir model is

qi/qsi

bici )

(7)

n

(1 -

qi/qsi)a ∑ i)1

πA ) RgT

i

and, for unary adsorption, eq 7 reduces to

bici )

qi/qsi

(8)

(1 - qi/qsi)ai

where ai and qsi are independent of temperature and bi has the same temperature dependence as eq 6. For ai ) 1, eq 7 reduces to the extended Langmuir model and eq 8 to the Langmuir model. Nitta et al.5 derived the multisite Langmuir model using a statistical thermodynamic method. Sircar,6 however, derived it simply from the classical Langmuirian kinetics. The mass action law was used to describe the reversible multicomponent adsorption process. The uptake rate of component i in a multicomponent system may be written as

d(qi/qsi) dt

ponent i at the same spreading pressure, π, and the same temperature as for the adsorbed mixture. The Gibbs adsorption isotherm approach yields the spreading pressure π as follows

∫0P

o

i

qi dPi ) constant Pi

(11)

where A is the surface area per unit mass of adsorbent, Pi ()Pyi) is the partial pressure of component i, qi is the measured adsorption isotherm of component i. When the single component isotherm is represented by the multisite Langmuir model, eq 8 is substituted in eq 11 together with the relation ci ) Pi/RgT to obtain

(

)

qi° πA ) Π ) (1 - ai)qi° - aiqsi ln 1 RgT qsi

(12)

where qi° is the adsorbed amount of component i at pressure Pi° and follows the multisite Langmuir isotherm

biPi° )

{RgTqi°}/{qsi} ({1 - qi°}/{qsi})ai

(13)

n

ad

) ki ci(1 -

qi/qsi)a - kideqi/qsi ∑ i)1 i

(9)

where kiad and kide are specific adsorption and desorption reaction rate constants, respectively. When the system reaches a dynamic equilibrium, namely, d(qi/qsi)/dt ) 0, the above equation reduces to eq 7 for multicomponent or to eq 8 for single component with bi ) kiad/kide. Ideal Adsorption Solution (IAS) Theory. IAS theory with the Langmuir model for unary data was described in detail and successfully applied to the prediction of binary isotherm of oxygen and nitrogen in our previous study.2 Since as discussed earlier, the Langmuir isotherm seemed inadequate for the present study, IAS theory with the multisite Langmuir model for unary data was used for predicting multicomponent equilibrium. IAS theory uses Raoult’s law to describe mutlicomponent adsorption equilibrium

On the adsorbed surface, the sum of mole fraction of each gas must be unity n

xi ) 1 ∑ i)1

(14)

Combining eqs 10, 13, and 14, we get n

∑ i)1

(

)

biPyi(1 - qi°/qsi)ai RgTqi°/qsi

)1

(15)

where xi and yi are mole fractions of component i in adsorbed phase and gas phase, respectively, P is the total pressure of gas phase, Pi° is the equilibrium gas-phase pressure corresponding to the adsorption of pure com-

Equations 12 and 15 have a set of (n + 1) nonlinear equations with (n + 1) unknowns, namely, qi° (i ) 1, 2, ..., n) and Π, which can be solved simultaneously using the nonlinear equation solver subroutine NEQNF in IMSL. NEQNF uses a modified Powell hybrid algorithm and a finite-difference approximation to calculate the Jacobian. After qi° is known, we can calculate Pi° from eq 13 and then xi from eq 10. The calculation of the total adsorbed amount, qT, is based on the fundamental definition of the ideal adsorbed solution7

(6) Sircar, S. Influence of adsorbate size and adsorbent heterogeneity on IAST. AIChE J. 1995, 41, 1135.

(7) Do, D. D. Adsorption analysis: Equilibria and kinetics; Imperial College Press: London, 1998.

Pyi ) Pioxi

(10)


1

n

)

qT

Langmuir, Vol. 19, No. 14, 2003 5725

xi

∑ i)1 q ° i

Therefore, the adsorbed amount for each component in a multicomponent mixture is obtained as

qi ) xiqT

[

]

∂q ji ∂cip ∂2cip 2 ∂cip + (1 - p) ) pDp + p ∂t ∂t ∂R2 R ∂R

)0

|

∂cip pDp ∂R

R)Rp

) kf(ci - cip|R)Rp)

}

yiPsys RgTsys

(22)

where ciim is an imaginary gas-phase concentration in equilibrium with the adsorbed phase concentration. It is imaginary because there is no gas inside the microparticle. This leads to the following material balance equation for component i in the micropores

{(

)}

qi ∂ciim ∂qi 1 ∂ 2 ) 2 r (Dco)i im ∂t ∂r r ∂r c i

(23)

In the above equation, the imaginary gas-phase concentration, ciim, is different in the single component, binary, and ternary systems and calculated from the appropriate form of the multisite Langmuir model. Boundary conditions for the microparticle balance are

|

∂ciim )0 ∂r r)0 3(Dco)i ∂ciim ∂ciim ) (kbo)i (q * - qi)|r)rc rc ∂r r)rc ∂qi i

|

}

(24)

In eqs 23 and 24, (Dco)i and (kbo)i are given by eqs 3 and 4, respectively. For a single component system, the boundary condition for component i at r ) rc is derived by substituting eq 8 in eq 24

)

r)rc

(18)

1 + (ai - 1)θi (kbo)i(θi* - θi)|r)rc bi(1 - θi)ai+1 (25)

In the case of a binary system, substituting the binary form of eq 7 in eq 24, the boundary conditions for components A and B at r ) rc are obtained

(19)

|

3(Dco)A ∂cAim rc ∂r

) r)rc

1 × bA(1 - θA - θB)aA+1

[kAA(θA* - θA) + kAB(θB* - θB)]|r)rc

In the above equation, ci is related to the operating pressure of the system by the following equation

ci )

qi ∂ciim Ji ) -(Dco)i im ∂r ci

|

Boundary conditions for macropore balance

R)0

The diffusion flux of component i to the microparticles is derived from the chemical potential theory by introducing an imaginary gas-phase concentration8

3(Dco)i ∂ciim rc ∂r

Mass balance in macropore for component i

|

(21)

(17)

Integral Uptake. The model equations for unary, binary, and ternary integral uptakes in adsorbent particles in DAB experiments are presented in this section. The multisite Langmuir model is used to represent the single component isotherms and predict the multicomponent isotherms. The following assumptions are used in the model development: (1) The ideal gas law applies and the system is isothermal. (2) In view of the high flow rate of adsorbate gas and negligible thickness of the bed in DAB measurements, it is reasonable to neglect concentration gradient and pressure drop between inlet and outlet of the adsorber. (3) The fluid and adsorbent solid phases are linked through an external film resistance. The external film resistance is relatively small and may be neglected. However, this approach is often numerically more advantageous than the alternative approach of applying the equilibrium boundary condition at the solid surface. (4) Molecular diffusion dominates in the macropores. (5) Both macro- and microparticles are spherical. (6) The transport behavior in the micropores is viewed as a series combination of barrier resistance confined at the micropore mouth followed by pore diffusional resistance in the interior of the micropores. (7) The chemical potential gradient is considered to be the driving force for diffusion across the micropore mouth and in the micropore interior. (8) The limiting micropore transport parameters are viewed as increasing functions of adsorbent loading according to eqs 3 and 4. Subject to these assumptions, the following equations constitute the model:

∂cip ∂R

∂q ji 3 ) - Ji|r)rc ∂t rc

(16)

|

3(Dco)B ∂cBim rc ∂r

(20)

where Psys and Tsys are operating pressure and temperature of the differential adsorption bed, respectively, yi is mole fraction of component i in the feed. q j i in eq 18 is the average adsorbed concentration of component i in the micropore. Mass balance at the micropore surface is given by

) r)rc

1 × bB(1 - θA - θB)aB+1

[kBA(θA* - θA) + kBB(θB* - θB)]|r)rc

}

(26)

where (8) Hu, X.; Do, D. D. Multicomponent adsorption kinetics of hydrocarbons onto activated carbon. Chem. Eng. Sci. 1993, 48, 1317.

5726


kAA ) (kbo)A[1 + (aA - 1)θA - θB] kAB ) (kbo)AaAθA kBA ) (kbo)BaBθB kBB ) (kbo)B[1 - θA + (aB - 1)θB]

Qinglin et al.

}

(27)

}

For a ternary system, substituting the ternary form of eq 7 in eq 24, the boundary conditions for components A, B, and C at r ) rc become

|

3(Dco)A ∂cAim rc ∂r

)

r)rc

1 × bA(1 - θA - θB - θC)aA+1

[kAA(θA* - θA) + kAB(θB* - θB) + kAC(θC* - θC)]|r)rc

|

3(Dco)B ∂cBim rc ∂r

)

r)rc

1 × bB(1 - θA - θB - θC)aB+1

[kBA(θA* - θA) + kBB(θB* - θB) + kBC(θC* - θC)]|r)rc

|

3(Dco)C ∂cCim rc ∂r

)

r)rc

1 × bC(1 - θA - θB - θC)aC+1

[kCA(θA* - θA) + kCB(θB* - θB) + kCC(θC* - θC)]|r)rc

(28)

where

kAA ) (kbo)A[1 + (aA - 1)θA - θB - θC] kAB ) (kbo)AaAθA kAC ) (kbo)AaAθA kBA ) (kbo)BaBθB kBB ) (kbo)B[1 - θA + (aB - 1)θB - θC] kBC ) (kbo)BaBθB kCA ) (kbo)CaCθC kCB ) (kbo)CaCθC kCC ) (kbo)C[1 - θA - θB + (aC - 1)θC]

} } }

∫01 3ciPζ2 dζ + (1 - P) ∫01 3qj iζ2 dζ qiP*

4

fmin ) (30)

(31)

(32)

where

q ji )

∫01 3qiη2 dη

qip* ) pcio + (1 - p)qi*

Results and Discussion Unary Equilibrium. Single component equilibrium parameters for the multisite Langmuir model were extracted by simultaneous nonlinear regression of oxygen, nitrogen, carbon dioxide, and methane data at all available temperatures with the constraint that aiqsi was equal for all four gases in a given adsorbent. The following nonlinear objective function was used for data regression

(29)

In the above equations, θi ) qi/qsi. θi* in eqs 25, 26, and 28 is the fractional coverage of the adsorption sites corresponding to qi*, which is the adsorbed amount in equilibrium with gas-phase concentration and can be calculated from the corresponding single component, binary, or ternary multisite Langmuir model. The nonlinear equation solver subroutine NEQNF in IMSL was used to calculate qi*. The fractional uptake was calculated from volume integration of the concentration profiles in the macropores and micropores

P qip ) qip*

gas concentration of component i and ζ ()R/RP) and η ()r/rc) are the dimensionless radial distances along macroparticle and microparticle, respectively. Model Solution. The model equations were written in the dimensionless forms and then discretized in space by the method of orthogonal collocation. This reduced the coupled partial differential equations to a set of coupled algebraic and ordinary differential equations. The ordinary differential equations were integrated in the time domain using Gear’s stiff variable step integration routine provided in the FORSIM package to obtain compositions in the gas and solid phases as functions of time at the collocation points along the macroparticle and microparticle radii. The nonlinear boundary equations at r ) rc were solved by using the nonlinear equation solver subroutine NEQNF given in the IMSL package. The number of collocation points is important in the accuracy of the numerical solution. The optimum numbers were determined by increasing collocation points until further change did not affect the solutions any more. In this study, 16 collocation points along both macroparticle and microparticle radii were used.

(33) (34)

In the above equations, qip and qip* are the adsorbed amount at a certain time and at equilibrium, respectively, for component i based on particle volume. cio is the feed

n

[(qjexp)k - (qjtheo)k]2 ∑ ∑ k)1 j)1

(35)

where k and j are the indexes for adsorbates and data points, respectively, superscripts exp and theo denote the experimental and theoretical equilibrium results. For each adsorbent, the minimization (eq 35) was carried out with respect to the experimental data for four gases measured over a wide pressure range at several temperatures by varying the 13 parameters, ai, bio, ∆Ui, and qs. The above nonlinear optimization problem was solved using the minimization subroutine DBCONF available in IMSL. DBCONF is a double precision subroutine and uses a quasi-Newton method and a finite-difference gradient to minimize a function of N variables subject to bounds on the variables. The theoretical equilibrium adsorbed amount, qjtheo, was calculated from eq 8 using the nonlinear equation solver subroutine DNEQNF in IMSL. DNEQNF is a double precision subroutine. During the minimization, initial guesses for the unknowns were systematically varied to ensure that the smallest possible residual was reached. The optimized parameters for the multisite Langmuir model for all four gases and the corresponding residuals are also given in Table 1. In addition, fits of the multisite Langmuir model to oxygen, nitrogen, carbon dioxide, and methane isotherm data on BF and Takeda CMS samples are shown in Figure 2. The isotherms of all four gases are represented very well by the multisite Langmuir model, which is also evident from consistent improvement in the residuals. The improvement in the residual is particularly remarkable for carbon dioxide. Compared to the Langmuir model fits in Figure 1 and those in our previous communication,1 significant improvement of multisite Langmuir model fit is obvious even from visual inspection in the high-pressure region. The number of adsorption sites (a values in Table 1)


Langmuir, Vol. 19, No. 14, 2003 5727

Figure 2. Measured equilibrium isotherms of pure gases on BF and Takeda CMS adsorbents and fits of the multisite Langmuir model: (a) oxygen (solid symbols) and nitrogen (open symbols) in BF CMS; (b) oxygen (solid symbols) and nitrogen (open symbols) in Takeda CMS; (c) carbon dioxide (solid symbols) and methane (open symbols) in BF CMS; (d) carbon dioxide (solid symbols) and methane (open symbols) in Takeda CMS. Solid and dotted lines are fits of the multisite Langmuir model for the solid and open symbols, respectively.

occupied by each adsorbate molecule in the two CMS samples is in the order oxygen e carbon dioxide < nitrogen < methane. It is also apparent from Table 1 that saturation capacities, qsi, of oxygen, nitrogen, and methane obtained from the multisite Langmuir model fit are significantly larger than those obtained from the Langmuir model fit. To make a quick check on the practicality of the very large saturation capacities obtained from the multisite Langmuir model, oxygen and nitrogen isotherms were measured at -35 °C in the pressure range of 10-19 bar. From the results shown in Figure 3, it is clear that the adsorbed amounts exceed the saturation capacities obtained from the Langmuir model fits. On the other hand, the predictions based on multisite Langmuir model parameters are indeed very good. Although the saturation capacities obtained from the multisite Langmuir model are still quite far from the data shown in Figure 3, we can at least conclude that the values may not be impractical. Prediction of Binary and Ternary Equilibrium. A logical start for this discussion would be to look at the improvement that a mulitsite Langmuir model can bring to the prediction of oxygen-nitrogen mixture equilibrium data presented in a previous communication,2 where IAS theory, with an individually fitted Langmuir model, was adjudged to be somewhat more accurate than the extended Langmuir isotherm. Representative comparison of the IAS theory results with multisite Langmuir model predictions for oxygen-nitrogen mixture adsorption on BF CMS is shown in Figure 4. Clearly, the oxygen data are much better predicted by the multisite Langmuir model. To a

Figure 3. Measured equilibrium adsorbed amount of oxygen and nitrogen in Takeda CMS at -35 °C are compared with predictions of multisite Langmuir model and Langmuir model using parameters in Table 1. Solid and dotted lines are predictions of multisite Langmuir and Langmuir models, respectively.

lesser extent, a similar conclusion is also apparent for nitrogen. Very small improvements were also noticed for the data on Takeda CMS, which is in line with the small improvement in the residuals of the single component data fits shown in Table 1. The binary experimental isotherm data of a methanecarbon dioxide mixture on both BF and Takeda CMS samples are compared with the following predictions in Figure 5: (i) multicomponent multisite Langmuir model, (ii) IAS theory using Langmuir model individually fitted to unary data, and (iii) IAS theory using multisite Langmuir model fitted to unary data. The single compo-

5728


Qinglin et al.

Figure 4. Binary isotherms for oxygen and nitrogen in BF CMS at 1.5 °C.

Figure 6. Ternary isotherms for carbon dioxide, methane and nitrogen in (a) BF CMS and (b) Takeda CMS at 30 °C. Predictions of multicomponent multisite Langmuir model and IAS theory using multisite Langmuir model fitted to unary data are completely indistinguishable.

Figure 5. Binary isotherms for carbon dioxide and methane in (a) BF CMS and (b) Takeda CMS at 30 °C. Predictions of multicomponent multisite Langmuir model and IAS theory using multisite Langmuir model fitted to unary data are completely indistinguishable.

nent isotherm parameters in Table 1 were directly used in the above predictions. The deviation between experimental results and predictions by IAS theory with Langmuir model is indeed very large. The reason is obviously the unsatisfactory individual fit of carbon dioxide isotherm data by the Langmuir model. This clearly suggests that IAS theory with Langmuir isotherm for the individual components is not suitable for the prediction of methane-carbon dioxide mixture equilibrium. However, both multicomponent multisite Langmuir model and IAS theory using multisite Langmuir model fitted to unary data are able to provide reasonable prediction of the binary experimental equilibrium results. In fact, the two predictions are indistinguishable. Sircar6 has shown that the prediction of multicomponent equilibrium using IAS theory based on multisite Langmuir model for the individual components may deviate significantly from the prediction of multicomponent multisite Langmuir model when the adsorbate sizes (ai) differ substantially. In his study, he used isotherms of CF2CL2-NH3 and NH3-CO2

mixtures on BPL carbon. The adsorbate size parameter (ai) values for CF2CL2, NH3, and CO2 were 6.5, 1.53, and 3.0, respectively. Since methane and carbon dioxide are much closer in size, such deviation was not observed in our study. The ternary experimental isotherm data of methanenitrogen-carbon dioxide mixture on both BF and Takeda CMS samples at 30 °C are compared with the following predictions in Figure 6: (i) multicomponent multisite Langmuir model and (ii) IAS theory using multisite Langmuir model fitted to unary data. The two model predictions are also indistinguishable in this case. It is interesting to note that in both the adsorbents, the predicted values for carbon dioxide are marginally lower, while those for methane are higher by similar extent. Since such trends were not clear for the binary results, it may not be appropriate to comment on these small deviations without further study. At this stage, it is reasonable to conclude that the two models can predict the ternary experimental equilibrium results reasonably well. From the above discussion, it is clear that both multicomponent multisite Langmuir model and IAS theory with multisite Langmuir model can satisfactorily predict the mixture equilibria examined in this study. However, since IAS theory adds additional complexity to calculations,2 only the multisite Langmuir model was used in the kinetic studies in this paper. Single Component Differential Uptake. The differential uptakes of CO2 and CH4 in the BF CMS sample were measured volumetrically over a wide range of adsorbent loadings. The uptake results and optimum fits of the pore, barrier, and dual models are shown in Figure 7. The barrier and the pore models are actually two extreme cases of the dual model. The dual model solution


Langmuir, Vol. 19, No. 14, 2003 5729

Figure 7. Single component differential uptakes of carbon dioxide and methane in BF CMS (a) at different temperatures at low coverage and (b) at various loadings at a given temperature.

reduces to that of the barrier model when a large value is assigned to the micropore diffusivity and vice versa. In our previous communication,1 the contribution from valve resistance in the constant volume apparatus was shown to be negligible. It is clear from Figure 7 that the barrier model fits the early part of uptake data but fails in the later part. On the other hand, the solution of the pore model is distinctly different from the shape of the experimental uptake in the early part but fits the later part of the data well. Undoubtedly, the dual model gives the best fit of the entire uptake data. The observations are entirely consistent with the oxygen and nitrogen uptake results in two of the Takeda and one of the BF CMS samples and those of methane and carbon dioxide in one of the Takeda CMS samples reported1 earlier from this laboratory. The limiting transport parameters of the dual model for CO2 and CH4 in BF CMS were extracted by minimizing the residual of model fit to the experimental data at low surface coverage (θ f 0). The limiting transport parameters and transport activation energies for O2, N2, CO2, and CH4 in BF and Takeda CMS are listed in Table 2. The transport parameters of the four gases in BF and Takeda CMS shown in Table 2 are apparently in the order oxygen > carbon dioxide > nitrogen > methane. However, the sizes of these molecules based on Lennard-Jones values are in the order oxygen (3.43 Å) < nitrogen (3.68 Å) < methane (3.82 Å) < carbon dioxide (4.00 Å). These values have been reported by Hirschfelder et al.9 Reid and Thomas10,11 have also pointed out this discrepancy between (9) Hirschfelder, J. O.; Curtiss, C. F.; Bird, R. B. Molecular theory of gases and liquids; Wiley: New York, 1954. (10) Reid, C. R.; Thomas, K. M. Adsorption rates of gases on a carbon molecular sieve used for air separation: linear adsorptives as probes for kinetic selectivity. Langmuir 1999, 15, 3206.

observed kinetic rates in CMS and corresponding LennardJones dimensions. It is important to note that the LennardJones diameters used by Reid and Thomas are somewhat different (oxygen (3.46 Å), nitrogen (3.64 Å), methane (3.80 Å), and carbon dioxide (3.30 Å)) and were taken from a different source.12 It can be seen that the values for carbon dioxide are significantly different, whereas the values are very close for the other gases. Reid et al.13 have also reported kinetic data for neon (2.75 Å), argon (3.40 Å), and krypton (3.60 Å) and concluded that even for these spherical molecules the Lennard-Jones diameters cannot explain the order of kinetic data (argon > neon > krypton). In a more recent publication, Reid and Thomas11 have reported molecular sizes in three dimensions obtained from zero integral neglect of differential overlap (ZINDO)14 calculation method for linear, planar, and tetrahedral molecules and have further confirmed that the observed kinetics cannot even be explained by these alternative molecular dimensions. In contrast, it is observed from this study that the transport coefficients of oxygen, carbon dioxide, nitrogen, and methane show an inverse trend with respect to the a values obtained from the combined regression of the multisite Langmuir model to their unary equilibrium data. (11) Reid, C. R.; Thomas, K. M. Adsorption kinetics and size exclusion properties of probe molecules for the selective porosity in a carbon molecular sieve used for air separation. J. Phys. Chem. B 2001, 105, 10619. (12) Armor, J. N. Carbon molecular sieves for air separation. In Separation Technology; Vansant, E. F., Ed.; Elsevier Science B.V.: Amsterdam, 1994; p 163. (13) Reid, C. R.; O’koye, I. P.; Thomas, K. M. Adsorption of gases on carbon molecular sieves used for air separation. Spherical adsorptives as probes for kinetic selectivity. Langmuir 1998, 14, 2415. (14) Webster, C. E.; Drago, R. S.; Zerner, M. C. Molecular dimensions for adsoptives. J. Am. Chem. Soc. 1998, 120, 5509.

5730


Qinglin et al.

Table 2. Transport Parameters of Gases for the Dual Model in BF CMS and Takeda CMS adsorbent BF CMS

adsorbate O2a N2a CO2 CH4

Takeda CMSa

O2 N2 CO2 CH4

a

temp (K)

Dco/rc2 × 103 (s-1)

302.15 267.65 253.15 302.15 283.15 275.15 303.15 283.15 263.15 343.15 323.15 303.15 267.65 253.15 302.15 283.15 273.15 303.15 283.15 263.15 333.15 303.15 283.15

6.84 3.27 1.85 0.43 0.18 0.13 1.59 1.01 0.66 0.018 0.0084 0.0025 4.78 2.80 0.58 0.22 0.13 2.10 1.06 0.42 0.0077 0.0026 0.0008

Dco′/rc2 (s-1)

Ed (kcal/mol)

βp

5.70

4.01

0

83.94

7.32

2.28

0.50

3.47

0

63.12

10.24

55.63

4.98

0

706.31

8.42

2.28

81.43

6.35

0

2.77

8.43

5.52

5.52

kbo × 102 (s-1)

kbo′ (s-1)

Eb (kcal/mol)

βb

24 7.57 4.57 1.06 0.51 0.44 2.64 1.14 0.43 0.065 0.021 0.0076 8.1 4.45 0.88 0.39 0.26 5.72 1.23 0.57 0.017 0.0033 0.0012

919.38

5.01

9.17

121.61

5.62

7.93

4180

7.21

7.81

7310

11.10

6.06

28.62

5.57

9.17

819.83

6.88

7.93

7.25

7.81

9.85

6.06

5607.6 468.74

Data taken from Huang et al.1

Concentration Dependence of Transport Parameters. When the adsorption equilibrium is represented by multisite Langmuir isotherm, the concentration dependence of micropore diffusivity according to the chemical potential gradient theory as the driving force takes the following form

Dc d ln c 1 + θ(a - 1) ) ) Dco d ln q 1-θ

(36)

The above equation reduces to Darken’s equation when a ) 1. In this derivation, the limiting diffusivity, Dco, is assumed to be independent of θ. Starting from the same principle, a similar equation can also be derived for the barrier coefficient

kb 1 + θ(a - 1) ) kbo 1-θ

(37)

In the above equation, the limiting barrier coefficient, kbo, is also assumed to be independent of θ. The dual transport parameters of oxygen, nitrogen, methane, and carbon dioxide in both BF and Takeda CMS samples obtained from uptake results measured at various levels of adsorbent loadings are plotted as a function of surface coverage (θ) in Figure 8. Plots of eqs 36 and 37 are also shown in Figure 8. The individual multisite Langmuir saturation capacity given in Table 1 was used to calculate the surface coverage. It is obvious that both the micropore transport parameters are increasing functions of adsorbed phase concentration. For carbon dioxide and oxygen, concentration dependence of micropore diffusivity approximately follows eq 36. In other cases, the dependence is stronger than that expected from eq 36 or eq 37. For oxygen, nitrogen, and carbon dioxide, concentration dependence of the barrier coefficient is stronger than that for micropore diffusivity. The two trends are, however, comparable for methane. These observations are very similar to those in our previous communication,1 in which single component equilibrium was represented by the Langmuir model. The similarity confirms that the very strong concentration dependence of transport parameters

in CMS samples is not an artifact of the theoretical isotherm model chosen. Combining the empirical equations (eqs 1 and 2) proposed to account for the concentration dependence of the thermodynamically corrected transport parameters with eqs 36 and 37 above, the following equations are obtained

Dc 1 + θ(a - 1) θ ) 1 + βp Dco* 1-θ 1-θ

(

kb 1 + θ(a - 1) θ ) 1 + βb kbo* 1-θ 1-θ

(

)

(38)

)

(39)

From the above equations, it is clear that the value of Dc/Dco* (or kb/kbo*) depends on both βp (or βb) and the equilibrium parameter (a). βp was extracted separately for each component by optimizing the fit of eq 38 to the combined experimental Dc/Dco* vs θ data from both adsorbents shown in Figure 8. Another fitting parameter, βb, was similarly obtained by optimizing the fit of eq 39 to the combined experimental kb/kbo* vs θ data from both adsorbents. The fits are also shown in Figure 8, and the extracted values are given in Table 2. In the following sections, two versions of the dual model will be used to analyze the integral uptake results. The dual model with βiP ) βib ) 0 will be referred to as Dual Model 1, while that with βiP and βib from Table 2 will be called Dual Model 2. It should be highlighted that Dual Model 1 assumes a chemical potential gradient as the driving force for diffusion with constant (pore size independent) limiting transport parameters. On the other hand, Dual Model 2 allows for pore size dependence of the thermodynamically corrected transport parameters according to eqs 3 and 4. Single component equilibrium parameters (multisite Langmuir isotherm) in Table 1 and transport parameters in Table 2 were directly used in model calculations. Single Component Integral Uptake. The integral uptakes of carbon dioxide and methane in BF and Takeda CMS were measured at 30 °C. In these measurements, the adsorbent particles were subjected to a step change


Langmuir, Vol. 19, No. 14, 2003 5731

Figure 8. Concentration dependence of the two transport parameters of the dual model in the two CMS samples. The dotted lines are the plots of eq 36 or eq 37. The continuous lines are the optimized fits of eq 38 and eq 39 with βp and βb as the fitting parameters, respectively. Thick and thin lines are for the solid and open symbols, respectively. In many cases, the thick and the thin lines completely overlap.

of 0-3 bar for carbon dioxide and 0-5 bar for methane. The results are compared with the dual model predictions in Figure 9. Although the pressure step change is not large, the impact of the stronger concentration dependence is still quite significant, which is evident from the difference between the two models. Dual Model 2 provides an

excellent prediction of the experimental results in both the CMS samples. Dual Model 1 is able to capture the trends of experimental results qualitatively, but the quantitative deviation is large. It is clear from Table 2 that βp values for all the four gases are smaller than the corresponding βb values. While this may suggest that at high loading the contribution of the barrier resistance to

5732


Qinglin et al.

Figure 9. Single component integral uptakes of carbon dioxide and methane in BF and Takeda CMS samples at 30 °C.

transport becomes negligible, it will nevertheless remain important in process calculations, as directly evident from the Dual Model 2 results for carbon dioxide for which βp ) 0. Binary Integral Uptake. The proposed multicomponent extensions (eqs 3 and 4) of the empirical equations to account for the concentration dependence of thermodynamically corrected micropore transport parameters have been validated with binary oxygen-nitrogen integral uptakes in the two CMS samples in our previous communication.2 In view of the improved fits of the oxygen and nitrogen unary equilibrium data and improvement in their binary equilibrium prediction by the mulitisite Langmuir model, its impact on the prediction of binary uptake of oxygen-nitrogen mixture was also investigated. Representative results are shown in Figure 10. The two isotherm models produce practically the same results for oxygen. For nitrogen, while the IAS theory does well in predicting the latter part of the uptake, the multisite Langmuir isotherm is able to strike a good balance in the entire range. To further investigate the effectiveness of eqs 3 and 4 in other systems, codiffusional uptakes of methanecarbon dioxide (70:30) mixture in BF and Takeda CMS samples were measured for a step change of 0-10 bar at 30 °C. Comparisons between the experimental results and model predictions are shown in Figure 11. The roll-up in carbon dioxide uptake is predicted by both the models. Carbon dioxide diffuses faster than methane into the micropores. Amount of carbon dioxide adsorbed initially exceeds the binary equilibrium limit corresponding to the composition of the bulk gas. The roll-up of carbon dioxide is this excess amount displaced by the slower diffusing methane. Since both the models use the same isotherm equations, it is not surprising that both of them predict the roll-up. The location and extent of roll-up are, however, determined by the counterexchange kinetics of the two species, and in this regard the Dual Model 2, other than

Figure 10. Binary integral uptakes of oxygen and nitrogen in BF CMS at -10 °C (mixture, O2:N2 ) 50:50; pressure step, 0-10 bar).

a curious deviation in the very early part, is much closer to the true behavior. In both the CMS samples, Dual Model 1 seems to qualitatively capture the experimental trends, but the quantitative departure is quite significant, particularly in the later parts of the uptakes. In comparison, the overall quantitative superiority of the Dual Model 2 is obvious. Of course the quantitative agreement is less satisfactory here compared to that obtained for oxygen-nitrogen mixture uptake. This is most likely due to the deviation in binary isotherm prediction. Overall, the concentration dependence of the thermodynamically corrected transport


Langmuir, Vol. 19, No. 14, 2003 5733

Figure 11. Binary integral uptakes of carbon dioxide and methane in BF and Takeda CMS samples at 30 °C (pressure step, 0-10 bar).

parameters seems to appreciably affect the binary uptakes of both gases, and the extended empirical models are able to capture these effects. It should be emphasized once again that only single component parameters from independent measurements were used in the binary predictions without any adjustment. Ternary Integral Uptake. Having verified the proposed empirical equations to account for the concentration dependence of the thermodynamically corrected transport parameters for binary uptake, these equations were further verified with ternary integral uptake results in both the adsorbents. The ternary uptake experiments, subjected to pressure perturbation of 0-10 bar in the external fluid phase, were also conducted at 30 °C. Mixtures of nitrogen, methane, and carbon dioxide in the proportion 45:39:16 and 72:20:8 were used in the DAB measurements in BF CMS and Takeda CMS, respectively. Since nitrogen is comparatively weakly adsorbed and carbon dioxide is very strongly adsorbed, it was necessary to maintain a high proportion of nitrogen in the mixture to ensure a balanced proportion of all gases in the adsorbed phase, which could be properly analyzed upon desorption. The experimental results and the model prediction are compared in Figure 12. As can be seen from Figure 12, there are two roll-ups in the ternary integral uptake, one for carbon dioxide and the other for nitrogen. Carbon dioxide is the fastest diffusing among the three components. More carbon dioxide than permitted by ternary equilibrium value is initially adsorbed. The first roll-up corresponds to the replacement of excess carbon dioxide by nitrogen and methane. Since nitrogen diffuses faster than methane, it is also initially adsorbed in excess of its ternary equilibrium limit. The second roll-up is this excess nitrogen replaced by methane. The slow approach of carbon dioxide and nitrogen toward equilibrium is controlled by the adsorption rate of methane. The two models are very close for carbon dioxide uptake in Takeda CMS, which is consistent with our intuition since carbon dioxide concentration in this run was very low. In all other cases, despite qualitative similarities, the vast quantitative superiority of Dual Model 2,

Figure 12. Ternary integral uptakes of carbon dioxide, methane, and nitrogen in BF and Takeda CMS samples at 30 °C (pressure step, 0-10 bar).

particularly for methane, is undeniable. Therefore, the impact of the concentration dependence of the thermodynamically corrected transport parameters is also very clear in the ternary uptakes and the effectiveness of the extended empirical models to account for these effects is further validated. Conclusions This study establishes the effectiveness of the multisite Langmuir isotherm as a reliable model for representing unary equilibrium data and predicting mixture equilibrium in carbon molecular sieves. This conclusion has been verified with oxygen, nitrogen, carbon dioxide, methane, and their binary and ternary mixture data on BF and Takeda CMS samples. Adsorption kinetics of carbon dioxide and methane in BF CMS reinforces our previous finding that transport of gases in the micropores of CMS is controlled by a dual resistance. A comprehensive evaluation of the transport of all four gases in the two CMS samples (Figure 8) also reinforces our previous observation that the thermodynamically corrected transport coefficients generally show strong dependence on the adsorbed phase. Our empirical approach to account for the concentration dependence of the thermodynamically corrected transport parameters has been extensively verified with unary, binary, and ternary integral uptake experiments involving large pressure steps. The only limitation of this approach is that we cannot yet predict β; it must be extracted individually for each component from unary experimental data. It is also observed that, although the saturation capacity of the same gas varies significantly between the two adsorbents, the difference in a values is practically insignificant (Table 1). It is also observed that the transport coefficients of the four gases show an inverse

5734


Qinglin et al.

trend with their a values. These observations seem to imply that the a values obtained from the regression of combined unary equilibrium data using multisite Langmuir model may be reliable relative size indicators of the gas molecules in the adsorbed state.

qs

Glossary

q°

a A b bo c co cp cT cim Dc Dco Dco* Dco′ Dp Eb Ed J kf kb kbo kbo* kbo′ kad kde P Pi Pi° Psys q qp

adsorption sites occupied by each adsorbate molecule in the adsorbed phase surface area per unit mass of adsorbent, cm2 Langmuir constant, cm3/(g mol) pre-exponential constant for temperature dependence of b, cm3/(g mol) gas-phase concentration, g mol/cm3 constant feed concentration in DAB measurements, g mol/cm3 gas concentration in the macropores, g mol/cm3 total concentration in the gas phase, g mol/cm3 imaginary gas-phase concentration in the micropore, g mol/cm3 micropore diffusivity, cm2/s limiting micropore diffusivity, cm2/s limiting micropore diffusivity in the smallest accessible pore (see eq 1), cm2/s pre-exponential constant for temperature dependence of diffusivity, cm2/s macropore diffusivity, cm2/s activation energy for diffusion across the barrier resistance at the pore mouth, kcal/mol activation energy for diffusion in micropore interior, kcal/mol diffusion flux, g mol cm-2 s-1 fluid phase mass transfer coefficient, cm/s barrier coefficient, s-1 limiting barrier coefficient, s-1 limiting barrier coefficient at the smallest accessible pore mouth (see eq 2), s-1 pre-exponential constant for temperature dependence of barrier coefficient, s-1 Langmuir adsorption rate constant, cm3 (g mol)-1 s-1 Langmuir desorption rate constant, s-1 pressure, atm partial pressure of component i, atm hypothetical pressure in the IAS theory that yields same spreading pressure for every component i in the mixture, atm pressure in the desorption system, atm adsorbed gas-phase concentration, g mol/cm3 adsorbed gas-phase concentration based on particle volume, g mol/cm3

qT

q* qp * q j r rc R Rg Rp t T Tsys ∆U x y

monolayer saturation capacity according to Langmuir isotherm or total adsorption site concentration in an adsorbent according to multisite Langmuir isotherm, g mol/cm3 total adsorbed amount defined by eq 16, g mol/ cm3 equilibrium adsorbed amount at pressure P°, g mol/cm3 equilibrium adsorbed amount based on microparticle volume, g mol/cm3 equilibrium adsorbed amount based on macroparticle volume, g mol/cm3 average adsorbate concentration in the micropore, g mol/cm3 radial distance coordinate of the microparticle, cm microparticle radius, cm radial distance coordinate of adsorbent pellet, cm universal gas constant ()82.05 cm3 atm (g mol)-1 K-1) adsorbent pellet radius, cm time, s temperature, K temperature of desorption system, K change of internal energy due to adsorption, kcal/ (g mol) mole fraction in adsorbed phase mole fraction in gas phase

Greek Letters β p , βb π Π p θ θ* ζ η

fitting parameters in the proposed empirical models defined by eqs 1 and 2 spreading pressure, atm constant defined in eq 12 particle void fraction fractional coverage of the adsorption sites ()q/ qs) fractional coverage of the adsorption sites corresponding to q* ()q*/qs) dimensionless parameters along the radius of macropore ()R/RP) dimensionless parameters along the radius of micropore ()r/rc)

Subscript i

ith component (A, B, and C for a ternary system)

Note Added after ASAP Posting. This article was released ASAP on 6/11/2003 with errors in Table 1. The correct version was posted on 6/17/2003. LA0270791

Binary and Ternary Adsorption Kinetics of Gases in Carbon Molecular

Recommend Documents