Exploring the Surface of Bioactive Glasses: Water ... - ACS Publications

Jul 16, 2008 - Surface Transformations of Bioglass 45S5 during Scaffold Synthesis for Bone ... Modeling the Water−Bioglass Interface by Ab Initio Mo...
0 downloads 0 Views 2MB Size
Download PDF
11936

J. Phys. Chem. C 2008, 112, 11936–11945

Exploring the Surface of Bioactive Glasses: Water Adsorption and Reactivity Antonio Tilocca*,† and Alastair N. Cormack‡ Department of Chemistry and Materials Simulation Laboratory, UniVersity College London, U.K., and New York State College of Ceramics, Alfred UniVersity ReceiVed: April 23, 2008; ReVised Manuscript ReceiVed: May 30, 2008

Car-Parrinello Molecular Dynamics simulations and structural optimizations have been carried out to investigate the surface of a highly bioactive phosphosilicate glass, containing Na and Ca modifier cations. The identification and characterization of the active surface sites was performed on the basis of their interaction with a water molecule. Spontaneous water dissociation was observed on the strongest surface sites, represented by three-coordinated Si atoms combined with a proton acceptor such as a nonbridging oxygen. Additional adsorption sites are represented by Na/Ca cation modifiers, which appear to be directly involved in the glass dissolution mechanism, by providing favorable pathways which allow water to penetrate under the surface. Small silica rings, 2- or 3-membered, are stable features of the surface, despite a favorable energetic balance for their opening upon water dissociation. The minimum energy paths and barriers for alternative 2 M ringopening mechanisms, obtained by means of String Method Car-Parrinello calculations, show that the opening of small rings is hindered by an energy barrier. Together with the availability of alternative, more favorable molecular adsorption sites where water can migrate before dissociation, this can result in some of the small rings not being open and hydroxylated upon immersion in a physiological environment, thus leaving them available to assist nucleation of dissolved Ca and P ions on the glass surface, as it had been previously proposed to interpret the initial stages of the bioactive fixation mechanism. 1. Introduction Silica-based glasses containing different amounts of Na, Ca and P ions have become a paradigm of second-generation surface-active biomaterials for medical applications. They combine high biocompatibility with the ability to induce the deposition of a layer of crystalline hydroxyapatite (HA) on their surface shortly after implant in vivo or in simulated body fluid (SBF); the HA layer, after incorporating collagen fibrils, produces a strong and stable bond between the glass and human hard and soft tissues, thus promoting fast integration of the implant with the host.1–5 These bioactiVe glasses also show high potential for third-generation tissue-engineering applications, where they can activate genes which stimulate healing of living tissues and cellular repair.6,7 The glass bioactivity, i.e., their ability to bond to bone and/or to induce tissue repair and regeneration, is usually assessed by measuring the rate of HA formation in vitro or in vivo;8–13 as this property is extremely sensitive to composition, particle morphology, surface texture, and thermal treatment of the glass, trial-and-error approaches are often employed to tailor the rate of HA formation toward a specific level, which might differ in different biomedical applications, or, analogously, to adjust the mechanical properties and the biodegradability of the glass toward a specified target.6,14 Although trial and error approaches, based on animals and clinical trials, are highly expensive and often ineffective, a more systematic and rational approach to bioactive glass optimization is hindered by the present lack of fundamental understanding of the effects leading to their different bioactive properties, which in turn depends on the lack of data concerning the atomistic structure of these materials. Following earlier attempts * Corresponding author: [email protected]. † University College London. ‡ Alfred University.

to identify structural features correlated to the bioactivity,15,16 only very recently have direct structural investigation on bioactive glasses, based on nuclear magnetic resonance (NMR)17 and Neutron/X-ray diffraction18 experiments, as well as Molecular Dynamics (MD) simulations19–25 started to shed a light on bulk features which can steer the bioactive behavior. For instance, our previous simulations suggested that the high bioactivity of the common 45S5 Bioglass reflects the high fragmentation of its phosphosilicate network: unlike less bioactive compositions, for 45S5 we recorded a high concentration of chain-like fragments, which can be released into the solution at a low energetic cost, thus accelerating the partial dissolution of the glass network which is required to start the bone-bonding mechanism.23,26 On the other hand, the tendence to form clusters and inhomogeneous Ca/P-rich regions in the bulk seems to mark the transition from bioactive to bioinactive glasses.19,25 The study of the structure of the bulk bioactive glasses, albeit crucial, represents only the starting point on the way to rationalize the bioactive properties of these materials; because the bone-bonding process starts and is eventually completed at the solid-liquid interface between the glass and the surrounding physiological environment, further studies need to focus directly on the properties of the glass surface. In fact, the first steps of the HA formation mechanism involve hydrolysis and partial dissolution of the silicate network, with release of modifier cations and silicate fragments in the physiological environment; this initial attack to the surface is crucial as it produces a hydroxylated SiO2-rich surface layer, which can then induce calcium phosphate deposition and HA crystallization.5,27,28 Surface analytic methods (scanning and transmission electron microscopy, vibrational and Auger spectroscopy), coupled with emission spectroscopy techniques to analyze the ions released in the physiological solution, can successfully be employed to probe the dynamical changes occurring at the glass-liquid

10.1021/jp803541j CCC: $40.75  2008 American Chemical Society Published on Web 07/16/2008

Exploring the Surface of Bioactive Glasses interface upon immersion.2,29 From a more fundamental perspective, a different task is the extraction of specific data on the structure of the surface of bioactive glasses before contact with the solution, as the active sites available on the dry surface will ultimately determine the rate at which the sequence of chemical processes involved in the bioactive fixation mechanism can proceed after immersion. In the present case, this kind of fundamental studies are highly challenging, due to the amorphous nature and multicomponent character of these biomaterials; the characterization of the surface sites can be performed effectively by measuring their strength of interaction with a gasphase polar molecular “probe”, using in situ experimental techniques.30,31 These methods provide an averaged picture of the adsorption sites, but a more detailed site-by-site characterization of the surface could be achieved using atomistic simulations, which in this work we apply to model the dry Bioglass surface and its interaction with an isolated water molecule: the ability of surface sites to adsorb and dissociate water is indeed a powerful indicator of their strength.32,33 In general, MD simulations represent one of the best computational approaches to obtain and investigate the bulk structure of modified silicate glasses: binary silicate glasses have successfully been modeled using classical34–38 and ab initio39–42 MD (AIMD) using an effective melt-quenching approach and periodic supercells, whereas simulations of multicomponent silicate glasses within this framework are less frequent. In fact, it is hard to incorporate all the complex interactions governing Na2O:CaO:SiO2:P2O5 amorphous phases into standard interatomic potentials, and our previous calculations highlighted the need to include ionic polarization, at least in the approximate form of a shell-model approach, in order to improve the representation of bulk structural properties of these materials in the medium-range.43 Classical MD simulations are still invaluable to tackle equilibrium processes involving relatively long time scales, such as diffusive dynamics, and large system sizes, such as the interface between the glass surface and bulk liquid water; however, standard empirical potentials cannot easily describe processes involving dynamical rearrangements to the electronic density, such as those brought about by bond breaking and forming upon water adsorption on the surfaces, and significant efforts are needed to develop advanced force fields which allow for water dissociation.44 Even though some structural insight can be achieved using classical MD,45,46 quantum-mechanical (QM) approaches are better suited to model the distorted coordination patterns and defects which characterize the surface of standard and bioactive silicate glasses,47–50 (as well as the bulk of the corresponding melts)24 and the local water reactivity.33 Among the various QM approaches, AIMD simulations using periodic slab models within Density Functional Theory (DFT), as in the Car-Parrinello (CP) method,51 represent an excellent combination of computational accuracy and efficiency, allowing us to probe the local structure and reactivity of specific surface sites, embedded in a realistic environment including long-range and many-body effects at an ab initio level.33,52 While AIMD has been often applied to model the dry53,54 and hydrated55–61 SiO2 surface, similar investigations of the surface of multicomponent silicate glasses are not as common,62 and, to our knowledge, no previous AIMD studies of the hydration of these surfaces are available. This again reflects the complexity of the problem: whereas the small number of active sites available on the surface of a crystalline solid can be mapped and studied in a relatively straightforward fashion,33 the surfaces of binary or ternary silicate glass expose a large number of potential adsorption sites, whose effective

J. Phys. Chem. C, Vol. 112, No. 31, 2008 11937 exploration using ab initio methods require significant computational resources. We have carried out extensive Car-Parrinello MD (CPMD) simulations and DFT structural optimizations of the 45S5 Bioglass surface; the hydrophilicity of several adsorption sites, including coordinative defects and small rings, has been probed by examining the dynamics and energetics of their interaction with a water molecule. The results provide the first ab initio description of this surface and of its hydration properties. The necessarily small size of the periodic system which can be modeled using ab initio techniques limits the scope of the present exploration to the interaction of an isolated water molecule with several adsorption sites identified on the glass surface. Even though the distribution of these sites in a small glass sample might not reflect their “true” relative proportion on a much larger surface area, the specific adsorption sites identified and investigated here, such as coordinative defects and small rings, are generally regarded as the main active features of these and other oxide surfaces,2,33,55,56 and therefore should provide an adequate representation of the bioglass surface. For the purpose of investigating the bioactivity of these sites, because the ab initio study of a realistic bioglass-SBF interface, including additional ions dissolved in a large sample of bulk liquid water, would require too large computational resources, a convenient alternative is modeling the local interaction of each site with an isolated water molecule. This represents an important step, providing the first details on the specific strength and reactivity of the various sites, which should be at least partially transferable to more complex environments and interfaces. Besides their direct impact in biomaterial science, these results are relevant also in the context of extending our atomistic understanding of the surface properties of silicate glasses, beyond pure silica compositions. For instance, new insight into the interaction of water with silicate glasses containing sodium and calcium modifiers can support a more fundamental understanding of the dissolution mechanism of these technologically important compositions. 2. Computational Methods The exact composition of the glass modeled in this work was 46.3SiO 2: 24.3Na 2 O: 26.8CaO: 2.4P 2 O 5 mol%, corresponding to the standard 45S5 Bioglass. Two different bulk samples, discussed and analyzed in refs22 and,24 were used to obtain different models of the glass surface: the first bulk structure had been generated by melt-and-quench classical MD followed by CP relaxation at room temperature,22 whereas the second one was obtained using a full ab initio scheme, where CPMD was employed in both the melt-and-quench phase and the final room temperature relaxation.24 Both bulk samples employed a periodic cubic cell with 11.63 Å sides, with the glass density at room temperature63 of 2.66 g cm -3. Details of the computational setup and procedures used to obtain and analyze the bulk glasses can be found in refs 22 and 24. It is important to remark that, even though the full-CP approach turned out to be essential to provide a very accurate picture of the melt phase and of the dynamical changes occurring during cooling, both final glass samples can be considered accurate: the polarizable force field used to obtain the classical sample had previously proven reliable to model modified silicate glasses,19,43 and the small structural difference between the two glasses are more likely to be related to the statistical error inherent to their small size. The CPMD simulations51 of the surfaces were carried out using the Car-Parrinello code of the Quantum-ESPRESSO package.64 The electronic structure was treated within the

11938 J. Phys. Chem. C, Vol. 112, No. 31, 2008 Generalized Gradient Approximation (GGA) to density functional theory (DFT), through the PBE exchange-correlation functional.65 Ultrasoft pseudopotentials were employed to represent core-valence electron interactions,66 including Na and Ca semicore shells. A plane-wave basis set was used, with cutoffs set to 30 and 200 Ry for the smooth part of the wave functions and the augmented charge, respectively; k-sampling was restricted to the Γ point. The MD time step δt and fictitious electronic mass µ were 7 and 700 au, respectively, using the deuterium mass for H atoms. The latter is a common way to increase the separation between electronic and ionic degrees of freedom in CPMD simulations, thus allowing for a longer time step. Similar computational approaches and methods have often been employed to model water adsorption on the surface of silica55,56,59 and of other oxides.33 A slab geometry was created from the bulk samples through the following steps: (i) the last configuration of the MD trajectory of the bulk glass was optimized with respect to the volume of the cubic cell, in order to obtain the optimal (theoretical) density at zero applied pressure; the latter turned out to be within 2% of the experimental density, thus confirming the good quality of the bulk samples. (ii) The optimized structure at the theoretical density was then cut to create a surface, by increasing the cell parameter along z by 14.5 Å, while keeping 3D periodic boundary conditions.67–69 Different unique surfaces, containing different numbers of undercoordinated Si defects, can be exposed by cutting the bulk unit cell at different z heights. The initial surfaces were then selected after checking that, after the cut, the top layers did not contain any two-coordinated Si, or isolated oxygen atoms; in other words, initial surfaces obtained by breaking two or more Si-O bonds of the same Si were considered much less stable and not allowed. (iii) H atoms were added to saturate the dangling Si-O bonds created in the bottom layer: it should be noted that a large fraction of nonbridging oxygens (NBOs) is already present in the bulk structure, to balance the charge excess created by the incorporation of modifier Na and Ca cations in the glass. These “structural” NBOs already present in the bottom layer were not saturated, unlike the NBOs created by the surface cut. (iv) Finally, the initial surface was relaxed in a 5 ps CPMD trajectory. In all the calculations, the atoms in the bottom 3 Å of the slab were kept frozen to their initial, bulk-like configuration. Water adsorption was investigated following an effective approach previously adopted to study water on a TiO2 surface:33,68,69 after releasing a water molecule near to a possible adsorption site, the system was left free to evolve in a room temperature CPMD trajectory of 4-5 ps, until the water molecule settled in the adsorption site or moved toward a more favorable one; this final configuration was then geometryoptimized, and the procedure repeated for another initial adsorption site. All the geometry optimizations following the initial MD relaxation were carried out by damped CPMD33,68 until the largest component in the ionic forces on movable atoms was less than 0.01 eV/Å. The energy barrier and minimum energy pathway (MEP) of water dissociative adsorption on a two-membered surface ring were calculated through the String Method of Car-Parrinello (SMCP) dynamics,70,71 as implemented in the QuantumESPRESSO package, version 3.2.3. In the SMCP approach, after identifying and optimizing the initial and final states of a process, a number of replicas of the system are distributed along a “string”, connecting the initial and final states at the two ends. Geometrical constraints are used to maintain the replicas evenly distributed along the string, while the perpendicular forces acting

Tilocca and Cormack

Figure 1. Time evolution of: (top panel) the instantaneous temperature of the systems; (bottom panel) the fraction of undercoordinated Si (nSi3c/ nSi) in the top 5 Å of surfaces A and B. t ) 0 corresponds to the initial unrelaxed surface.

on each replica are simultaneously optimized through damped CP dynamics. Ten replicas were used in the present SMCP calculations. 3. Results 3.1. Relaxation Dynamics and Structure of the Dry Surface. As noted before, both glass surfaces were obtained by CPMD relaxation of suitable slabs cut from the bulk supercell, and they differ for the method used to obtain the initial bulk glass. In the following, “A” shall label the surface whose corresponding bulk glass had previously been obtained by classical MD melt-and-quench, whereas the bulk glass for surface “B” had been generated using a full CPMD procedure. The surface relaxation can be monitored by plotting the fraction of undercoordinated Si3c and NBO species as a function of the simulation time (Figure 1), where t ) 0 corresponds to the start of the CP dynamics for the as-created, unrelaxed surface. The figure shows that Si3c defects are rapidly healed, with their initial fraction roughly halved in 0.5 ps; the relaxation of the surface can be considered complete in ∼1 ps for surface B and ∼1.5 ps for surface A. The slightly longer relaxation time observed for surface A, probably related to the higher initial fraction of defects, prompted us to extend the corresponding MD trajectory up to ∼4.2 ps, to ensure that a stable surface configuration was obtained and no further rearrangements would occur. The progressive self-passivation of the surfaces is accompanied by a corresponding increase in the system temperature (no temperature control was imposed during this phase), which levels off at ∼ 570 K for surface A and 470 K for surface B, with a trend closely following the evolution of the number of defects. The observed short time scale for surface relaxation is in agreement with experimental and MD studies of pure silica glass surfaces:72,73 due to the rapid relaxation, it can be safely assumed that the surface relaxation is complete before water adsorption takes place, and therefore water will normally interact with a fully relaxed dry surface. The atomic composition of the two surface samples is highlighted in the z-profiles of Figure 2, showing the number

Exploring the Surface of Bioactive Glasses

Figure 2. z-profiles of the number of O, Si, P, Na and Ca atoms in the two surface samples A and B (see the text for their definition): for each sample, the slab was divided in slices 2.25 Å thick along z, and the numbers of atoms found in each slice was averaged over the corresponding MD trajectory. The horizontal dashed lines mark the corresponding bulk values.

of each atomic species within slices 2.25 Å thick, averaged over the last 2 ps of the corresponding CPMD trajectory. The top layers (z > 0 Å) correspond to the relaxed surface; from the abrupt decrease in the atomic population above z ) +4 Å one can roughly mark the surface region as the uppermost 3-4 Å of the slabs. As a further consistency check, in order to confirm that the overall thickness of the slab is adequate, Figure 2 also shows that the compositions in the central portion (-2 < z < +2 Å) of the slabs are close to the corresponding bulk values (represented by dashed lines). The structure of the surface can be further characterized through the z-profiles of Si and O species in different coordination, shown in Figure 3. Besides NBOs, which as mentioned before characterize the bulk of these modified glasses as well, a few residual three-coordinated Si atoms (Si3c) are present in both surfaces after relaxation, whereas no three-coordinated O atoms are formed. These defects are visible in the side views of the surfaces in Figure 4. The figures also highlight the presence of small rings on the surfaces: a three-membered (3M) ring is found in surface A, and a twomembered (2M) ring is in surface B. While the two Si atoms of the 2M-ring are fully coordinated, the 3M-ring exposes two undercoordinated Si3c: we have recently observed that undercoordinated Si defects are often found associated to small rings in the melt precursor of bioactive glasses;24 even though both Si3c and small rings in the melt are effectively healed during the cooling phase and therefore not commonly observed in the bulk glasses at room temperature,24 their occurrence in the surface of this composition reflects their metastable nature highlighted in the melt. Their formation is related to the low availability of network formers in the present bioactive composition (containing only ∼45% SiO2): when several Si-O bonds are simultaneously broken, either due to the exposure of the surface or to the high temperature of the melt, partial structural relaxation can be effectively achieved by saturating some of the dangling bonds enclosing them in small rings. Metastable disiloxane and trisiloxane rings have been identified

J. Phys. Chem. C, Vol. 112, No. 31, 2008 11939

Figure 3. z-profiles of the number of 1-,2- and 3-cooordinated O and of 2-, 3-, and 4-coordinated Si atoms in the two surface samples A and B (see the text for their definition): each slab was divided in slices 2.25 Å thick along z, and the numbers of atoms found in each slice was averaged over the corresponding MD trajectory.

Figure 4. Side view of surface A (top) and B (bottom). Si, O, Na, Ca and P atoms are represented as blue, red, gray, green and yellow spheres, respectively.

experimentally on silica surfaces;76–78 unlike most siloxane Si-O-Si bonds, which are chemically inert, the internal strain makes the siloxane bonds in 2- and 3-membered rings much more reactive toward water and other molecules;53,55,57,59,61,79,80 these sites are also deemed to induce calcium phosphate nucleation on bioactive surfaces.2,49,74,75 3.2. Water Adsorption and Dissociation. Despite the relatively small size of the slabs, the present computational

11940 J. Phys. Chem. C, Vol. 112, No. 31, 2008 approach, involving periodically repeated supercells and ab initio representation of many-body forces, provides a realistic picture of the local environment typical of as-formed bioactive glass surfaces. Several potentially active sites are exposed on the dry surfaces obtained by the procedure described in the previous section; we can investigate the reactivity of these sites by using water as an effective probe.33 Our approach involves using standard CP molecular dynamics to effectively explore the complex energy hypersurface of the exposed glass: a water molecule is released close to a possible adsorption site, and its dynamics followed for ∼3-5 ps at room temperature; this time is normally enough for the water molecule to settle in a stable adsorption minimum, which is not necessarily close to the initial site. The finite temperature dynamics allows the molecule to move spontaneously toward regions of lower (free) potential energy, which are then identified rapidly and effectively by direct inspection of the MD trajectory. After roughly identifying the local minimum corresponding to (closer to) the initial arrangement of water on the surface, the exact minimum geometry and energy are obtained by structural optimizations. This approach turns out to be especially effective to examine the highly heterogeneous surface of multicomponent glasses, exhibiting many surface sites with a very flexible local structure, whose individual study would be too computer time-consuming. As in other systematic approaches, such as accurate mapping the local electrostatic potential,81 the purpose is to optimize the identification of the adsorption sites of glass surfaces: in the present work, this optimization is achieved using finitetemperature ab initio MD to drive the water probe toward the most stable sites, screening out less favorable regions of the surface. 3.2.1. Water on Surface A. Four water adsorption modes, as identified on surface A, are shown in Figure 5. In structure (A1), the water molecule has significantly penetrated inside the surface, with its oxygen OW found 2.25 Å below the uppermost NBO. In this configuration, both water hydrogens are involved in hydrogen bonds (Hbs) with surface oxygens, namely a BO belonging to the 3M-ring and a NBO on the opposite side. The molecule is further encaged within two Na ions (at ∼2.3 Å) and a Ca ion at 3.6 Å. In structure (A2), the water molecule is adsorbed on the top of a siloxane bond, and it donates a Hb to an NBO protruding above the surface; as in A1, a modifier cation (Ca in this case) is found at a short distance (2.37 Å) from the water oxygen. These Na(Ca)-OW distances are close to the typical modifier cations-O distances found in bulk bioactive glasses.18,19 In structure (A3) water has transferred a hydrogen to another surface NBO protruding upward from the surface, although there is still a strong Hb linking OW and the H: during the MD trajectory one actually observes fluctuations where the proton is transferred along the Hb from the OW to the NBO and back, as if the two structures were dynamically equivalent. Once again, a modifier Ca cation is located at short distance (2.28 Å) from the OW. As discussed before, a prominent feature exposed in surface A is a 3M-ring, incorporating a pair of three-coordinated Si atoms. This site turns out to be very reactive: about 0.2 ps after releasing a water molecule 2.5 Å above one of the two undercoordinated Si (Figure 6a), a hydrogen is transferred from the water directly to the second Si3c in the ring, leading to a Si-OH silanol linked to a Si-H moiety (Figure 6b). The newly formed SiH bond is highly vibrationally excited, but it remains stable for the remaining 3 ps of the trajectory, where the release of its excess vibrational energy to the other degrees of freedom of the slab leads to an average system temperature around 500

Tilocca and Cormack

Figure 5. Water adsorption modes on surface A. (A1): water has significantly penetrated within the surface, forming hydrogen bonds with a bridging oxygen belonging to the 3M-ring and with an adjacent nonbridging oxygen. (A2): water is coordinated near to a Ca cation and donates a Hb to a surface NBO. (A3): one of the two water O-H bonds is transferred to a surface NBO, leaving a metastable hydroxyl close to a Ca ion. (A4): optimized structure after water adsorption and dissociation on the 3M-ring. Si, O, Na, Ca and P atoms are represented as blue, red, gray, green and yellow spheres, respectively.

K. Hydrosilane Si-H groups have been observed on the surface of activated silicas, and give rise to characteristic IR bands.82 The optimized final structure is labeled A4 in the bottom of Figure 5. Water adsorption and dissociation leads in this case to the saturation of both Si3c defects: the corresponding strong driving force apparently makes this dissociation mechanism more favorable than the hypothetical alternative dissociative mechanism involving opening of the 3M-ring. Indeed, unlike the spontaneous (barrierless) process observed here at room temperature, water dissociation through opening of the 3M-ring would not restore full tetrahedral coordination to the ringforming Si, and therefore it is less favorable, even taking into account the strain energy released by opening the ring. As will be discussed below, the opening of a fully saturated small ring in these glass compositions, although thermodynamically favorable, involves a significant kinetic barrier, at least when only a single water molecule is involved. Therefore the present MD trajectories suggest that, in the presence of coordination defects,

Exploring the Surface of Bioactive Glasses

J. Phys. Chem. C, Vol. 112, No. 31, 2008 11941

Figure 6. Four snapshots extracted from the MD trajectories on surfaces A and B. Top panels: initial and final state for water dissociation on 3-fold coordinated Si embedded in a 3-membered ring on surface A. Bottom panels: initial and final state for water dissociation on 3-fold coordinated Si on surface B. The molecular adsorption initial states (a) and (c) are not stable minima and are rapidly transformed into the corresponding dissociative products (b) and (d), whose optimized structures are A4 in in Figure 5 and B1 in Figure 7, respectively.

Figure 7. Structure-optimized water adsorption modes on surface B. (B1): after water adsorption and dissociation on the 3-fold-coordinated Si. (B2): after molecular adsorption on a Si atom of the 2M-ring. (B3): molecular adsorption and formation of a H-bonded complex bridging between the 2M-ring and an adjacent NBO. (B4): dissociative adsorption and opening of the 2M-ring. Si, O, Na, Ca, P, and H atoms are represented in blue, red, gray, green, yellow, and white colors, respectively. Only atoms close (within a 5 Å radius) to the adsorption sites have been extracted from the periodic structure and plotted.

alternative water dissociation pathways may prevail, leading to saturation of these defects while at the same time keeping the ring closed. 3.2.2. Water on Surface B. Water dissociation takes place immediately after releasing the molecule above a Si3c defect in surface B (Figure 6c,d); the process is favored by the presence of a nearby NBO bonded to an adjacent Si: this NBO accepts a proton from the water molecule through a transient 6M-ring intermediate, and the resulting hydroxylated structure is stable over a subsequent MD run of 2 ps at 300K. The final optimized structure with the two vicinal silanols is shown in Figure 7B1. On the other hand, when the water molecule is placed above the 2M-ring, no stable adsorption is observed in the subsequent MD trajectory, regardless of the initial distance/orientation of

the molecule with respect to the ring. For instance, we have attempted to start the dynamics with the molecule on either 4-fold-coordinated Si of the ring, and in both following MD trajectories the molecule moved away from the 2M-ring and adsorbed on a different modifier cation site nearby, similar to those discussed in the previous section (A1-A3 structures in Figure 5). As it will be shown below, this is a consequence of the metastable nature of the ring adsorption state involving a 5-fold-coordinate Si (shown in Figure 7B2), and the presence of a significant energy barrier separating this state from the more stable state (see below) with the 2M-ring opened through water dissociation. In other words, even though water chemisorption through 2M-ring opening is thermodynamically favorable, a single water molecule approaching the 2M-ring may not immediately dissociate, especially if several other favorable adsorption sites are available nearby, as in this case. Another possible state with water molecularly adsorbed on the 2M-ring emerged in a short MD trajectory started with the molecule located above and donating a H-bond to the uppermost BO in the ring. During the trajectory, a second Hb was formed between the other water hydrogen and another adjacent NBO; the optimized structure of this H-bonded complex is shown in Figure 7B3. In order to obtain the hydroxylated structure produced by ring opening, the topmost bridging oxygen in the 2M-ring was removed, and an OH group was attached 1.85 Å above each of the two Si atoms in the ring. This initial configuration was then used to start a RT trajectory of ∼2 ps, which showed that the opened-ring structure with two vicinal silanols is stable: the Si-OH groups, albeit very mobile and occasionally interacting through a H-bond during the trajectory, did not condense back into molecular water. The optimized structure with two Hbonded vicinal silanols is shown in Figure 7B4. 3.2.3. Opening of 3M-Ring. In order to study the opening of a 3M-ring, rather than on the 3M-ring containing undercoordinated Si in surface A, we focused on a defect-free 3M-ring which we were able to identify on another slab obtained from a larger bulk sample.83 The optimized structure of the 3M-ring is shown in Figure 8C1: the plane formed by the three Si atoms is roughly normal to the surface plane, thus resulting in a Si-O bond considerably more exposed and therefore prone to attack by water than the others. Therefore, in order to model the effect of ring opening by water dissociative adsorption, an OH and a

11942 J. Phys. Chem. C, Vol. 112, No. 31, 2008

Tilocca and Cormack

Figure 9. Potential energy profiles along the minimum-energy paths for the two ring-opening mechanisms discussed in the text. Dashed line: mechanism 1, starting from the H2O-Si(ring) complex. Full line: mechanism 2, starting with the H2O-BO(ring) complex. Letters 1a1d and 2a-2d label selected replicas whose structure is shown in Figures 10 and 11, respectively. A cubic polynomial was used for interpolation between the adjacent energy points.

Figure 8. (C1) Optimized structure of an exposed fully saturated 3Mring; (C2) corresponding hydroxylated structure obtained upon water adsorption and dissociation with ring opening. Si, O, Na, Ca, P, and H atoms are represented in blue, red, gray, green, yellow, and white colors, respectively. Only atoms close to the 3M-ring have been extracted from the periodic structure and plotted.

TABLE 1: Water Adsorption Energies ∆Eads (eV); ∆Eads ) -[Eint - Eslab - Ewat].a

a

system

∆Eads (eV)

A1 A2 A3 A4 B1 B2 B3 B4 C2

0.75 1.31 0.59 3.13 2.45 0.04 0.44 1.16 0.89

The different systems are labeled as in Figures 5, 7, and 8.

proton were attached to the Si and O of the exposed Si-O bond, respectively, and the structure was optimized, leading to the configuration shown in Figure 8C2. A short MD run at room temperature followed, during which the hydroxylated structure was stable. 3.2.4. Adsorption Energies. After structural optimizations, the water adsorption energy was calculated as: ∆Eads ) -[Eint - Eslab - Ewat], where Eint is the total energy of the interacting slab/water system, Eslab is the energy of the bare slab, and Ewat is the energy of an isolated water molecule calculated in the same supercell used for the surface.68 The calculated adsorption energies are tabulated in Table 1. The table shows that water dissociative adsorption at Si3c defects is the most favorable process, either on the 3M-ring of surface A (∆Eads )3.13 eV) or on the Si3c site in surface B (∆Eads )2.45 eV). The higher ∆Eads on surface A reflects the higher stabilization produced by saturating two dangling Si instead of

just one as for surface B. If water dissociation is coupled to the opening of a 2M-ring (as in structure B4), the resulting adsorption energy of 1.16 eV shows that the release of strain energy in the ring makes this process highly favorable, albeit significantly less than water dissociation at undercoordinated Si defects; the opening of a fully saturated 3M-ring (structure C2) also appears energetically favorable, even though 0.26 eV less than the more strained 2M-ring. Considerable stabilization (1.3 eV) results from molecular adsorption modes in close proximity of modifier Ca cations, as in structure A2; molecular adsorption on Na, as in A1, is less favorable, but still rather stable. On the other hand, the initial molecular adsorption on the siloxane fragments of the small rings is not as favorable (0.04-0.44 eV for structures B2 and B3, respectively), which explains the low residence times of water physisorbed on these sites during the room temperature MD. Nevertheless, the analysis in the next section confirms that metastable configurations where water is directly physisorbed on Si-BO-Si siloxane linkages must necessarily be involved in the opening of 2Mrings. 3.2.5. Determination of Energy Barriers for 2M-Ring opening. Unlike defects, coordinatively saturated rings do not spontaneously dissociate water during the room temperature dynamics, despite the higher stability of the final hydroxylated state, thus suggesting that a kinetic barrier to water dissociation and ring opening may slow down the process. This is especially important for water adsorption on the Bioglass surface, because of the availability of many alternative sites, involving modifier cations, where the molecule can be attracted and trapped. In order to quantitatively confirm this possibility, the energy barrier for opening the 2M-ring on the glass surface B was evaluated using the SMCP approach. Two reaction paths were considered: in the first, the initial state was the one in Figure 7B2, with molecular water adsorbed on a Si atom in the ring, whereas the initial state of the second path was the structure in Figure 7B3, with water H-bonded to the exposed BO in the ring. The common final state for the two paths was the one shown in Figure 7B4, with two silanols replacing the 2M-ring. The potential energy profiles along the minimum energy path (MEP), obtained by SMCP minimization, are shown in Figure 9; the corresponding reaction mechanisms are summarized in Figures 10 and 11, by means of the optimized structures of selected replicas taken from the corresponding MEP. The energy profiles clearly show that in both cases an energy barrier is present, as it had been already inferred from the MD simulations. The barrier to ring-opening mechanism 1, involving a (metastable)

Exploring the Surface of Bioactive Glasses

J. Phys. Chem. C, Vol. 112, No. 31, 2008 11943

Figure 10. Minimum Energy Paths (MEP) for 2M-ring opening mechanism 1, starting from H2O coordinated to Si in the ring. Structures 1a-1d are optimized replicas along the MEP, corresponding to the points labeled in Figure 9.

Figure 11. Minimum Energy Paths (MEP) for 2M-ring opening mechanism 2, starting from H2O coordinated to the exposed BO in the ring. Structures 2a-2d are optimized replicas along the MEP, corresponding to the points labeled in Figure 9.

5-coordinated Si intermediate, is smaller than the barrier to mechanism 2, where the initial interaction of the dissociated water with the ring does not involve a Si, but a H-bond to a BO. The relatively high energy barrier (1.35 eV) makes mechanism B less likely, despite the more stable initial state; therefore, even though an initial H-bond complex can be initially formed between water and a BO, water will eventually move apart and attach to a Si, forming the transient intermediate with a 5-fold-coordinated Si, from which water dissociation and ring opening can proceed with a lower barrier, of 0.62 eV. The opening of a 2 M ring on a pure amorphous SiO 2 surface had been modeled using constrained MD in ref 55, where the authors determined a barrier of 0.32 eV going from an initial state with water physisorbed on the silicon site of the ring to the final hydroxylated state. The lower barrier, compared to the corresponding one obtained in this work,85 might denote that 2M-rings on the Bioglass surface can be less reactive than on amorphous SiO2; however, another CPMD calculation led to an estimated activation energy of 0.9 eV,59 and actually our result turns out to be closer to the barrier calculated by Walsh et al. using cluster models,86 which also identify a similar transition state with a (metastable) 5-coordinated Si. 4. Discussion The ab initio simulations of the dry Bioglass surface and its interaction with a water probe provide a very detailed picture of the surface sites available before the glass is immersed in the physiological environment. While most experimental investigations study the bioactivity of these materials by looking at the physicochemical transformation occurring at their surface after immersion, the current study allows us to focus on the structure of the surface and to directly compare the activity of the different sites accessible to a water molecule. The emerging picture of the dry surface is a rather complex one: due to the low silica content and the high network fragmentation of the bulk glass, silicate groups close to the surface tend to associate in small (2- and 3-membered) rings in order to partially saturate the exposed dangling bonds, although not all coordinatively unsaturated Si are eliminated by surface relaxation, and highly reactive Si3c defects are left exposed. On the relaxed surface, exposed silicate fragments are alternated with regions populated by a higher concentration of modifier Na and Ca cations, which provide stable adsorption sites for molecular water.

The strongest adsorption sites are Si3c defects: at room temperature, water immediately dissociates upon contact with an Si3c, although it appears that the close proximity of a strong basic site, such as an NBO, is essential to help the dissociation by accepting the water proton.56,58,81,86 Even though the Si3c defects examined in our surface samples satisfied this requirement, there is the possibility that some isolated Si3c, not accompanied by an NBO, might be also present on the dry surface. On these sites, water will adsorb molecularly, despite a presumably high adsorption energy, unless additional water molecules are available to assist the dissociation process by accepting and transferring the proton.60,69 The calculated dissociative adsorption energy at defects (2.45 eV) could be compared with the experimental values of the initial enthalpy of adsorption on silica obtained from calorimetric methods: on the strongest sites characteristic of quartz dusts, likely associated to defects, water dissociation releases 2.1-2.5 eV.32 Other strong adsorption sites, not present on pure silica surfaces, are provided by modifier Na and (especially) Ca cations: previous MD simulations had shown that the introduction of Na modifier cations in silica glass introduces additional hydrophilic adsorption sites, of intermediate strength between network defects and sites with lower affinity for water, in agreement with our results.73 On the bioactive glass surface, modifier cations appear to behave as Lewis acids, coordinating the oxygen of a water molecule, which is adsorbed molecularly, although in some cases the adsorption can lead to a pseudodissociated state, where a water proton is transferred to a neighboring NBO, along a H-bond (Figure 5A3). Even though the metastable nature of the latter state seems to indicate that the stabilization of the resulting hydroxyl by the cation is not as ideal as it would be on a Si atom, these adsorption sites represent an important feature of the surface of bioactive glasses. In fact, the MD runs often highlighted a preference of water for binding to these cation sites, rather than to the exposed siloxane chains, provided that nearby NBOs were also available for a hydrogen bond. Because these Na/Ca-rich surface regions also appear relatively open, they could provide thermodynamically favorable channels for water to directly penetrate under the surface and inside the glass network, without first having to break surface siloxane bridges in order to open its way to the inner regions. An indication of this possibility is given in one of the MD trajectories, where water initially released above

11944 J. Phys. Chem. C, Vol. 112, No. 31, 2008

Tilocca and Cormack are associated with Si3c defects: this association is not unlikely, and it had been already noticed in the corresponding melt.24 5. Final Remarks

Figure 12. Distributions of O-Si-O and Si-O-Si angles calculated from the MD runs; (black) total distribution; (red) calculated for 2Mrings only; (green) calculated for 3M-rings only.

the surface rapidly finds a favorable entrance to the glass and after a few picoseconds is found in the configuration shown in Figure 5A1. This effect could have important consequences for the mechanism of dissolution of bioactive glasses in an aqueous environment, such as in the human host. Although saturated siloxane chains are exposed in our sample surfaces, they do not appear as stable adsorption sites as those involving NBOs: the hydrophobic character of dehydrated, defect-free silica surfaces, dominated by unstrained Si-O-Si linkages, is indeed well established.32,80,87,88 Furthermore, the present MD runs and static optimizations highlighted that the initial water physisorption on small 2 or 3 M rings is rather weak, in qualitative agreement with similar calculations of water interacting with 2M-ring in pure silica surfaces.56,57 These rings, especially the disiloxane ones, are usually considered highly reactive, on the basis of their considerable internal strain.79 The small rings formed on the present surface samples are indeed highly strained: Figure 12 shows that the O-Si-O and Si-O-Si angle distributions of our slabs include features at 90° associated to the internal angles in 2M-rings, and the angles associated to the 3M-ring are also significantly smaller than the unstrained O-Si-O and Si-O-Si ones, centered at 109° and 134°, respectively. In fact, the hydroxylation energy for 2Mring opening (1.16 eV) confirms that the release of strain energy makes the opening of 2M-rings thorugh water dissociation thermodynamically favorable also on the surface of bioactive glasses. The calculated hydroxylation energy is less favorable than the 1.7 eV value estimated for the pure silica glass surface in previous ab initio CPMD simulations,55 but it compares well with another recent DFT estimate (1.06 eV),57 and experimental energies around 1-1.5 eV were obtained for the dissociative adsorption enthalpy of water on highly dehydrated amorphous silicas.32 Notwithstanding the favorable energetic balance, a significant kinetic barrier slows down the reaction: the minimum energy path for the opening of a 2M-ring, determined using high-level SM optimizations, indicates that preadsorption of water on a ring Si is required, forming a 5-fold-coordinated intermediate, which then transforms into the hydroxylated structure overcoming a barrier of around 0.6 eV. The presence of this barrier might result in some of these small rings remaining closed upon contact with water at room temperature, especially due to the presence of more favorable physisorption sites, such as the modifier cations discussed before, where molecular water might migrate before being able to overcome the barrier to dissociation/ring opening. Our MD simulations also show that the barrier to opening of small rings on the surface can be lowered or completely removed when these rings

The present simulations represent the first atomistic exploration of the dry surface of bioactive glasses, performed using realistic slab models and state-of-the-art ab initio methods employing periodic supercells. Several unique features of these biomaterials have been highlighted by modeling the interaction of a water probe with the surface, resulting in a significant step toward a more fundamental understanding of their bioactive behavior. Due to the low silica content, the surface of bioactive silicate glasses is characterized by significant roughness and heterogeneity: one of the most significant features emerged here is the presence of favorable adsorption sites mostly populated by Na/Ca cations, which can attract small polar molecules such as water, allowing them to penetrate under the surface. This seems to denote a specific role of these sites in the partial dissolution of the glass surface, which is a key step in the bioactive fixation mechanism. Additional surface features are the small, strained 2- and 3-membered silicate rings; the simulations show that, despite a thermodynamic driving force toward the opening of these rings through water dissociation, a water molecule approaching the surface will initially interact weakly with these rings, due to the presence of a significant energy barrier separating the intermediate physisorbed state from the stable hydroxylated state. The presence of more favorable molecular adsorption sites, such as the modifier cations, will in some cases result in the water molecule being pulled away from the rings before it can be dissociatively chemisorbed. This is especially relevant in the context of the bioactive mechanism: because small rings can directly favor the nucleation and deposition of a Ca/P film from the physiological solution,2 their availability on the surface of bioactive glass, even after the material is immersed in an aqueous environment, could be crucial for the development of the bone-bonding interface. Finally, the role of three-coordinated Si defects was highlighted: spontaneous water dissociation occurs at these sites, provided that they are associated with a proton acceptor such as a nonbridging oxygen: it is likely that water will be initially be driven toward these strong sites, where dissociation will produce vicinal surface silanols, and in some cases also hydrosilane groups; when Si3c are associated with small rings, they can lower or remove the kinetic barrier to ring opening. After these strong chemisorption sites are fully hydroxylated, additional water molecules are likely to be adsorbed molecularly, although cooperative interactions between several water molecules can also induce further dissociation and surface hydroxylation, and mixed molecular/dissociated water adsorption might eventually be produced.69 Although the investigation of these effects is beyond the scope of this paper, which is aimed at exploring the features of the as-created, dry surface of bioactive glasses, a further essential step in any fundamental study of these biomaterials will necessary involve the interface with bulk liquid water, whose ab initio modeling requires much larger computational resources. Acknowledgment. A.T. thanks the UK’s Royal Society for financial support (University Research Fellowship). Computer resources on the HPCx service were provided via the UK’s HPC Materials Chemistry Consortium and funded by EPSRC (portfolio grant EP/D504872). References and Notes (1) Hench, L. L.; Andersson, O. H. in An Introduction to Bioceramics; Hench, L. L.; Wilson, J., Eds,: Singapore; World Scientific, 1993.

Exploring the Surface of Bioactive Glasses (2) Cerruti, M.; Sahai, N. ReV. Min. Geochem. 2006, 64, 283. (3) Vallet-Regi, M.; Ragel, C. V.; Salinas, A. J. Eur. J. Inorg. Chem. 2003, 1029. (4) Hench, L. L. Science 1984, 226, 630. (5) Hench, L. L. J. Am. Ceram. Soc. 1998, 81, 1705. (6) Boccaccini, A. R.; Chen, W.; Lefebvre, L.; Gremillard, L.; Chevalier, J. Faraday Discuss. 2007, 136, 27. (7) Hench, L. L.; Polak, J. M. Science 2002, 295, 1014. (8) Ebisawa, Y.; Kokubo, T.; Ohura, K.; Yamamuro, T. J. Mater. Sci. Mat. Med. 1990, 1, 239. (9) Cerruti, M.; Bianchi, C. L.; Bonino, F.; Damin, A.; Perardi, A.; Morterra, C. J. Phys. Chem. B 2005, 109, 14496. (10) Vogel, M.; Voigt, C.; Gross, U. M.; Muller-Mai, C. M. Biomaterials 2001, 22, 357. (11) Fujibayashi, S.; Neo, M.; Kim, H.-M.; Kokubo, T.; Nakamura, T. Biomaterials 2003, 24, 1349. (12) Lebecq, I.; Desanglois, F.; Leriche, A.; Follet-Houttemane, C. J. Biomed. Mater. Res. A 2007, 83, 156. (13) Xin, R. L.; Leng, Y.; Chen, J. Y.; Zhang, Q. Y. Biomaterials 2005, 26, 6477. (14) Chen, Q. Z.; Thompson, I. D.; Boccaccini, A. R. Biomaterials 2006, 27, 2414. (15) Strnad, Z. Biomaterials 1992, 13, 317. (16) Hill, R. G. J. Mater. Sci. Lett. 1996, 15, 1122. (17) Leonova, E.; Izquierdo-Barba, I.; Arcos, D.; Lopez-Noriega, A.; Hedin, N.; Vallet-Regi, M.; Eden, M. J. Phys. Chem. C 2008, 112, 5552. (18) FitzGerald, V.; Pickup, D. M.; Greenspan, D.; Sarkar, G.; Fitzgerald, J. J.; Wetherall, K. M.; Moss, R. M.; Jones, J. R.; Newport, R. J. AdV. Func. Mater. 2007, 17, 3746. (19) Tilocca, A.; Cormack, A. N.; de Leeuw, N. H. Chem. Mater. 2007, 19, 95. (20) Lusvardi, G.; Malavasi, G.; Menabue, L.; Menziani, M. C.; Pedone, A.; Segre, U. J. Phys. Chem. B 2005, 109, 21586. (21) Linati, G.; Lusvardi, G.; Malavasi, G.; Menabue, L.; Menziani, M. C.; Mustarelli, P.; Segre, U. J. Phys. Chem. B 2005, 109, 4989. (22) Tilocca, A.; de Leeuw, N. H. J. Phys. Chem. B 2006, 110, 25810. (23) Tilocca, A.; Cormack, A. N.; de Leeuw, N. H. Faraday Discuss. 2007, 136, 45. (24) Tilocca, A. Phys. ReV. B 2007, 76, 224202. (25) Tilocca, A.; Cormack, A. N. J. Phys. Chem. B 2007, 111, 14256. (26) Clark, A. E.; Hench, L. L. J. Biomed. Mater. Res. 1976, 10, 161. (27) Kokubo, T.; Kim, H.; Kawashita, M. Biomaterials 2003, 24, 2161. (28) Cho, S. B.; Miyaji, F.; Kokubo, T.; Nakanishi, K.; Soga, N.; Nakamura, T. J. Mater. Sci. Mat. Med. 1998, 9, 279. (29) Ramila, A.; Vallet-Regi, M. Biomaterials 2001, 22, 2301. (30) Cerruti, M.; Magnacca, G.; Bolis, V.; Morterra, C. J. Mater. Chem. 2003, 13, 1279. (31) Cerruti, M.; Bolis, V.; Magnacca, G.; Morterra, C. Phys. Chem. Chem. Phys. 2004, 6, 2468. (32) Bolis, V.; Fubini, B.; Marchese, L.; Martra, G.; Costa, D. J. Chem. Soc. Faraday Trans. 1991, 87, 497. (33) Tilocca, A.; Selloni, A. J. Chem. Phys. 2003, 119, 7445. (34) Huang, C.; Cormack, A. N. J. Chem. Phys. 1991, 95, 3634. (35) Smith, W.; Greaves, G. N.; Gillan, M. J. J. Chem. Phys. 1995, 103, 3091. (36) Horbach, J.; Kob, W.; Binder, K. J. Phys. Chem. B 1998, 103, 4104. (37) Lammert, H.; Heuer, A. Phys. ReV. B 2005, 72, 214202. (38) Mead, R. N.; Mountjoy, G. J. Phys. Chem. B 2006, 110, 14273. (39) Tilocca, A.; de Leeuw, N. H. J. Mater. Chem. 2006, 16, 1950. (40) Donadio, D.; Bernasconi, M.; Tassone, F. Phys. ReV. B 2004, 70, 214205. (41) Charpentier, T.; Ispas, S.; Profeta, M.; Mauri, F.; Pickard, C. J. J. Phys. Chem. B 2004, 108, 4147. (42) Du, J.; Corrales, L. R. J. Chem. Phys. 2006, 125, 114702. (43) Tilocca, A.; de Leeuw, N. H.; Cormack, A. N. Phys. ReV. B 2006, 73, 104209. (44) Mahadevan, T. S.; Garofalini, S. H. J. Phys. Chem. B 2007, 111, 8919. (45) Levine, S. M.; Garofalini, S. H. J. Chem. Phys. 1987, 86, 2997. (46) Zeitler, T. R.; Cormack, A. N. J. Cryst. Growth 2006, 294, 96. (47) Giordano, L.; Sushko, P. V.; Pacchioni, G.; Shluger, A. L. Phys. ReV. B 2007, 75, 024109. (48) Donadio, D.; Bernasconi, M.; Boero, M. Phys. ReV. Lett. 2001, 87, 195504. (49) Sahai, N.; Tossell, J. A. J. Phys. Chem. B 2000, 104, 4322.

J. Phys. Chem. C, Vol. 112, No. 31, 2008 11945 (50) Cerruti, M.; Morterra, C.; Ugliengo, P. J. Mater. Chem. 2004, 14, 3364. (51) Car, R.; Parrinello, M. Phys. ReV. Lett. 1985, 55, 2471. (52) Tilocca, A.; Di Valentin, C.; Selloni, A. J. Phys. Chem. B 2005, 109, 20963. (53) Ceresoli, D.; Bernasconi, M.; Iarlori, S.; Parrinello, M.; Tosatti, E. Phys. ReV. Lett. 2000, 84, 3887. (54) Rignanese, G.-M.; De Vita, A.; Charlier, J.-C.; Gonze, X.; Car, R. Phys. ReV. B 2000, 61, 13250. (55) Masini, P.; Bernasconi, M. Phys. ReV. B 2002, 14, 4133. (56) Rignanese, G.-M.; Charlier, J.-C.; Gonze, X. Phys. Chem. Chem. Phys. 2004, 6, 1920. (57) Du, M.-H.; Kolchin, A.; Cheng, H.-P. J. Chem. Phys. 2004, 120, 1044. (58) Du, J.; Cormack, A. N. J. Am. Ceram. Soc. 2005, 88, 2532. (59) Mischler, C.; Horbach, J.; Kob, W.; Binder, K. J. Phys.: Condens. Matter 2005, 17, 4005. (60) Ma, Y.; Foster, A. S.; Nieminen, R. M. J. Chem. Phys. 2005, 122, 144709. (61) Yang, Y.; Wang, E. G. Phys. ReV. B 2006, 73, 035406. (62) Ganster, P.; Benoit, M.; Delaye, J.-M.; Kob, W. Surf. Sci. 2008, 602, 114. (63) Handbook of Glass Data O. V. Mazurin.; M. V. Streltsina.; T. P. Shvaiko-Shvaikovskaya., Eds: Elsevier, 1993. (64) P. Giannozzi et al., http://www.quantum-espresso.org. (65) J.P., Perdew.; K., Burke.; M., Ernzerhof. Phys. ReV. Lett. 1996, 77, 3865. (66) Vanderbilt, D. Phys. ReV. B 1990, 41, 7892. (67) With the 2% volume expansion of the bulk cell and the 14.5 Å elongation of the c parameter, the dimensions of the slab were 11.7 Å × 11.7 Å × 26.2 Å; taking into account the subsequent introduction of H atoms in the bottom layer, this effectively corresponds to a vacuum gap of at least 10 Å between a water molecule adsorbed on the top layer and the periodic image of the H-saturated bottom layer; we have checked that no significant differences resulted when a larger vacuum separation was used. (68) Tilocca, A.; Selloni, A. Langmuir 2004, 20, 8379. (69) Tilocca, A.; Selloni, A. J. Phys. Chem. B 2004, 8, 4743. (70) Kanai, Y.; Tilocca, A.; Selloni, A.; Car, R. J. Chem. Phys. 2004, 121, 3359. (71) Tilocca, A.; Selloni, A. ChemPhysChem 2005, 8, 1911. (72) D’Souza, A. S.; Pantano, C. G. J. Am. Ceram. Soc. 1999, 82, 1289. (73) Leed, A. E.; Pantano, C. G. J. Non-Cryst. Solids 2003, 325, 48. (74) Hench, L. L.; West, J. K. Annu. ReV. Mater. Sci. 1995, 25, 37. (75) Sahai, N.; Anseau, M. Biomaterials 2005, 26, 5763. (76) Bunker, B. C.; Haaland, D. M.; Michalske, T. A.; Smith, W. L. Surf. Sci. 1989, 222, 95. (77) Wallace, S.; West, J. K.; Hench, L. L. J. Non-Cryst. Solids 1993, 152, 101. (78) Chiang, C. M.; Zegarski, B. R.; Dubois, L. H. J. Phys. Chem. 1993, 97, 6948. (79) B.A., Morrow.; B., A.; Cody, I. A.; Lee, L. S. M. J. Phys. Chem. 1976, (80), 2761. (80) Grabbe, A.; Michalske, T. A.; Smith, W. L. J. Phys. Chem. 1995, 99, 4648. (81) Leed, A. E.; Sofo, J. O.; Pantano, C. G. Phys. ReV. B 2005, 72, 155427. (82) Pacchioni, G.; Vitiello, M. Phys. ReV. B 1998, 58, 7745. (83) Tilocca, A.; Cormack, A. N. manuscript in preparation (84) Feuston, B. P.; Garofalini, S. H. J. Appl. Phys. 1990, 68, 4830. (85) Even though the authors investigate the effect of using two different reaction coordinates, their initial state in both paths is analogous to the one in Figure 7B2 of the present paper, and therefore both their estimated activation energies should be compared with the corresponding value of 0.62 eV obtained for path A in this work. It should be remarked that the method used in the present paper to locate the MEP is independent on the choice of an arbitrary reaction coordinate, and it allows us to better discriminate between alternative reaction paths, with different initial or final states. (86) Walsh, T. R.; Wilson, M.; Sutton, A. P. J. Chem. Phys. 2000, 113, 9191. (87) Iler, R. K., The Chemistry of Silica; New York: Wiley Interscience, 1979. (88) Bakaev, V. A.; Steele, W. A. J. Chem. Phys. 1999, 111, 9803.

JP803541J