Transport Behavior of Oxygen and Nitrogen through Organasilicon

phase Optimized Molecular Potentials for Atomistic Simulation Studies (COMPASS) force field was used to construct the polymers. Diffusion coefficients...
1 downloads 0 Views 580KB Size
J. Phys. Chem. B 2006, 110, 17375-17382

17375

Transport Behavior of Oxygen and Nitrogen through Organasilicon-Containing Polystyrenes by Molecular Simulation Qing Lin Liu* and Yu Huang Department of Chemical and Biochemical Engineering, Collage of Chemistry & Chemical Engineering, Xiamen UniVersity, Xiamen, 361005, China ReceiVed: May 24, 2006; In Final Form: July 17, 2006

Molecular dynamics (MD) simulations have been used to study the transport properties of oxygen and nitrogen in the para-substituted polystyrenes which possess one to four Si atoms in each substituent. The Condensedphase Optimized Molecular Potentials for Atomistic Simulation Studies (COMPASS) force field was used to construct the polymers. Diffusion coefficients were obtained from molecular dynamics (NVT ensemble) with up to 3 ns simulation times. After molecular dynamic simulation, the trajectories of the small molecules in the polymer matrix were obtained. Then diffusion coefficients have been calculated from the Einstein relationship revealing a considerable agreement between the simulated and the experimental data. And solubility coefficients have been calculated by the Grand Canonical Monte Carlo (GCMC) method. The solubility of oxygen increased with increasing Si content in the polymer membrane. The para-substituted polystyrenes with a branched substituent at the R-position showed higher permeability than those of the nonbranched ones. The higher the glass transition temperature (Tg) of the membrane, the larger the diffusion coefficients of oxygen and nitrogen obtained.

Introduction The rapidly increasing speed of computing over the past decade, combined with its equally rapidly decreasing cost, is revolutionizing science and engineering, so much so that computational science is increasingly regarded by many as the equal partner of the two traditional sciences, theoretical and experimental.1 In recent years, the growth of computational power and the development of theoretical models have made the computer simulation undergo a remarkable development. Gas separation through a membrane is one of the most attractive industrial processes because of the low energy consumption and ease of operation.2 One of the phenomena that, until now, is still not completely understood is the transport of gas molecules through membranes because experimental techniques cannot easily give details at the molecular level.3 Because of the static structure on an atomistic scale and the dynamic behavior on a time scale ranging from picoseconds (ps) to nanoseconds (ns) of the separation system (membrane plus penetrants) to determine the small molecule membrane separation mechanisms, it is impossible to get any direct experimental data about what is happening in these dimensions during the separation processes.4 To rectify this situation, the method of atomistic simulation has been widely used for understanding and predicting the structure-property relationships of materials. This method may provide a guideline for designing new polymers with high mechanical strength, high permselectivity, and high permeability. Atomistic simulation techniques have proven to be a useful tool for the understanding of structureproperty relationships of materials and MD can be used particularly for detailed descriptions of the complex morphologies and transport mechanisms. A number of MD and GCMC simulations of diffusion and sorption of small molecules in polymer membrane have been * Address corrspondence to this author. E-mail: [email protected].

reported in the literatures. Most of these simulations have focused on the investigation of the transport properties of small molecules in amorphous membranes such as polypropylene, polyimide, poly(dimethylsiloxane), and so on in resf 5-23. In this study, the main goal of MD simulation is essentially to make predictions of the physical properties of the polymer and give better insight into the mechanism of diffusion. Molecular dynamics simulation (NPT ensemble) has been used to relax the initial packing modeling. The trajectories of the small molecules in the polymer matrix were obtained from molecular dynamics (NVT ensemble) with use of up to 3 ns simulation times. Diffusion and solubility coefficients have been calculated by the MD and GCMC method, respectively. And the diffusion and sorption coefficients were calculated to discuss the effect of Si content on the permeability. Simulation Details Simulation of Polymer Microstructures. The simulations were performed with Material Studio (MS) software of Accelrys Inc. For the construction of the polymer chains and the amorphous cells, Materials Visualizer and Amorphous Cell modules were used. The calculation was carried out on the SGI workstations. The chemical structures of the repeat units of the seven Si-containing polystyrenes used in this study are shown in Figure 1. And their physical properties are listed in Table 1. An “all atom” model was used in the present study while the “united atom” model leads to an overestimate of diffusion coefficients of small penetrant molecules in polymers.13 The structure of Si-containing polystyrene monomers was build by the Materials Visualizer module. The potential types of atoms for these repeat units indicated in Figure 1 were elucidated according to the COMPASS force field.25 The monomers were initially subjected to energy minimization to seek the minimum of the potential energy surface of microstructure, respectively.

10.1021/jp063174x CCC: $33.50 © 2006 American Chemical Society Published on Web 08/15/2006

17376 J. Phys. Chem. B, Vol. 110, No. 35, 2006

Liu and Huang

Figure 1. Repeat units of Si-containing polystyrenes.

TABLE 1: Characteristics of the Organosilicon-Containing Polystyrene Samples Used in This Study24 polymer

10-5Mwa

Tgb (°C)

dc (g/cm3)

Si content (%)

FFVd

PS1 PS2 PS3 PS4 PS5 PS6 PS7

0.94 4.4 3.9 2.6 0.67 2.0 5.5

85 148 82 86 78 102 124

0.938 0.923 0.935 0.932

14.7 21.4 25.1 27.6 13.0 19.4 18.4

0.179 0.167 0.141 0.134

a Mw was determined by GPC, using polystyrene standard samples. Tg was determined by DSC (scanning rate, 20 deg/min). c Density was measured by the sink and float test. d FFV was calculated by the following equation of FFV ) (V - V0)/V, where, V and V0 denote the molar specific volume at 298 and 0 K, respectively. V0 was estimated by van der Waals volume, using the equation of V0 ) 1.3Vw.

b

The general procedure used during all the energy minimization stages in this study starts with steepest-descent and changes to conjugate-gradient as the energy derivative on any atom decreases to 1000 kcal/mol. As this value reaches 10 kcal/mol, it switches to the Newton-Raphson, and then stops when it is less than 0.001 kcal/mol or when at least 50 000 steps are completed. After the energy minimization of monomers, the repeat units were then used to construct a polymer chain containing 30-60 monomer units and were subsequently subjected to energy minimization. To reduce the chain end effect on the simulations only a single long chain was utilized.3 It should be mentioned that the restriction of the possible number of atoms for membrane module leads to a characteristic side length of the packing model of just a few nanometers. Therefore, to avoid severe surface effects it is necessary to introduce periodic boundary conditions.4 The initial polymer periodic boundary condition structures are grown with the Amorphous Cell modules, which implement a modification of the rotational isomeric state (RIS) method of Theodorous and Suter.26,27 The generation of an initial polymer structure is based on the “selfavoiding” random walk method of Theodouous and Suter.26 According to this method, three backbone-chain atoms of a polymer together with their bonds are initially “placed” in a cubic unit cell in some random orientation. Next, the length of the backbone chain is increased stepwise by adding one atom at a time to the growing polymer chain. During this process, the allowed rotation states (orientations) of successive bonds between adjacent atoms are determined from probability functions by energy consideration by using the standard Monte Carlo method. To avoid packing algorithm related catenations and spearings of aromatic units, it was necessary to start the packing at a very low density of 0.1 g/cm3.4 And for the purpose of avoiding the generation of high energy configurations, the lookahead feature of the Amorphous Cell module was used with a value of 5, with all 5 added at each step, and the maximum number of look-ahead configurations was set to 50. All the

Figure 2. Typical polymer packing model: (a) initial polymer packing model and (b) polymer model after energy refinement.

microstructures were assumed to contain one polymer chain, three oxygen molecules, and three nitrogen molecules in each cell. The reason for using three small molecules is to get the mean squared displacement (MSD) by calculating from the trajectories of the three penetrant molecules of the same species. That is to say we used three small molecules of each kind in these simulations to obtain the average diffusion coefficient. Following the cell construction at the initial density, the Discover module of the MS software was used for calculation

Transport Behavior of O and N through Polystyrenes

J. Phys. Chem. B, Vol. 110, No. 35, 2006 17377

TABLE 2: Details of the Model at 298 K

polymer

no. of repeat units

no. of atoms in the model

no. of inserted penetrant molecules

initial density before refinement (g‚cm-3)

density after refinement (g‚cm-3)

lateral dimensions before refinement (Å)

lateral dimensions after refinement (Å)

PS1 PS2 PS3 PS4 PS5 PS6 PS7

60 45 35 30 60 45 40

1862 1937 1927 2012 1922 1982 2030

6 6 6 6 6 6 6

0.10 0.10 0.10 0.10 0.10 0.10 0.10

0.927 0.915 0.934 0.933 0.926 0.928 0.934

57.45 58.10 57.94 58.74 58.64 58.98 58.21

27.35 27.78 27.51 27.90 27.92 28.29 27.75

TABLE 3: Transport Properties of Nitrogen in Organosilicon-Containing Polystyrenes at 298 K

d

polymer

Dcala

Dexpb

Scalc

Sexpb

Pcald

Pexpb

PS1 PS2 PS3 PS4 PS5 PS6 PS7

5.24 6.67 1.07 4.01 2.99 3.8 2.34

4.82 6.09 0.34 2.7

3.11 3.49 5.13 3.73 4.54 5.49 3.67

1.53 1.65 3.3 2.1

17.8 32.7 4.33 15.0 13.7 20.9 8.22

6.6 10.1 1.2 5.4 4.8 7.8 3.99

2.4 1.99

3.4 2.01

a In 10-7cm2 s-1. b Reference 24. c In 10-3 cm3 (STP) cm-3 cm Hg-1. In Barrer.

TABLE 4: Transport Properties of Oxygen in Organosilicon-Containing Polystyrenes at 298 K

d

Figure 3. (a) Mean squared displacement s(t) of the oxygen in a PS1 matrix as a function of time and (b) log s(t) vs log t plot for the diffusion of oxygen in a PS1 matrix.

and minimization of the potential energy of the penetrant/ polymer system. The total potential energy of the initial polymer is usually much higher than that of the practical polymer it simulates. This high energy is due to the nonbonded interactions between the atoms of the polymer chains, which arise from the overlap of atoms in the simulated initial microstructure.13 So all initial packings were subject to extensive equilibration procedures. First, several cycles of NPT (constant particle number, temperature, and pressure) dynamics of the Discover module were performed at pressures of thousands of bars and 298 K to scale the cells to the specified experimental density, which are given in Table 1. Then the Basic Refine Protocol of Amorphous Cell module was used. A brief NVT (constant particle, volume, and temperature) dynamics run at 298 K and a final energy minimization were performed to minimize the first energy minimization. Subsequently, annealing was performed by the NPT dynamics procedure of Temperature Cycle Protocols at 1 bar. During annealing, the cells were first heated

polymer

Dcala

Dexpb

Scalc

Sexpb

Pcald

Pexpb

PS1 PS2 PS3 PS4 PS5 PS6 PS7

11.56 13.79 2.53 7.23 9.32 10.78 8.37

11.1 13.0 1.0 4.8

5.23 6.17 6.92 6.67 4.96 6.31 5.78

2.33 2.9 3.9 3.8

75.5 134.7 17.5 48.2 46.2 68.0 48.4

25.9 37.7 4.0 18 17 26 18.1

7.4 6.51

3.4 2.78

a In 10-7 cm2 s-1. b Reference 24. c In 10-3 cm3 (STP) cm-3 cm Hg-1. In Barrer.

by 50 K increments from 298 to 800 K which is well above the Tg of the polymer. The duration of the NPT dynamics simulations at each temperature was 150 ps. And then a NVT dynamics was used to further relax the polymer structure. Figure 2 shows the cubic cell structures of the system of PS1 before and after the energy refinement. In all runs, Andersen’s method and Berendsen’s method were used for temperature and pressure control, respectively.28,29 During the simulations, a large charge group-based cutoff of 13 Å, with the spline and buffer widths 2 and 1.5 Å, respectively, was applied to evaluate the nonbonded interactions. The time step was 1 fs. Diffusion, Solubility, and Permeability. Diffusion of oxygen and nitrogen molecules in the equilibrated cells was determined by NVT dynamics at 298 K. The cells were relaxed by NPT dynamics at 298 K for 1000 fs. Then a NVT dynamics was subsequently performed at the density of the model polymer systems for 3 ns to study the motion details of molecules in the cell. From the trajectories recorded at 1 ps intervals, the consecutive positions of the penetrant molecules diffusing in the polymer matrix are computed as a function of time. The diffusion coefficients for the penetrant can be calculated by means of the Einstein relationship:

〈|r(t) - r(0)|2〉 tf∞ 6t

D ) lim

(1)

where r(0) is the initial positional coordinate of the penetrant molecule in the polymer matrix, and r(t) is the positional

17378 J. Phys. Chem. B, Vol. 110, No. 35, 2006

Liu and Huang

Figure 4. Displacement of oxygen (a) and nitrogen (b) in PS1 matrix as a function of the simulation time.

coordinate of this molecule after a time t, and |r(t) - r(0)| represents the displacement of the penetrant molecule during time t. The MSD (s(t) ) 〈|r(t) - r(0)|2〉) means the given penetrant molecule averaged over all possible time origins t ) 0 and all simulated trajectories of the penetrant. Because eq 1 relies on the assumption of a random walk for penetrant through the polymer matrix, one has to consider the possible effect of anomalous diffusion first reported by Mu¨ller-Plathe et al.30 The criterion for a model system reaching the normal diffusive regime is that the slope of log s(t) vs log t approaches unity. As shown in Figure 3b, the slope of log s(t) vs log t is 0.57 suggesting an anomalous diffusion at shorter times and increases to approximately unity (1.06) at about 300 ps. So eq 1 can be used to determine the D-value from MD simulation in this study, as indicated in Figure 3a. The sorption of the penetrant molecule through the polymeric matrix was calculated by the Sorption module of MS Software. Both van der Waals and Coulombic forces were included in the sorption simulation. In this simulation, gas molecules were inserted into the interface of the equilibrated cells. The gas concentration probability is based on the calculation of the energy change between the new configuration and the previous one. Solubility coefficients are obtained from the GCMC method that uses a Metropolis algorithm for accepting or rejecting configurational moves (rotation and translation of the sorbate molecule) as well as for sorbate insertion and deletion.31 The probability of gas creation and destruction was expressed as

[ (

P ) min 1; exp

)]

NikT -∆E ( ln kT fi V

(2)

Figure 5. Trajectories of oxygen (a) and nitrogen (b) in PS1 matrix: the size of the cubes corresponds to the volume of the simulated basic volume element.

where ∆E is calculated from the sum of nonbonded potential terms, k is the Boltzmann’s constant, V is the cell volume, fi is the fugacity of component i in the gas phase, and Ni is current number of component i in the membrane cell. The probability of translation movement in the cell is given by

[

P ) min 1; exp

(-∆E kT )]

(3)

A random sorbate molecule was chosen in the cell for the gas rotation move. The rotation axis was chosen randomly, and the molecule rotated by a random amount within the cell. The new gaseous configuration, based on the energy change, was obtained with the same probability applied to the translation movement as mentioned above.32 By performing the simulation over a

Transport Behavior of O and N through Polystyrenes

J. Phys. Chem. B, Vol. 110, No. 35, 2006 17379

Figure 6. Plots of the oxygen (a) and nitrogen (b) solubility coefficient vs Si content in the polystyrene membranes.

Figure 7. Plots of the oxygen (a) and nitrogen (b) diffusion coefficient vs Si content in the polystyrene membranes.

range of pressures (or fugacities for nonideal gases) from 0.05 to 1 atm, a sorption isotherm can be obtained in the form of a plot of the concentration of sorbed gas, C, as a function of pressure (or fugacity) at constant temperature. At each pressure, 1 000 000 steps of GCMC calculations were performed with use of the initial equilibration period of 100 000 steps. The solubility coefficient, S, is then obtained from the limiting slope of the sorption isotherm at zero pressure as22

nitrogen transport in the organosilicon-containing polystyrenes. Comparison with the corresponding experimental values is also given; the estimated and the experimental values did not differ significantly. Most of the ratios of the estimated results over the experimental data are not considerably greater than a factor of 3. Judging these results one can notice that the possible errors of experimental solubility and diffusivity data can be quite high.33 This is due to the difficulties in obtaining really amorphous polymeric materials. On the other hand, the accuracy of calculated values may be affected by several factors such as the force field used, the overall simulation time, the cutoff, and other assumptions and approximations embodied in the MD and MC simulations.33 It is, therefore, generally accepted that a coincidence between measured and simulated diffusivity and solubility values with a factor of 3 to 5 is acceptable.30 That confirms, in turn, that these models can be considered to be sufficiently equilibrated. The simplest way of studying the diffusion of an individual penetrant is to inspect its path through the spaces. To get a more detailed insight into different models of penetrant molecules, it can also make sense to investigate displacement |R(t)| ) {|r(t) - r(0)|2}1/2 for an individual penetrant.3 Figure 4 shows the displacement of oxygen and nitrogen diffusion in the Poly1. We can obtain from the picture that the “jump” motions of the oxygen are greater than those of nitrogen. That may make the diffusion coefficient of oxygen much larger than that of nitrogen. This is because the diffusion of small penetrant molecules in a polymer can be visualized as occurring by the hopping mechanism.13 We can also obtain the information of the small

S ) lim(C/p) pf0

(4)

where C is in units of cm3 (STP)/cm3 polymer, and p is the pressure. Once the diffusion and solubility coefficients are obtained, the permeability coefficient, P, can be calculated by the relationship

P ) DS

(5)

Results and Discussion Table 2 lists the detail of models. Compared to the corresponding experimental values of density in Table 1, the estimated densities are in close agreement with the experimental data. The estimated data show a deviation from the experimental data between -1.1% for PS1 and 0.1% for PS4. Therefore, in the given case, the investigated models can be considered to be reasonably equilibrated. Tables 3 and 4 list the values of diffusion, solubility, and permeability coefficients at 298 K estimated form MD and GCMC simulations for oxygen and

17380 J. Phys. Chem. B, Vol. 110, No. 35, 2006

Liu and Huang

Figure 9. Plots of the oxygen (a) and nitrogen (b) diffusion coefficient vs Tg in the polystyrene membranes.

Figure 8. The probability distribution of oxygen (a) and nitrogen (b) calculated in the PS1 matrix.

molecules’ transfer behaviors from Figure 5. The trajectory of the oxygen is much lager than that of the nitrogen, which can also show us the diffusivity of oxygen and nitrogen. Figure 6 shows the impact of silicon content on the solubility of oxygen and nitrogen. However, the tendency for the solubility coefficients of oxygen changing with the silicon content is very different from that for nitrogen. As shown in Figure 6, the solubility coefficients of oxygen increased with increasing silicon content, while that of nitrogen did not have the same tendency. Figure 7 shows the variation of diffusion coefficients of oxygen and nitrogen with the silicon content in the organosilicon-containing polystyrenes. It suggests that the variation of these for both oxygen and nitrogen has a similar trend. It is found that the silicon atom has an empty d-orbital. In the siloxane linkage, the electrons of the oxygen’s lone pair interact with the empty d-orbital to increase the charge density on the Si atom.34 On the other hand, because there is no lone pair in the C-Si linkage in an organosilicon moiety, the empty d-orbital

has the opportunity to interact with the oxygen molecule in the membrane.35 That may explain why only the solubility coefficient of oxygen increases with the silicon content and also why the solubility coefficient of oxygen is higher than that of nitrogen, as indicated in Figure 8. From Figure 9, we can observe the relationship between the Tg of the polymer and the diffusion of the oxygen and nitrogen in the membrane. Yampol’skii et al.36 reported independently that the free volume of the glassy polymers increased with increasing Tg. That may explain the relationship between Tg and diffusion coefficient of oxygen and nitrogen. In our case of organosilicon-containing polymers, introducing the branching structure at the para-position of the benzene ring increased the Tg (PS1 < PS2, PS3 < PS4, and PS5 < PS6). In other words, the polymer having methylene linkages in the substituents showed fairly low Tg due to the flexibility of the substituents (PS1, PS3, and PS5), as indicated in Figure 11. According to the relationship between free volume and Tg in the glass polymer mentioned previously, the polymers having a flexibile substituents (PS1, PS3, and PS5) have a comparatively lower free volume. That may explain the tendency of the diffusion coefficient (PS2 > PS1, PS4 > PS3, and PS6 > PS5). Table 5 lists the values of solubility selectivity, diffusion selectivity, and permselectivity of oxygen vs nitrogen. The selectivity of oxygen vs nitrogen for solubility and diffusion coefficients has no obvious relationship. We can also see from Figure 10 that there is no clear relationship between the selectivity of oxygen vs nitrogen to Si content. However, the reason may be due to the complex interaction of solubility with

Transport Behavior of O and N through Polystyrenes

J. Phys. Chem. B, Vol. 110, No. 35, 2006 17381

Figure 11. MSD of substituents of the polymers as a function of time.

TABLE 5: Selectivity Values for Oxygen/Nitrogen Separation on the Seven Polystyrenes polymer

SO2/SN2

DO2/DN2

PO2/PN2

PS1 PS2 PS3 PS4 PS5 PS6 PS7

1.68 1.77 1.35 1.79 1.10 1.15 1.57

2.21 2.07 2.36 1.80 3.12 2.84 3.58

3.72 3.66 3.19 3.22 3.37 3.26 5.62

influenced by the side chain mobility of the polystyrene, the slower the mobility of the side chain, the higher Tg and larger free volume will be. To synthesize a high oxygen and nitrogen permeability and permselectivity membrane, we may choose high Si content and lower mobility organosilicon-containing polystyrenes. Acknowledgment. The support of National Nature Science Foundation of China Grant no. 50573063 and the research fund for the Doctoral Program of Higher Education (no. 2005038401) in preparation of this article is gratefully acknowledged. References and Notes

Figure 10. Plots of the solubility selectivity (a), diffusion selectivity (b), and permselectivity (c) vs Si content in the polystyrene membranes.

diffusion. The result shows that higher Si content may increase the oxygen’s but not the nitrogen’s solubility coefficient and the lower mobility of the side chain can increase the diffusion coefficient. This shows us clearly that to synthesize a high oxygen and nitrogen permeability and permselectivity membrane we must use the organosilicon-containing polystyrenes. Conclusion The solubility coefficient of oxygen increases with increasing Si content because of the strong affinity of the C-Si linkage toward oxygen molecules. The diffusion coefficient is strongly

(1) Peter, T. C. Molecular Dynamics Simulation of Realistic Systems. Fluid Phase Equilib. 1996, 116, 237-248. (2) Vansant, E. F.; Dewolfs, R. Gas Separation Technology: Proceedings of the International Symposium on Gas Separation Technology; Elsevier Publishing Company: Amsterdam, The Netherlands, 1990. (3) Tocci, E.; Bellacchio, E.; Russo, N.; Diroli, E. Diffusion of gases in PEEKs membranes: molecular dynamics simulations. J. Membr. Sci. 2002, 206, 389-398. (4) Hofmann, D.; Fritz, L.; Ulbrich, J.; Schepers, C.; Bohning, M. Detailed-atomistic molecular modeling of small molecule diffusion and solution processes in polymeric membrane materials. Macromol. Theory Simul. 2000, 9, 293-327. (5) Kucukpinar, E.; Doruker, P. Molecular simulations of small gas diffusion and solubility in copolymers of styrene. Polymer 2003, 44, 36073620. (6) Mu¨ller-Plathe, F. Diffusion of penetrants in amorphous polymers: A molecular dynamics study. J. Chem. Phys. 1991, 94, 3192-3199. (7) Mu¨ller-Plathe, F. Molecular dynamics simulation of gas transport in amorphous polypropylene. J. Chem. Phys. 1992, 96, 3200-3205. (8) Cuthbert, T. A.; Wagner, N. J.; Paulaitis, M. E.; Murgia, G. M. D Simulations of Penetrant Diffusion in Amorphous Polypropylene: Diffusion Mechanisms and Simulation Size Effects. Macromolecules 1999, 32, 50175028. (9) Boyd, R. H.; Pant, P. V. K. Molecular packing and diffusion in polyisobutylene. Macromolecules 1991, 24, 6325-6331. (10) Mu¨ller-Plathe, F.; Rogers, S. C.; van Gunsteren, W. F. Diffusion coefficients of penetrant gases in polyisobutylene can be calculated correctly by molecular-dynamics simulations. Macormolecules 1992, 25, 6722-6724.

17382 J. Phys. Chem. B, Vol. 110, No. 35, 2006 (11) Tamai, Y.; Tanaka, H.; Nakanishi, K. Molecular Simulation of Permeation of Small Penetrants through Membranes. 1. Diffusion Coefficients. Macromolecules 1994, 27, 4498-4508. (12) Sok, R. M.; Berendsen, H. J. C.; van Gunsteren, W. F. Molecular dynamics simulation of the transport of small molecules across a polymer membrane. J. Chem. Phys. 1992, 96, 4699-4704. (13) Charati, S. G.; Stern, S. A. Diffusion of gases in silicone polymers: Molecular dynamics simulations. Macromolecules 1998, 31, 5529-5535. (14) Gee, R. H.; Boyd, R. H. Small penetrant diffusion in polybutadiene: a molecular dynamics simulation study. Polymer 1995, 36, 14351440. (15) Han, J.; Boyd, R. H. Molecular packing and small-penetrant diffusion in polystyrene: a molecular dynamics simulation study. Polymer 1996, 37, 1797-1804. (16) Hofmann, D.; Ulbrich, J.; Fritsch, D.; Paul, D. Molecular modelling simulation of gas transport in amorphous polyimide and poly(amide imide) membrane materials. Polymer 1996, 37, 4773-4785. (17) Zhang, R.; Mattice, W. L. Atomistic modeling of the diffusion of small penetrant molecules in the bulk amorphous polyimide of 3,3′,4,4′benzophenonetetracarboxylic dianhydride and 2,2-dimethyl-1,3-(4-aminophenoxy)propane. J. Membr. Sci. 1995, 108, 15-23. (18) Mu¨ller-Plathe, F. Diffusion of water in swollen poly(vinyl alcohol) membranes studied by molecular dynamics simulation. J. Membr. Sci. 1998, 141, 147-154. (19) Tamai, Y.; Tanaka, H. Permeation of small penetrants in hydrogels. Fluid Phase Equilib. 1998, 144, 441-448. (20) Kim, W. K.; Mattice, W. L. Static and Dynamic Behavior of H2O and O2 Penetrants in a Polybenzoxazine. Macromolecules 1998, 31, 93379344. (21) Fried, J. R.; Sadat-Akhavi, M.; Mark, J. E. Molecular simulation of gas permeability: poly(2,6-dimethyl-1,4-phenylene oxide). J. Membr. Sci. 1998, 149, 115-126. (22) Fried, J. R.; Ren, P. The atomistic simulation of the gas permeability of poly(organophosphazenes). Part 1. Poly(dibutoxyphosphazenes). Comput. Theor. Polym. Sci. 2000, 10, 447-463.

Liu and Huang (23) Fried, J. R.; Goyal, D. K. Molecular Simulation of Gas Transport in Poly[1-(trimethylsilyl)-1-propyne]. J. Polym. Sci.: Part B: Polym. Phys. 1998, 36, 519-536. (24) Nagasaki, Y.; Hashimot, Y.; Kato, M.; Kimijima, T. Gas permeation properties of organosilicon-containing polystyrenes. J. Membr. Sci. 1996, 110, 91-97. (25) COMPASS-Condensed Phase Optimized Molecular Potentials for Atomistic Simulation Studies; Distributed by Accelrys Inc., 1999. (26) Theodorou, D. N.; Suter, U. W. Detailed molecular structure of a vinyl polymer glass. Macromolecule 1985, 18, 1467-1478. (27) Theodorou, D. N.; Suter, U. W. Atomistic modeling of mechanical properties of polymeric glasses. Macromolecule 1986, 19, 139-154. (28) Andersen, H. C. Molecular dynamics simulations at constant pressure and/or temperature. J. Chem. Phys. 1980, 72, 2384-2393. (29) Berendsen, H. J. C.; Postama, J. P. M.; Van Gunsteren, W. F.; DiNola, A.; Haak J. R. Molecular dynamics with coupling to an external bath. J. Chem. Phys. 1984, 81, 3684-3690. (30) Mu¨ller-Plathe, F. Permeation of polymers-a computational approach. Acta Polym. 1994, 45, 259-2932. (31) Metropolis, N.; Rosenbluth, A. W.; Rosenbluth, M. N.; Teller, A. H.; Teller, E. Equation of State Calculations by Fast Computing Machines. J. Chem. Phys. 1953, 21, 1087-109. (32) Tung, K. L.; Lu, K. T.; Ruaan, R. C.; Lai, J. Y. MD and MC simulation analyses on the effect of solvent types on accessible free volume and gas sorption in PMMA membranes. Desalination 2006, 192, 391400. (33) Tocci, E.; Hofmann, D.; Paul, D.; Russo, N.; Diroli, E. A molecular simulation study on gas diffusion in a dense poly(ether-ether-ketone) membrane. Polymer 2001, 42, 521-533. (34) Rochow, E. G. Silicon and Silicones; Springer-Verlag: Berlin, Germany, 1990. (35) Elschenbroich, C. H.; Salzer, A. Organometallics; VCH: Weinheim, Germany, 1992. (36) Yu, P.; Yampol’akii, Y. P.; Volkov, V. V. Studies in gas permeability and membrane gas separation in the Soviet Union. J. Membr. Sci. 1991, 64, 191-228.