First-principles simulation, meaning density-functional theory calculations with plane waves and pseudopotentials, has become a prized technique in condensed-matter theory. Here I look at the basics of the suject, give a brief review of the theory, examining the strengths and weaknesses of its implementation, and illustrating some of the ways simulators approach problems through a small case study. I also discuss why and how modern software design methods have been used in writing a completely new modular version of the CASTEP code.

Recent developments in density functional theory (DFT) methods applicable to studies of large periodic systems are outlined. During the past three decades, DFT has become an essential part of computational materials science, addressing problems in materials design and processing. The theory allows us to interpret experimental data and to generate property data (such as binding energies of molecules on surfaces) for known materials, and also serves as an aid in the search for and design of novel materials and processes. A number of algorithmic implementations are currently being used, including ultrasoft pseudopotentials, efficient iterative schemes for solving the one-electron DFT equations, and computationally efficient codes for massively parallel computers. The first part of this article provides an overview of plane-wave pseudopotential DFT methods. Their capabilities are subsequently illustrated by examples including the prediction of crystal structures, the study of the compressibility of minerals, and applications to pressure-induced phase transitions. Future theoretical and computational developments are expected to lead to improved accuracy and to treatment of larger systems with a higher computational efficiency.

Iterative diagonalization of the Hamiltonian matrix is required to solve very large electronicstructure problems. Present algorithms are limited in their convergence rates at low wave numbers by stability problems associated with large changes in the Hartree potential, and at high wave numbers with large changes in the kinetic energy. A new method is described which includes the e6'ect of density changes on the potentials and properly scales the changes in kinetic energy. The use of this method has increased the rate of convergence by over an order of magnitude for large problems. I. INTRQDUCTIQNThe motivation for this work has its foundations in attempts to model periodic ce11s containing many atoms while using a plane-wave basis. The number of plane waves required is moderate when modeling a material such as silicon, but is much larger when some of the atoms are first-row or transition-metal atoms. The nonlocal pseudopotentials associated with these atoms are quite deep and require high-wave-number plane waves to describe the resulting wave functions. The dual requirements of simultaneously describing a large cell volume and a high kinetic energy lead to an explosive growth in the number of basis functions. The need for a new algorithm became apparent when simulations had to be abandoned due to computational difhculties related to the 1arge number of plane waves, as well as the charge instabilities associated with large systems.We outline an adaptation to systems which have orthonormality constraints of an iterative conjugate-gradient method known to be eScient in very large minimization problems. The practical use of this method has been limited by the difhculty of determining the exact linear combinations of vectors which minimize the total energy. We describe a simple new technique for minimizing the total energy in a self-consistent manner. II. BACKGROUNDThe method of self-consistent fields commonly used in quantum-mechanical simulation is gracefully nonlinear.In the usual formulation, the wave function of the system to be studied is expressed as a product of single-particle eigenstates. The energy of the system is a function of a set of coeKcients of basis functions. While there exist many stationary points in the energy, for unpolarized electrons there are no false minima. As long as the ground state is not orthogonal to the starting point, the exact ground state is reachable from any starting set of coefticients by following a path of decreasing energy. Since the energy is at least a quartic function of the coe%cients, no direct method of solution exists and approximate solutions must be iterated until they no longer change. To begin the standard procedure, an initial electron-electron potential is determined by some sort of guess and a Hamiltonian matrix is generated. The Harniltonian matrix is diagonalized and the lowest eigenvectors are occupied. An electron density and corresponding electron-electron potential are generated from the solutions, and the process begins again. If one simply uses the outp...

The spin-polarized plane-wave pseudopotential method, based on density-functional theory, has been used to calculate the electronic band structures and the optical absorption spectra of nitrogen-doped and oxygen-deficient anatase TiO2. The calculated results are in good agreement with our experimental measurements. These ab initio calculations reveal that the optical absorption of nitrogen-doped TiO2 in the visible light region is primarily located between 400 and 500 nm, while that of oxygen-deficient TiO2 is mainly above 500 nm. These results have important implications for the understanding and further development of photocatalytic materials that are active under visible light.

We describe and test a novel molecular dynamics method which combines quantum-mechanical embedding and classical force model optimization into a unified scheme free of the boundary region, and the transferability problems which these techniques, taken separately, involve. The scheme is based on the idea of augmenting a unique, simple parametrized force model by incorporating in it, at run time, the quantum-mechanical information necessary to ensure accurate trajectories. The scheme is tested on a number of silicon systems composed of up to approximately 200 000 atoms.

A new method is presented for performing first-principles molecular-dynamics simulations of systems with variable occupancies. We adopt a matrix representation for the one-particle statistical operator Γ, to introduce a "projected" free energy functional G that depends on the Kohn-Sham orbitals only and that is invariant under their unitary transformations. The Liouville equation [Γ,Ĥ] = 0 is always satisfied, guaranteeing a very efficient and stable variational minimization algorithm that can be extended to non-conventional entropic formulations or fictitious thermal distributions.In recent years, the range of problems that can be studied with quantitative accuracy using the methods of computational solid state physics has expanded dramatically. It is now possible to calculate many materials properties with an accuracy that is often comparable to that of experiments. This degree of confidence is based on the fundamental quantum-mechanical treatment offered by density-functional theory (DFT) [1], coupled with the availability of increasingly powerful computers and with the development of algorithms tuned towards optimal performance [2].The application of these methods and techniques to metallic systems has nonetheless encountered several difficulties, that have made progress slower than for the case of semiconductors and insulators. The discontinuous variation of the orbital occupancies across the Brillouin zone (BZ) makes the occupation numbers rather ill-conditioned variables, and the self-consistent solution of the screening problem can suffer from several instabilities. The absence of a gap in the energy spectrum and the requirement of an exact diagonalization for the Hamiltonian matrix everywhere in the BZ (in order to assign the occupation numbers) introduce "slow frequencies" in the evolution of the orbitals towards the ground state and preclude the straightforward extension to metals of algorithms which performed well for insulators. Smearing the Fermi surface with a finite electronic temperature [3] allows for an improved BZ sampling, but only partially alleviates the problems alluded to above.In this Letter, we introduce a new approach which solves many of these problems in a natural way, and which provides a general and efficient framework for obtaining the ground state of a Kohn-Sham Hamiltonian at a finite electronic temperature. The typical context is the Mermin formulation for the Fermi-Dirac statistics [4], but the method also applies when generalized entropic functionals are introduced [5], as it is often the case for metallic systems. Other applications include DFT studies of insulators or semiconductors with thermally excited states [6], and fractional quantum Hall states [7]. The language of ensemble-DFT [8] is used, and an orbitalbased variational algorithm for the minimization to the ground state is developed and implemented. Dramatic improvements are obtained in the convergence of the energies and especially of the Hellmann-Feynman forces.Within ensemble-DFT, the Helmholtz free energy fun...

Total-energy pseudopotential calculations are used to study the imaging process in noncontact atomic-force microscopy on Si͑111͒ surfaces. At the distance of closest approach between the tip and the surface, there is an onset of covalent chemical bonding between the dangling bonds of the tip and the surface. Displacement curves and lateral scans on the surface show that this interaction energy and force are comparable to the macroscopic Van der Waals interaction. However, the covalent interaction completely dominates the force gradients probed in the experiments. Hence, this covalent interaction is responsible for the atomic resolution obtained on reactive surfaces and it should play a role in improving the resolution in other systems. Our results provide a clear understanding of a number of issues such as ͑i͒ the experimental difficulty in achieving stable operation, ͑ii͒ the quality of the images obtained in different experiments and the role of tip preparations and ͑iii͒ recently observed discontinuities in the force gradient curves. ͓S0163-1829͑98͒08336-2͔

When a brittle material is loaded to the limit of its strength, it fails by the nucleation and propagation of a crack(1). The conditions for crack propagation are created by stress concentration in the region of the crack tip and depend on macroscopic parameters such as the geometry and dimensions of the specimen(2). The way the crack propagates, however, is entirely determined by atomic- scale phenomena, because brittle crack tips are atomically sharp and propagate by breaking the variously oriented interatomic bonds, one at a time, at each point of the moving crack front(1,3). The physical interplay of multiple length scales makes brittle fracture a complex 'multi-scale' phenomenon. Several intermediate scales may arise in more complex situations, for example in the presence of microdefects or grain boundaries. The occurrence of various instabilities in crack propagation at very high speeds is well known(1), and significant advances have been made recently in understanding their origin(4,5). Here we investigate low-speed propagation instabilities in silicon using quantum-mechanical hybrid, multi-scale modelling and single-crystal fracture experiments. Our simulations predict a crack- tip reconstruction that makes low-speed crack propagation unstable on the ( 111) cleavage plane, which is conventionally thought of as the most stable cleavage plane. We perform experiments in which this instability is observed at a range of low speeds, using an experimental technique designed for the investigation of fracture under low tensile loads. Further simulations explain why, conversely, at moderately high speeds crack propagation on the (110) cleavage plane becomes unstable and deflects onto ( 111) planes, as previously observed experimentally(6,7)