Journal Club Theme of February 2009: Finite Element Methods in Quantum Mechanics
Welcome to the February 2009 issue. In this issue, we will discuss the use of finite elements (FEs) in quantum mechanics, with specific focus on the quantum-mechanical problem that arises in crystalline solids. We will consider the electronic structure theory based on the Kohn-Sham equations of density functional theory (KS-DFT): in real-space, Schrödinger and Poisson equations are solved in a parallelepiped unit cell with Bloch-periodic and periodic boundary conditions, respectively. The planewave pseudopotential approach is the method of choice in such quantum-mechanical simulations, but there has been growing interest in recent years on the use of various real-space mesh approaches, of which finite elements are gaining increasing prominence. Most of us are very well-aware of the use of finite elements in solid and structural mechanics applications, but might be less familiar with its place or potential in quantum-mechanical calculations. To bridge this gap, I provide a brief introduction to the topic and sketch the main ingredients, which is supplemented by links to three review articles for further details. Parts of these journal articles would be very accessible to readers of iMechanica who are familiar with FE and other grid-based methods. In the September 15, 2008 Journal Club issue, the computational challenges in electronic-structure calculations have been outlined by Vikram Gavini; please take a look at the nice overview of DFT that he has written.
The solution of the equations of DFT provide a means to determine material properties completely from quantum-mechanical first principles (ab initio), without any experimental input or tunable parameters. This facilitates fundamental understanding and robust predictions for a wide range of properties across the gamut of material systems. The quantum-mechanical problem in a crystalline solid consists of determining the electronic charge density (corresponding wavefunction and effective potential) for a system consisting of positively charged nuclei, surrounded by negatively charged electrons. In the all-electron problem, the Coulomb potential Z/r diverges (solution has cusps and oscillates rapidly near nuclear positions), and hence is not readily amenable to accurate numerical calculations for even moderate system sizes. For first principles computations in crystalline solids, the pseudopotential approach is widely adopted in most quantum molecular dynamics codes. The electrons in the inner shell (core electrons) are tightly bound and do not contribute to any valence bonding. In the pseudopotential approach, the core electrons are frozen in their atomic state, and the divergent Coulomb potential is replaced by a modified potential (pseudopotential) such that the valence states close to the core are less oscillatory, but do not change outside the core region. Only pseudoatomic wavefunctions for the valence (outermost shells) states are solved for in the Schrödinger equation, which makes it numerically tractable via planewaves or with real-space mesh techniques.
On using the pseudopotential approach, the KS-DFT equations consist of the solution of single-particle like Schrödinger equations, which are coupled to a Poisson equation. The single-particle Schrödinger equation (atomic units used) for the ith state is:
subject to boundary conditions consistent with Bloch's theorem. In the above equation, and are the pseudoatomic wavefunction and energy eigenvalue, respectively. The total effective potential now consists of a local ionic part, a non-local ionic part, the electronic (Hartree) potential and the exchange-correlation potential (due to many-body interactions and Pauli's exclusion principle). Typically, to enable predictive capabilities, energy eigenvalues need to be computed within 1 mHa (chemical) accuracy (1 Hartree ~ 27.21 eV). The single-particle pseudowavefunctions are squared sum to form the charge density, which is used in the Poisson equation to solve for the Coulomb potential (local ionic and electronic contributions). Once again is formed and the process is repeated until self-consistency is attained (charge density and effective potential do not change). The total energy, forces, etc., can now be computed to enable quantum molecular dynamics simulations.
The solution of the Poisson equation scales linearly with the number of degrees of freedom N, but the Schrödinger solution varies as the cube of N (see this plot). For self-consistency, since repeated Poisson and Schrödinger equations are solved for and in excess of 1000 eigenfunctions may be required for large systems (~100 atoms or more), the solution of the Schrödinger equation is the limiting step in the solution of the equations of DFT. Unlike linear equations that are relatively easy to solve, accurate eigensolutions for thousands of eigenpairs are computationally demanding and algorithmically challenging, since the eigenfunctions also have to meet the orthogonality constraint.
From Planewaves to Real-Space Mesh Techniques
As the preceding discussion indicates, the efficient solution of the Schrödinger equation in DFT computations is of paramount importance. This need is especially pronounced for metallic systems with heavier atoms and under extreme conditions (high-temperature and/or -pressure) whose pseudopotentials are deep and sharply localized. In such instances, planewaves (Fourier bases) are less than ideal since they have the same resolution everywhere, and hence real-space approaches with their ability to have variable resolution in space (adaptive) become more attractive. More importantly, for large-scale electronic-structure calculations, basis-sets that are compactly-supported such as finite elements or wavelets lead to structured sparse system matrices and hence are more amenable to iterative solution schemes and to implementation on massively parallel computational platforms. As a variational method, finite elements are systematically improvable and convergence of the energy eigenvalues is from above (min-max theorem). Furthermore, boundary conditions (periodic, Bloch, Dirichlet or a combination of these) are easily incorporated within the weak formulation, which make FEs attractive for modeling crystalline solids, molecules, clusters, or surfaces. These attributes make finite elements particularly promising for density functional calculations, and hence the optimism that in the upcoming years there will emerge wider interest in exploring finite element methods in quantum-mechanical computations. I close with a reading list, and a short summary of each review article.
- T. L. Beck (2000), "Real-Space Mesh Techniques in Density-Functional Theory," Reviews of Modern Physics, Vol. 72, Number 4, pp. 1041–1080. [arXiv] [Journal]
Beck's review discusses finite-difference and finite element formulations in DFT. Emphasis is placed on Poisson and nonlinear Poisson-Boltzmann equations in electrostatics and on solutions of Hartree-Fock and KS-equations (eigenvalue problems) of DFT. First, Beck provides a context for real-space mesh techniques by describing some of the prominent developments so far on electronic structure methods (planewaves, Gaussian, LCAO, etc.). Then, the theory behind KS-DFT and classical electrostatics is described. Second-order FD is briefly discussed, and then a higher-order finite-difference method is used to solve the Poisson equation with a singular charge density. The finite element formulation for the Poisson equation is presented. Multigrid solvers are attractive for real-space methods and the main feature of a multigrid method are discussed. Beck presents energy minimization techniques to solve the Schrödinger eigenproblem, and provides a historical perspective on the developments in finite-difference and finite element methods for self-consistent calculations. Lastly, attention is given to time-dependent DFT calculations.
- J. E. Pask and P. A. Sterne (2005), "Finite Element Methods in Ab Initio Electronic Structure Calculations," Modelling and Simulation in Materials Science and Engineering, Vol. 13, pp. R71–R96. [Journal]
A PDF of this paper is uploaded (courtesy of Pask). Pask and Sterne review finite element bases and their use in the self-consistent solution of the Kohn-Sham equations of DFT. The solution of the Schrödinger and Poisson equations is discussed, with particular attention to the imposition of the required Bloch-periodic and periodic boundary conditions, respectively. The use of these solutions in the self-consistent solution of the Kohn-Sham equations and computation of the DFT total energy is then discussed, and applications are given. To impose Bloch-periodic boundary conditions, the wavefunction is rewritten in terms of a periodic function u(x), and the variational (weak) form for the Schrödinger equation in terms of u(x)is derived. The weak form within the unit cell and expressions for the contributions of the non-local part of the pseudopotential are presented. In general, the trial and test functions can now be complex-valued functions, and the overlap matrix S and the local part of the Hamiltonian H are Hermitian. Band structure for a Si pseudopotential is presented and the optimal sextic convergence in energy for cubic finite elements is demonstrated. In the context of crystalline calculations, particular attention is given to the handling of the long-range Coulomb interaction via use of neutral charge densities in the Poisson solution. Self-consistent finite element solutions for the band structure of GaAs are shown, and uniform convergence with mesh refinement to the exact self-consistent solution is demonstrated.
- T. Torsti et al. (2006), "Three Real-Space Discretization Techniques in Electronic Structure Calculations,"Physica Status Solidi. B, Basic Research, Vol. 243, Number 5, pp. 1016–1053. [arXiv] [Journal]
In this paper, Torsti and co-workers compare and contrast the performance of finite-differences, finite elements, and wavelets in electronic-structure calculations. The computational problems under consideration to be solved are the single-particle Schrödinger equation for the pseudowavefunction, and the Poisson equation for the Coulomb potential. Basic theory on FD, FE, and wavelets is first presented; hierarchical finite element bases are also touched upon. Solution approaches for linear equation solvers and for generalized eigenproblem solvers are discussed in significant detail. Applications of finite-differences to quantum dots, surface nanostructures and positron calculations are presented, whereas finite element solutions for all-electron calculations for molecules are performed. Finally, the similarities and differences between the different methods of discretization are nicely summarized.