Computer simulation on supercritical carbon dioxide fluid a potential model for the benzene-carbon dioxide system from ab initio calculations

EOgltlgglA ELSEVIER

Fluid Phase Equilibria 104 (1995) 375-390

C o m p u t e r S i m u l a t i o n on Supercritical C a r b o n Dioxide Fluid A potential model for the benzene-carbon dioxide system from ab initio calculations Jun-Wei Shen, K.B. Domafiski, Osamu Kitao, Koichiro Nakanishi

Division of Molecular Engineering, Faculty of Engineering, Kyoto University, Sakyo-ku, Kyoto 606-01, Japan Keywords: potential model, ab initio calculations, Monte Carlo, benzene, carbon dioxide Received 18 May 1994; accepted in final form I August 1994

ABSTRACT Ab initio q u a n t u m chemical calculations were done for a system of a benzene molecule and a carbon dioxide molecule for a lot of relative configurations. Since intermolecular interactions in this system were considered to be small at some configurations, the calculations were carried out with 6-31G* basis set and the 2nd order M¢ller-Plesset p e r t u r b a t i o n method. The resultant potential surface was explained in detail using four representative configurations. At the short intermolecular distance, the interaction showed strong dependence on the relative configurations. On the basis of these data, an analytical pair potential function was proposed for the b e n z e n e - C O 2 system. Using this potential model, we have done Monte Carlo calculations for the CO2 fluid containing one benzene molecule at several supercritical conditions. Our model was confirmed to reflect the shape of a plate-like molecule on the solvent structures.

INTRODUCTION As one of chemical engineering processes, supercritical fluid extractions, nowadays, have wide applications to separation and purification of materials which are particularly sensitive to temperature and inadequate to be treated by conventional distillation procedures. For volatilization and separation of relatively heavy organic compounds, the supercritical processes are expected to significantly reduce the energy requirements (Kohn and Savage, 1979). Concerning the supercritical extraction techniques, the optimum design and efficient operations need knowledge and prediction on the characteristic features of solvents around solutes. Although a large number of researches on supercritical fluid extractions have been already reported at even the several production stages, there is little quantitative fundamental understanding on the fluid structures in both experimental and theoretical studies (Bruno and

0378-3812/95/$09.50 © 1995 - Elsevier Science B.V. All rights reserved SSDI 0378-3812 (94) 02662-9

376

J.-W. Shen et al. / Fluid Phase Equilibria 104 (1995) 375-390

Ely, 1991). The major difficulties in the experimental studies come from high compressibility of the mixture around the critical point. Near the point, the experimental data give only drastic changes in the phase behavior of mixtures but no detailed information. On the other hand, high non-ideal properties of the mixtures cause severe problems in the theoretical treatments. Usually, the solutes are heavy organic compounds and the solvents are relatively small molecules; this situation gives no satisfactory mixing rule. As a result, conventional empirical potential models give poor results for such fluids (Shing et al., 1982, 1983). Therefore, computer simulations with a set of reliable potential models are expected to give new aspects on the supercritical extractions. As a supercritical extraction solvent, CO2 (carbon dioxide) is now commonly used in the practical applications due to the tractable critical conditions (31°C and 73.8 bar). Using ab initio quantum chemical calculations, we have already constructed our potential model for the CO2 fluid and confirmed that the model well reproduce the experimental phase diagram (Domafiski et al., 1994). Although the CO2 fluid is known to have good properties as a supercritical solvent, the extraction mechanism has not been explained in detail at the molecular level. In order to understand how the CO2 fluid dissolves several kinds of solutes and evaporates from solutions, the detailed changes in the solution environments must be investigated with an appropriate solute along the changes in temperature and pressure of the system. A benzene molecule is adopted as the solute in this study. Our interest is to elucidate effects of the plate-like molecule on the CO2 solvent structures. In attempting to interpret behaviors of the supercritical fluid around the solute molecule, a key quantity is an interaction energy between the solute and solvent molecules. For this purpose, a reliable potential function is required for the benzene-CO2 system. There is, however, no extensive potential function. Therefore, we started our research from constructing a potential model for the benzene-CO2 system. This paper is organized as follows. First, we perform ab initio quantum chemical calculations at extremely large numbers of relative configurations for the benzene-CO2 system. The potential energy surface obtained for the system is discussed in detail using four representative configurations. Next, a potential model is proposed for the solute-solvent system. Using iso-energy contour maps by this potential model, we summarize several information to consider future simulation results. Finally, we report primary results of MC (Monte Carlo) simulations for the CO2 fluid containing one benzene molecule at several supercritical conditions.

J.-W. Shen et al. I Fluid Phase Equilibria 104 (1995) 375-390

377

C A L C U L A T I O N A L DETAILS

Ab initio quantum chemical calculations A system constituted by a benzene molecule and a CO2 molecule is considered to have a small interaction due to the dispersion force at some configurations. Moreover, to make an appropriate potential model, our potential surface should smoothly continue to the infinite distance limit of this system. Therefore, we adopted a relatively flexible basis set, 6-31G*, and MP2 (the second-order perturbation m e t h o d by M¢ller-Plesset theory). First, we optimized the structure of a benzene molecule with the same basis set and wave function level; the results are Rc_c=1.395 /~ and Rc_H=l.087 ~ , respectively. In the same t r e a t m e n t for a CO2 molecule, R c - o is 1.180/~ (Domafiski et al., 1994). All calculations in this paper used these data. Z w

/,,,

i i

~-X

/ Y

Fig.1 Definition of relative configurations between a benzene molecule and a CO2 molecule. Our c o m p u t e r simulations need an analytical function which describes an interaction energy of the b e n z e n e - C O 2 system with a wide range of configurations. In order to get the sufficient information, the complex potential surface requires ab initio calculations at extremely large numbers of configurations. Five parameters describe the relative configurations as depicted in Fig.l, in which a benzene molecule is located at the origin of coordinates and "r" is a distance between the center of the benzene molecule and t h a t of a CO2 molecule. The other four parameters are defined as follows:

378

J.- W. Shen et al. / Fluid Phase Equilibria 104 (1995) 375-390

1. 0 is a rotation angle of the O'Z' axis around the OZ axis, 2. ¢ is an angle between the OZ axis and the intercenter distance "r" 3. a is a rotation angle between the CO2 molecule and the O'Z' axis, 4.13 is a rotation angle of the CO2 molecule around the O'Z'. For one fixed distance "r", all angles change from 0 ° to 360 ° with an increment of 45 ° (¢, a , and/3) or 30 ° (0). This procedure generates 2,208 configurations of the b e n z e n e - C O 2 system, which are reduced to 45 unique configurations due to the s y m m e t r y of the benzene molecule. For each configuration, we calculated 10 points at different intermolecular distances "r" from 2.5 to 10.0/~. The resultant 450 unique points correspond to 22,080 points on the potential surface. All the calculations were carried out using Gaussian 90 package program (Frisch et al., 1990).

A potential function On the basis of above ab initio calculations, we proposed a potential model for the b e n z e n e - C O 2 system. There are m a n y possibilities to select an analytical potential function for this system. Here we adopted a site-site LennardJones model with Coulomb term such as

Vab initio(r)= ~ ~

i=l j = l

4eij

rij)

-

rij

J

+

rij J

,

(1)

where the i's and f s are atomic sites on the benzene and CO2 molecules, respectively. The rij is the intersite distance; eij, aij, and qi,j are energy, size, and Coulomb parameters, respectively. The constraints, qc = -2qo for the CO2 molecule and qH = --qc for the benzene molecule, are imposed on the Coulomb parameters.

Monte Carlo simulations We carried out MC calculations with a constant-NVTensemble by Metropolis procedure (Metropolis et al., 1953). The system consists of one benzene molecule and 256 CO2 molecules with a periodic b o u n d a r y condition. The benzene molecule is a rigid body of Dsh s y m m e t r y and each CO2 molecule is also a linear rigid body. These structures are due to our optimization. The interaction energy between the benzene and CO2 molecules was described by our new potential function and t h a t between CO2 molecules was described by

J.-W. Shen et al. /Fluid Phase Equilibria 104 (1995) 375-390

379

our previous results (Domafiski et al., 1994). We performed the simulations for several different temperatures from 310K to 400K at a density of 0.5749 g.cm -3. The temperature region is slightly higher than the experimental critical temperature of the CO2 fluid (304.15K). After confirming an equilibrium of the system, we used 100,000 MC steps for several subsequent analyses. Here one MC step generates randomly new 256 configurations in which every COs molecule is moved once, while the benzene molecule is kept at the initial point and configuration. RESULTS AND DISCUSSION

Potential surface from ab initio calculations Present ab initio data gave a lot of interesting information for the benzeneCO2 system. To summarize characteristic points of the potential surface, we use four representative configurations, CF1 -~ CF4. Those are depicted in Figs.2(a) ,-~ (d) with the corresponding potential energy curves. In these figures, dotted lines are potential energy curves due to HF (Hartree-Fock) wavefunction and solid lines are those due to MP2 wavefunction. The HF describes an electrostatic interaction and the MP2 adds an electron correlation effect to the HF. The electron correlation modifies the electrostatic interaction at HF level by deformation of an electron cloud and explains a dispersion interaction. Therefore, from comparisons between the dotted and solid lines, we can investigate what kind of interaction works at the specific configuration. CF1 is the most stable configuration among all our present calculations. There is a minimum o f - 4 . 6 kJ.mo1-1 at 3.7 /~ with the HF level and an electron correlation effect makes the minimum deeper to -12.1 kJ.mo1-1 and inner to 3.3 /~. Both electrostatic and dispersion interactions work well for this configuration. When the COs molecule approaches the benzene molecule in CF1, two different kinds of electrostatic interactions occur; one is an attractive interaction between the C-COs (carbon atom of the CO2 molecule) (6+) and the EC-B (electron cloud of the benzene molecule) ( 6 - ) and another is a repulsive interaction between the O-COs (oxygen atoms of the CO2 molecule) ( 6 - ) and the EC-B (6-). Since the former is more significant than the latter, the HF wavefunction gives a sufficient deep minimum. When the COs molecule rotates around the C8 axis of the benzene molecule, little difference exists in the interaction energy. This confirms that molecular partial charge does not localize at the carbon atom positions of the benzene molecule. Therefore, addition of the electron correlation easily deforms the EC-B and reduces the electrostatic repulsion between the O-COs and the EC-B. More-

J.-W. Shen et al. / Fluid Phase Equilibria 104 (1995) 375-390

380

over, a dispersion force causes the additional attractive interaction between the CO2 molecule and the benzene molecule. As a result, a deep potential m i n i m u m occurs in this configuration. The flexible EC-B plays a significant role in this CF1. 15 O

E

I

10 _(a)

I

15

I

HF . . . . MP2 - -

' l

0

0

-~

-5

-5

/ ~

-10

-10

.CFU½

r-

-15 15

I

I

,

,

,

I-

lo |(c)

,

I

|i

i

/',

I

-15

i

15

. PHF 2 ....

1

0

-5

CF3

.............................

"

i

{d)

I

I

I

I

I

I

I

I

HF . . . . MP2

-

CF4 .......................

-10 i

2

....

-5

~5 -10 0

.....

........

n

-15

MP2

10 5

a~

'

5

eUJ

~,

'

10

I

5

<°>

4 6 Distance [A]

I

8

< ~ ~

-15 0

i 2

I

i

4 6 Distance [A]

I

8

10

Fig.2 Ab initio potential energy curves and the corresponding configurations. Dotted lines are due to HF and solid lines are due to MP2: (a) CF1, (b) CF2, (c) CF3, and (d) CF4. CF2 changes a shape of the potential energy curve by addition of an electron correlation effect. The HF wavefunction gives a simple repulsive curve, but the MP2 wavefunction has a shallow m i n i m u m of-2.5 k J . m o l - l a t 4.5 •. W h e n the CO2 molecule approaches the benzene molecule in CF2, a strong electron repulsive interaction occurs between the O-CO2 ( 5 - ) and the EC-B ( 5 - ) . W i t h a similar mechanism to the case in CF1, the electron correlation effect reduces the electron repulsive interaction and causes a dispersion interaction between the benzene and CO2 molecules. Since a distance in Fig.2 is defined as a length between the CBM (center of the benzene molecule) and

J.-W. Shen et al. / Fluid Phase Equilibria 104 (1995) 375-390

381

the C-CO2, at the minimum position in CF2, the lower O-CO2 is located at 3.3/~ above the benzene plate. This distance is almost the same as the position where a deep minimum occurs in CF1. This coincidence is reasonable because the dispersion force finally determines the minimum position for both CF1 and CF2. Addition of an electron correlation effect also changes a shape of the potential energy curve in CF3. A simple repulsive curve at the HF wavefunction changes to an attractive curve with a shallow minimum of-0.9 kJ.mo1-1 at 5.5/~ in the MP2 wavefunction. In this CF3, the C-CO2 approaches the H-B (hydrogen atom of the benzene molecule) and causes two kinds of electrostatic repulsive interactions; one is due to the C-CO2 (6+) and the H-B (6+) and another is due to the O-CO2 ( 6 - ) and the EC-B (6-). Although the electron correlation effect reduces the electron repulsive interaction and causes a dispersion interaction, the quantity is not so large as that of CF2. The small change is explained by the fact that CF3 cannot use flexibility of the EC-B so much as CF2. However, there is also an interesting coincidence of distance between this CF3 and previous two configurations. When the CO2 molecule is captured at the minimum at 5.5 /~, the distance between the C-CO2 and the H-B is 3.0/~. Namely, in both sides of the benzene plate, the dispersion force makes the O-CO2 and the EC-B interact almost at the same distance as CF1 and CF2. CF4 is explained by a different mechanism from the above three configurations. The HF wavefunction gives a shallow minimum of-1.5 kJ.mo1-1 at 6.3 /~, and the MP2 wavefunction has a deeper minimum of-4.2 kJ.mo1-1 at 6.2 /~. In this CF4, one of the O-CO2 ( 6 - ) approaches the H-B (6+) and causes an electrostatic attractive interaction with the H-B. Although an electron correlation effect makes the minimum deeper than the HF case, this configuration is described mainly by the electrostatic interaction. A dispersion force does not work well in this CF4. As summarized in the previous paper (Domafiski et al., 1992), a CO2 molecule is characterized by the partial charges at the carbon atom (6+) and the oxygen atoms (6-). On the other hand, a benzene molecule is described by the flexible electron cloud on the plate ( 6 - ) and the hydrogen atoms (6+). The potential energy curve profile strongly depends on the relative configuration of the benzene-CO2 system, and the contents are explained by the combinations of those same or different partial charges between both molecules. Moreover, the dispersion force works well and determines the minimum position of the potential energy curves at some configurations where the O-CO2 can interact with the EC-B.

J.-W. Shen et al. / Fluid Phase Equilibria 104 (1995) 375-390

382

A potential model

We fitted the present ab initio data to an analytical potential function using site-site LJ + Coulomb interaction terms as mentioned above. As a total, 450 ab initio data for the benzene-CO2 system have been used for optimization of coefficients in equation (1). The optimized results are given in Table 1 and summarized in Fig.3. In this figure, calculated energy values by the present potential function are plotted against original ab initio data for each configuration. Table 1 The potential parameters of the benzene-CO2 system.

i,j

C-C C-O H-C H-O

a,j (/~)

Eli (kJ.mo1-1) 1.22 x 10-4 8.97x10 -1 3.95x10 -23 6.14x10 -1

5.9838 3.0207 177.68 2.3314

0

q, qj (.~ . kJ.mo1-1)

-76.065 38.032 76.065 -38.032

2O

E

•

•

C

.o_ "5

10

t-

I.t_ "0 N

"~

0

Q..

0 "6 ~ -10 ¢-

UJ

-10

0

10

20

ab initio Energy [kJ/mol]

Fig.3 Comparison between energies by the optimized function and the ab initio data.

383

J.-W. Shen et al. / Fluid Phase Equilibria 104 (1995) 375-390

(a) I

I

I

I

I

I

81

I

I

I

I

(b) I

I

I

\

6

"~ 4 N

_

8

-11.0 ..."" J"" -4.2

I

I

I

I

:-d;,., I

I

8

I

I

,,

I

I

I

I

I

I

,< N

4

2 -

-11.0

..-""

@ i,/" 0

"2

"

4 X [A]

6

8

0

"2

"

4 X [A]

6

8

Fig.4 Iso-energy contour maps of a CO2 molecule around a benzene molecule. Top: perpendicular to the benzene plate, bisecting C-C bond; middle: perpendicular to the benzene plate, including C-H bond; bottom: on the benzene plate. (a) A CO2 molecule approaches the benzene parallel to the plate. (b) A CO2 molecule approaches the benzene perpendicular to the plate. The detailed orientations of the CO2 molecule are explained in the text. Numbers in these figures indicate the depth of potential wells in kJ.mo1-1. The dotted diagonal lines on the top and middle maps show the division at an angle of 30° from the benzene plate. The dotted lines on the bottom maps correspond to the lines with 30 ° from C-H bond of the benzene molecule.

384

J.-W. Shen et al. /Fluid Phase Equilibria 104 (1995) 375-390

To check how the potential model reflects relative configurations of the benzene-CO~, system, iso-energy contour maps were prepared in Figs.4(a) and (b), where a benzene molecule is placed at the center of these maps. Fig.4(a) explains the case when a CO2 molecule approaches the benzene molecule parallel to the plate. In those three figures, the CO2 molecule is included on the map. In the bottom of Fig.4(a), the CO2 molecule is on the line which connects the CBM and the C-CO2. These cases described in Fig.4(a) include CF1 and CF4 as the special configurations . On the other hand, Fig.4(b) explains the case when a CO2 approaches the benzene molecule perpendicular to the plate. In the top and middle of Fig.4(b), the CO2 molecule is included on the map. In the bottom of Fig.4(b), the C-CO2 is on the map and the CO2 molecule itself is perpendicular to the map. These cases described in Fig.4(b) include CF2 and CF3 as the special configurations . Characteristic points in Fig.4(a) are a deep potential well above the benzene plate and a shallow potential well beside the H-B. In the top of Fig.4(a), the former has a minimum of-11.0 kJ-mo1-1 at (0.0/~, 3.3/~) and the latter has a minimum of-4.2 kJ.mo1-1 at (5.5/~, 0.0/~). In the middle of Fig.4(a), the shallow well shifts to (6.1 /~, 0.0/~) due to the existence of a C-H bond. The minimum is -3.5 kJ.mo1-1 on this map. In the bottom of Fig.4(a), there is an interesting minimum of-4.3 kJ.mo1-1 at (5.0/~, 2.3/~), which is at 5.5/~ from the CBM. In this configuration, one of the O-CO2 faces the benzene molecule and makes an electrostatic attractive interaction with two H-B's. We call the configuration CF5, which is depicted in Fig.5(a) with the potential curve by the optimized function. Fig.4(b) is characterized by a shallow potential well in the intermediate configuration between CF2 and CF3; we call it CF6, which is depicted in Fig.5(b). In this CF6, there is a potential well which has a minimum of-4.8 kJ.mo1-1 at (4.0 A, 2.7/~) in the top of Fig.4(b). This minimum is at 4.8/~ from the CBM. As summarized in the previous section, CF2 and CF3 show only shallow minimum due to an electrostatic repulsive interaction. However, at this CF6, one of the O-CO2 can make an electrostatic attractive interaction with the H-B. Moreover, the C-CO2 can also make an electrostatic attractive interaction with the EC-B. In the middle of Fig.4(b), the well slightly shifts to (3.8/~, 3.3/~), which is at 5.0/~ from the CBM with a minimum of-4.6 kJ.mo1-1. This shift is explained by reducing the numbers of hydrogen atoms which effectively interact with the O-CO2. In the bottom of Fig.4(b), there is also an interesting minimum of -2.2 kJ.mo1-1 at (4.1 /~, 2.4/~), which is at 4.8 /~ from the CBM. In this configuration, at both sides of the benzene plate, the O-CO2 make a dispersion interaction with the EC-B. We call the configuration CF7, which is depicted in Fig.5(c). Since the minimum position

J.-W. Shen et al. / Fluid Phase Equilibria 104 (1995) 375-390

385

is 3,5 /~ from the benzene ring, the situation due to the dispersion force is almost same as those of CF1, CF2, and CF3. 15

15

E

10

10

2

s

5

o

I

I

I

I

I

I

0

III

"~

(b)

-5

-5

~ -10

-10

t'-

El_

I

-15

I

I

I

I

-15 0

15 o

E

2

4 6 Distance [A]

I

8

10

10 5

~

e..-

0

LLI

-~

-5

N -10 -15 0

I 2

I I 4 6 Distance [,&]

I 8

Fig.5 Potential energy curves by the optimized function: (a) CFb, (b) CF6, and (c) CF7. These results indicate that our model can reflect fairly several changes in the potential energy surface with the relative configuration of the benzene-CO2 system. Following observations in Fig.4 are extremely important for investigating characteristics of solvent structures around a solute. 1. At the short distance region from the CBM, there are noticeable potential wells. These wells capture CO2 molecules around the benzene molecule. Since locations of wells are different between Fig.4(a) and (b), these wells are expected to control relative orientations of CO2 molecules to the benzene molecule. Therefore, in the fluid structure, several characteristic features are expected to appear around these regions.

386

J.-W. Shen et al. / Fluid Phase Equilibria 104 (1995) 375-390

2. Beyond 7 /~ from the CBM, there is no significant character on these contour maps. Therefore, in the fluid structure, no difference is expected between these regions and the pure CO2 fluid.

Results from MC simulations In this section, we present our preliminary results of MC simulations and discuss orientations of CO2 molecules around a benzene molecule. The fluid structures are summarized by an RDF (radial distribution function) from the CBM. Since a CO2 molecule has a symmetry of Dooh, we can define two kinds of RDFs; one is between the CBM and the C-CO2 (g(r)center-c) and another is between the CBM and the O-CO2 (g(r)center-o). In order to examine how a plate-like molecule causes heterogeneity in distribution of the CO2 fluid, we divided the space around a benzene molecule into two regions, namely, in-plane part and out-of-plane part. The definition is based on whether the C-CO2 in question is in the benzene plate or in the normal direction to the benzene plate. Fig.6 describes these two regions. Therefore, four kinds of RDFs are defined: (a) out-of-plane part of g(r)center-c, (b) out-of-plane part of g(r)center-o, (c) in-plane part of g(r)center-c, and (d) in-plane part of g(r)center-o. In Figs.7(a)~(d), each RDF is shown as a function of temperature from 310K to 400K.

normal direction

"•-or-plane

l

~- e

Fig.6 Definition of in-plane and out-of-plane parts of the CO2 fluid around a benzene molecule. In Fig.7(a), the first peak is located at 3-6 /~. From the top and middle of Fig.4(a) and (b), two configurations, CF1 and CF6, can contribute to this peak. W h e n a CO2 molecule approaches a benzene molecule parallel to the

387

J.-W. Shen et al. / Fluid Phase Equilibria 104 (1995) 375-390

plate as CF1, the molecule is captured by a deep well at (0.0 /~, 3.3 /~) in the top and middle of Fig.4(a). This potential well has an extremely large interaction energy of-11.0 kJ.mol-land firmly controls an orientation of the CO2 molecule to the benzene plate at the lower temperature. The higher temperature of the system reduces the control and awkwardly makes the RDF lower. On the other hand, a COs molecule approaches a benzene molecule perpendicular to the plate as CF6, the CO2 molecule can feel a shallow well at 4.8-5.0/~ from the CBM. The outer part of the first peak can be partially explained by the approach to the edge of benzene plate. 3

i

i

i

(a) ,o

i

3

(b)

310K

2

~

320K . . . .

-~,":.:,,\ .':-

350K . . . . .

?

2

I

I

1

A

....

//"t,~

_

I

350K . . . . .

E o

"c"

0 0

2

3

I

4

r[Al

I

6

8

0

0

0

2

4

6

I

(c)

3

I

(d)

31 OK

I

.

'310K " I 320K . . . . I 350K . . . . . |

:::: ,o

8

r[A]

2

oIL

G)

2

o

0

...... I....... i 2

4 r [A]

I

I~'~~

6

8

...... ,...... J- .... i....... i ..... 0

0

2

4

r [,h,]

6

8

Fig.7 Temperature dependence on RDFs: (a)out-of-plane part between the CBM and the C-CO2, (b)out-of-plane part between the CBM and the O-CO2, (c)in-plane part between the CBM and the C-COs, and (d)in-plane part between the CBM and the O-COs. In Fig.7(b), there is a sharp first peak at 3-5/~. Since this RDF, g(r)center-

0

388

J.-W. Shen et al. / Fluid Phase Equilibria 104 (1995) 375-390

o, measures both O-CO2, the clear peak explains t h a t the oxygen atoms located at the almost same distance from the CBM. If a CO2 molecule approaches a benzene molecule parallel to the plate, the CF1 is consistent with the peak. On the other hand, the approach perpendicular to the plate in CF6 should cause a broad distribution to the peak. Therefore, in the first peak of the out of plain space, the CO2 molecule is expected to approach to the benzene molecule parallel to the plate. In Fig.7(c), the first peak starts from 4 ~ and the m a x i m u m is located at 5.5 /~. Since a shallow potential well in the middle of Fig.4(a) is at 6.1 /~, the CF4 does not contribute to the first peak. The location of the first peak indicates the contribution of CF5 and CFT. The former has a m i n i m u m at 5.5 /~ from the CBM and the latter has t h a t at 4.8/~, respectively. Therefore, at the first peak of the in-plane space, a CO2 molecule is expected to approach a benzene molecule from the configurations which avoid the H-B. The R D F in Fig.7(d) has a different structure from other R D F s in Fig.7(a) ,-~ (c). In the last case, the height of the first peak around 5 /~ is relatively low even at the lower t e m p e r a t u r e and a big shoulder is seen from 5/~ to 7 /~. These results reflect t h a t the position from the CBM is relatively different between two O-CO2's. Fig.7(c) suggests the contributions of CF5 and CF7. If a CO2 molecule approaches a benzene molecule in CFT, a shallow well in the b o t t o m of Fig.4(b) should cause a sharp peak at 5 /~ in the RDF. The first peak of Fig.7(d) is explained partially by this CFT. If a CO2 molecule approaches a benzene molecule in the CF5 and is captured by another shallow well in the b o t t o m of Fig.4(a), the outer O-CO2 is located at 6 - 7 / ~ from the CBM. The vague peak in Fig.7(d) indicates a possibility that a CO2 molecule approaches a benzene molecule parallel to the plate in the in-plane space, too. Further investigations have been in progress for the b e n z e n e - C O 2 system by focusing on how the benzene molecule determines the structure of the neighboring CO~. fluid and how the fluid structure depends on the temperature and density of the system. CONCLUDING REMARKS Our final purpose of this work is to investigate the CO2 fluid structures around a benzene molecule in the supercritical region of the system. For this purpose, we tried to construct a solute-solvent intermolecular potential function which is usable to computer simulations of the system. First, ab initio q u a n t u m chemical calculations were carried out for a syst e m of a benzene molecule and a CO2 molecule with large numbers of relative configurations. We have adopted 6-31G* basis set and included electron cor-

J.-W. Shen et al. / Fluid Phase Equilibria 104 (1995) 375-390

389

relation effects by MP2 theory because some configurations are considered to be characterized by only the dispersion interaction. The detailed features of these ab initio calculations are explained by using four representative configurations. A CO2 molecule has partial charges at the carbon a t o m (6+) and the oxygen atoms ( 6 - ) . On the other hand, a benzene molecule is characterized by the flexible electron cloud on the plate ( 6 - ) and the hydrogen atoms (6+). The potential curve profile is mainly described with an electrostatic interaction caused by the combinations of those same or different partial charges between b o t h molecules. Moreover, the dispersion force plays an i m p o r t a n t role in determining the m i n i m u m position of the potential curve at some configurations. Next, the ab initio d a t a were fitted to an analytical potential function usable in molecular simulations. We confirmed t h a t our model sufficiently reproduced the ab initio d a t a of each configuration. Moreover, iso-energy contour maps were prepared to obtain an overview of the potential model. These maps were extremely insightful to consider the solvent structures around the solute benzene. Finally, using the potential models, we have done preliminary MC calculations for the CO2 fluid containing one benzene molecule at the supercritical region. In order to examine the effect of the plate-like shape of the benzene molecule, the fluid structures around a solute were divided into the in-plane and out-of-plane parts. The analyses of RDFs from the CBM indicated that our model considerably reflected a shape of the plate-like molecule on the solvent structures. ACKNOWLEDGEMENT The present calculations have been carried out at Supercomputer Laboratory at Institute for Chemical Research, Kyoto University. This study has been supported by a Grant-in-Aid for Scientific Research from the Ministry of Education, Science, and Culture, Japan. REFERENCES Bruno, T.J., Ely, J.F., 1991, Supercritical Fluid Technology: Reviews in Modern Theory and Applications. CRC Press, Inc., Boca Raton Ann Arbor Boston London. Domafiski, K.B., Kitao, O., Nakanishi, K., 1992. The peculiar potential surface of the carbon dioxide dimer. Chem. Phys. Left., 199: 525. Domafiski, K.B., Kitao, O., Nakanishi, K., 1994. A new potential model for carbon dioxide from ab initio calculations. Molecular Simulation, 12: 343. Frisch, M.J., Head-Gordon, M., Trucks, G.W., Foresman, J.B., Schlegel, H.B., Raghavachari, K., Robb, M.A., Binkley, J.S., Gonzalez, C., DeFrees, D.J., Fox, D.J., Whiteside, R.A.,

390

J.-W. Shen et al. / Fluid Phase Equilibria 104 (1995) 375-390

Seeger, R., Melius, C.F., Baker, J., Martin, R.L., Kahn, L.R., Stewart, J.J.P., Topiol, S., Pople, J.A., 1990. Gaussian 90. Gaussian Inc., Pittsburgh PA. Kohn, P.M., Savage, P.R., 1979, Supercritical fluids try for CPI applications. Chem. Eng., 86: 41. Metropolis, N., Rosenbluth, A.W., Rosenbluth, M.N. Teller, A.H., Teller, E., 1953. Equation of state calculations by fast computing machines. J. Chem. Phys., 21: 1087. Shing, K.S., Gubbins, K.E., 1982, The chemical potential in dense fluids and fluid mixtures via computer simulation. Mol. Phys., 46: 1109. Shing, K.S., Gubbins, K.E., 1983, The chemical potential in non-ideal liquid mixtures. Computer simulation and theory. Mol. Phys., 49: 1121.