Energy & Fuels 1992,6, 771-778
77 1
Computer Simulation of the Molecular Structure of Bituminous Coal G.A. Carlson Fuel Science Department, Sandia National Laboratories, Albuquerque, New Mexico 87185 Received June 4,1992. Revised Manuscript Received August 26, 1992
Molecular modeling techniques have been used to study the three-dimensional structure and the energetics of the Given, Wiser, Solomon, and Shinn molecular models of bituminous coal. These studies demonstrate the importance of van der Waals (vdW) interactions and hydrogen bonding in the formation and stabilization of coal macromolecular structure. VdW interactions are responsible for most of the stabilization in these models, and a tendency for the relative importance of vdW interactions to increase with increasing rank is indicated. The strengths of vdW interactions for a series of aromatic and saturated ring molecules are determined and related to those in the coal structures. Physical densities and microporosities calculated for the simulated coal structures are in good agreement with those observed experimentally.
Introduction The structure of coal is inherently complex and varies widely depending on the origin, history, and age or rank of the particular coal examined. Nonetheless, because of the relationship of coal's structure to its reactivity in combustion, pyrolysis and liquefaction processes, there have been many studies to define its molecular (chemical) and conformational (physical) structure and properties. A number of these studies have focused on the molecular details of bituminous coal.14 In these studies, structures were derived using data from a variety of sources, including coal atomic composition; analysis of products from chemical reactions or from liquefaction or pyrolysis treatments of coal samples;GPC and VPO; and spectroscopic analyses using NMR, IR, and FTIR. The molecular models developed over the past 30 years are typically based on clusters. Clusters are molecular substructures consisting of, on the average, three or four fused aromatic rings that contain appropriate numbers of 0, N, and S atoms and have attached short-chain aliphatic groups. The clusters are connected together by hydroaromatic, etheric, or short aliphatic linkages. Because coal is a heterogeneous material, these models are understood to be representations of average coal structures and not to be quantitatively accurate descriptions of a particular coal structure. These molecular structures, while providing much helpful information on the chemical nature of coal, do not provide data on its three-dimensional structure or on the intercluster interactions that provide the basis for many of the physical properties of coal. For example, the glassy nature of coal, the glass-to-rubber transition that coals go through when heated, and the nature of coal-solvent interactions are not explained adequately by a knowledge of coal molecular structure alone. To obtain a more complete picture of the physical characteristics of coal (1) Given, P. H. Fuel 1960,39,147-153. (2) W k r , W. H. NATO ASI Series C 1983,124,325350, (3) Solomon, P. R. New Approaches in Coal Chemistry; ACS Sympoeium Series No. 169; American Chemical Society: Washington, DC, 1981; pp 61-71. (4) Shinn, J. H. Fuel 1984,63, 1187-1196.
related to these molecular structures, Spiro constructed three-dimensional representations of the Given,l Wiser,2 and Solomon3 structures using space-filling physical models.5 Because of steric problems, most of the molecular structures had to be altered somewhat before they could be constructed. Spiro's models provided insight into the intercluster interactions and the degree of ring alignment (stacking) that might occur in bituminous coal. These observations led to a plausible explanation for coal plasticity based on sliding of ring systems, facilitated by small aliphatic fragments released in the early stages of pyrolysis reactions. A limitation of these space-fillingmodel studies was that energeticsof the various structures and structural conformations could not be determined. Thus, the relative probability of the various three-dimensional structures was not established. Recently, with the development of molecular modeling software: it has become possible not only to visualize molecular structures in three dimensions (on the computer screen), but also to calculate energetically favorable structural conformations using molecular mechanics and molecular dynamics methods. Molecular modeling techniques are being used widely today to provide insight into the structure, properties and interactions of biomolecules (e.g., enzymes, proteins, inhibitors)' in order to guide the design and synthesis of pharmaceuticals. In the area of fuel science, variations of these techniques have been applied to visualize asphaltene structures (but not to calculate minimum-energy structures)Band to construct average molecular representations of kerogen macromole c u l e ~and, ~ recently, bituminous coalslO based on ana~
(5) Spiro, C. L. Fuel 1981,60, 1121-1126. (6) Fruhbeis, H.; Klein, R.; Wallmeir, H. Angew. Chem., Znt. Ed. Engl. 1987,26,403-418. (7) Bajorath, J.; Kraut, J.; Li,Z.; Kitaon, D. H.; Hagler, A. T. Proc. Natl. Acad. Sci. 1991,8%,6423-6426. (8) Robinson, K. K. Proceedings of the Electric Power Research Institute Conference on Coal Structure, Palo Alto, CA; EPRI, 1987; pp 3-1-3- 11. (9) Faulon, J. L.; Vandenbroucke, M.; Drappier, J. M.; Behar, F.; Romero, M. Adu. Org. Geochem. 1990,16,981-993. (10) Faulon, J. L.; Hatcher, P. G.; Wenrel, K. A. Prepr. Pap.-Am. Chem. SOC.,Diu. Fuel Chem. 1992, 37(2), 900-907.
0887-0624/92/2506-0771$03.00/00 1992 American Chemical Society
772 Energy & Fuels, Vol. 6, No. 6, 1992
lytical data. In a preliminary report, we constructed and energy-minimized four previously proposed molecular models of bituminous coal.ll Here, we present more detailed information on these computer-generated molecular models, including the results of molecular dynamics calculations and the importance of van der Waals and hydrogen-bonding interactions.
Experimental Section Molecular modeling studies were performed using BIOGRAF (Molecular Simulations, Sunnyvale, CA), a software program allowing construction, visualization, and energy calculations. Studies were performed on a Silicon Graphics 4D/320GTXB workstation. The modeling program allows the user to construct molecules on a computer screen and to visualize them in three dimensions (usingrotational and translational motion and depthcueing of the computer images) in a variety of formats (stick; dot surface; space-filling representations). The user can also carry out energy minimization and molecular dynamics calculations using a force field approximation. The force field methodology employed in this software (DREIDING, DREIDING/A)12 is described in greater detail below. The DREIDING force field used in the majority of calculations in the current study is a very general force field that can be used for a large number of atom types. It allows calculation of the total energy of a molecular structure as a sum of bonding interactions (E, = stretch; E b = bend; Et = torsion; and Ei = inversion) and nonbonding interactions (E, = van der Waals; E, = electrostatic; and Eh = hydrogen bonds).
E = (E, + E, + E, + E,) + (E,+ E, + E,) The DREIDING force field does not include higher-order terms such as stretch-bend interactions, atomic dipoles, etc. In common with other force field methods, it does not permit bonds to be broken or made during calculations. Because van der Waals (vdW) and hydrogen bond (HB) interactions are considered in most detail in this paper, their method of calculation will be discussed more fully than the other interactions. VdW interactions affect all atom pairs in a structure except for those that are already coupled through bond, angle, or torsion interactions. VdW interaction energies are calculated in DREIDING using a Lennard-Jones 12-6 potential, with a different vdW bond strength (well depth) DO,scaled distance p= Rm/Ro, and vdW bond length Ro assigned for each different A-B atom pair.12
E , = Do{pE- 2p;}S A switching function S is employed to gradually "cut off" the vdW calculations for atom pairs more than 8 A apart. Beyond this distance, the interactions become relatively unimportant because of the l / F dependence. For the larger structures in the current study, use of this cutoff distance reduced the calculational time by roughly &fold,while affectingthe calculatedvdW energies by less than 3 % . Hydrogen bond interaction energies (&) are calculated using a Lennard-Jones 12-10 potential based on the scaled distance (PDA = RDdRh, with RDA the distance between donor and acceptor and Rh the optimal hydrogen bond distance) and an additional dependence on the angle e D n A between donor atom, hydrogen atom, and acceptor atom.12
In these calculations, all oxygen, nitrogen, and sulfur atoms are (11)Carleon, G. A.; Granoff, B. Prepr. Pap.-Am. Chem. SOC.,Diu. Fuel Chem. 1989,34(3), 780-786. (12) Mayo, S. L.;Olafaon,B. D.; Goddard, W. A. J . Phys. Chem. 1990, 94.8897-8909.
Carlson considered to be potential hydrogen bond acceptors and donors. For all three atom types, the equilibrium distance Ro between donor and acceptor is defined to be 2.75 A, and the equilibrium hydrogen bond strength DOis fixed at 9.00 kcaVmol (DREIDING) or 9.50 kcal/mol (DREIDING/A). After a molecular model is built on the computer, the minimumenergy structure or conformation can be identified using energy minimization or molecular dynamics techniques.6J2 In energy minimization, the positions of atoms in a structure are perturbed in a stepwise process, using either conjugate gradient or steepest descents techniques, until a minimum-energy conformation is found. This method is most successfulwhen applied to relatively simple structures, for which a single low-energy conformation exists. For more complicated structures, which often include a number of low-energyconformationsseparated by energy barriers, molecular dynamics techniques are employed. In molecular dynamics calculations, the individual atoms in a structure are allowed to move with velocitiesrepresentative of vibrational and rotational motion at a chosen "temperature". The interactions of the atoms are defined as before by the force-field equation. Newton's equations of motion are then solved to determine the motion of each atom as a function of time. The motion of the a t o m allows them to surmount potential energy barriers between different structures and thus not to remain trapped in a local energy minimum. In these calculations,the energiesare evaluated and velocities are reassigned frequently (at roughly 1-fa intervals; 1 fs represents a small fraction of a vibrational frequency for a molecule), causing the structures to twist and fold in ways that generally tend to optimize the nonbonding interactions (van der Waals, electrostatic, and hydrogen bonds), while maintaining appropriate bond lengths and angles. Molecular dynamics calculations typically require many thousands to millions of evaluation steps, representing the equivalent of picoseconds to nanosecondsof molecular motion, to sample a reasonable number of the possible minimum-energy conformations. In the current studies, the "temperature" in the dynamics calculations was initially held at 300 K for 30 pa, and then ramped downward from 300 to0K over a 3 0 - p period (using1-fa dynamics steps). The initial 30-pa period permitted the coal structures to move from the initial structure through a number of different conformations, generally becoming more folded (Le., threedimensional) and lower in energy with time. During the last 30-pa period, the structure relaxed into one minimum energy conformation as the temperature dropped. This conformation may be the global energy minimum or a local energy minimum, because there are likely to be a number of similarly low-energy conformations for a molecule as complex as coal. Our experience is that most of these low-energy conformations have similar properties (e.g., energy, density, and porosity). A 60 OOO-step, 60-ps dynamics run for a 1311-atom structure, with over 100 OOO separate calculations performed per step, required about 25 h of computation. Most calculations in this work were performed using DREIDING, the most recent (1990) version of the force field. Only the vdW calculations of model ring compounds (and the associated calculations of coal vdW interactions) were done with the earlier DREIDING/A force field. Results from the two force fields were similar in most regards (the absolute energies from DREIDING/A tended to be about 10% lower, for undetermined reasons). Electrostatic interactions were found not to contribute significantly to the energy minimization in early DREIDING calculations for the coal models, and are not included in the calculations reported here. After minimum-energy structures were obtained, a separate computer calculation was utilized to determine their physical densities. In this calculation, a rectangular box just large enough to contain the structure was fist defined and divided into uniform 1A3 volume elementa. Then all volume elements were assigned, by comparison with the coordinates of all atoms in the structure, to either structural volume, void volume, or external volume. Initially, volume elements 'within" the structure were defined as those included within a surface stretched across all outer atomic
Molecular Structure of Bituminous Coal
Energy &Fuels, Vol. 6, No. 6,1992 773
''1L
Figure 1. Molecular models of bituminous coal from the literature: ref3; (d) Shinn model, ref 4. volume elements less than 4 8, apart. Void volume was then defined as all volume elements that were within the structure in a pore or opening a t least 4 8, in minimum dimension. This minimum pore size was chosen to approximate the size of aprobe atom or molecule that might be used in a density determination. The structural volume was defined as all remaining volume elements that were within the structure. All volume elements not within the structure were assigned to external volume. Once the volume elements were all assigned, the density of the structure was calculated from the molecular weight of the structure, multiplied by Avogadro's number, and divided by the appropriate volume. "True" densities, as might be determined experimentally using helium displacement, were calculated from the structural volume. Microporosity was calculated from the void volume and expressed as a percentage of the total volume within the structure. "Particle" densities, as would be measured by displacement of mercury, could not be calculated, because our models are much too smallto allow a determination of mesoporosity or macroporosity.
Modeling Results and Discussion Three-dimensional representations based on four postulated bituminous coal molecular models, those of Given, Wiser, Solomon, and Shinn,1-4 were created, energyminimized, and analyzed using molecular modeling techniques. These four models, which range in size from about 200 to 1300 atoms, were chosen as representative of those that have appeared in the literature over the past 30 years.
6H
OH
Given model, ref 1;(b) Wiser model, ref 2; (c) Solomon model,
The models from the literature are shown in Figure 1. Two of the four models were modified somewhat prior to computer modeling. To satisfy steric requirements, the Given' model was modified slightly by replacing a strained quaternary carbon bond at the center of the structure with a ternary bond, as suggested by S p i r ~ .A~later version2 of the Wiser coal model with the disulfide bridge replaced by an ethylene linkage was used. Solomon's model3 was significantly modified by building three replicates of the smallest molecular fragment and connecting one to each of the three indicated open bonds on the two larger fragments. This increased the size of the model and hence the potential for intercluster interactions. The Shinn model4 was used without change. For each computer model, any peripheral open bonds indicated in the original model were capped with methyl groups. Figure 2 shows the computer modeling representations of the Shinn molecular structure. The model was built on the computer in a quasi-two-dimensional conformation, shown from the top view in Figure 2a. Following a total of 60 ps molecular dynamics and subsequent energy minimization (Figure 2b), the computer model has obtained a more three-dimensional shape and is much more compact. The Shinn model is illustrated because it represents the largest and most complex of the molecular models studied. Computer-modeling results for the Wiser and modified Solomon models showed structural folding
Carlson
774 Energy & Fuels, Vol. 6, No. 6, 1992
a
b
-
Figure 2. Computer modeling of 1311-atom Shinn model. (a) Computer model as originally created, top view, shown both as stick figure and as solid-filled representation. (b) Energy-minimized computer model after 60-ps molecular dynamics calculation, top view, shown both as stick figure and as solid-filled representation. Atoms are color-coded turquoise or gray (carbon), white (hydrogen), red (oxygen), blue (nitrogen), and yellow (sulfur). 3000 2500
2 0
2000
2 cd
8 B
3
1500 1000
5oo
E, , !,
E;;Y,
I , ,
, ,
,,,,,,,,
Covalent Energy Non-Covalent Energy
,
,I
0 0
10
20
30
40
50
60
Time (ps) Figure 3. Total, covalent, and noncovalent energy terms for the Shinn computer as a function of time during a 60-ps molecular dynamics calculation. The points represent the energies of the structure at 1-ps intervals, following an additional 50 steps of energy minimization of the structure. For the first 30 ps, the “temperature” of the dynamics run is 300 K. For the last 30 ps, the temperature ramps down from 300 to 0 K a t a rate of 1K/O.l PS.
similar to Shinn. The Given model, because of the stiff linkages between molecular clusters, was unable to fold into a compact structure. Figure 3 shows the total energy, the energy of bonding interactions, and the energy of nonbonding interactions for the Shinn structure during the 60-ps molecular
dynamics run. The energiesshown represent the dynamics results at 1-ps intervals, following 50-step minimization of the structure after each dynamics interval. The first 30 ps of dynamics was run at a constant 300 K temperature. The last 30 ps of dynamics was run with the temperature lowered from 300 to 0 K, at 1K/0.1 ps. For about the first 20 ps, a rapid decline in total energy with time was observed as the structure rapidly folded into a more compact form. Following this, the energy remained relatively constant, as the structure moved between more or less equivalent conformations. The energy breakdown into bonding and nonbonding interactions shows that most of the reduction in energy with time is due to reductions in nonbonding energy (hydrogen-bonding and van der Waals interactions), while bonding energy (stretches, bends, torsions, etc.) remains relatively constant after the first few picoseconds. Table I shows the molecular parameters for all four models, along with the energies calculated for the minimum energy conformations. The energies are shown normalized for easier comparison. Three of the models show similar (normalized)minimum energies. The higherenergy exception is the Given model, which did not fold significantly and thus had limited non-bonding-energy stabilization. Figure 4, illustrating a portion of the energy-minimized Wiser structure, demonstrates graphically the importance of nonbonding interactions in structural folding. A
Molecular Structure of Bituminous Coal
Energy & Fuels, Vol. 6, No. 6, 1992 775
Table I. Molecular Parameters and Calculated Energies for Bituminous Coal Models parameter Given Wiser modified Solomon Shinn “bituminous” type of coal “low rank” PSOC 170 vitrinite-rich high-vol bituminous 192 393 no. of atoms 396 1311 1492 2967 mol wt 3020 9956 0.70 0.74 0.66 0.71 CAr/CTot 0.21 0.28 0.34 0.40 HAJ(HA~ + HAI) wt fraction C 0.820 0.782 0.823 0.789 H 0.053 0.059 0.056 0.057 0.113 0 0.107 0.090 0.119 0.019 0.014 N 0.009 0.015 0.032 S 0.021 0.019 normalized formula C100H7709.8N2.0 C100H90010.9N1.6S1.6 ~100H8lO8.2Nl.0Sl.0 C100H87011.3N1.7S0.9 normalized energy of minimized 2.07 1.78 1.75 1.65 struct, (kcal/atom)
E
0
c
0
I
I‘ Figure 4. A portion of the minimum-energy conformation of the Wiser structure, demonstrating the importance of vdW and hydrogen bonding in stabilizing the three-dimensionalstructure. The parallel stacking of two adjacent ring structures (hydroaromatic near the top, aromatic near the bottom, shown side-on) illustratesthe favorable vdW interactionsof fused ring structures. An instance of hydrogen bonding is shown by the dashed line between a phenolic -OH species on the lower aromatic structure and a carbonyl oxygen on the upper hydroaromatic structure. Atom color coding is same as in Figure 2.
hydrogen bond between a phenolic -OH (lower center) and a carbonyl oxygen (right center) is seen in the figure. Attractive van der Waals (vdW)interactions are illustrated by the stacking of two adjacent ring structures, one near the top and one near the bottom of the figure. Because the nonbonding interactions appeared to provide the main driving force for the development of threedimensional structure and minimum-energy conformations, these were explored in more depth. Table I1 details the hydrogen bond interactions in three of the models, again normalized for ease of comparison. The Given model is not included, because it did not undergo folding or hydrogen bonding to any significant extent. Each model detailed in Table I1 contains roughly the same relative number of hydrogen bond acceptors and donors. However, the relative number of hydrogen bonds formed in the energy-minimized modified Solomon structure was very much lower than in the other two models. Possible reasons for this observation are discussed below. Intercluster hydrogen bonds (Le., those between donors and acceptors on different clusters) are differentiated in the table from total hydrogen bonds, ’because only intercluster bonds
Table 11. Van der Waals and Hydrogen Bond Interactions for Bituminous Coal Models modified parameter Wiser Solomon Shinn H-bond donors/100 atoms 4.1 3.6 3.0 H-bond acceptors/ 100 atoms 6.6 6.9 5.1 4.3 total H bonds/ 100 atoms 6.1 0.8 intercluster H bonds/ 100 atoms 3.8 2.3 0.3 mol wt/intercluster H bond 201 332 2967 -123 total intercluster H-bond energy (kcal) -38 -6.8 intercluster H-bond energy/ -2.6 -4.1 -6.8 intercluster H bond (kcal) total vdW interactions 23947 103015 21532 total AvdW energy (kcal) -61 -226 -80 AvdW energy/vdW interaction (kcal) -0.003 -0.004 -0.002
stabilize folded structures. Intracluster hydrogen bonds (i.e., those between donors and acceptors on the same cluster) do not provide a driving force for folding of the structure, because they can be formed with or without structural folding. The relative numbers of intercluster hydrogen bonds for the Wiser and Shinn structures are in agreement with solvent swelling data,13 which suggest roughly one hydrogen bond is present for each 200-300 molecular weight (-25-40 atoms) of Illinois No. 6 and Bruceton coals. Experimentally, hydrogen bond strengths vary widely (typically, from 2 to 10 kcal), depending on the donors and acceptors involved and their associated molecular structures.14 The DREIDING force field treats all hydrogen bonds as equivalent, with a maximum hydrogen bond energy of 9.0 kcal/mol for optimallyoriented interactions, regardless of atom type. However, because of nonoptimal angles or distances, the average energies per hydrogen bond calculated for the different models considered here ranged from 2.6 to 6.8 kcal/mol. This range of hydrogen bond energies is generally in agreement with the values expected for hydrogen bonds in ~ 0 a l . l ~ Table I1 also shows calculated vdW interaction data for the same three models. The total number of vdW interactions calculated includes all pairs of atoms (except for bonded atom pairs and their next-nearest-neighbors) that are separated by less than 9 A. The energiesof specific vdW interactions vary greatly, depending on the distance of separation; however, it is clear that a very large number of vdW interactions, distributed over the entire structure, are involved in creating this important force for energy minimization. (13) Larsen, J. W.; Green, T. K.; Kovac, J. J. Org. Chem. 1985, 50, 4729-4735. (14) Pimentel, G. C.; McClellan, A. L. The Hydrogen Bond; W . H. Freeman: San Francisco, CA, 1960.
Carlson
776 Energy & Fuels, Vol. 6, No. 6, 1992
Table 111. Comparison of van der Waals and Hydrogen-Bonding Interactions in Bituminous Coal Models modified parameter Wiser Solomon Shinn AH-bond energyi100 atoms (kcal) -9.7 -1.8 -9.4 AvdW energyi100 atoms (kcal) -15.5 -20.2 -17.2 AH-bond energyiAvdW energy 0.63 0.09 0.55
Relative energy changes due to vdW and hydrogenbonding interactions are compared in Table 111. Although there are quantitative differences in the data for the three models, vdW interactions appear to provide a stronger driving force for energy minimization in each case than do hydrogen bonds. This is an interesting conclusion, since hydrogen bonding seems generally to have been regarded as the more important nonbonding interaction in bituminous c0al.'~J6 On the other hand, other noncovalent forces including vdW interactions have also been shown to have importance in bituminous coal s t r u ~ t u r e s . On ~~J~ the basis of the analysis of boiling point data, White and SchmidtIg concluded that vdW interactions were the dominant intermolecular forces in the liquefaction products from bituminous coal. Although there are some differences in molecular composition expected for bituminous coal and its liquefaction products (liquefaction products would have less oxygen, hence less hydrogen bonding),their results nonetheless support the contention that vdW interactions should also be important interactions in coal. The markedly lower stabilization energy attributable to hydrogen bonding in the modified Solomon structure may be related to the somewhat (-25%) lower number of hydrogen bond donors and acceptors in this model (Table 11) coupled with a higher degree of aromaticity (Table I). Because hydrogen bonds are pairwise interactions, lowering the number of both donor and acceptor atoms could have a more than linear effect on the probability of hydrogen bond formation. Also, aromatic structures have a higher vdW stabilization energy than nonaromatic structures, as shown later in this section. Finally, in a coal structure containing large fused ring systems (the clusters in the modified Solomon model contained on the average five or six fused rings), vdW and hydrogen bonding will be competitive effects. That is, a structural conformation maximizing hydrogen bonding will in all probability be quite different from one maximizing vdW interactions. Thus, in the modified Solomon model, the combination of less species able to hydrogen bond, coupled with a higher driving force for vdW interactions due to high aromaticity, appears strongly to favor vdW interactions. If these findings can be generalized, they suggest a strong effect of coal rank on the relative importance of vdW vs hydrogen bonding. Experimental evidence supports a shift with rank from a strong hydrogen bond dependence for low-rank coals20,21 to a dominance by vdW interactions between polynucleararomatic systems for high-rank ~ o a l s . ~ ~ ~ ~ ~ ~~
(15) Liotta, R.; Rose, K.; Hippo, E. J. Org. Chem. 1981, 46, 277-283. (16) Stenberg, V. I.; Baltisberger, R. J.; Patal, K. M.; Raman, K.; Woolsey, N. F. In Coal Science; Gorbaty, M. L., Larsen, J. W., Wender, I., Eds., Academic Press: New York, 1983; Vol. 2, pp 125-171. (17) Brenner, D. Fuel 1985, 64, 167-173. (18) Nishioka, M.; Larsen, J. W. Energy Fuels 1990, 4 , 100-106. (19) White, C. M.; Schmidt, C. E. Fuel 1987, 66, 1030-1035. (20) Larsen, J. W.; Baskar, A. J. Energy Fuels 1987, 1 , 230-232. (21) Padrich, T. D.; Lockwood, S.J. Prepr. Pap-Am. Chem. SOC., Diu.Fuel Chem. 1984, 29(1), 153-160. (22) Quinga,E. M. Y.; Larsen, J. W. Energy Fuels 1987, 1 , 300-304.
b
= Aromatics Saturated Rin'gs = Tetralin = Paracvclophanes 1 = Coal Structur es
0
0=
I
M
.e
I ,
250 0
$
50.0
\
I -1.25-1
-0.75-0.50-0.25
'
1
0
1
0.25 0.50 0.75 1
1.25
vdW Energy (kcal/Carbon Atom) Figure 5. Van der Waals stabilization energies calculated for the interactions of pairs of molecules from several families of fused ring hydrocarbons using the DREIDING/A force field. The stabilization energy is the difference in energy for the molecular pair a t optimum spacing, relative to the energy for the pair at large distance. For ease in comparing the data from different families, the energies are normalized to the number of carbon atoms per molecule for aromatics, hydroaromatics, and tetralin; to half the number of carbon atoms for the paracyclophanes; or to the number of non-hydrogen (heavy) atoms per coal cluster. Likewise, for paracyclophanes the molecular weight is taken as half the molecular weight of the parent molecule, and for coal the average cluster molecular weight is used.
To cast more light on the role of vdW interactions in bituminous coal and its liquefaction products, a variety of ring compounds covering the range of cluster sizes believed to be present in bituminous coal was also investigated. Both aromatic and fully-saturated ring compounds, ranging in size from benzene and cyclohexane (CS)to coronene and perhydrocoronene (c24),were studied. One partially saturated molecule (tetralin) was also included for comparison. Calculations of vdW interactions were made for energy-minimized pairs of molecules, first well-separated and then at the optimum position for vdW interaction. The difference in energies was taken to be the vdW stabilization energy. The minimum-energy conformations for the molecules were generally with the rings parallel and overlapping, usually with some rotation of the rings relative to one another. VdW calculations were also performed on [2.2lparacyclophane, [3.3lparacyclophane, and [4.4lparacyclophane, a series of highly-strained molecules consisting of pairs of paraxylenes having the methyl groups at both ends linked together by zero, one, or twomethylene groups, respecti~ely.~~ The paracyclophanes were intended to represent the extreme case of coal clusters coupled together in a highly-strained fashion. For comparison with the vdW calculations for the ring compound pairs, each paracyclophane molecule was viewed as if it were a speciallyconstrained (linked) pair of para-substituted benzene molecules. Thus, an energy calculation for each paracyclophane was compared with that for the two substituted benzene molecules resulting from breaking the paracyclophane end linkages and adding hydrogens to cap the resulting radicals. The difference in energy of the wellseparated product molecules relative to the paracyclophane molecule was taken to be the vdW stabilization energy. (23) Mallya, N.; Stock, L. M. Fuel 1986, 65, 736-738. (24) Cram, D. J.; Hornby, R. B.; Truesdale, E. A.; Reich, H. J.; Delton, M. H.; Cram, J. M. Tetrahedron 1974, 30, 1757-68.
Molecular Structure of Bituminous Coal
Figure 5 shows the vdW stabilization energy data, normalized by the number of carbon atoms in the molecule (normalized by half the number of carbon atoms in the paracyclophanes). The data are plotted against molecular weight for each ring compound pair and against half the molecular weight for each paracyclophane. Also included in Figure 5 is the vdW stabilization energy for each coal model, plotted against the average molecular weight per molecular cluster in the structure. Here, the vdW data are normalized to the number of non-hydrogen atoms. The average molecular weight per cluster was determined from the total molecular weight of the model divided by the number of clusters in the model. The number of clusters was determined by inspection. The order of vdW stabilization (from most to least) is the Shinn, the modified Solomon, the Wiser, and the Given models. Three of the models show relatively similar stabilization energies. The exception is the Given structure, which folded very little during minimization and thus showed very little change in vdW energy. As seen in Figure 5, the vdW stabilization energies of different types of ring systems are distinctly different. For the aromatic ring compounds, there is a continuous increase in normalized vdW stabilization energy with ring size. This is perhaps not surprising, because for each atom there is a larger number of vdW interactions with the adjacent molecule as the ring size increases. For the fullysaturated ring compounds, on the other hand, there is a lesser degree of vdW stabilization compared with aromatics, and the increase in stabilization with molecular weight is not significant. Tetralin is found to be intermediate vdW stabilization. It appears that the out-ofplane hydrogens in the saturated ring compounds prevent the molecules from approaching to optimal distances for vdW interactions. The distances between equivalent carbons in the energy-minimized saturated ring pairs are from 4.7 to 5.1 A, significantly different from the C-C vdW distance (3.9A) or the distance between carbons in the aromatic pairs (3.5-3.7 A). The most striking results are for the paracyclophanes, which actually show a vdW repulsion rather than a stabilization for the energyminimized structures. The [2.2lparacyclophane is particularly strained, due to the relatively short -C2H4linkages to each benzene ring. The [3.31- and [4.41paracyclophanes, with longer linkages, are much less strained, although the vdW interactions are still modestly repulsive. By comparison with the various molecular pairs, three of the coal structures show vdW stabilizations that are intermediate in magnitude. This appears reasonable, since coal clusters contain both aromatic and hydroaromatic molecular species. Further, the clusters are coupled together (although not nearly so rigidly as the paracyclophanes), somewhat reducing their ability to undergo rearrangement to give optimum vdW interactions. The comparison of the coal structure data with that of the molecular pairs may not be completely justified, because the coal structures evaluated have multiple intercluster interactions, rather than the single pairwise interactions considered for the ring compounds. However, separate calculations have shown that the normalized vdW interactions for a given molecular pair are nearly the same (within a few percent) as those for energy-minimized groupings of three, four, or more molecules. Thus, the closest pairwise interaction is by far the dominant one in
Energy & Fuels, Vol. 6, No. 6, 1992 777 Table IV. True Density and Microporosity of Bituminous Coal Models parameter
true density (g/cm3) microporosity ( % )
Wiser 1.27 0.7
modified Solomon 1.31 0.0
Shinn
1.23 0.2
these systems. It also should be mentioned that a fullyoverlapping stacking of aromatic ring compounds is not generally the minimum-energy configuration observed experimentally; a herringbone stacking motif (with only partially-overlapping rings) is often favored, because of quadrupolequadrupole interaction^.^^ The current calculations, using a relatively simple force-field approximation, do not include quadrupole interactions. Finally, based on the minimum-energy conformations for the three models that folded effectively,"true" densities and microporosities were calculated. The results are given in Table IV. The calculated true densities averaged 1.27 f 0.04g/cm3,effectively the same as the experimentallydetermined true densities for bituminous coal (1.28-1.33 g / ~ m ~ )and , ~ 6for vitrinite macerals, the primary constituents of bituminous coal (1.25-1.30 g / ~ m ~ The ) . ~mi~ croporosities calculated for the structures were less than 1% in every case. This agrees with recent experimental data, which suggest less than 1% microporosity (pores < 1.5 nm diameter) in Pittsburgh No. 8 bituminous coal.28 Previous experimental data indicated values of 4-6 % for microporosity of hvC Becuase the quantitative measurement and evaluation of microporosity are at best very difficult, the correct value for microporosity of bituminous coal is uncertain. The modeling results are at least not incompatible with the available experimental data. It might be argued that the size of the structures in the current simulations is too small to accurately determine the microporosity that is present. It is certainly true that we are unable to simulate larger pores, which may be responsible for the largest part of the pore volume observed experimentally in most coaL28 .
Summary and Conclusions Four molecular models of bituminous coal have been studied using molecular modeling techniques. The threedimensional minimum-energy conformations obtained after extensive molecular dynamics calculations for three of the molecular models have provided insights into the types of bonding that provide the rigid structure observed for coals. In particular, both vdW interactions and hydrogen bonding appear to be powerful driving forces for the formation and maintenance of three-dimensional structure and rigidity, although vdW interactions appear to provide a somewhat greater degree of energy stabilization. The relative importance of vdW and hydrogen bonding may be very dependent on coal rank. Our simulations support a very limited degree of microporosity in the bituminous coal structure. Because our calculated structures are limited in size and the experimental data are not unambiguous, this subject requires more study. (25) Miller, J. H.; Mallard, W. G.; Smyth, K. C. J. Phys. Chem. 1984, 88,4963-4970. (26) Mahajan, 0. P. Coal Porosity In Coal Structure; Meyere, R. A., Ed.; Academic Press: New York, 1982; pp 51-59. (27) Crelling, J. C. Proc. 1987 Int. Conf. Cool Sci. 1987, 119-122. (28) White, W. E.; Bartholomew, C. H.; Hecker, W. C.; Smith, D. M. Adsorpt. Sci. Technol. 1990, 7, 1tN-209. (29) Gan, H.; Nandi, S.P.; Walker, P. L., Jr. Fuel 1972,51,272-277.
778 Energy & Fuels, Vol. 6, No.6,1992
In this paper, we have not considered the concept of network (cross-linked) coal structures. Computer simulation of network coal structures is in progress and will be the subject of a later publication. These network studies in particular may cast new light on the subject of microporoeity in coal structure. Additionally, we have not directly treated the subject of solvent interactions with coal structures. This is an obvious extension of the methodology presented here, because solvent interactions are based on the same nonbonded interactions considered
Carlson here. Coal-solvent interactions of cross-linked structures are probably most relevant. Modeling of solvent interactions and solvent swelling is planned for future studies.
Acknowledgment. This work was supported by the US. Department of Energy and Sandia National Laboratories under Contract DE-AC04-76DP00789. We thank Dr. Jean-Loup Faulon for contributions to the program for density calculations, and Dr. Malvina Farcasiu and Professor John Larsen for many helpful discussions.
