ERNEST ORLANDO LAWRENCE BERKELEY NATIONAL LABORATORY

Molecular Simulations of Biogenic Nanoparticles

Current Research Team

Research Goals

Because of their ubiquitous presence in natural materials and their strong surface reactions, nanoparticles produced by bacteria and fungi in natural settings figure importantly in a broad range of environmental phenomena, from global climate change to contaminant remediation. Our research is designed to provide molecular-scale information about the structure and surface chemistry of these nanoparticles (in particular, layer type minerals), based on state-of-the-art computer simulation using tested codes and realistic models of layer type minerals. Computer simulations allow us to explore molecular details that spectroscopy cannot probe or fully elucidate. In particular, they often bring out features of nanoparticle electronic structure that are important to materials science applications.
 

Methods

We have been developing an approach to the surface geochemistry of smectitie, Mn oxide, and Fe sulfide based on molecular-level simulation. Our research, supported by the DOE-BES program at Lawrence Berkeley National Laboratory, is aimed at understanding the mechanisms by which these minerals interact with water molecules, cations and organic molecules in aqueous subsurface environments. Our principal techniques, molecular dynamics (MD) and density functional theory (DFT) simulations, are well-established, essential components of research in theoretical physical chemistry.

The underlying approach in MD simulations is to construct potential functions that model all of the known interactions in a system of ions, atoms, and molecules, then devise a strategy for sampling the phase space of the interacting system in order to compute its properties. In a MD simulation, the phase space of the system is sampled through numerical integration (Beeman algorithm) of the Newton-Euler equations of motion for each molecular species, which is performed consistently with the suite of potential functions assumed.

In recent years the interpretation of both diffraction and spectroscopic data on clay minerals has been facilitated by a class of independent, first-principles (sometimes termed "ab initio") quantum-mechanical simulations, which now have achieved sufficient accuracy to predict crystallographic properties of clay minerals without recourse to empirical parameterization. The ingredients of this scheme are the atomic nuclei and electrons whose interactions are described by the density functional theory (DFT) formulation of quantum mechanics. DFT-based methods are able to describe structure and bonding properties to a high degree of accuracy, including structural trends. Indeed, these methods represent the only practical quantum-mechanical approach for studying complex materials such as clay minerals and metal oxides.

Our MD calculations utilize the code MOLDY, developed by K. Refson. Our ab initio calculations are performed using a parallel version of the code CASTEP, developed by the UKCP (United Kingdom Car-Parrinello) consortium. Numerical calculations have been carried out under DOE support on Cray XT4, Cray XE6, and IBM iDataPlex supercomputers at the NERSC (National Energy Research Scientific Computing Center).

top

Research Projects

Magnetic ordering in mackinawite

Mackinawite (FeS) is a layer type mineral (P4/nmm space group) composed of sheets of edge-sharing FeS₄ tetrahedra stacked along the c-axis and held together through van der Waals forces (Fig. 1a). In nature, this mineral is typically the first solid phase to be precipitated during iron sulfide mineralization mediated by sulfate-reducing bacteria. Mackinawite is known to influence the mobility of environmentally-important trace elements through sorption and redox processes, and it figures importantly in the well-known iron-sulfur hypothesis about the origin of life. Its Mössbauer spectrum exhibits one singlet without quadrupole and hyperfine coupling, indicating that there is no magnetic ordering from room temperature down to 4.2 K. Our DFT calculations, however, show a magnetic moment on the Fe ions and that the ground state of mackinawite should be both metallic and magnetic (Fig. 1b). A similar conflict between Mössbauer spectra and DFT calculations has been reported for Fe-based high-temperature superconductors, such as FeSe (T_c = 12 K), LaO_0.89F_0.11FeAs (T_c = 26 K), and CeO_0.85FeAs (T_c = 46.5 K), which are isostructural with mackinawite and for which it is known that magnetism is suppressed by strong magnetic moment (spin) fluctuations.

If Fe ions in metallic systems like mackinawite have a magnetic moment, multiplet splitting in the Fe 3s core-level photoemission spectrum occurs due to the exchange interaction of core 3s electrons with the net spin of the 3d electrons. The Fe 3s core-level photoemission spectrum of mackinawite, obtained at the Advanced Light Source (Fig. 1c), indeed showed exchange-energy splitting (ΔE_3s) consistent with a magnetic moment on the Fe ions. All of these computational and experimental results indicate the existence of strong itinerant spin fluctuations in mackinawite. Because these spin fluctuations are believed to the main mediators of electron-pairing in the Fe-based superconductors, we conjecture that mackinawite may in fact be a superconductor.

Figure 1. (a) Mackinawite structure (Purple = Fe; Orange = S, a = b = 3.674 Å, c = 5.033Å). (b) DFT-calculated total energy of mackinawite per formula unit vs. the lattice parameter in different magnetic orderings (FM: ferromagnetic; NM: non-magnetic; AFM: antiferromagnetic order). In the AFM ordering, three types were considered. The vertical line represents the experimental lattice parameter. (c) Experimental Fe-3s core-level photoemission spectrum of mackinawite. ΔE_3s is the multiplet splitting energy, equal to 2.9 eV.
 

Mn-vacancy-induced photoconductivity of layer type manganese oxides

Manganese(IV) oxides are known to impact a broad range of biogeochemical processes, mainly through their high capacity for metal sorption and their facile oxidation of organic and inorganic compounds. Mn oxides found in weathering environments and natural waters are produced mainly by bacteria and take on layer structures, i.e., stacks of sheets of edge-sharing MnO₆ octahedra. An important structural characteristic of these oxides is the presence of Mn(IV) cation vacancies whose charge deficit is typically compensated by metal cations or protons. These vacancies have long been identified as strong adsorption sites for metals, but they may also play an important role in Mn redox biogeochemistry, particularly photo-induced redox reactions. Because detailed electronic band structure is the key to understanding photo-induced redox reactions, Mn(IV) oxides both vacancy-free and containing cation vacancies charge-compensated by protons were investigated using ab initio quantum mechanics simulations as realized in the code, CASTEP.

The ab initio study showed that a Mn(IV) vacancy reduces the band gap energy (Fig. 2a) and separates photo-excited electrons and holes (Fig. 2b). A reduction in band-gap energy generates more pairs of electrons and holes upon illumination, and the distinct separation of charge carriers enhances their transfer before loss by recombination. Therefore, Mn(IV) vacancies enhance effectively the photoconductivity of layer type MnO₂, facilitating photo-redox reactions between the mineral and inorganic or organic compounds.

Recent studies in materials science indicate that synthetic layer type Mn(IV) oxide nanoparticles with cation vacancies like those found in the biogenic MnO₂ minerals are semiconductors that produce photocurrent under visible light stimulation, thus making them very attractive for applications in energy storage, solar cell fabrication, and catalysis. Our prediction of effective band gap energy reduction by cation vacancies indicates that photocurrent production by these layer type MnO₂ nanoparticles can be optimized by the control of vacancy concentrations during laboratory synthesis.

Figure 2. Calculated electronic properties of vacancy-free MnO₂ (left) and MnO₂ containing 12.5 or 3.3 % Mn(IV) vacancies charge-compensated by protons (right). (a) Band structures, showing reduction of the band gap energy (denoted by a red double arrow) from 1.3 to 0.3 eV. (Solid lines denote bands for spin-up electrons while dashed lines denote bands for spin-down electrons.) (b) Charge distributions for valence-band maximum hole states (orange) and conduction-band minimum electron states (blue), showing that a cation vacancy compensated by protons (H) effectively separates photo-induced charge carriers against recombination.
 

Surface complexes of zinc in birnessite interlayers

Biogeochemical cycling of zinc (Zn) is strongly influenced by sorption on hexagonal birnessite. Spectroscopic studies have identified that Zn forms both tetrahedral (Zn^IV) and octahedral (Zn^VI) triple-corner-sharing surface complexes (TCS) at Mn(IV) vacancy sites in hexagonal birnessite (Fig. 3). The octahedral complex is expected because it is similar to that of Zn in the Mn oxide mineral, chalcophanite (ZnMn₃O₇·3H₂O). However, the reason for the occurrence of the four-coordinate Zn surface species remains unclear. We examined this issue using spin-polarized DFT.

Total energy, magnetic moments of vacancy O, and electron-overlap populations between Zn and vacancy O all have proved that Zn^IV-TCS is stable in birnessite without a need for Mn(III) substitution in the octahedral sheet and that it is more effective in reducing undersaturation of surface O at a Mn vacancy than is Zn^VI-TCS. Geometry-optimization of a hypothetical monohydrate mineral, ZnMn₃O₇·H₂O, which has a similar structure to chalcophanite but contains only tetrahedral Zn, and comparison with chaclophanite strcuture suggest that occurrence of tetrahedal Zn along with octahedral Zn in hexagonal birnessite is more related to the Zn-birnessite stability mediated by interactions of Zn ions (e.g., through H-bonds) with neighboring Zn ions or with adjacent Mn octahedral sheets and stacking modes of the sheets, rather than to the existence of Mn(III) substitution or local Zn bond strength difference.

(a)                                                 (b)                                                  (c)

Figure 3. Red: O; White: H; Teal: Zn; Purple octahedra: Mn(IV)O₆ octahedra. Structure of Zn complexes in (a) chalcophanite (ZnMn₃O₇·3H₂O), viewed perpendicularly to the c-axis. Isolated surface complexes of Zn at a Mn(IV) vacancy of birnessite interlayers: (b) tetrahedral Zn triple-corner-sharing (TCS) complex, with one H₂O directly coordinated to Zn and two extra H₂O H-bonding to the complex; (c) octahedral Zn TCS complex, with three H₂O coordinated to Zn.
 

Nickel partitioning at Mn(IV) vacancies in birnessite

Besides forming TCS surface complexes, Ni is also able to enter a cation vacancy, in effect substituting for the missing Mn (Ni-Sub). The ratio of Ni-Sub/Ni-TCS is observed to increase with pH and time of equilibration, which can be attributed to a pH-dependent activation energy for the TCS → Sub transformation, but this hypothesis is challenging to test by experiment. Our DFT study of the reaction path (Fig.4) found that: (1) an energy barrier exists between Ni-TCS and Ni-Sub and (2) the height of the barrier increases with the number of protons bound at the vacancy site (i.e., it is pH-dependent). The activation energy explains the predominance of Ni-TCS species over Ni-Sub in short-term experiments. Note also in Fig. 4 the increasing depth of the Ni-Sub energy minimum relative to Ni-TCS as pH increases.

Figure 4. Reaction path in the Ni-TCS to Ni-Sub transformation at (a) low pH and (b) high pH. The reaction coordinate marked 0.0 corresponds to Ni-TCS while the reaction coordinate marked 1.0 corresponds to Ni-Sub. Low pH and high pH conditions were represented by differing numbers of protons bound at the vacancy site. The energy (ΔE) is expressed relative to the total energy of Ni-TCS (reaction coordinate 0.0).
 

Structure and dynamics of water and solutes near the surfaces of smectite clay minerals

The plane of oxygen atoms on the cleavage surface of a 2:1 clay mineral is called a siloxane surface. This plane is characterized by hexagonal symmetry among its constituent oxygen ions. Associated with the siloxane surface is a roughly hexagonal cavity, formed by the bases (triads of oxygen ions) of six corner-sharing silica tetrahedra. This cavity has a diameter of about 0.26 nm. If there is no structural charge localized near a cavity, it can bind water molecules attracted to the proton in the structural OH group nestled inside it, or it can bind hydrophobic molecules, such as methane, that are attracted to its oxygen atoms in preference to those in more polar water molecules. If there is structural charge localized near a cavity, then interlayer cations and water molecules both are attracted to this charge and may compete. In aqueous solutions that contain several types of cations, competition of cations (such as Na⁺, Ca²⁺, and CaCl⁺ ion pairs in Fig. 5) for the structural charge sites results in the well-known cation exchange capacity of clay minerals. The first few statistical monolayers of water on siloxane surfaces (Fig. 5) are structured by their propensity to solvate adsorbed cations while forming hydrogen bonds with surface O atoms.

Figure 5. Density map of Na⁺ (dark blue), Ca²⁺ (light blue) and Cl^- (yellow) ions during a 200 ps interval of a molecular dynamics simulation of a 1.6 M Na-Ca-Cl brine in a 6 nm pore between two smectite clay particles (only one side of clay particle is shown).  A snapshot of the location of water O atoms near the surface (red spheres) indicates that Na⁺ and Ca²⁺adsorb as outer-sphere surface complexes.
 

Diffusion of solutes in water

The transport of solutes through porous media and their reactivity at nanoparticle surfaces can be kinetically controlled by molecular diffusion in bulk water or near hydrated mineral surfaces. Molecular diffusion in water can also result in significant isotopic fractionation, allowing the inference of geochemical fluxes from observed isotopic distributions. This emerging area of experimental research is being developped by our collaborator, Frank Richter (University of Chicago). Essential to the interpretation of these distributions is molecular-scale information about the magnitude and mechanisms of solute diffusion and diffusive fractionation, which can be obtained from molecular dynamics simulations.

Isotopes of noble gases dissolved in groundwater are widely used as geochemical indicators of hydrologic transport processes and climatic change. Despite the reliance of these applications on a knowlege of the diffusion coefficients of noble gas isotopes in liquid water, very few measurements of these critical parameters have been reported, primarily because of significant analytical difficulties. Molecular dynamics simulations, which are not limited in this way, can then be very useful tools for examining the isotopic mass dependence of noble gas diffusion coefficients.

CO₂-aqueous brine interfaces

Properties of supercritical CO₂-aqueous brine interfaces (such as interfacial tension or interfacial mass transfer rates and equilibria) strongly influence the behavior of CO₂ plumes in geologic media and, more broadly, may determine the success of carbon capture and storage (CCS) as a CO₂ abatement technology.  Molecular dynamics simulations can complement experimental techniques by probing these interfacial properties and their molecular-scale basis over broad ranges of temperature, pressure, and brine chemistry (Fig. 9).

Figure 9. Molecular dynamics simulation of a 2 M Na-Cl aqueous brine in contact with supercritical CO₂.  The figure shows (1) a snapshot of all CO₂ molecules (gray) and of the water molecules near the left CO2-water interface, (2) a density map of Na⁺ (blue), Cl^- (yellow) and water O (red) in the liquid water phase during a 20 ps time interval, and (3) the trajectories of two water molecules, one diffusing into the CO₂ phase (left side of the figure) and one diffusing along the interface (right side of the figure).
 

CO₂-water interfacial tension

Intererfacial tension (IFT) between supercritical CO₂ and the aqueous phase in a reservoir is the principal physical property that determines the storage capacity of the reservoir and the storage pressure above which CO₂ may flow through the sealing cap rock overlying the storage unit. Very limited high-quality experimental IFT data are available, so alternative methods of predicting IFT over the wide range of storage unit pressures (5-45 MPa) and temperatures (30-110 °C) are needed. Thus we are using MD simulation to investigate the interfacial properties of CO₂-H₂O and CO₂-brine mixtures. Simulations of the CO₂-water system have shown that the strength of the short-range part of the CO₂-H₂O interaction controls the adsorption of CO₂ to the fluid-fluid interface, which in turn controls the IFT. Adsorption is apparent in Figure 10, where carbon dioxide density is enhanced adjacent to the water phase. Our simulations also reproduce the experimentally-observed increase in IFT with increasing salinity. 

Figure 10.  (Top) Snapshot of a molecular dynamics simulation cell containing CO₂ (EPM2 model) and H₂O (SPC/E model) molecules equilibrated at 70 °C and 20 MPa. Oxygen, carbon, and hydrogen atoms are represented by red, cyan, and gray spheres respectively. (Bottom) Time-averaged supercritical CO₂ and liquid water density of the same simulation as a function of distance along the cell. The CO₂ density is significantly enhanced near the interface.

 top

Publications

Kwon, K.D., Refson, K., Bone, S., Qiao, R., Yang, W., Liu, Z., Sposito, G., 2011. Magnetic ordering in mackinawite (tetragonal FeS): evidence for strong itinerant spin fluctuations, Phys. Rev. B  83, 064402.

Kwon, K.D., Refson, K., Sposito, G., 2010. Surface complexation of Pb(II) by hexagonal birnessite nanoparticles. Geochim. Cosmochim. Acta 74, 6731-6740.

Kwon, K.D., Sposito, G., 2010. Reactivity of biogenic manganese oxide for metal sequestration and photochemistry: computational solid state physics study. J. Miner. Soc. Korea 23, 161-170.

Peña, J., Kwon, K.D., Refson, K., Bargar, J.R., Sposito, G., 2010. Mechanisms of nickel sorption by a bacteriogenic birnessite, Geochimica et Cosmochimica Acta 74, 3076-3089.

Bourg, I.C., Sposito, G., 2010.  Connecting the molecular scale to the continuum scale for diffusion processes in smectite-rich porous media. Environ. Sci. Technol. 44, 2085-2091.

Bourg, I.C., Richter, F.M., Christensen, J.N., Sposito, G., 2010.  Isotopic mass-dependence of alkali metal cation diffusion coefficients in water. Geochim. Cosmochim. Acta 74, 2249-2256.

Kwon, K.D., Refson, K., Sposito, G., 2009. On the role of Mn(IV) vacancies in the photoreductive dissolution of hexagonal birnessite. Geochim. Cosmochim. Acta 73, 4142-4150.

Kwon, K.D., Refson, K., Sposito, G., 2009. Zinc surface complexes on birnessite: A density functional theory study. Geochim. Cosmochim. Acta 73, 1273-1284.

Kwon, K.D., Refson, K., Sposito, G., 2008. Defect-induced photoconductivity in layered manganese oxides: A density functional theory study. Phys. Rev. Lett. 100, 146601.

Bourg, I.C., Sposito, G., 2008. Isotopic fractionation of noble gases by diffusion in liquid water: Molecular dynamics simulations and hydrologic applications. Geochim. Cosmochim. Acta 72, 2237-2247.

Bourg, I.C., Sposito, G., Bourg, A.C.M., 2008. Modeling the diffusion of Na+ in compacted water-saturated Na-bentonite as a function of pore water ionic strength. Appl. Geochem. 23, 3635-3641.

Bourg, I.C., Sposito, G., Bourg, A.C.M., 2007. Modeling cation diffusion in compacted water-saturated sodium bentonite at low ionic strength. Environ. Sci. Technol. 41, 8118-8122.

Bourg, I.C., Sposito, G., Bourg, A.C.M., 2007. Modeling the acid-base surface chemistry of montmorillonite. J. Colloid Interface Sci. 312, 297-310.

Bourg, I.C., Sposito, G., 2007. Molecular dynamics simulations of kinetic isotope fractionation during the diffusion of ionic species in liquid water. Geochim. Cosmochim. Acta 71, 5583-5589.

Sutton, R., Sposito, G., 2006. Molecular simulation of humic substance-Ca-montmorillonite complexes. Geochim. Cosmochim. Acta 70, 3566-3581.

Bourg, I.C., Sposito, G., Bourg, A.C.M., 2006. Tracer diffusion in compacted water-saturated bentonite. Clays Clay Miner. 54, 363-374.

Sutton, R., Sposito, G., 2005. Molecular structure in soil humic substances: The new view. Environ. Sci. Technol. 39, 9009-9015.

Sutton, R., Sposito, G., Diallo, M.S., Schulten, H.R., 2005. Molecular simulation of a model of dissolved organic matter. Environ. Toxicol. Chem. 24, 1902-1911.

Sposito, G., Park, S.-H., Refson, K., 2005. NERSC supercomputers are being used to predict the atomic structures of environmental nanoparticles that play important roles in global geochemical cycles. DOE Greenbook: Needs and Directions in High Performance Computing for the Office of Science, U.S. Department of Energy, Office of Science, June 2005 (PPPL-4090), pp. 30-31.

Park, S.-H., Sposito, G., 2004. Molecular modeling of clay structure and surface chemistry. In S.A. Auerbach, K.A. Carrado, and P.K. Dutta (eds.), Handbook of Layered Materials Science & Technology, Chap. 2. Marcel Dekker, New York.

Refson, K., Park, S.-H., Sposito, G., 2003. Ab initio computational crystallography of 2:1 clay minerals: 1. Pyrophyllite-1Tc. J. Phys. Chem. B 107, 13376-13383.

Park, S.H., Sposito, G., 2003. Do montmorillonite surfaces promote methane hydrate formation?: Monte Carlo and molecular dynamics simulations. J. Phys. Chem. B 107, 2281-2290.

Sutton, R., Sposito, G., 2002. Animated molecular dynamics simulations of hydrated caesium-smectite interlayers. Geochem. Trans. 3(9), 73-80.

Park, S.-H., Sposito, G., 2002. Structure of water adsorbed on a mica surface. Phys. Rev. Lett. 89, 85501.

Park, S.-H., Sposito, G., Sutton, R., Greathouse J., 2001. Density functional theory calculations on the structures of 2:1 clay materials. Earth Sciences Division 2000-2001 Annual Report, Lawrence Berkeley National Laboratory, p. 22.

Sutton, R., Sposito, G., 2001. Molecular simulation of interlayer structure and dynamics in 12.4 Å Cs-smectite hydrates. J. Colloid Interface Sci. 237, 174-184.

Greathouse, J.A., Refson, K., Sposito, G., 2000. Molecular dynamics simulation of water mobility in magnesium-smectite hydrates. J. Am. Chem. Soc. 122, 11459-11464.

Park, S.-H., Sposito, G., 2000. Monte Carlo simulation of total radial distribution functions for interlayer water in Li-, Na-, and K-montmorillonite hydrates. J. Phys. Chem. B 104, 4642-4648.

top