PaperlayerView exact revision on arXiv

Source-linked AI summary

SIESTA: recent developments and applications

Alberto García, Nick Papior, Arsalan Akhtar, Emilio Artacho, Volker Blum, Emanuele Bosoni, Pedro Brandimarte, Mads Brandbyge, J. I. Cerdá, Fabiano Corsetti, Ramón Cuadrado, Vladimir Dikan, Jaime Ferrer, Julian Gale, Pablo García-Fernández, V. M. García-Suárez, Sandra García, Georg Huhs, Sergio Illera, Richard Korytár, Peter Koval, Irina Lebedeva, Lin Lin, Pablo López-Tarifa, Sara G. Mayo, Stephan Mohr, Pablo Ordejón, Andrei Postnikov, Yann Pouillon, Miguel Pruneda, Roberto Robles, Daniel Sánchez-Portal, Jose M. Soler, Rafi Ullah, Victor Wen-zhe Yu, Javier Junquera

arXiv:2006.01270v1physics.comp-phcond-mat.mtrl-sci

TL;DR

SIESTA is reviewed to document its current status and recent improvements in capabilities, performance, ease of use, and interoperability. The review describes core and newly implemented features, illustrates applications including electronic stopping, and notes that localized atomic orbitals enable reduced-scaling methods while basis-set convergence remains non-uniform.

  • Problem

    The paper addresses the need to assess SIESTA’s present status and document steps improving its capabilities, performance, ease of use, and visibility.

  • Method

    The review synthesizes SIESTA’s localized atomic-orbital approach, recent electronic-structure and performance enhancements, interoperability efforts, and application examples.

  • Results

    The review reports applications spanning diverse systems and capabilities, including successful first-principles electronic-stopping simulations for protons and antiprotons in LiF.

  • Takeaways & Limitations

    Strictly localized orbitals reduce the number of basis functions and produce sparsity, enabling high performance and reduced-scaling methods in SIESTA.

  • Takeaways & Limitations

    Atomic-orbital basis sets lack a unique convergence procedure, although practical accuracy comparable to well-converged plane-wave calculations can be attained.

Abstract

from arXiv · show

A review of the present status, recent enhancements, and applicability of the SIESTA program is presented. Since its debut in the mid-nineties, SIESTA's flexibility, efficiency and free distribution has given advanced materials simulation capabilities to many groups worldwide. The core methodological scheme of SIESTA combines finite-support pseudo-atomic orbitals as basis sets, norm-conserving pseudopotentials, and a real-space grid for the representation of charge density and potentials and the computation of their associated matrix elements. Here we describe the more recent implementations on top of that core scheme, which include: full spin-orbit interaction, non-repeated and multiple-contact ballistic electron transport, DFT+U and hybrid functionals, time-dependent DFT, novel reduced-scaling solvers, density-functional perturbation theory, efficient Van der Waals non-local density functionals, and enhanced molecular-dynamics options. In addition, a substantial effort has been made in enhancing interoperability and interfacing with other codes and utilities, such as Wannier90 and the second-principles modelling it can be used for, an AiiDA plugin for workflow automatization, interface to Lua for steering SIESTA runs, and various postprocessing utilities. SIESTA has also been engaged in the Electronic Structure Library effort from its inception, which has allowed the sharing of various low level libraries, as well as data standards and support for them, in particular the PSML definition and library for transferable pseudopotentials, and the interface to the ELSI library of solvers. Code sharing is made easier by the new open-source licensing model of the program. This review also presents examples of application of the capabilities of the code, as well as a view of on-going and future developments.

I. INTRODUCTION.

SIESTA combines strictly localized atomic-orbital bases, norm-conserving pseudopotentials, real-space grids, and selectable solvers to provide flexible electronic-structure simulations. This review surveys recent methodological, performance, interoperability, and licensing improvements while documenting the method’s strengths and boundaries.

  • A. Theory background and context: Its localized basis enables reduced-scaling algorithms and tunable accuracy and cost, from exploratory calculations to simulations comparable with plane-wave methods.The method offers solver choices ranging from cube-scaling diagonalization to reduced-scaling alternatives.
  • I. INTRODUCTION.: The review presents new core features including DFT+U, spin-orbit interaction, hybrid functionals, improved TD-DFT, DFPT, and new analysis capabilities.It also reports substantial performance improvements in TRANSIESTA and interfaces to new electronic-structure solvers.
  • I. INTRODUCTION.: Recent work expands interoperability through pseudopotential sharing, Wannierization, multiscale interfaces, external libraries, and the Electronic Structure Library ecosystem.The PSML pseudopotential standard, ELSI solver interface, and new open-source licensing model support code and data sharing.
  • A. Theory background and context: SIESTA uses finite-support atomic-like orbitals, norm-conserving pseudopotentials, and a real-space grid to construct sparse Hamiltonian and overlap matrices.The grid represents charge densities and potentials and computes Hamiltonian terms that are not handled as two-center integrals.
  • B. Overview of Siesta capabilities: Typical SIESTA calculations use 10-20 basis orbitals per atom, while empty vacuum regions are essentially free of basis-function cost.This supports systems ranging from dozens to hundreds of atoms on modest hardware, including calculations using cubic-scaling diagonalization.
  • B. Overview of Siesta capabilities: SIESTA retains important scope boundaries: atomic-orbital convergence is non-uniform, pseudopotentials can be disadvantageous for core-electron effects, and periodic Poisson treatment requires repeated supercells for low-dimensional systems.Users remain responsible for choosing basis sets, density functionals, and pseudopotentials, although optimization tools and curated databases are provided.

III. RECENT DEVELOPMENTS IN SIESTA

Recent SIESTA development combines a more open distribution model with interoperability improvements and expanded electronic-structure capabilities. These changes include standardized pseudopotential exchange, new solver interfaces, and methods addressing correlated systems.

  • New distribution model and development infrastructure: SIESTA adopted GPL open-source licensing and more transparent Git-based development, enabling broader contribution and easier code access.
  • Interoperability: PSML formalizes norm-conserving pseudopotential data with metadata and provenance, while libPSML enables interoperability between SIESTA and ABINIT.
  • DFT+U for correlated systems: The LDA+U implementation treats localized d or f states with a Hubbard model while retaining approximate DFT treatment for other valence electrons.
  • DFT+U for correlated systems: 3.08 eV optical band gap versus 1.08 eV with bare GGA-PBE brings NiO close to its 3.10 eV experimental absorption onset.
  • DFT+U for correlated systems: The Ni magnetic moment reaches 1.67 µB, within the experimental range and above the 1.39 µB bare GGA-PBE result.

D. Van der Waals functionals

SIESTA extends its localized-orbital framework with efficient van der Waals functionals and hybrid-functional support. The section also describes Gaussian-based integral evaluation and related electronic-structure implementations.

  • Van der Waals functionals: Van der Waals computation scales as O(N logN) rather than O(N²) through polynomial and Fourier expansions of the nonlocal interaction kernel.
  • Van der Waals functionals: The extended van der Waals kernel depends on densities, density gradients, and interpoint distance, and has been applied to nanotubes, hydrogen adsorption, and liquid water.
  • Hybrid functionals: HSE06 adds screened nonlocal Hartree-Fock exact exchange to semilocal density functionals, with Gaussian-fitted numerical atomic orbitals used to evaluate four-center integrals.
  • Hybrid functionals: Non-negligible electron-repulsion integrals are calculated in the first SCF cycle and stored, preserving sparsity in the Hartree-Fock-exchange Hamiltonian.
  • Hybrid-functional results: HSE06 opens the band gaps of bulk Si and BaTiO3 relative to the semilocal functional.
  • Spin-orbit coupling: Spin-orbit coupling is implemented non-perturbatively with full spinors and can support magnetic crystalline anisotropy calculations.

G. New electronic-structure solvers

SIESTA has broadened its solver ecosystem from conventional diagonalization to reduced-scaling methods through a unified ELSI interface. This integration expands solver access, parallelization, and performance options.

  • SIESTA performance is usually dominated by solving the generalized eigenvalue problem, while finite-support bases make Hamiltonian and overlap construction comparatively lightweight.
  • ELPA provides a drop-in alternative to SCALAPACK, using an additional tridiagonal-conversion step for better scalability and significant speedups.
  • FOE obtains the density matrix from a polynomial expansion of the Fermi–Dirac function without diagonalization.
  • PEXSI uses pole expansion and selected inversion to reduce complexity to at most O(N^2) for dense systems and O(N) for quasi-one-dimensional systems.
  • ELSI unifies access to eigensolvers and density-matrix solvers, including ELPA, OMM, PEXSI, EigenExa, MAGMA, SLEPc, and NTPoly.
  • The SIESTA–ELSI interface adds parallelization over k-points and spins, while newer PEXSI versions use fewer poles and parallelize over chemical-potential trial points.

H. Time dependent DFT

SIESTA’s real-time TD-DFT implementation propagates occupied states in a non-orthogonal, possibly moving basis. A two-stage scheme separates electronic propagation from basis changes and supports parallel execution.

  • Real-time TD-DFT was initially implemented with a Crank–Nicolson propagator for wave-function coefficients.
  • The finite-step propagation uses overlap and Hamiltonian matrices in SIESTA’s non-orthogonal basis, with simplified forms available for sufficiently small time steps and smoothly varying Hamiltonians.
  • For moving atoms, the implementation accounts for changes in both the basis set and the Hilbert space through a two-step procedure involving Löwdin orthonormalization.
  • The Crank–Nicolson algorithm is unitary by construction, preserving orthonormality regardless of the time-step size.
  • The implementation propagates fully occupied states and supports parallelization over k-points and collinear spin.

1. Parallelization

TD-DFT propagation uses dense-matrix infrastructure and inherits scaling similar to conventional diagonalization-based SIESTA runs. In moving-atom simulations, the Löwdin basis-change step is the main computational bottleneck.

  • The two-stage propagation uses matrix operations, inversion, and overlap-matrix diagonalization, with occupied-state coefficients stored in a rectangular matrix.
  • MatrixSwitch abstracts parallel matrix manipulation while invoking BLACS and SCALAPACK dense-matrix libraries.
  • SIESTA’s native compressed-sparse-column matrices can be converted to dense one-dimensional block-cyclic form without MPI communication, whereas two-dimensional redistribution requires communication.
  • The propagation scales similarly to conventional diagonalization-based DFT because both rely on analogous routines in the same dense-matrix library.
  • 18% of runtime was spent in Crank–Nicolson on 30 processors and 25% on 316, while the Löwdin step remained the most expensive operation.
  • For atomic motion, the Löwdin step limits TD-DFT performance; future work targets alternative basis changes, iterative inversion, and improved algorithms.

3. Improved real-time propagators

SIESTA offers progressively more elaborate real-time propagators for systems with strongly varying electronic densities. A two-step half-step Hamiltonian scheme improves energy conservation at an approximately 35% CPU-time increase.

  • The basic propagator is especially suited to small density perturbations, whereas heavily driven systems require schemes with better time-reversal preservation.
  • Propagation accuracy can be increased by dividing each time step into n sub-steps for the underlying first-order expansion.
  • The newest scheme propagates to half-step, explicitly evaluates the half-step Hamiltonian, then completes the full step using that Hamiltonian.
  • ~35% higher CPU time enables better energy conservation when the Kohn–Sham potential varies strongly in time.
  • In double-ionized uracil, the standard implementation fails to handle the highly excited system properly, while extrapolation approaches CPMD-level energy conservation and the two-step scheme improves it further.
  • For smaller time steps, the two-step scheme provides even clearer energy-conservation improvements.

4. Electronic stopping of atomic projectiles

SIESTA’s TD-DFT enabled successful first-principles simulations of electronic stopping for atomic projectiles, while recent DFPT developments extend phonon calculations and improve scaling for sufficiently large systems.

  • Electronic stopping: SIESTA’s TD-DFT enabled the first explicit first-principles simulation of electronic stopping for protons and antiprotons in insulating LiF.The proton–antiproton sign difference produced a stopping-power difference beyond linear response, known as the Barkas effect.
  • Electronic stopping: Electronic stopping simulations were later extended with improved TD-DFT versions, including deformation-density analysis for a proton moving through bulk Ge.Figure 11 shows electron deformation density around a proton traveling leftwards at 1 a.u.
  • Density-functional perturbation theory: The DFPT implementation was merged into SIESTA and computes Γ-point phonons with LDA or GGA at costs comparable to ground-state calculations.Perturbed Hamiltonian and overlap matrix elements follow the ground-state methodology.
  • Density-functional perturbation theory: The alternative DFPT approach targets systems with a gap and changes density with O(N_b^2 M) work after an initialization step scaling as O(N_b^3).The supplied passages describe the initialization as producing perturbation-independent quantities used across atomic displacements.
  • Density-functional perturbation theory: A reduced-scaling DFPT algorithm reuses perturbation-independent quantities to extract the whole dynamical matrix more efficiently for large systems.For carbon fullerenes, a preliminary serial calculation placed the crossover near 650 atoms; parallel initialization can lower this threshold.

J. TranSiesta

TRANSIESTA extends SIESTA to biased open-boundary transport with multiple electrodes, thermoelectric effects, phonon transport, and alternative matrix-inversion algorithms. Its BTD implementation improves scaling and can substantially outperform dense methods for large systems.

  • Transport capabilities: TRANSIESTA uses nonequilibrium Green functions with bulk electrodes to perform biased open-boundary calculations and support Ne ≥ 1 electrodes.The rewritten code also supports thermoelectric calculations, real-space calculations without k-points, and phonon transport.
  • Transport capabilities: Different electrode temperatures and chemical potentials allow self-consistent calculations of thermoelectric effects from separate reservoirs.The spectral function carries electrons from each electrode, while electrode Fermi functions encode reservoir conditions.
  • Algorithms and performance: TRANSIESTA offers BTD, MUMPS sparse, and LAPACK dense inversion algorithms, with BTD providing linear scaling for constant-width systems.The BTD method can outperform dense inversion by a factor of 100 and reduces memory requirements.
  • Electrode representations: Real-space self-energies describe semi-infinite electrodes around defects or line structures without image couplings or k-point sampling.They are presented as superior to Brillouin-zone-integrated quantities for these settings.
  • Transport extensions: TBTRANS supports user-defined tight-binding models, while PH-TRANS and INELASTICA extend the framework to phonon and inelastic electron–phonon transport.INELASTICA treats phonon-excitation effects perturbatively in postprocessing.
  • Algorithms and performance: For the largest pristine-graphene system, BTD achieved more than 40 times the speed of LAPACK in both EGF and NEGF performance characterization.MUMPS gained speed after 5,000 orbitals.

K. Wannierization

SIESTA’s WANNIER90 interface connects Bloch-state calculations to maximally localized Wannier functions and enables efficient electronic-property and multiscale workflows. Wannier-based second-principles methods scale to very large systems but assume invariant bond topology and have limited explicit electronic applications so far.

  • Wannierization: SIESTA can call WANNIER90 as a library or use it for postprocessing to generate maximally localized Wannier functions.The interface makes the unitary transformations between Bloch and Wannier states directly accessible.
  • Wannierization: Wannier functions provide an exact tight-binding representation of Bloch-band dispersion and supply parameters for automatic multiscale simulations.Their construction minimizes a localization functional based on quadratic spreads around Wannier centers.
  • Second-principles multiscale methods: Second-principles DFT uses first-principles models to simulate systems containing tens or hundreds of thousands of atoms through the independent SCALE-UP code.SIESTA supplies model ingredients through MODELMAKER, including force fields, Wannier Hamiltonians, electron–lattice terms, and electron–electron interactions.
  • Second-principles multiscale methods: SPDFT expands the DFT energy around a reference density using a small deformation density, with corrections calculated efficiently as an approximation to the full DFT energy.The zeroth-order term equals the full DFT energy for the reference density, while first- and second-order terms represent band and interaction changes.
  • Second-principles multiscale methods: SPDFT is restricted to processes with invariant bond topology because its deformation density must remain small when bonds are not created or destroyed.This scope boundary is stated for the stationary second-principles problem.
  • Applications and scope: Explicit electronic degrees of freedom have appeared in few SPDFT publications, whereas lattice applications include thermal conductivity and ferroelectric phenomena in perovskites.Reported examples include SrTiO3 and PbTiO3, including strain-dependent domain structures.

M. Scripting and integration in external frameworks

SIESTA exposes its functionality through Lua, AiiDA workflows, and postprocessing utilities, enabling rapid algorithm prototyping, automated simulations, interoperability, and richer analysis of electronic-structure results.

  • Lua interface: Embedded Lua scripts can control SIESTA execution and data without recompiling the host code.Lua is lightweight and interfaces directly with the program’s data structures and routines.
  • Lua interface: The FLOS Lua library implements custom molecular-dynamics modes, geometry relaxations, and optimization schemes that users can reuse or extend.Shared user-level scripts accelerate new-functionality development relative to integrating changes into the core code base.
  • Lua interface: Lua enabled non-intrusive implementations of standard NEB, DNEB, and variable-cell NEB algorithms for transition-state searches.Earlier SIESTA implementations required significant code changes and did not reach the mainstream version.
  • AiiDA integration: The open-source aiida-siesta plugin prepares inputs, runs SIESTA calculations, and supports workflows with provenance, data handling, and sharing.AiiDA workflows include standard materials properties such as band structures.
  • AiiDA integration: An AiiDA workflow can generate a synthetic STM image from a structure, with its execution represented as an automatically generated graph.The passage identifies this as a workflow example rather than a reported performance result.
  • Analysis utilities: SIESTA utilities provide fat bands, spin-texture analysis, non-collinear and spin–orbit wave-function analysis, band unfolding, and real-space visualization.Full unfolding is available even for non-periodic systems such as liquids modeled in large simulation cells.
  • Analysis utilities: SISL postprocesses SIESTA outputs, converts data formats, prepares density-matrix inputs, and builds reduced tight-binding models for large-scale simulations.Hamiltonian manipulation can retain selected band-structure features and calculate far-field currents with little accuracy loss.
  • External interoperability: SISL interoperates with ASE by converting between their geometry objects.SISL was initially developed for SIESTA and TRANSIESTA output but now supports other DFT codes.

3. Other post-processing and visualization utilities

SIESTA's post-processing ecosystem provides visualization and analysis tools for structures, electronic properties, lattice dynamics, phonon spectra, thermodynamic quantities, and optical response. These utilities extend core calculations into interpretable physical observables and experimental comparisons.

  • Visualization: Visualization converters transform SIESTA structure and grid data for XCrySDen and VESTA, including interpolated spatial functions and graphical analyses.The tools support clipping, translations, isosurfaces, contour plots, and Fermi-surface visualization.
  • Lattice dynamics: Lattice-dynamics utilities animate vibration eigenvectors and analyze large-supercell phonons through densities of modes and symmetry or q-vector projections.Projections can distinguish longitudinal and transverse trends, reveal broadened phonon dispersions, and emphasize modes relevant to infrared or Raman spectra.
  • Lattice dynamics: The phdos and vibent tools extract atom- or polarization-resolved mode densities and calculate temperature-dependent vibrational free-energy and entropy contributions.Examples include Be1/3Zn2/3Se, impurity-substituted CZTS, and a Ni4 molecular magnet.
  • Lattice dynamics: The velcf tool obtains phonon frequencies from velocity autocorrelation functions computed from molecular-dynamics histories.A 1000-step simulation at 600 K produced spectra similar to frozen-phonon results for the example shown.
  • Optical properties: Linear-response TDDFT connects induced density to external perturbations through interacting and non-interacting response functions, enabling optical-property calculations from SIESTA orbitals.The PySCF-NAO implementation has been used for optical properties of metallic clusters and plasmon size-dependence.

5. Thermal transport by the AEMD method

SIESTA includes post-processing and pseudopotential-based methods for extracting chemical information from core-level shifts. These methods offer initial- and final-state approximations with different treatments of core-hole relaxation.

  • Core-level shifts: Two schemes are implemented in SIESTA for calculating core-level shifts within a pseudopotential approach.Core-level shifts help analyze changes in the local and chemical environments of atoms and complement experimental information.
  • Initial-state approximation: The initial-state approximation estimates environmental shifts from changes in crystal-potential expectation values using atomic core-state wavefunctions.It neglects electronic relaxation in the presence of the core hole and averages splittings caused by lost spherical symmetry.
  • Final-state approximation: The final-state approximation includes core-hole relaxation by comparing total energies from separate calculations using a core-hole pseudopotential at different sites.A special pseudopotential with a missing core electron is required for each site.
  • Applications: The implemented methodology has been applied to shifts induced by hydrogen bonding in organic molecules.

O. Software-engineering advances and partnerships

SIESTA's software-engineering work addresses maintainability, performance, interoperability, and community development while supporting demanding application areas. The review highlights both practical advances and the complexity introduced by a broader software ecosystem.

  • Software engineering: Reference-counted data structures simplify bookkeeping for richer molecular-mechanics and self-consistent-field iterations.
  • Performance: A mixed MPI/OpenMP programming model improves arithmetic-communication balance and has produced significant speedups for large TranSiesta systems.
  • Partnerships: Participation in MaX and the Electronic Structure Library supports modular, reusable, and exascale-oriented development focused on domain-specific performance.
  • Applications: SIESTA's multi-terminal TranSiesta implementation enabled transport simulations that matched experimental transconductance spectra and assigned observed resonances.The application concerned atomic-scale in-plane transport on reconstructed Ge(001), with collection distances as short as 30 nm.
  • Applications: SIESTA simulations of PbTiO3/SrTiO3 superlattices and low-dimensional materials reproduced vortex structures and identified multiple closely spaced charge-density-wave structures.For 2H-NbSe2, six structures within a few meV were compatible with the experimental 3×3 modulation.

D. Siesta in biology: pilin proteins as conductors

SIESTA has been applied to biological systems, including electrostatics and transport-relevant behavior in pilin proteins under aqueous conditions. The example combines molecular mechanics with first-principles calculations for a large solvated system.

  • Broader scope: The review situates this application within SIESTA's broader use in molecular biology and continuing development toward modular, maintainable, and interoperable workflows.
  • Biological motivation: Pilin's singly oriented alpha helix aligns peptide-bond dipoles along the protein axis, creating a polarization with electrostatic consequences.
  • Method: The study combined long molecular-mechanics simulations in 0.1 M NaCl at pH 7 with SIESTA calculations validating the molecular-mechanics field.The wet system contained 4580 atoms, while a 944-atom vacuum relaxation supported validation.
  • Results: Aqueous surroundings replace pilin's homogeneous vacuum potential ramp with significant long-wavelength fluctuations while preserving a sizeable gap.Coherent transport is not likely, but frontier orbitals evolve suggestively toward enhanced diffusive electron transport.
Loading 2006.01270v1…