Hostname: page-component-848d4c4894-v5vhk Total loading time: 0 Render date: 2024-06-25T13:19:39.825Z Has data issue: false hasContentIssue false

Molecular Orbital Calculations of Water-Ice Surface Interactions*

Published online by Cambridge University Press:  30 January 2017

Patricia L. M. Plummer*
Affiliation:
Department of Physics and Graduate Center for Cloud Physics, University of Missouri-Rolla, Rolla, Missouri 65401, U.S.A.
Rights & Permissions [Opens in a new window]

Abstract

Detailed knowledge of the energetics of the interaction of a water monomer with an ice surface is important in understanding the factors which control both the nucleation and the growth of ice. In this work we model this interaction by using a small section of an ice surface and semi-empirical quantum-mechanical calculations of the CNDO/INDO type to investigate the energetics of a gas-phase molecule approaching the surface. In the attempt to simplify the system as much as possible while retaining the essential features of the "real" physical situation, 13 water molecules arranged in the lattice positions for ice Ih were used to represent the surface. Electron-density plots for both basal and prism surfaces were generated. Potential-energy curves were calculated for a monomer approaching a number of different sites on each type of surface. Monomer-surface binding energies varied from approximately 17 to 33 kJ/mol for the surface sites examined. Ion pairs of the H30+-OH- type were included in the surface and their effect on the energy of interaction with an incoming monomer calculated. For these studies no distortion of the ice surface was allowed. The energies found follow the same trends as observed for a free monomer-ion interaction.

Résumé

Résumé

La connaissance détaillée de l'aspect énergétique de l'interaction entre un monomère d'eau et une surface de glace est importante dans la compréhension des facteurs qui contrôlent àla fois la nucléation et la croissance de la glace. Dam ce travail, nous avons construit un modèle d'interaction en utilisant une petite section de la surface de glace et des calculs semi-empiriques de mécanique quantique du type CNDO/INDO utilisés pour calculer les caractéristiques énergétiques d'une molécule en phase gazeuse approchant la surface. De manière à simplifier le plus possible le système tout en retenant les facteurs essentiels de la situation physique "réelle", 13 molécules d'eau placées dans les positions du réseau de glace Ih sont utilisées pour représenter la surface. Des courbes de density électronique pour les plans de base et prismatiques ont été obtenus. Des courbes d'énergie potentielle ont été calculées pour un monomère approchant un nombre de sites différents sur chaque type de surface. Les énergies de liaison monomère-surface varient approximativement entre 17 et 33 kJ/mol pour les têtes dc surface étudiés. Des paires d'ions de type OH--H3O+ ont été inclus à la surface et leur effet sur I'énergie de l'interaction avec un monomére à proximité de I'interface est calculé. Pour ces études aucune distorsion de la surface de glace n'a été envisagée. Les énergies trouvées suivent les mêmes tendances que pour une interaction ion-monomère libre.

Zusammenfassung

Zusammenfassung

Genaue Kenntnis der Energieverhältnisse der Wechselwirkung eines Wassermonomers mit einer Eisoberfläche ist wichtig für das Verständnis der Faktoren, welche sowohl die Keimbildung als auch das Wachstum von Eis beherrschen. In dieser Arbeit bilden wir diese Wechselwirkung im Modell nach. Wir nehmen einen kleinen Ausschnitt einer Eisoberfläche und halbempirische quantenmechanische Berechnungen des CNDO/INDO Typs, urn die Energieverhältnisse eines Moleküls der Gasphasc beim Annähern an die Oberfläche zu untersuchen. Beim Versuch, das System soweit wie möglich zu vereinfachen und zugleich die wesentlichen Kennzeichen des "wirklichen" physikalischen Zustandes zu bewahren, wurden dreizehn Wassermoleküle benützt, die auf den Gitterplätzen von Eis Ih angeordnet waren, urn die Oberfläche darzustellcn. Elektronendichtediagramme wurden for Basis- und Prismenoberfläche aufgestellt. Es wurden Kurven der potentiellen Energie für ein Monomer berechnet, das sich verschiedenen Stellen auf jeder der beiden Oberflächenarten annähert. Die Bindungsenergien Monomer-Oberfläche für die untersuchten Stellen reichten von annähernd 17 bis 33 kJ/mol. Ionenpaare des Typs H3O+-0H- wurden in die Oberfläche mit eingeschlossen und ihre Wirkung auf die Wechselwirkungsenergie mit einem auftreffenden Monomer berechnet. Bei diesen Untersuchungen wurde keine Verzerrung der Eisoberfläche zugelassen. Die gefundenen Energien folgen dem selben Trend wie die beobachtetc freie Monomer-Ion-Wechselwirkung.

Type
Research Article
Copyright
Copyright © International Glaciological Society 1978

Introduction

The theoretical investigation of intermolecular interactions, particularly those involving hydrogen bonds, is of fundamental importance to chemistry, biology, and solid state and surface physics. Because of its importance in many chemical and biological processes, the study of the properties of water in all of its phases has been an area of intensely active research for the past decade. The earlier theoretical studies of ice (Weissmann and Cohan, 1965) were confined to a model of the lattice in which there were only two or three interacting molecules. Often, only the electrons involved in the hydrogen bond were included in the calculations. However, since 1972, several authors (Crowe and Santry, 1973; Santry, 1973; Cotterill and others, 1973; Scott, 1975; Plummer and Stein, 1976) have shown the power of semi-empirical molecular-orbital calculations for studying larger systems of ice and ice-like lattices. These calculations concentrated on determinations of structural parameters, molecular electronic charge distributions, and stabilization energy for idealized infinite ice lattices or small ice-like aggregates. It was found that these semi-empirical quantum-mechanical calculations were in surprisingly close agreement with those properties which could be compared with experimental values. Approximate molecular-orbital calculations of the CNDO/INDO type (Pople and Beveridge, 1970) give good agreement (Plummer, 1972 and unpublished results) with the results of ab initio calculations of a small basis set of small clusters of water molecules, and even predict the non-additive effects of hydrogen bonding found near the Hartree-Fock limit of the extended basis calculation (Hankins and others, 1970), thus giving confidence that the use of these techniques is justified and does yield physically meaningful results.

The nature of the surface of ice itself has long been a subject of controversy (Fletcher, 1970, 1973; Jellinek, 1972; Nason and Fletcher, 1975). The only semi-quantitative theoretical treatment of the ice surface has been developed by Fletcher (1970, 1973) but this addressed the question of the presence or absence of a liquid-like layer on the surface and not that of the stability of the different planes of ice nor the energetics of the surface-monomer interaction. Since the interaction of the surface with an incoming or adsorbed water monomer is of major importance in understanding the processes of nucleation and crystal growth, this report will concentrate on these aspects of the ice problem. The monomer-surface interaction is modelled by using a small section of an ice surface and then applying semi-empirical quantum-mechanical calculations of the CNDO/INDO type to investigate the energetics of a gas-phase molecule approaching the surface. Preliminary investigations of the water-ice surface using a ten-molecule cluster to represent the ice surface have been reported (Plummer and Stein, 1976) and compared with the results from an empirical pair-potential model (Plummer and others, 1976). The latter model is discussed in more detail in Kiefer and Hale (1977). Additional comparisons between the quantum-mechanical calculations and the pair potential will be made in a future report (Plummer, in press).

In an attempt to simplify the system as much as possible while retaining the essential features of the "real" physical situation, an idealized section of the ice surface was represented by thirteen water molecules arranged in the lattice positions of ice Ih. Questions into the effects of surface relaxation will be dealt with in subsequent work in which molecular dynamics will be used. The problem of the nature of "real" ice surface is not specifically addressed, but these studies together with the molecular-dynamics studies can be very instructive and provide valuable data on a molecular level.

Potential-energy curves were calculated for a monomer approaching different sites on both the basal and prism planes. Three-dimensional electron-density plots of the ice surfaces were generated. The specific geometries and surface structures considered are discussed. Then, details of the calculations are discussed together with a brief description of the CNDO technique, and we discuss the results of these calculations. Conclusions and implications for future work are included in the last section.

When a small group of molecules are used to represent the ice surface, we adopt the assumptions that the oxygen sublattice of ice Ih extends to the surface without re-arrangement, this appears to be experimentally justified for sufficiently low temperatures, (Firment and Somorjai, 1975),* and that a molecule approaching the surface is influenced primarily by

first and second nearest neighbors of the approach site. This latter point is amplified and justified later.

The Model

For the calculation of the interaction energies of a water monomer with the basal plane, a group of thirteen water molecules was arranged on a two-level ice Ih lattice as illustrated in Figure I. The positions of the oxygens are indicated by large circles while the hydrogen positions for a specific basal surface are shown as small filled circles. The c-axis is normal to the plane and the oxygens represented by filled circles lie 0.9 Å (1 Å = 10-10 m) below the plane. The oxygen-oxygen distances were fixed at 2.76 Å. The intramolecular O-H separation was 1.0 Å, slightly greater than in the vapor (0.957 Å), and no distinction was made in O-H separation between "free" hydrogens and those participating in hydrogen bonds. The nearest-neighbor distance between molecules in the same plane (dashed lines) was 4.5 Å, the lattice constant a 0. All H-O-H angles were opened up to 109.50°, the tetrahedral angle. The area of this section of the surface is approximately 80 ÅZ. For the basal face, four distinct proton configurations have been examined: (I) a configuration in which all the molecules in the first layer (open circles) had protons pointing out of the plane parallel to the c-axis, (2) a configuration in which no protons pointed out of the plane, (3) a configuration in which the protons on the central site pointed away from the c-axis while the remaining protons alternated up and down (referred to as BA) and (4) a second alternating configuration of surface protons in which the remaining protons were arranged so as to give the maximum surface stability (Plummer, 1972 and unpublished results). The calculations reported here were for the BA surface as it appeared to be the best representation of an "average" basal surface of ice. Additionally, entropy considerations would make a geometry like BA, with more proton disorder than the optimal surface, physically more probable. It has also been suggested (Pruppacher and Pflaum, 1975) that such a surface would be more conducive to the nucleation of new planes of ice than would a more polar surface.

Fig. 1. Section of the basal surface of ice Ih showing the proton configurations for an alternating (BA) surface. The primary bonding site is indicated by 1. The dashed lines connect molecules in the same plane and illustrate the hexagonal symmetry of the lattice.

For the basal plane, the interaction energies for the monomer approaching four specific types of surface sites were calculated. The site labelled 1 in Figure I was considered as the primary bonding position and the potential surface was then calculated for a variety of rotational orientations of the approaching monomers. The remaining three "sites" were not epitaxial positions but were of interest in "mapping out" the energy contours of the surface. These non-epitaxial positions should play a role in surface diffusion and in any liquid-like layer on the surface. One such position was in the center of one of the six-membered rings. A second position was directly above a second-layer oxygen. In these two cases as well as for the primary site, the incoming monomer was oriented so as to donate one of its hydrogens for bond formation with the surface. Thus the dipole of the approaching monomer was oriented approximately 52° from the surface. The fourth surface site was located half-way between site 1 and a first-layer nearest-neighbor molecule. The monomer approaching this site was oriented so both hydrogens pointed toward the surface.

For evaluating the interaction energy of a monomer with a prism face of ice Ih, the molecules represented the surface were arranged as shown in Figure 2. The c-axis is in the plane of the molecules and the plane of the oxygens represented by the closed circles lies 1.38 Å below those represented by open circles. Only one configuration of surface protons was considered for the prism face. Again the primary prism site is designated as 1 (Fig. 2).

Fig. 2. Section of the prism surface of ice Ih. The primary bonding site is indicated by 1.

The only other prism site examined had the monomer approaching with both hydrogens available for bonding to the surface. It should be noted that, since the hydrogen bonds which extend from the surface are not perpendicular to the surface for the prism face, nearest-neighbor surface molecules are less suitable for bifurcated bonding with an approaching monomer.

In addition to the several neutral basal and prism surface geometries, ion-pair defects were introduced into the surface and the interaction energies of these determined with an approaching monomer. The ion-pair defects were produced by moving a proton from one end of the hydrogen bond to the other (see Fig. 3). The central site was always one molecule of the ion pair. In one such defect pair, the surface molecule (first layer) had an H3O+ configuration while its companion in the second layer was OH-. In a second configuration the surface molecule was negative OH- while the H3O+ was in the second layer. A third ion-pair configuration has the ion pair separated by a neutral water and both ions were in the first layer of molecules (see Fig. 4). Since such defects are present in concentrations of 1 in 106 (Hobbs, 1974, chapters 2 and 6) their effects on the interaction should be studied. If the effects are sufficiently small, the ion defects could be neglected in any statistical study of the ice surface- vapor- surface cluster system

Fig. 3. Sections of the basal surface illustrating two configurations of ion pair defects.

Fig. 4. Basal surface illustrating separated ion defects.

The Calculations

For these studies of the monomer-ice surface interaction all valence electrons, both of the monomer and of the molecules used to represent the surface, were included. The atomic orbital basis for these calculations was of the conventional Slater valence-set type, that is, is orbitals on all the hydrogens and 2S, 2p x , 2p y orbitals centered on each oxygen. The problem then consisted of calculating the expectation value of the energy operator for a variety of monomer positions using product wave functions describing the 42 nuclei and 112 electrons. In order to make the problem tractable (and not too expensive) a number of the multicenter integrals were neglected (neglect of differential overlap) and the two electron repulsion integrals were parameterized from experimental atomic-spectra data. The specific technique used has been widely applied to molecular systems and systems of molecular aggregates. These approximate molecular-orbital theories are known by the acronyms CNDO, INDO, or NDDO depending on whether the neglect is "complete", "intermediate", or limited to "neglect of diatomic differential overlap" (Pople and Beveridge, 1970).

We have used both the CNDO and INDO methods and have found the qualitative predictions to be the same from both. The results reported here are CNDO/2 calculations (the 2 identifies the specific parametrization used for the integrals) which could be more readily compared with and calibrated against previous results for water systems (Hoyland and Kier, 1969; Dobosh, 1974; Scott, 1975).

The calculations were performed on an IBM 370/168 computer and required approximately 1.5 minutes of CPU time per monomer-surface geometry. The calculation of the electron 'density for a specific surface configuration required approximately 1.2 min of CPU time, with an additional 0.15 min to generate a three-dimensional plot.

The stability of a given surface geometry was determined by subtracting 13 times the energy of an individual surface monomer from the energy for the 13-molecule aggregate. The binding energy of the monomer E b with the surface was determined by subtracting the energy of the surface E s plus the energy of the monomer E n , from the total energy for the 14-molecule aggregate E t

()

The net atomic charge on any atom I a was equal to the core charge on the atom ?a minus the gross electron population on that atom P a

()

For all of the calculations, the geometry of the incoming monomer was fixed at an O-H separation of 1.0 Å and a H-O-H angle of 105°, and these were not allowed to change upon bond formation. Additionally, no relaxation of the surface was allowed, either for a change in surface geometry or upon interaction with the approaching monomer.

For each basal or prism site the interaction energies were calculated for 12 to 20 monomer surface distances surrounding the minimum separation. A least-squares fit to the resulting potential curves allowed the calculation of vibrational frequencies for a monomer-surface stretch. For the principal bonding sites on each face, the energy as a function of rotation about the monomer-surface bond was calculated and barriers to rotation estimated.

Results and Discussion

These calculations have used a small cluster of 13 water molecules to model the surface of ice Ih. The use of clusters of atoms to model the solid substrate in chemisorption investigations (Johnson and Messmer, 1974 ; Messmer, 1977) produced results in reasonable agreement with experiment. Since we were attempting to model the interaction of a water monomer with an infinite ice surface, we needed to determine whether the properties sought could be attributed to localized properties of the extended solid and whether the number of molecules used to represent the surface was sufficient to reproduce the short-range effects of the monomer-surface interaction.

It has been found in simulations of ice (Santry, 1973) that intermolecular electron exchange and delocalization contributions are quite short range, dominated by nearest-neighbor interactions and convergent by second-nearest-neighbor inclusion. We found the same trends for the monomer-surface interactions. For the basal plane, 94% of the interaction energy could be attributed to nearest-neighbor interactions with an additional 5.4% from second-nearest-neighbor interactions. For the prism face, convergence was not quite so rapid, and third-nearest-neighbor interactions were needed to provide greater than 99% of the interaction energy. With the cluster employed here the monomer interaction included fourthnearest neighbors for the principal bonding sites on each face and third-nearest neighbors elsewhere. Thus, it appears reasonable to assume that the energetics of the molecule close to the surface could be modelled using this small section of the surface. The principal limitation imposed by the use of the small cluster is that it cannot be expected faithfully to mimic relaxation of the extended crystal when defects are incorporated. Thus no relaxation was allowed in these structures, and our attention was directed to the perfect ice surface.

Monomer-basal-surface interactions

To visualize more easily how the surface would appear to an incoming water molecule, the electron densities at a number of points on the surface were calculated. Then the density in any plane could be displayed as a three-dimensional plot in which the relative magnitude of the density is represented by the height of the peaks in the plot. One such plot for the basal (BA) face is shown in Figure 5. This plot illustrates the electron density perpendicular to the c-axis in the plane of the upper layer of molecules (refer to Fig. I). The area of the plot is 180 Å2 and the distance between grid points is 0.25 Å. The viewer is oriented 30° above the plane and 45° clockwise from the orientation depicted in Figure 1. The slight asymmetries in the plot are due to the disorder in the proton configurations. As expected, the lone pair of the oxygen atoms produced the peak densities, and these together with the effect of molecular dipoles served to orient the incoming monomer. However, the magnitude of these effects was quite small beyond 4 or 5 Å. To illustrate the result that nearest-neigh b or interactions dominate and how short-range many of the surface features are, Figure 6 shows the electrondensity plot I Å above the surface.

Fig. 5. Electron-density map of the basal (BA) surface of ice Ih. The magnitude of the density in the plane of the upper level of molecules is represented by the height of the peaks.

Fig. 6. Electron-density map 1 Å above the basal (BA) ice Ih surface.

The change in electron density with distance above the plane has been examined more closely in the vicinity of the approach sites. In Figure 7 the density is plotted against distance above the plane for points along the line connecting the primary bonding site and the site for bifurcated bonding. The center peak is located at the primary site; bifurcated bonding occurs midway between the center and right-most peaks, with the two hydrogens in the plane and pointing toward the peaks. Figure 8 illustrates the electron density for an approach trajectoryover a second-layer molecule. The approach site is at the center peak; the other points lie on a line parallel to that of Figure 7. The scale of peak heights has been increased in Figure 8 to show more detail.

The most stable configuration of basal (BA) and monomer had the monomer 2.55 Å above the primary site with the non-bonded monomer hydrogen oriented 180° from the surface molecule H-O-H bisector. The potential curve for this interaction is shown in Figure 9 and has an energy minimum of -33.3 kJ/mol (34.7 x 10-2 eV). The barrier to rotation about the monomer- surface bond was 2. 1 kJ/mol (2.17 x 10-2 eV) which is approximately kT. The stretching frequency for the bond, calculated via a least-squares fit to the potential curve, was 280 cm-1.

Fig. 7. Electron-density map perpendicular to the basal (BA) surface along the trajectories approaching the primary and bifurcated sites.

Fig. 8. Electron-density map for the approach to the second layer bonding site. The axes are parallel to those in Figure 7.

The second most stable monomer-surface geometry was the one in which the monomer was in a bifurcated position and had both hydrogens participating in bonds to the surface. The distance of the monomer oxygen from the surface at the energy minimum was 1.7 Å and the energy was -28.4 kJ/mol (0.30 eV). The potential curve for the monomer approaching this site together with that for the monomer over a second-layer molecule are shown in Figure 10. Figure 11 shows the potential barrier for rotation between the primary site with one hydrogen bond and the bifurcated site with two monomer-surface bonds. The barrier height was 11.3 kJ/mol (0.117 eV) above the primary site. Other orientations of the monomer dipole above the bifurcated site were investigated but all had weaker bonding to the surface.

The remaining two basal sites for which the monomer-surface interaction was calculated were not true bonding configurations since no monomer-surface hydrogen bond was formed. However these are both local attractive minima in the potential surface. The attraction interaction is approximately one-half that of the primary bonding site. The magnitude and exact location of these local minima are highly dependent on the orientation of both the surface protons and those of the incoming water monomer.

Fig. 9. Potential energy versus monomer-surface distance for a monomer approaching the primary bonding site on the basal surface.

Fig. 10. Potential energy versus monomer-surface distance for a monomer approaching (A) the birfurcated site and (B) the "second" layer site on the basal surface.

Fig. 11. Potential energy barrier for transition between the primary and bifurcated basal sites.

Monomer-prism-surface interactions

For the prism face, electron-density plots of the surface were again generated in order to assess in a better way the effects of the surface on an incoming monomer. The plot obtained for the prism plane is shown in Figure 12 where the density is plotted parallel to the c-axis and through the upper layer of oxygen atoms (see Fig. 2). The magnitude of the interaction of the monomer with the primary site is found to be slightly less than for the basal face (-33 kJ/mol) although the difference is at the limit of the accuracy of the calculations. However, bonding to the bifurcated site was substantially less for the prism face amounting to an energy of 16.3 kJ/mol. The barrier for rotation about the monomer surface bond at the primary site was about the same as for the basal face: 2.09 kJ/mol. The vibrational frequency calculated for elongation of this bond (266 em-1) is slightly less than for the basal plane.

Fig. 12. Electron-density map for prism face of ice Ih.

Monomer-ion-pair-defect interactions

Ordinarily the Bernal-Fowler rules apply in ice, that is, each water molecule is surrounded by two near hydrogens and two hydrogens at approximately twice the distance of the near hydrogens. These conditions guarantee the neutrality of each molecule as well as assuring the presence of one and only one hydrogen per bond. However, if this were always the case it would be hard to explain the observed electrical properties of ice. Bjerrum (1951) suggested that certain types of point defects, which violated the Bernal-Fowler rules, could be present and that their migration could provide an explanation for many of the electrical properties of ice (Camplin and Glen, 1973; Glen and Paren, 1975). The types of point defects suggested were ion-pair defects and orientational defects. In the former, a proton is transferred from one water molecule to its neighbor producing an H30+ and an OH- pair. If the negative ion then captures a proton from an adjacent molecule before recombination occurs, the ionic defect has been somewhat stabilized and we have the first step in charge migration. In the case of the orientational defects, rotation of a single molecule can put two hydrogens on the same bond (a D-defect) and simultaneously create a bond having no hydrogens (an L-defect). In this report we will discuss only the ionic defects. The interactions of the monomer with a surface containing orientational defects as well as isolated ions will be reported later (Plummer, in press).

As described earlier, a study of ion-pair defects on the basal surface was performed. An upper bound for the energy required to produce the defects is given by

()

where E d is the defect energy, E s n is the energy of a neutral surface, and E s i is the energy associated with an ionic surface. This number will be an over-estimate, since it includes only the energy required for ionization and charge separation but no compensating relaxation of the surface. The values obtained ranged from 600 to 850 kJ/mol (6.3 to 8.8 eV), about the same magnitude as the proton affinity of water, 670 kJ/mol (7 eV).

Fig. 13. Potential energy versus monomer-surface distance for a monomer approaching a positive surface site top) and a negative surface site (bottom) of basal ice Ih having ion-pair deject sites.

The primary basal site for monomer interaction was always one of the ion pairs and was the only site of monomer approach included in these studies. In the first case, a proton was moved from a second-layer molecule producing an H30+ ion at the surface having all protons pointing into the surface. This orientation was not conducive to hydrogen-bond formation with an incoming monomer but nevertheless exhibited an attractive interaction in excess of 20.9 kJ/mol. When the proton was shifted from the primary site to a second-layer molecule producing an OH- primary site, a very strong bond with the incoming monomer resulted. The potential energy as a function of distance for these interactions is illustrated in Figure 13. The negative primary site also exhibited a higher barrier to rotation about the monomer-ion bonds than found with the neutral site.

A third ion-pair configuration was produced by shifting the proton from the second-layer H30+ to an adjacent upper-layer molecule (see Fig. 4). This had the effect of separating the ion pair and having both ion defects in the upper layer of molecules. The energy required to produce the separated defect was 210 kJ/mol (2.2 eV). When the ions were separated in this fashion, the strength of the interaction with an incoming monomer was essentially the same as for a monomer with a "free" OH- ion and greater than in the second type of configuration.

Conclusions

In this study we have modelled the interaction of a water monomer with an extended, perfect ice surface represented by a cluster of 13 water molecules whose oxygens were arranged at the lattice sites for hexagonal ice. Summaries of the results appear in Tables I and II. We have shown that the interaction energies are dominated by nearest-neighbor effects and hence this model appears quite reasonable.

The attraction energies were calculated for a variety of possible attachment sites on both prism and basal faces. It was found that epitaxial sites, those which allowed a hydrogen bond to be formed between the monomer and the surface, produced the largest interaction energies. These energies were approximately the same (-33.4 kJ/mol) for the two ice crystal surfaces.

Table I. Summary of Results for Basal Surface-Energy Calculations

Table II. Summary of Results for Energy Calculation for Defect Surfaces

The second-most-stable configuration has two "bonds" formed between the monomer and the basal surface. This bifurcated structure was separated from the most-stable configuration by an energy barrier of 6.3 kJ/mol. This orientation of the monomer is probably most frequently observed as a metastable or intermediate structure for a molecule in transition between two of the epitaxial bonding sites. This position can be regarded as a local minimum in the potential for surface diffusion. Its existence and the moderate barrier between it and a more stable epitaxial structure, makes reasonable the consideration of a "bipedal walk" process for the surface diffusion of water on the basal surface of ice. The situation is quite different on the prism face where a similar bifurcated geometry was found to be less stable than "nonbonding" configurations on the basal ice surface.

The other two approach sites considered for the basal surface were not true bonding configurations, but they proved to be attractive sites with interaction energies in the range of - (17-21) kJ/mol. The magnitude of this attraction will be highly dependent on the orientation of both the monomer and the surface hydrogens and will be investigated more extensively in future work. Since these are local minima in the surface potential, they can attract and physically adsorb a monomer if its incoming energy is not too great. The monomer could be quasi-bound for some time before being re-evaporated or diffusing across the surface to be incorporated into an epitaxial position. These sites can also provide added stabilization for the addition of a second molecule to a monomer already bound to a primary site thereby enhancing the initial step in cluster formation.

From these calculations we can estimate a "static" diffusion rate and residence time for a monomer on the basal surface. For these evaluations, we use the expressions of Burton and others (1951).

()

and

()

where χs is the diffusion time, α is the distance between nearest neighbors in the same plane, v is the stretching frequency, v' is the diffusion frequency, E A and E D are the average binding and diffusion energies, and τs is the residence time. For the average binding energy E A we used that associated with a primary site, 32 kJ/mol. The diffusion barrier was taken as to kJ/mol and v' and v were 450 cm-1 and 28o cm-1 respectively. For a temperature of 260 K, we find the mean path length to be approximately 270α, where α is 4.5 Å, the distance between nearest neighbors in the same plane. The residence time is estimated as 7 x 10-7 s. It should be noted that these are estimates of dynamic quantities from static calculations. Additionally, they should be regarded as an upper limit to these quantities as we have taken the bottom of the potential well as our value for E A so that zero-point vibrational energies have not been excluded. Probably vibrational excitation should be included as well. However, these refinements did not seem justified in the light of the other approximations inherent in Equations (4) and (5). The estimate of χs = 270α may be compared with the estimate of Burton and others of χs = 400α for the (111) face of an f c.c. crystal.

We also examined the effect on the monomer-surface interaction of incorporating ion-pair defects in the surface. As one might expect, we found the interaction energies to be greater, ranging from three to five times as strong as with the neutral site. The potential wells were also steeper, yielding higher stretching frequencies. The energies required to produce the defects compare well with the experimental proton affinity of water and also with calculated barriers to proton transfer (Weissmann and Cohan, 1965). Future studies will address the question of lattice relaxation and defect mobility. It would be expected that the presence of ion-defect pairs or isolated ions in the surface would significantly alter diffusion rates near those sites. Indeed, because of the strength of the interaction, the results of these calculations predict no diffusion from the ionic sites.

In summary, I emphasize the following points: First, while the strengths of the bondings between the basal and prism planes are roughly the same, the studies reported here suggest that the diffusion rates and residence times may be substantially different for the two crystal faces. Therefore, factors such as supersaturation and flux of monomers to the surface, as well as the temperature of the surface and vapor, will play a large role in determining the ice-crystal growth habit. Such factors should be studied in a dynamical fashion and such calculations are under way. Secondly, the ice surface should be viewed not as having isolated strong-bonding regions separated by repulsive regions but rather as consisting of attractive regions separated by low energy barriers. The non-bonding local minima provide the energy mechanism whereby a molecule approaching a surface does not have to attach initially to a primary bonding site but can be attracted and remain in contact with the surface long enough to either move to a more favorable bonding position on the surface or become attached to a molecule or cluster already bonded to the surface. This, in turn, again suggests that surface diffusion would be a primary mechanism for surface cluster formation and growth. Third, the presence of these non-bonding sites together with the essentially free rotation about monomer-surface bonds provide motions not present in the crystal, and hence may be the molecular origin of a liquid-like layer on the ice surface.

While these calculations are only a model for the real interaction, they represent a first attempt to model this interaction from a molecular point of view using quantum-mechanical calculations. The results reported do provide considerable information about the monomer - surface interaction and within the limits of the model allow us to draw preliminary conclusions about processes on the surface of ice.

Acknowledgements

The author would like to acknowledge the help of E. Stein and K. Whisnant for their computational assistance and work with the ten-molecule surfaces. She wishes to thank the Computer Facilities at the University of Missouri for much of the computer time needed for these studies. She also thanks Dr James L. Kassner, Jr for his continued encouragement of this work and Dr 0. R. Plummer, for many useful discussions and suggestions.

Discussion

W. B. KAMB: What you call a "bifurcated site" does not involve bifurcated hydrogen bonds as that term is conventionally used; you could perhaps better call it a "doubly-bonded site".

Also, the particular interaction between the H30+ ion and monomer that you described is one where no H-bond is formed. If this is definitely more stable than the H-bonded configuration, it might suggest that, in the ice crystal, H30+ ions and L-defects should be permanently associated.

P. L. M. PLUMMER: Perhaps it would be better to refer to the single-proton-down configuration as a "single donor" and the two-proton-down configuration as a "double donor" configuration, then an oxygen-down orientation could be a "single" or "double" acceptor depending on its attachment to the surface.

I cannot address your second question as to whether the H30+ and L-defects should be permanently associated in the bulk, but for the surface-monomer case, where the monomer is free to rotate so as to optimize the interaction when all the protons on H3O+ are directed into the solid, the L-defect configuration is more stable. It should be recalled, however, that the all-proton-down orientation is only one of four for H30+, the others are orientated for hydrogen-bond formation with an incoming monomer which would, of course, be a much stronger interaction.

J. W. GLEN: Whilst agreeing with Dr Kamb's second remark, I disagree with the first as I believe the "bifurcated" approach was to two molecules both in the top layer and thus tending to a four-oxygen ring. As the apparent failure of the monomer to give hydrogen bonding in one orientation of the H30+ ion on the surface is potentially important, as Kamb has pointed out, I am interested to know whether you tried OH- ions in both orientations (proton out of the surface as well as proton in the surface).

PLUMMER: No, I have not as yet, but that calculation will be done. However, for the OH-, either orientation is favourable for hydrogen-bond formation with the monomer.

J. E. BERTIE: Your calculations will certainly be used by those interested in the formation and growth of ice crystals from the vapour. It is particularly important, therefore, that an estimate of the effect of the approximations in your treatment be given. I note that you used the CNDO and INDO approximations to the ab initio molecular-orbital treatment. These are very approximate in principle and are nowadays only used in cases which cannot be treated by the full molecular-orbital method. Since the energies that you calculated are small (≲40 kJ/mol), it is particularly necessary to give an estimate of the effect of your approximations and their values. Can you give this estimate please?

PLUMMER: It is certainly true that the CNDO/INDO techniques are semi-empirical in nature and must be used with great care. The results are much more reliable for comparing different structures within the same types of systems where you are looking at energy differences, and are not nearly as good for either absolute total energies or for comparisons between systems (say water systems against organics). It is my feeling that the values for bonding energies are probably good to 10 or 15% and that differences greater than 4 kJ/mol are real. Those less than 1.2 kJ/mol may well not be. We have compared the CNDO with other semi-empirical methods, empirical potentials, and with ab initio calculations, when possible, and these are the ranges of accuracy which seem reasonable from these comparisons. The separations are surely about to% too small.

Footnotes

*

Supported in part by Atmospheric Research Section, National Science Foundation and the National Aeronautics and Space Administration.

* The model is probably valid at substantially higher temperatures than those studied experimentaly ≈ 155 K). This is being investigated using molecular dynamics simulations between temperatures of 160 and 265 K.

References

Bjerrum, N. 1951. Structure and properties of ice. Kongelige Danske Videnskabernes Selskab. Matematiskbsiske Medddelser, Bd. 27, Nr. I.Google Scholar
Burton, W. K., and others. 1951. The growth of crystals and the equilibrium structure of their surfaces, by W. K. Burton, N. Cabrera and F. C. Frank. Philosophical Transactions of the Royal Society of London, Ser. A, Vol. 243, No. 866, p. 299-358.Google Scholar
Camplin, G. C., and Glen, J. W. 1973. The dielectric properties of HF-doped single crystals of ice. (In Whalley, E., and others, ed. Physics and chemistry of ice: papers presented at the Symposium on the Physics and Chemistry of Ice, held in Ottawa, Canada, 14-18 August 072. Edited by E. Whalley, S. J. Jones, L. W. Gold. Ottawa, Royal Society of Canada, p. 256-61.)Google Scholar
Cotterill, R. M. J., and others. 1973. Computer studies of perfect-crystal properties and defect structures in ice I, [by] R. M. J. Cotterill, J. W. Martin, 0. V. Nielson and 0. B. Pedersen. (In Whalley, E., and others, ed. Physics and chemistry of ice : papers presented at the Symposium on the Physics and Chemistry of Ice, held in Ottawa, Canada, 14-18 August 1972. Edited by E. Whalley, S. J. Jones, L. W. Gold. Ottawa, Royal Society of Canada, p. 23-27.)Google Scholar
Crowe, R. W., and Santry, D. P. 1973. Molecular orbital theory for infinite molecular aggregates: application to idealozed hexagonal and cubic ices. Chemical Physics, Vol. 2, No. 3, p. 304-20.Google Scholar
Dobosh, P. A. 1974. Quantum Chemistry Program Exchange. Program No. 141. [Available from QCPE, Indiana University, Bloomington, Indiana 47401.]Google Scholar
Firment, L. E., and Sornorjai, G. A. 1975. Low-energy electron diffraction studies of molecular crystals. The surface structure of vapor-grown ice and naphthalene. Journal of Chemical Physics, Vol. 63, No. 2, p. 1037-38. Fletcher, N. H. 1970. The chemical physics of ice. Cambridge, University Press. (Cambridge Monographs on Physics.)Google Scholar
Fletcher, N. H. 1973. The surface of ice. (In Whalley, E., and others, ed. Physics and chemistry of ice: papers presented at the Symposium on the Physics and Chemistry of Ice, held in Ottawa, Canada, 14-18 August 1972. Edited by E. Whalley, S. J. Jones, L. W. Gold. Ottawa, Royal Society of Canada, p. 132-36.)Google Scholar
Glen, J. W., and Paren, J. G. 1975. The electrical properties of snow and ice. Journal of Glaciology, Vol. 15, No. 73, p. 15-38.Google Scholar
Hankins, D., and others. 1970. Water molecule interactions, [by] D. Hankins and J. W. Moskowitz and F. H. Stillinger. Journal of Chemical Physics, Vol. 53, No. 12, P. 4544-54.Google Scholar
Hobbs, P. V. 1974. Ice physics. Oxford, Clarendon Press.Google Scholar
Hoyland, J. R., and Kier, L. B. 1969. Semiempirical quantum mechanical studies on some hydrogen bonded systems. Theoretica Chimica Acta, Vol. 15, Fasc. 1, p. 1-11.Google Scholar
Jellinek, H. H. 1972. The ice interface. (In Home, R. A., ed. Water and aqueous solutions; structure, thermodynamics, and transport processes. New York, etc., Wiley-Interscience, p. 65-107.)Google Scholar
Johnson, K. H., and Messmer, R. P. 1974. Clusters, chemisorption, and catalysis. Journal of Vacuum Science and Technology, Vol. 11, No. 1, p. 236-42.Google Scholar
Kiefer, J., and Hale, B. N. 1977. The water monomer on the basal plane of ice Is: an effective pair, central force potential model of the static interaction. Journal of Chemical Physics, Vol. 67, No. 7, p. 3206-13.Google Scholar
Messmer, R. P. 1977. Molecular cluster approach to some solid-state problems. (In Segal, C. A., ed. Modern theoretical chemistry. Vol. 8. New York, Plenum Press, p. 215-46.)Google Scholar
Nason, D., and Fletcher, N. H. 1975. Photoemission from ice and water surfaces. Quasiliquid layer effect. Journal of Chemical Physics, Vol. 62, No. 11, p. 4444-49.Google Scholar
Plummer, P. L. M. 1972. Cooperative effects of bonding in small water clusters. Bulletin of the American Physical Society, Ser. 2. Vol. 17, No. 5, p. 696. [Abstract.]Google Scholar
Plummer, P L M. In press. A molecular orbital study of water-ice interactions. II. Effects of proton order and orientational defects. Journal of Chemical Physics. Google Scholar
Plummer, P. L. M., and Stein, E. M. 1976. Energetics of water-ice surface interactions: a quantum mechanical study. (In International Cloud Physics Conference, July 26-30, 1976, Boulder, Colorado, U.S.A. Boston, American Meteorological Society, p. 130-35.)Google Scholar
Plummer, P. L. M., and others. 1976. Models for water and ice clusters on simple surfaces, by P. L. M. Plummer, B. N. Hale, J. Kiefer, E. M. Stein. (In Kerker, M., ed. Colloid and interface science. Vol. 2. Aerosols, emulsions, and surfactants. Proceedings of the International Conference on Colloids and Surfaces-50th Colloid and Surface Science Symposium, held in San Juan, Puerto Rico, on June 21-25, 1976. New York, etc., Academic Press Inc., P. 45-62.)Google Scholar
Pople, J. A., and Beveridge, D. L. 1970. Approximate molecular orbital theory. New York, etc., McGraw-Hill Book Co., Inc. (McGraw-Hill Series in Advanced Chemistry.)Google Scholar
Pruppacher, H. R., and Pflaum, J. C. 1975. Some characteristics of ice nucleation active sites derived from experiments with a ferroelectric substrate. Journal of Colloid and Interface Science, Vol. 52, No. 3, p. 543-52.Google Scholar
Santry, D. P. 1973. Molecular-orbital studies on ice Ih. (In Whalley, E., and others, ed. Physics and chemistry of ice: papers presented at the Symposium on the Physics and Chemistry of Ice, held in Ottawa, Canada, 14-18 August 1972. Edited by E. Whalley, S. J. Jones, L. W. Gold. Ottawa, Royal Society of Canada, p. 19-22.)Google Scholar
Scott, B. F. 1975. Molecular orbital studies on ice-like structures. Theoretica Chimica Acta, Vol. 38, Fasc. 2, p. 85-100.Google Scholar
Weissmann, M., and Cohan, N. V. 1965. Molecular orbital study of the hydrogen bond in ice. Journal of Chemical Physics, Vol. 43, No. 1, p. 119-23.Google Scholar
Figure 0

Fig. 1. Section of the basal surface of ice Ih showing the proton configurations for an alternating (BA) surface. The primary bonding site is indicated by 1. The dashed lines connect molecules in the same plane and illustrate the hexagonal symmetry of the lattice.

Figure 1

Fig. 2. Section of the prism surface of ice Ih. The primary bonding site is indicated by 1.

Figure 2

Fig. 3. Sections of the basal surface illustrating two configurations of ion pair defects.

Figure 3

Fig. 4. Basal surface illustrating separated ion defects.

Figure 4

Fig. 5. Electron-density map of the basal (BA) surface of ice Ih. The magnitude of the density in the plane of the upper level of molecules is represented by the height of the peaks.

Figure 5

Fig. 6. Electron-density map 1 Å above the basal (BA) ice Ih surface.

Figure 6

Fig. 7. Electron-density map perpendicular to the basal (BA) surface along the trajectories approaching the primary and bifurcated sites.

Figure 7

Fig. 8. Electron-density map for the approach to the second layer bonding site. The axes are parallel to those in Figure 7.

Figure 8

Fig. 9. Potential energy versus monomer-surface distance for a monomer approaching the primary bonding site on the basal surface.

Figure 9

Fig. 10. Potential energy versus monomer-surface distance for a monomer approaching (A) the birfurcated site and (B) the "second" layer site on the basal surface.

Figure 10

Fig. 11. Potential energy barrier for transition between the primary and bifurcated basal sites.

Figure 11

Fig. 12. Electron-density map for prism face of ice Ih.

Figure 12

Fig. 13. Potential energy versus monomer-surface distance for a monomer approaching a positive surface site top) and a negative surface site (bottom) of basal ice Ih having ion-pair deject sites.

Figure 13

Table I. Summary of Results for Basal Surface-Energy Calculations

Figure 14

Table II. Summary of Results for Energy Calculation for Defect Surfaces