By using a combination of constant pressure ab-initio molecular dynamics simulations and density-functional linear response theory we analyze the high pressure phases of molecular solid hydrogen. We use a gradient corrected LDA, and, in the case of ab-initio molecular dynamics, a freshly implemented efficient technique for Brillouin zone sampling. An extremely good k-point sampling turns out to be crucial for obtaining the correct ground state. This approach allows us to optimize simultaneously the orientational degrees of freedom, the lattice constants, and the space group. This can be done either by a local optimization technique, or by running molecular dynamics (MD) trajectories. The MD allows for the system to undergo structural transformations spontaneously. In the lower pressure, namely for the broken symmetry phase (BSP or phase II), we find that the best candidates are the "quadrupolar" orthorhombic structures of Pca2_1 and P21/c symmetry. Equation of state, lattice parameters and vibronic frequencies, are in very good agreement with experimental observations (see J. Kohanoff et al.).
We are presently extending our simulations to higher pressures, with the aim of understanding the nature of phase III (H-A) and of the long-sought hypothetical metallic phase.
We have used constant-pressure ab-initio molecular dynamics to simulate the high-pressure induced transformation of graphite into diamond (see S. Scandolo et al.). The two simulation cells adopted contain 48 and 64 carbon atoms, described by a state-of-the-art LDA pseudopotential scheme. Within a few ps simulation time, and under increasing pressure, a spontaneous conversion from the initial hexagonal graphite to a final diamond structure is obtained. We find that the conversion path proceeds trough a preliminary sliding of the graphite planes into an orthorhombic stacking wherefrom an abrupt collapse and buckling of the planes leads to both cubic and hexagonal (Lonsdaleite) forms of diamond in comparable proportions. This result confirms earlier speculations by Wheeler and Lewis [Mat. Res. Bull. 10, 687 (1975)]. The mutual orientation of the initial and final phases - a crucial indicator of the actual conversion path - is identical to that observed in shock-wave experiments.
A new metallic phase of carbon, metastable at terapascal pressure, is also discovered at higher pressures (see S. Scandolo et al.). The simulation shows that, upon fast compression, diamond survives in a metastable state up to about 3 TPa, where it collapses into a six-fold coordinated structure, named SC4, metallic, never considered so far. Although SC4 is not the lowest energy phase at these pressures, it is argued that large activation barriers may hinder the transition of diamond to other phases at lower or comparable pressure.
In order to understand the behavior of diamond in the so-called Diamond Anvil Cell (DAC), widely used in high-pressure research, we are currently studying the more realistic case of a combined hydrostatic and uniaxial compression of diamond (see S. Scandolo et al.). In particular, we have calculated the nonlinear stress-vs-strain curves for tetragonal (001) distorsion of the diamond lattice (see Zhao Ji-Jun et al.). These curves are of fundamental importance in the modeling of diamond-anvil cell devices.
The behavior of solid oxygen in the pressure range between 5-116 GPa is studied by ab-initio simulations, showing a spontaneous phase transformation from the antiferromagnetic insulating delta-O2 phase to a non-magnetic, metallic molecular phase ( S. Serra et al.). The calculated static structure factor of this phase is in excellent agreement with X-ray diffraction data in the metallic zeta-O2 phase above 96 GPa [Y. Akahama et al, Phys. Rev. Lett. 74, 4690 (1995)]. We thus propose that zeta-O2 should be base centered monoclinic with space group C2/m and 4 molecules per cell, suggesting a re-indexing of experimental diffraction peaks. Physical constraints on the intermediate-pressure epsilon-O2 phase are also obtained.
We have simulated the polymerization of crystalline acetylene into polyacetylene and its subsequent transformation into into a-C:H under pressure shock (see M. Bernasconi et al.). We find that crystalline acetylene polymerizes at 15 GPa, with formation of both cis and trans chains. Subsequent compression shows that polyacetylene develops interchain links, causing a gradual saturation of C-C bonds, and ending up at 50 GPa with a-C:H containing about 80 % of spł carbons. The sp˛ to spł conversion is irreversible and is not undone reverting back to zero pressure. The final a-C:H is a wide gap insulator and, unlike the conventionally generated a-C:H, is highly anisotropic keeping some memory of the original polyacetylene chains axis.
Knowledge of the equation of state, phase diagram and other physical properties such as the electric conductivity of Methane, Water and Ammonia at high pressures and temperatures is of great interest in planetary physics. It is in fact known that the interiors of planets like Uranus and Neptune are mainly composed of such molecules, in solar proportions. Pressures and temperatures inside these planets are thought to range from some GPa and about 1000 K, to several hundred GPa and several thousand kelvin. Using our constant-pressure molecular dynamics tool we have first investigated the behavior of methane in these conditions, and found that, contrary to the current interpretation of shock-wave experiments, methane dissociates below 100 GPa into a mixture of hydrocarbons, and it fully separates into hydrogen and carbon only above 300 GPa. This may explain the anomalous abundance of hydrocarbons in the atmosphere of Neptune ( F. Ancilotto et al.).
Simulations on water and ammonia along the planetary isentrope show that they instead behave as fully dissociated ionic, electronically insulating fluid phases, which turn metallic for temperatures exceeding 7000 K for water and 5500 K for ammonia. At lower temperatures, the phase diagrams of water and ammonia exhibit a superprotonic solid phase between the solid and the ionic liquid ( C. Cavazzoni et al.).
By performing constant-pressure deformable-cell ab-initio molecular dynamics simulations we have studied the pressure-induced chemical instability of CO above 5 GPa [Katz et al, J. Phys. Chem. 88, 3176 (1984)]. The simulation shows that, contrary to previous speculations, polymerization proceeds without CO bond dissociation. The resulting polymer consists of a disordered network of small polycarbonyl (CO)_n chains connecting fivefold C4O cycles ( S. Bernard et al.). The computed vibrational spectra and electronic gap agree very well with (and shed light on) very recent experimental data obtained at Livermore.
The behavior of Iron at earth's core conditions is far from being fully understood [see e.g. O.L. Anderson, Science 278, 821 (1997)]. Open issues include the melting temperature at the inner core boundary, whose estimates range from 5000 to 9000 K, and the hypothetical presence of a new solid phase above 200 GPa and 4000 K, whose presence would reconcile most of the experimental data. We are currently undertaking extensive ab-initio simulations of both liquid and solid iron, in order to clarify such issues.
Gray tin (alpha-Sn) surfaces are not very widely studied, perhaps because good single crystals have not been historically available. In recent years, alpha-Sn(111) has been successfully grown on InSb(111), up to a thickness of about 30 monolayers (roughly 100 A) [T. Osaka, et al, Phys. Rev. B 50, 7567 (1994)]. Good alpha-Sn(111) surfaces of films more than four monolayers thick were found to display a 3x3 reconstruction when grown at lower temperatures. However they irreversibly switched to a 2x2 reconstruction upon annealing at high temperatures.
We have carried out (see Z.Y. Lu et al. ('96)) an ab-initio study of alpha-Sn (111), with the aim of predicting and understanding its structure, reconstructions, and electronic states. We consider a variety of structural possibilities, and optimize them by moving atoms according to Hellmann-Feynman forces. Our results indicate that the unreconstructed surface is highly unstable, while a variety of reconstructions compete for the true ground state. Extrapolated trends from diamond to Si to Ge are well borne out, with a 2x1 pi-bonded chain reconstruction prevailing in the absence of adatoms, and a c(4x2) or (2x2) basic adatom-restatom reconstruction otherwise. We thus propose that the observed Sn(111) 2x2 reconstruction could be explained by the adatom-restatom model.
For the observed 3x3 reconstruction, the scenario appears less obvious. We considered two different models: a dimer-adatom-stacking fault (DAS) model and a distorted sqrt{3} x sqrt{3}, inspired by existing Charge-Density-Wave(CDW) distorted surfaces. Both models yield a higher surface energy when compared to the 2x2 adatom-restatom model, the DAS only slightly worse than the other ( Z.Y. Lu et al. ('98b)). Thus the nature of the metastable 3x3, possibly a result of growth kinetics, remains an open issue.
We have also considered the (001) surface of alpha-Sn. RHEED data [W. T. Yuen, Semic. Sci. Technol. 5, 373 (1990)] on alpha-Sn(100) exhibit a variety of reconstructions with periodicities 2x1, p(2x2), and c(4x4), attributed to possible ordering of dimers, in analogy to Si and Ge (100) surfaces. Our calculations ( Z.Y. Lu et al. ('98a)) show that surface dimers indeed form, accompanied by a large energy gain of 0.618 eV/(surface atom) with respect to the ideal surface. As in Si and Ge, the dimer is buckled, but in alpha-Sn the amount of buckling is surprisingly large, 1.0 A, to be compared with 0.4 A (Si) and 0.74 A (Ge). A frozen phonon calculation predicts a corresponding surface dimer rocking mode at 4.8 THz. The surface core level shift was found to be 0.6 eV for the up dimer atom. In the ground state of alpha-Sn(100) we find that dimers tend to order ``antiferromagnetically''. Calculations show that the most favored states with asymmetric buckled dimers are the c(4x2) and p(2x2) antiferro-reconstructions, found to be nearly degenerate.
We have carried out ab-initio electronic structure calculations on the surface of sqrt{3} adsorbates. In particular, we have addressed the issue of metal-insulator instabilities, charge-density-waves (CDWs) or spin-density-waves (SDWs), driven by partly filled surface states and their 2D Fermi surface, and/or by the onset of magnetic instabilities. The focus is both on the newly discovered commensurate CDW transitions in the Pb/Ge(111) and Sn/Ge(111) structures, and on the puzzling semiconducting behavior of the Pb/Ge(111), K/Si(111):B and SiC(0001) surfaces. In all cases, the main factor driving the instability appears to be an extremely narrow surface state band. We have carried out so far preliminary calculations for the Si/Si(111) surface, chosen as our model system, within the gradient corrected local density (LDA+GC) and local spin density (LSD+GC) approximations, with the aim of understanding the possible interplay between 2D Fermi surface and electron correlations in the surface + adsorbate system. Our spin-unrestricted results show that the sqrt{3} paramagnetic surface is unstable towards a commensurate density wave with periodicity 3x3 and magnetization 1/3 ( S. Scandolo et al. ) . We are currently extending the calculations to the Sn/Ge(111) surface.
The role of correlations has been investigated also with model many-body calculations based on the Hubbard-Holstein model , and Refs. G. Santoro et al. ('98) and G. Santoro et al. ('99)
This project starts in 1994, when mature parallel architectures and operating systems, with a good cost/benefit ratio, become available and the High Performance Computing (HPC) community was moving towards parallel programming. The first release (late 1995) of our Car-Parrinello parallel code has been developed on the SISSA IBM SP2 machine. Communications among processor were implemented using the pvm library. The speed-up was good for a low number of processors ( less than 16 ). At the beginning of the 1997 a second release was written and ported on the CRAY T3E parallel architecture. Key features of this version were a communication layer, introduced to isolate the CP routines from the communication library, and an improved scalability obtained with a better data distribution. The introduction of the communication layer has improved the portability of the code, as a matter of fact, the code has been ported with little effort on many architectures. Supported architectures are : IBM SP2, CRAY T3E, FUJITSU VPP700, SGI ORIGIN2000, SUN ULTRA HPC 10000. During the 1998 the data distribution and communications have been revised and optimized. In particular we have reduced the syncronizations and rewritten the fft driver for a new distribution of the data. Now we are revriting the code using FORTRAN 90 to change the programming paradigm, from a procedural programming to a modular programming. This change is necessary to introduce data hiding and data abstraction, techniques that make the code more easy to maintain and upgrade. Moreover a modular programming paradigm is essentials if more people are working on it. From the point of view of the efficiency modular programming enhance the locality of the data being used (small strides), and therefore optimize the use of caches. Performace and Optimizations of the code.