What are the contributions mentioned in the paper "Octopus: a tool for the application of time-dependent density functional theory" ?

In this paper, the European Theoretical Spectroscopy Facility ( ETSF ) has been used for the analysis of the properties of nano-structures and bio-Systems at Surfaces.

(Open Access) Octopus: a tool for the application of time-dependent density functional theory (2006) | Andrea Castro

phys. stat. sol. (b) 243, No. 11, 2465–2488 (2006) / DOI 10.1002/pssb.200642067

Feature

Article

Feature Article

octopus: a tool for the application of time-dependent density

functional theory

Alberto Castro

*, 1, 2

, Heiko Appel

1, 2

, Micael Oliveira

2, 3

, Carlo A. Rozzi

2, 4

, Xavier Andrade

2, 5

Florian Lorenzen

, M. A. L. Marques

2, 3

, E. K. U. Gross

1, 2

, and Angel Rubio

2, 6

Institut für Theoretische Physik, Fachbereich Physik, Freie Universität Berlin, 14195 Berlin, Germany

European Theoretical Spectroscopy Facility (ETSF)

Departamento de Física, Universidade de Coimbra, Rua Larga, 3004-516 Coimbra, Portugal

CNR-INFM National Research Center on nanoStructures and bioSystems at Surfaces (S3),

Via Campi 213A, 41100 Modena, Italy

Laboratoire des Solides Irradiés, École Polytechnique, 91128 Palaiseau, France

Departamento de Física de Materiales, Facultad de Ciencias Químicas (UPV/EHU),

Centro Mixto CSIC-UPV/EHU and Donostia International Physics Center, San Sebastián 20018, Spain

Received 9 February 2006, revised 30 March 2006, accepted 4 April 2006

Published online 9 June 2006

PACS 71.15.–m, 71.15.Mb, 73.21.La, 73.22.–f, 73.50.Fq, 78.20.Bh

We report on the background, current status, and current lines of development of the

octopus project. This

program materializes the main equations of density-functional theory in the ground state, and of time-

dependent density-functional theory for dynamical effects. The focus is nowadays placed on the optical

(i.e. electronic) linear response properties of nanostructures and biomolecules, and on the non-linear re-

sponse to high-intensity fields of finite systems, with particular attention to the coupled ionic-electronic

motion (i.e. photo-chemical processes). In addition, we are currently extending the code to the treatment

of periodic systems (both to one-dimensional chains, two-dimensional slabs, or fully periodic solids),

magnetic properties (ground state properties and excitations), and to the field of quantum-mechanical

transport or “molecular electronics.” In this communication, we concentrate on the development of the

methodology: we review the essential numerical schemes used in the code, and report on the most recent

implementations, with special attention to the introduction of adaptive coordinates, to the extension of our

real-space technique to tackle periodic systems, and on large-scale parallelization. More information on

the code, as well as the code itself, can be found at

http://www.tddft.org/programs/octopus/.

1 Introduction

Both density-functional theory (DFT) [1, 2], and time-dependent density-functional theory (TDDFT)

[3, 4] have enjoyed a steady increase of their popularity ever since they were born, in the sixties and

eighties respectively. The reason is that both theories achieve, for many problems, an unparalleled bal-

Corresponding author: e-mail: alberto@physik.fu-berlin.de, Phone: +49 030 838 53028, Fax: +49 030838 55258

The European Theoretical Spectroscopy Facility (ETSF) is an initiative of the NANOQUANTA EC 6th Framework Network

of Excellence. It is aimed to build a “knowledge center” to deal with problems of spectroscopy and the properties of electronic

excited states in matter, particularly nanostructures, as well as nanoelectronics and the energetics of atomic motion on the na-

nometre scale.

http://www.cmt.york.ac.uk/nanoquanta/, http://www.cmt.york.ac.uk/etsf/etsf_home.htm

2466 A. Castro et al.: octopus: a tool for the application of TDDFT

ance between accuracy and computational cost. Although the scope of applicability of traditional Quan-

tum Chemistry techniques, or of Quantum Monte-Carlo procedures, have also increased in recent years

[5, 6], DFT/TDDFT is still the method of choice for large systems (e.g., molecular systems of biological

interest) undergoing complex processes.

Correspondingly, numerous software packages that solve DFT/TDDFT equations are available [7].

Among them,

octopus [8] is one with special focus on TDDFT. In the present Article, we describe the

current status of this project, and the aims of its developing process. In brief, some of the key aspects that

describe

octopus are:

Target problems:

(i) Linear optical (i.e. electronic) response of molecules or clusters.

(ii) Non-linear response to classical high-intensity electromagnetic fields, taking into account both the

ionic and electronic degrees of freedom.

(iii) Ground-state and excited state electronic properties of systems with lower dimensionality, such as

quantum dots.

(iv) Photo-induced reactions of molecules (e.g., photo-dissociation, photo-isomerization, etc.).

(v) In the immediate future, extension of these procedures to systems that are infinite and periodic in

one or more dimensions (polymers, slabs, nanotubes, solids), and to electronic transport.

Theoretical base:

(i) The underlying theories are DFT and TDDFT. Also, the code may perform dynamics by consider-

ing the classical (i.e. point-particle) approximation for the nuclei. These dynamics may be non-adiabatic,

since the system evolves following the Ehrenfest path.

(ii) Regarding TDDFT, we have implemented two different approaches: On the one hand, the “stan-

dard” TDDFT-based linear-response theory, which provides us with excitation energies and oscillator

strengths for ground-state to excited-state transitions. On the other hand, we have also implemented the

explicit time-propagation of the TDDFT equations, which allows for the use of large external potentials,

well beyond the range of validity of perturbation theory.

Methodology:

(i) As numerical representation, we have chosen to work without a basis set, relying on numerical

meshes. Nevertheless, auxiliary basis sets (plane waves, atomic orbitals) are used when necessary.

Recently, we have added the possibility of working with non-uniform grids, which adapt to the inho-

mogeneity of the problem, and of making use of multigrid techniques to accelerate the calculations. The

adaptive coordinates implementation will be discussed in some detail in Section 4.

(ii) For most calculations, the code relies on the use of pseudopotentials [9]. We currently allow for

two types: Troullier–Martins [10], and Hartwigsen–Goedecker–Hutter [11].

(iii) In addition to being able to treat systems in the standard 3 dimensions, 2D and 1D modes are also

available. These are useful for studying, e.g., the two-dimensional electron gas that characterizes a wide

class of quantum dots.

Technical aspects:

(i) The code has been designed with emphasis on parallel scalability. In consequence, it allows for

multiple task divisions. We will comment on this aspect in Section 5.

(ii) The language of most of the code is Fortran 90 (almost 50.000 lines at present). Other languages,

such as C or Perl, are also used.

(iii) We have struggled to employ only standard and portable tools. The resulting code may run on

virtually any Unix-like platform.

(iv) The package is licensed under the GNU General Public License (GPL). In consequence, it is

available for use, inspection, and modification for anyone, at

http://www.tddft.programs/octopus/.

phys. stat. sol. (b) 243, No. 11 (2006) 2467

Feature

Article

Section 2 summarizes the algorithms used for ground-state calculations, while the next section handles

response properties. Sections 4–6 report some of the recent additions to the package: adaptive coordinates

to numerically represent the problem, support for parallel calculations (a lot of effort is being put onto the

scalability of the computations to large number of processors), and the treatment of periodic systems.

2 Ground-state DFT calculations

2.1 Kohn–Sham equations

Kohn–Sham (KS) DFT [1, 2] provides the ground state one-particle density

of a system of

elec-

trons exposed to an external potential

()v r

, by identifying it with the density of a non-interacting system

of electrons subject to the so-called KS potential,

()v r

. This system, being non-interacting, may be solved

through a set of one-particle equations, the KS equations [1] (atomic units will be used throughout):

() () ( 1... )

iii

hiNϕεϕ==,,,rr

(1)

() | ()|

n ϕ

(2)

The functions

and the real numbers

are the KS orbitals and KS eigenvalues, respectively. The KS

state is the single Slater determinant built from those orbitals. The KS Hamiltonian is given by

KS KS

[] ()hn v=- — + .r

(3)

The first term, the kinetic operator, is approximated in a real-space formulation by a finite difference

formula – details about this will be given in Section 4. The KS potential is usually separated as follows:

KS Hartree xc

() () []() []()vvvnvn=+ + .rr r r

(4)

2.1.1 External potential

The external potential

()v r

is typically the sum of the Coulomb potential generated by each of the nuclei.

In a pseudopotential formulation, this includes both local and non-local components. For an atom

positioned at

, the pseudopotential

is the sum of a local operator

local

and a set of non-local projec-

tors described by atom-centered functions

local

| | () () | ()vv

κκ

αα αα

ϕϕξϕξ·Ò= +·Ò.

rrr r (5)

Note that these projectors are typically well localized in real space, so their action is computationally

feasible and faster than in a plane wave formulation.

The code also allows for other “user-defined” external potentials. For example, one can attempt to

model the solvent environment of a given system with the electrostatic potential generated by a set of

point charges and/or dipoles (e.g. to model a chromophore in its protein environment [12]). This is the

basic principle of the so-called QM/MM techniques [13]. Also, the user may define a model potential

describing a two-dimensional quantum dot, and can specify it simply by writing down its mathematical

function in the input file.

2.1.2 Hartree potential

The second term of (4),

Hartree

[]()vnr

is more time-consuming. There are various ways of obtaining this

potential numerically, and we have investigated and implemented some of them [14]. These are ex-

plained in more detail in Section 2.2.

2468 A. Castro et al.: octopus: a tool for the application of TDDFT

2.1.3 Exchange-correlation potential

The exchange-correlation (xc) potential,

xc 0

[]vn

, is an unknown functional of the density, and has to be

approximated. It is the first functional derivative of the exchange-correlation energy functional:

[]()

δ()

=.r

(6)

We have incorporated in

octopus a wide variety of possible functionals, ranging from the standard local

density approximation (LDA) [15] and generalized gradient approximations [16], to the state-of-the-art

orbital-dependent functionals [17–19]. “Traditional” LDAs and GGAs are easy, since they are explicit

functionals of the density and its gradient. The more recent orbital-dependent functionals, however, are

explicit functionals of the KS orbitals (and so they are implicit functionals of the density through the

orbitals) and require the use of the optimized effective potential method (OEP) [17, 18]. We have im-

plemented these functionals in

octopus. Both the Krieger, Li and Iafrate (KLI) [19] approximation and

the full solution of the OEP equation [20] (still in experimental phase) are available.

We are now extending the set of functionals to cope with current-density functionals [21]. At this

point, we emphasize that along with the

octopus distribution we provide a standard “exchange and corre-

lation library,” written in C. All (TD)DFT codes require an equivalent piece of software, and in our opin-

ion, it would be mutually beneficial to share an open, reliable library. We expect that this may be a first

step towards this goal.

2.1.4 Eigensolvers

Once we know how to construct the real-space representation of the Hamiltonian for a “trial” density

(or, in fact, for a trial set of KS orbitals

from which the density is generated), we are faced with the

problem of solving the Kohn–Sham Eq. (1) for the

lowest lying eigenpairs of this Hamiltonian opera-

tor. In real space this amounts to the solution of an eigenproblem for large sparse matrices. The literature

in this field is abundant [22], and we have tried several schemes in

octopus. The following are available

in the current version of the code: conjugate-gradients based schemes [23], Lanczos-based algo-

rithms [24] and the Jacobi–Davidson procedure [25].

2.1.5 Mixing

We are left with the mixing of the density, which is essential for the convergence of the self-consistent

procedure. For that purpose, we employ some standard techniques. Essentially, one has to build recur-

sively a series of densities

()i

that converges to the solution density

. Each new density is generated

through a prescription of the form:

(1) (1) () (1) ( )

[...]

iiiiis

nGnnn n



++--

=,,,,,

(7)

where

(1)i



is the density obtained from Eq. (2) using the Kohn–Sham orbitals of step

1i +

. The simplest

example of such a prescription is the so-called “linear mixing” [26], for which Eq. (7) takes the form:

(1)

(1) ()

(1 )

nnnαα



=- + .

However, octopus allows for more sophisticated procedures – we refer the

reader to the original references: the generalized Broyden algorithm of Johnson [27], and the “guaranteed

reduction” Pulay algorithm [28].

2.1.6 Spin

All the previous equations were written considering no spin polarization. However,

octopus is also able

to perform calculations using spin-density functional theory, either considering complete spin alignment

throughout the system or not. This latter case requires the use of the generalized local spin-density the-

ory [29]. The wave functions are then described as two-component spinors

(

)

() () ()Φϕϕ=,rrr

where

the components are complex wave functions.

phys. stat. sol. (b) 243, No. 11 (2006) 2469

Feature

Article

Finally we also mention the possibility to perform calculations including external magnetic fields. As

noted above, current-density functionals are being implemented, but it is already possible to perform

calculations including a static uniform magnetic field and using the standard density functionals.

2.2 Hartree potential

In three dimensions, the Hartree potential may be represented in two equivalent forms: as the integral:

Hartree

()

[]() d

vn r

(8)

or as the solution of Poisson’s equation:

Hartree

[]() 4π()vn n—=-.rr

(9)

There are various ways in which these equations may be solved, and we have investigated and imple-

mented some of them [14].

2.2.7 Conjugate gradients

This amounts to solving Eq. (9) via a conjugate gradients algorithm. This poses the problem of the

boundary conditions for

. The standard solution is to obtain the boundary conditions by calculating the

value of

at points around the simulation box by making use of a multipole expansion representation of

the density

: For points outside, the potential is given by

Hartree

(1)

4π 1

() ()

lm lm

lml

vYrQ

•

==-

ÂÂ

d()()

lm lm

QrrYn=,

where

are spherical harmonics. octopus now offers an alternative: we subtract from

a sum of densi-

ties

lm lm

, where

are the multipoles of

, and where

are auxiliary known charge distributions

whose

()lm

-moment is unity, and whose other moments are zero:

lm lm

lml

nn Qn

==-

=- .

ÂÂ

(11)

For a sufficiently large integer

has negligible boundary conditions, so that

Hartree

[]vn

may be calcu-

lated with the usual Laplacian with zero boundary conditions. Since Poisson’s equation is linear,

Hartree Hartree Hartree

[] [] [ ],

lm lm

lml

vnvn Qvn

==-

ÂÂ

(12)

the functions

Hartree

[]

can be obtained exactly (see Ref. [14] for explicit analytical expressions for

and

Hartree

[]

2.2.8 Multigrids

Still in real-space, as a recent addition,

octopus now also allows for the use of the multigrid method [30,

31]. Multigrid is a linear scaling iterative method to solve elliptic problems. The base of this scheme is to

use a group of different grids that have less points than the original grid where the problem is discretized.

(10)

Octopus: a tool for the application of time-dependent density functional theory

Figures

Citations

ABINIT: First-principles approach to material and nanosystem properties

Density functional theory for transition metals and transition metal chemistry

BerkeleyGW: A massively parallel computer package for the calculation of the quasiparticle and optical properties of materials and nanostructures

Progress in Time-Dependent Density-Functional Theory

Time-dependent density-functional theory for molecules and molecular solids

References

Generalized Gradient Approximation Made Simple

Self-Consistent Equations Including Exchange and Correlation Effects

Inhomogeneous Electron Gas

The inhomogeneous electron gas.

Self-interaction correction to density-functional approximations for many-electron systems

Related Papers (5)

Density-Functional Theory for Time-Dependent Systems

Efficient pseudopotentials for plane-wave calculations

Generalized Gradient Approximation Made Simple

Self-Consistent Equations Including Exchange and Correlation Effects

Self-interaction correction to density-functional approximations for many-electron systems

Frequently Asked Questions (1)

Q1. What are the contributions mentioned in the paper "Octopus: a tool for the application of time-dependent density functional theory" ?