Molecular dynamics (MD) method is used to investigate structural transformation and the loss of intermediate range order in SiO2 glass at very large positive pressures and the modification of SiO2 glass network at very large negative pressures. The nature of molecular vibrations in solid C60 has been studied with tight binding molecular dynamics (TBMD) method. Implementations of simulation algorithms on parallel computers are also discussed.

In a-SiO2 at high pressures, the height of the first sharp diffraction peak in S(q) is considerably diminished and its position shifts to larger wave vectors. At twice the normal density, Si-O bond length increases, Si-O coordination changes from 4 to 6, and O-Si-O band-angle changes from 109° to 90°. This is clearly a tetrahedral to octahedral transformation, which is observed recently by Meade, Hemley, and Mao in their diffraction experiments using synchrotron radiation.

MD simulations of porous silica are carried out in the density range 2.2 - 0.1 g/cm3 Internal surface area, pore surface-to-volume ratio, gyration radius, and fractal dimension are studied as a function of density. Simulations are in good agreement with the experimental results obtained by x-ray scattering. The results reveal a crossover of the structural correlations between short- to intermediate-range (< 8 Å) and fractal- to large-scale-regime (10 ~ 100 Å).

Dispersion and density of states (DOS) of inter- and intra-molecular phonons are calculated for orientationally ordered and disordered solid C60 using the TBMD method. Inter-molecular phonon DOS extends up to 7.6 meV and shows libron peaks at 2.4 meV and 3.7 meV in the orientationally ordered phase. Orientational disorder softens libron modes. Intra-molecular phonons below 70 meV also show significant dispersion. Our results are in good agreement with the recent inelastic neutron scattering experiments.

MD is a numerical approach which involves the solution of Newton's equations for particles in the system. The multiple-time-step (MTS) approach reduces this computation significantly by exploiting the different time scales for short-range and intermediate-range interactions. Using the linked-list scheme, parallel algorithms are designed to implement on the in-house 8-node iPSC/860, a MIMD (multiple instruction multiple data) machine. Our group has also developed a quantum dynamical simulation scheme for Computational Nanoelectronics based on the quantum molecular dynamics (QMD) method. The QMD algorithm has been implemented on the in-house 8,192-node MasPar, a SIMD (single instruction multiple data) architecture.

A scheme of developing ab initio many body potentials based on total energy calculations within density functional theory (DFT) is presented and demonstrated for transition metal alloys. An ab initio interatomic potential for Ni/Cr alloys is constructed with no input from experimental data. Molecular dynamics simulations have been performed to study thermal expansions. The coefficient of thermal expansion (CTE) has been calculated over a wide range of temperature, and good agreement is obtained between theory and experiment.

A new method is presented for the generation of valence-only local hamiltonians, or pseudohamiltonians, within the DFT-LDA framework. At the moment these promise to be very useful tools to overcome the problem of eliminating the core electrons from many QMC calculations.

We have investigated total energy expressions that consist of a term describing the eigenvalue sum and a pair potential term. Such expressions can be used to fit the results of first principles total energy calculations at given structures, and then obtain the total energy of another configuration of atoms avoiding the complexity of further ab initio calculations.

In this work we present a method of fitting APW total energy results to a non-orthogonal tight-binding Hamiltonian from which the sum of the eigenvalues is derived and to a pair potential represented by a 5th order polynomial. We fit total energies for the fee and bec structures and we then obtain the elastic constants Cij in good agreement with both full potential LAPW calculations and experiment. We present results of this method for Ir.

We propose a new scheme to extract “optimal” interatomic potentials starting from a large number of atomic configurations (and their forces) obtained from first-principles calculations. The method appears to be able to overcome the difficulties encountered by traditional fitting approaches when using rich and complex analytical forms, and constitute a step forward towards large-scale simulations of condensed matter systems with a degree of accuracy comparable to that obtained by ab initio methods. A first exploratory application to aluminum is presented.

Parallel computing offers new capabilities for using molecular dynamics (MD) to simulate larger numbers of atoms and longer time scales. In this paper we discuss two methods we have used to implement the embedded atom method (EAM) formalism for molecular dynamics on multiple-instruction/multiple-data (MIMD) parallel computers. The first method (atom-decomposition) is simple and suitable for small numbers of atoms. The second method (force-decomposition) is new and is particularly appropriate for the EAM because all the computations are between pairs of atoms. Both methods have the advantage of not requiring any geometric information about the physical domain being simulated. We present timing results for the two parallel methods on a benchmark EAM problem and briefly indicate how the methods can be used in other kinds of materials MD simulations.

We present a method for the self-consistent determination of inter-atomic force-constants. Using non-local ab-initio pseudopotentials to represent the ion-electron interaction and linear response theory to calculate the self-consistent change in the electron density, we are able to calculate the dynamical matrices at arbitrary points in the Brillouin zone. Diagonalization of the dynamical matrix yields phonon eigenvectors and eigenvalues, and fourier inversion yields the real-space interatomic force-constants. We present numerical results for the phonon-dispersion of a variety of metals.

The reconstruction of high-angle twist grain boundaries on the four densest atomic planes in gold are investigated using the recently developed method of grand-canonical simulated quenching. It is found that the grain boundaries on the two densest planes, (111) and (100), do not reconstruct, while those on the (110) and (113) planes do.

The melting of crystalline silicon and the cooling of liquid silicon are investigated using Molecular Dynamics. Both the Stillinger-Weber (SW) potential and the Tight-Binding Bond Model are used to calculate the forces. The electrical properties are investigated using an empirical pseudopotential method with a plane wave basis. The melting point of the solid is found to be about 2300K. The dependency of this temperature with cell size is investigated. On cooling, there are changes in some of the properties of the liquid: the energy per particle decreases, the diffusion constant decreases, and the low frequency electrical conductivity decreases slightly as the temperature decreases. Between 1180K and 980K the liquid undergoes a transition to a glassy phase. There are large changes in the pair correlation function, the SW three-body energy distribution, the diffusion constant, the density of electron single particle states and the electrical conductivity. All of these changes are consistent with increased tetrahedral bonding.

Within a tight-binding molecular dynamics scheme we investigate pure amorphous silicon (a-Si) obtained by direct quenching from the melt. Using different rates for the cooling process, we demonstrate that both structural and electronic properties of a-Si depend on the sample preparation. Possible size-effects are also investigated using 64- and 216-atom supercells. Finally, we discuss the reliability and transferability of the present scheme for large scale simulations of covalent materials.

We calculate the trends in the Jahn-Teller (JT) parameters of the Si vacancy as a function of the breathing mode coordinate. Using ab-initio local-density-functional theory and Hellman--Feynman forces we find the couplings between the JT modes and the breathing mode to be significant. Since the breathing mode coordinate is very sensitive to pressure we predict a pressure-sensitive JT contribution to the effective electronic correlation energy U. Furthermore we find that the tetragonal JT energy increasingly dominates over the trigonal energy as the breathing coordinate decreases due to a softening of the tetragonal spring constant.

A quantum dynamical simulation method is developed to investigate coupled many electron-phonon systems. Both electron and phonon wave functions are numerically propagated in time. The method is applied to the study of resonant tunneling of an electron through double quantum dots. Phonon-induced electron localization is observed. The space- splitting Schrödinger solver and dynamical-simulated-annealing Poisson solver are implemented on an 8,192-node MP-1 computer from MasPar.

We have employed a two-dimensional computer model to study the effect of volume fraction of second phase constituents on load transfer (stiffness) and strength in brittle short-fiber composites, i.e. composites containing a random distribution of aligned fibers, and brittle particulate composites. We find that the efficiency of load transfer to the second phase consituent increases with volume fraction in particulate composites, while it decreases for short-fiber composites. The strength of brittle particulate composites is found to decrease, while the strength of brittle short-fiber composites marginally increases only at fiber volume fractions equal or greater than 0.25.

Flexible boundary codes for the atomistic simulation of dislocations and other defects have been developed in the past mainly by Sinclair [1], Gehlen et al.[2], and Sinclair et al.[3]. These codes permitted the use of smaller atomic arrays than rigid boundary codes, gave descriptions of core non-linear effects and allowed fair assessments of the Peierls stress for dislocation motion. Green functions (continuum or discrete) or surface traction forces were used to relax the boundary atoms.

A much simpler approach is followed here. Core and mobility effects at the boundary are accounted for by a dipole tensor centered at the dislocation line, whose components constitute six more parameters of the minimization process. Results are presented for [100] dislocations in NiAl. It is shown that, within the limitations of the technique, reliable values of the Peierls stress are obtained.

Large scale molecular dynamics simulations on a massively parallel computer are performed to investigate the mechanical behavior of 2-dimensional materials. A model embedded atom many- body potential is examined, corresponding to “ductile” materials. A parallel MD algorithm is developed to exploit the architecture of the Connection Machine, enabling simulations of > 106atoms. A model spallation experiment is performed on a 2-D triagonal crystal with a well-defined nanocrystalline defect on the spall plane. The process of spallation is modelled as a uniform adiabatic expansion. The spall strength is shown to be proportional to the logarithm of the applied strain rate and a dislocation dynamics model is used to explain the results. Good predictions for the onset of spallation in the computer experiments is found from the simple model. The nanocrystal defect affects the propagation of the shock front and failure is enhanced along the grain boundary.

We report on local and global processes during the mechanical deformation of packings of spheres. We demonstrate that solid friction imposes a counter-rotation interaction to two adjacent grains during the deformation in order to minimize losses. This is equivalent to introduce an anti-ferromagnetic spin-like interaction in the structure; the search for the deformation process is quite similar to that one of finding the groundstate of an anti- ferromagnetic disordered structure with special boundary conditions. We demonstrate experimentally that boundaries may produce leading frustration processes; in this case, they govern most of the deformation mode and cancel mostly the influence of bulk mechanical properties (at least in the case when global “antiferromagnetic” order of the bulk is allowed by the absence of frustration). Frustration due to triangular local structure is described and existence of “anti-ferromagnetic-like” domains is proved; these domains are induced by the adequation between the boundary conditions and the deformation process and the domain wall is a stress-localization band as named in classical soil mechanics.

A method for determining the intrinsic resistance to dislocation emission and glide in atomic models is presented and applied to EAM models of nickel containing a single-ended crack tip. The method is an adaptation of the Peierls approach in which the work done on the glide plane is computed as a function of dislocation position and its derivative is the intrinsic resistance to dislocation motion. The results indicate that there are two parts to the resistance dislocation glide from the crack tip: a distinct barrier to injection of the dislocation from the crack tip and a periodic resistance associated with the lattice. When a hydrogen interstitial is placed on the slip plane near the crack tip the barrier height is reduced but its width is increased.

The cohesion of a grain boundary (GB) is believed to be the controlling factor limiting theductility of high-strength metallic alloys, and particularly W. Intergranular embrittlement isusually associated with segregation of impurities at the GBs. Impurities present in ppmconcentrations can result in a dramatic decrease in plasticity. This paper reviews recent results onboth semi-empirical and first-principles modelling of the energetics and the electronic structuresof impurities on a Σ3 (111) GB in W. Our calculations have shown that impurities, such as N, O,P, S, and Si weaken the intergranular cohesion resulting in “loosening” the GB. The presence ofB and C on the contrary, enhances the interatomic interaction across the GB. The so-called site-competition effect should play an important role affecting impurity distribution in W GBs.Among the impurities analyzed, B in the GB has the lowest energy, and thus would tend todisplace other impurity atoms from the GB. Microalloying with 10-50 ppm B may be an effectiveway of improving tungsten's ductility. These results are important for understanding thefundamental physics of intergranular embrittlement.

First principles theories of alloy phase equilibrium have been successfully used in recent years to compute temperature-composition phase diagrams for solid state phases. One particular approach, originating with the successful phenomenological Ising models to describe the alloy Hamiltonian, uses a cluster expansion of the configurational energy in terms of short-ranged pair and many-body interactions. The approach is deeply rooted in our ability to compute accurate total energies of relatively complex compounds, using density functional theory in the local approximation, from which effective interactions may be obtained. Fundamental aspects of the method which include convergence of the cluster expansion, treatment of the configurational entropy and description of vibrational modes are reviewed. Applications of the theory are given for binary alloys in the Ru-Zr-Nb system using the Linear Muffin Tin Orbital method for the total energy calculations, the Cluster Variation method for the description of the configurational entropy, and the Debye-Gruneisen approximation for the vibrational modes. The results are used to compute the equilibrium phase diagram for the Zr-Nb system and to assess current experimental data on phase stability in the Ru-Nb system. In the latter case, the calculations indicate that the DO19 structure is a likely candidate structure for the experimentally observed hexagonal-based compound Ru3Nb. Investigation of the energies and interactions of tetragonal structures as a function of the c/a ratio suggest the L10 structure as a likely candidate for the observed tetragonal phase near 1:1 stoichiometry.

The Monte Carlo technique was used to calculate the phase diagram of the fee superstructures in the Pd-V system using up to fourth nearest neighbor, concentration independent, pair and multiplet interactions. The interactions were computed by the method of Direct Configurational Averaging using a Linearized Muffin-Tin Orbital Hamiltonian cast into the tight binding form. The phase diagram was computed with a fast Monte Carlo simulation technique using environment sampling. The two fee ground states in experimental diagram are predicted to be stable. The computed transition temperatures are higher than those found experimentally.