SIESTA is a density-functional code able to perform efficient electronic structure calculations and ab initio molecular dynamics  simulations of molecules and solids. SIESTA’s efficiency stems from the use of strictly localised basis sets and from the implementation of linear-scaling algorithms which can be applied to suitable systems.

A very important feature of the code is that its accuracy and cost can be tuned in a wide range, from quick exploratory calculations to highly accurate simulations matching the quality of other approaches, such as plane-wave methods. The possibility of treating large systems with some first-principles electronic-structure methods has opened up new opportunities in many disciplines. SIESTA incorporates also a subsystem (TranSIESTA) for obtaining the electronic structure under nonequilibrium conditions, as appropriate for electronic transport.

During the MaX project, a number of new functionalities were added (Full Spin-Orbit interaction, TDDFT, Linear Response, N-terminal TranSIESTA), and new performance improvements achieved (through the use of new solver libraries and by the implementation of a hybrid MPI-OpenMP parallelization scheme).
The program’s flexibility was enhanced by the implementation of a scripting framework based in the Lua language, and by the encapsulation of several libraries (notably GridXC, which includes the efficient implementation of van der Waals functionals) that have been contributed to the Electronic Structure Library, an initiative supported by the ECAM CoE.

SIESTA is now distributed under the terms of the GPL license and has become quite popular, being increasingly used by researchers in geosciences, biology, and engineering (apart from those in its natural field of materials physics and chemistry). Currently there are several thousands users all over the world. The mailinglist reaches over a thousand subscribers, and the paper describing the method1 has received more than 6700 citations so far. Several SIESTA schools have been organised over the years.

One can distinguish two different sections in the code: the electronic-structure solver (basically obtaining the density-matrix from the Hamiltonian and overlap matrices), and the Hamiltonian setup.

  1. The electronic-structure solver was traditionally based on direct diagonalisation with the ScaLaPack library, with cubic operation count and moderate parallel scaling. Recently, we have implemented a new algorithm based on a pole expansion of the Fermi function plus selected inversion of the sparse matrices involved (PEXSI)2. This method allows massively parallel scaling and features a much smaller memory footprint and a reduced complexity: at most quadratic for three-dimensional systems, and simply linear with system size for quasi one-dimensional ones. Within the MaX project we have also integrated within SIESTA the CheSS linear-scaling solver.
  2. The Hamiltonian setup is almost perfectly linear-scaling by construction, due to the finite-range of the basis set. Its parallel efficiency depends on the pattern of distribution of the orbital and real-space-grid data. Thanks to work done in the past few years at BSC on near-optimal distributions on the grid, and within the MaX project on hybrid OpenMP-MPI operation, the parallel scaling of the setup of H is quite good, and this section of the code represents only a small fraction of the computational time for large systems.

Fig. 1 Algorithm scaling of the PEXSI solver in SIESTA compared with the ScaLaPack-based solver, for systems of different effective dimensionality: a DNA strand, a graphene-BN bilayer, and liquid water.

Fig. 2 Parallel strong scaling of SIESTA-PEXSI and the (Scalapack) diagonalization approach when they are applied to a DNA chain and a Graphene-Boron Nitride stack, a prototypes of quasi one-dimensional and two-dimensional systems. “ppp” stands for the number of MPI processes used in each pole computation, and the nearly straight lines represent increasing levels of parallelization over poles. Due to the large memory requirements of the (dense) diagonalization approach, a substantial minimum number of processors is needed, hence the limited domain of the diagonalization curves.

1 J. M. Soler, E. Artacho, J. D. Gale, A. Garcia, J. Junquera, P. Ordejon, and D. Sanchez-Portal, The SIESTA method for ab initio order-n materials simulation,
J. Phys.: Condens. Matter 14, 2745 (2002).

2 L. Lin, A. Garcia, G. Huhs, and C. Yang, Siesta-PEXSI: Massively parallel method for efficient and accurate ab initio materials simulation without matrix diagonalization, J. Phys.: Condens. Matter 26, 305503 (2014).