Which quantity can Octopus compute

Octopus is a software package for density-functional theory (DFT), and time-dependent density functional theory (TDDFT). As such, it can compute many quantities. These quantities can be integrated quantities, scalar fields, vector fields. Some are accessible for equilibrium (ground-state) calculations, whereas some are only computed for time-dependent calculations.

Below, we list the quantities that can be computed by Octopus and provide relevant equations and references as well as link to details on how to compute these quantities using Octopus.

Ionic output


Input example:

Outputs file containing the coordinates of the atoms treated within quantum mechanics. If OutputFormat = xyz, the file is called geometry.xyz; a file crystal.xyz is written with a supercell geometry if the system is periodic; if point charges were defined in the PDB file (see PDBCoordinates), they will be output in the file geometry_classical.xyz. If OutputFormat = xcrysden, a file called geometry.xsf is written.


Input example:

Outputs file forces.xsf containing structure and forces on the atoms as a vector associated with each atom, which can be visualized with XCrySDen.

Equilibrium electronic output

This series of output quantities are obtained by specifying the variable block Output in the input file. Other ones can be obtained for time-dependent calculations using the variable block TDOutput, see below.


Input example:

outputs Kohn-Sham potential, separated by parts. File names are:

For $v_{s}$ and $v_{xc}$, a suffix for spin is added in the spin-polarized case.

It is also possible to output the gradient of the potential using Output = potential_gradient

Electronic density

Input example:

Outputs the electronic density. The output file is called density-, or lr_density- in linear response.

In the polarized case, this is given by $$ n_\sigma(\vec{r}) = \sum_{i=1}^{N_\sigma} \left| \phi_{i\sigma}(\vec{r}) \right|^2 ,$$

whereas in the noncollinear case, it is $$ n_{\sigma\sigma'}(\vec{r}) = \sum_{i=1}^{N_\sigma} \phi_{i\sigma}(\vec{r}) \phi^*_{i\sigma'}(\vec{r}) .$$

In the spinor case, four real components are printed, labelled sp-1 to sp-4 , corresponding to $n_{\uparrow\uparrow}$, $n_{\downarrow\downarrow}$, $\Re{n_{\uparrow\downarrow}}$, and $\Im{n_{\uparrow\downarrow}}$.

Kinetic energy density

Input example:

Outputs kinetic-energy density, defined as:

$$ \tau_\sigma(\vec{r}) = \sum_{i=1}^{N_\sigma} \sum_{\vec{k}} w_{\vec{k}} \left| \vec{\nabla} \phi_{i\sigma}(\vec{r}) \right|^2 .$$

The index $\sigma$ is the spin index for the spin-polarized case, or if you are using spinors. For spin-unpolarized calculations, you get the total kinetic-energy density. The previous expression assumes full or null occupations. If fractional occupation numbers, each term in the sum is weighted by the occupation. Also, if we are working with an infinite system, all k-points are summed up, with their corresponding weights. The files will be called tau-sp1 and tau-sp2 , if the spin-resolved kinetic energy density is produced (runs in spin-polarized and spinors mode), or only tau if the run is in spin-unpolarized mode.

Current density

Input example:

Outputs the total current density, defined as

$$ \vec{j}(\vec{r}) = \sum_{i=1}^{N} \sum_{\vec{k}} w_{\vec{k}} \Big( \Im \left[ \langle i,\vec{k} | \vec{\nabla} + \frac{|e|}{c} \vec{A}(t) | i,\vec{k} \rangle \right] + \Re\left[ \langle i,\vec{k} | e^{-i\frac{|e|}{c}\vec{A}(t)\cdot\vec{\hat{r}}} [\vec{\hat{r}}, \hat{V}_{\rm NL}] e^{i\frac{|e|}{c}\vec{A}(t)\cdot\vec{\hat{r}}} | i,\vec{k} \rangle \right] \Big) ,$$

where $\hat{V}_{\rm NL}$ is the nonlocal part of the Hamiltonian (pseudopotential, Fock operator, etc.), and $\vec{A}(t)$ is a potential vector potential used during the simulation. In order to control if the nonlocal part of the pseudopotential (and other nonlocal operator contributions) are included or not, please refer to the variable CurrentDensity. Note that other contributions are added for mGGA and hybrid functionals.

The output file is called current- . For linear response, the filename is lr_current- .

Energy density

Input example:

Outputs the total energy density to a file called energy_density. This contains contribution from the kinetic energy density, the external potential, the Hartree, and the $xc$ part. Note that the energy density is only correct for LDA and GGA functionals. For more complex functionals, some contributions might be missing, especially for potential-only functionals which are not derived from an energy functionals.


Octopus offers different ways to output the wavefunction. This can be on real space grid, in Fourier space, or looking at the density of each individual orbitals.

Input example:

outputs the wavefunctions on the real-space grid. Which wavefunctions are to be printed is specified by the variable OutputWfsNumber. The output file is called wf- , or lr_wf- in linear response.

Input example:

outputs the modulus squared of the wavefunctions. The output file is called sqm-wf- . For linear response, the filename is sqm_lr_wf- .

Input example:

outputs wavefunctions in Fourier space. This output is experimental! This is only implemented for the ETSF file format output. The file will be called wfs-pw-etsf.nc .

Electron localization function (ELF)

Octopus can compute the electron localization function (ELF). This can be computed both in the static case and in the time-dependent case. The definition of the time-dependent case can be found in Ref.[ Burnus, T. and Marques, M. A. L. and Gross, E. K. U., Time-dependent electron localization function, Phys. Rev. A 71 010501 (2005); ]. Note that the images and videos from this paper can be found here.

Input example:

In the general case, the ELF is defined as:

$$ f_{\rm ELF}(\vec{r},t) = \frac{1}{1+[D_\sigma(\vec{r},t)/D_\sigma^0(\vec{r},t)]^2} ,$$


$$D_\sigma^0(\vec{r},t) = \tau_\sigma^{\rm HEG}[n_\sigma(\vec{r},t)] = \frac{3}{5}(6\pi^2)^{2/3}n_\sigma^{5/3}(\vec{r},t),$$

is the kinetic energy density of a homogeneous electron gas of density $n_\sigma$. Here $\sigma$ refers to the spin index.

The function $D_\sigma$ is given in the general case by

$$ D_\sigma(\vec{r},t) = \tau_\sigma(\vec{r},t) - \frac{1}{4}\frac{[\nabla n_\sigma(\vec{r},t)]^2}{n_\sigma(\vec{r},t)} - \frac{j_\sigma^2(\vec{r},t)}{n_\sigma(\vec{r},t)} ,$$

where the last term is the current term of the time-dependent ELF, that is also present for complex-valued wavefunctions, like in current-carrying ground states. In order to control the use (or not) of the current term in the ELF, please refer to the variable ELFWithCurrentTerm.

The output file is called elf- , or lr_elf- in linear response, in which case the associated function D is also written, as lr_elf_D- . Only in 2D and 3D.

Octopus can also output the basins of attraction of the ELF.

Input example:
The output file is called elf_rs_basins.info . Only in 2D and 3D.

Bader population analysis

Input example:

Outputs Laplacian of the density which shows lone pairs, bonded charge concentrations and regions subject to electrophilic or nucleophilic attack. See RF Bader, Atoms in Molecules: A Quantum Theory (Oxford Univ. Press, Oxford, 1990). The corresponding basins are also outputed as bader_basins .

Electronic pressure

Input example:

It is possible to output electronic pressure, as defined [ Tao, Jianmin and Vignale, Giovanni and Tokatly, I. V., Quantum Stress Focusing in Descriptive Chemistry, Phys. Rev. Lett. 100 206405 (2008); ].

To be more precise, what the code outputs is the function

$$L(\vec{r}) = \frac{1}{2}\Big(1+\frac{\tilde{p}}{\sqrt{1+\tilde{p}^2}} \Big) ,$$

where we defined here $\tilde{p} = (p^{\rm S}+p^{\rm xc})/p^{\rm TF}$, with the independent particle pressure $p^{\rm S}$ given by

$$p^{\rm S}(\vec{r}) = \frac{1}{3}\tau(\vec{r}) - \frac{1}{4} \nabla^2 n(\vec{r}) ,$$

and the exchange-correlation pressure $p^{\rm xc}$ given by

$$p^{\rm xc}(\vec{r}) = n(\vec{r})v_{\rm xc}(\vec{r}) - e_{\rm xc}(\vec{r}),$$

where $e_{\rm xc}(\vec{r})$ is the exchange-correlation energy density.

Please note that the original paper does not include the exchange-correlation pressure in the definition of $L(\vec{r})$. Moreover, this quantity is only valid for spherically-symmetric atoms. At the moment, this is only available for spin-unpolarized calculations and for the Kohn-Sham DFT theory level . Also, if the functional does not have an energy, this won’t work.

Matrix elements

Outputs can output a series of matrix elements from the electronic states. What is output can be controlled by the OutputMatrixElements variable, which needs to be specified in the input file.

Input example:

At the moment, the code support the following options for OutputMatrixElements:

Octopus can also output transition-potential approximation (TPA) matrix elements, using $\vec{q}$-vector specified by MomentumTransfer.

Input example:

Density of states

Input example:
Outputs density of states. See DOSEnergyMax, DOSEnergyMin, DOSEnergyPoints, and DOSGamma.

Dipole-moment density and Polarizability density

Outputs dipole-moment density dipole_density- , or polarizability density alpha_density- in linear response. If ResponseMethod = finite_differences, the hyperpolarizability density beta_density- is also printed.

Input example:

Grid-point coordinates

Outputs values of the coordinates over the grid. Files will be called mesh_r- followed by the direction.

Input example:

Exchange-correlation density

Outputs the xc density, which is defined as the charge density that generates the XC potential, i.e., $-1/4\pi$ times the Laplacian of the XC potential. The files are called nxc .

Input example:

Heat current density

Outputs the total heat current density. The output file is called heat_current- .

Input example:

Electron-photon correlation function

Outputs the electron-photon correlation function. The output file is called photon_correlator .

Input example:


todo: document J_flow option!

Exchange correlation torque

Input example:

Outputs the exchange-correlation torque. This output is only possible for the spinor case and in the 3D case. $$ \vec{t}(\vec{r}) = \vec{m}(\vec{r})\times \vec{B}_{\rm xc}(\vec{r}) ,$$

where $\vec{m}$ is the magnetization vector, and $\mathbf{B}_{\rm xc}$ is the exchange-correlation magnetic field.

Brillouin-zone resolved outputs

It is also possible to output certain quantities resolved in momentum space, for periodic systems. The quantities are obtained on the k-point grid used to sample the Brillouin zone. At the moment, only specific output formats are possible.


Outputs the total current density resolved in momentum space. The output file is called current_kpt- .


Outputs the electronic density resolved in momentum space.


Outputs the eigenvalues resolved in momentum space, with one file for each band.

Time-dependent outputs

This series of output quantities are obtained by specifying the variable block TDOutput in the input files.


Outputs the (electric) multipole moments of the density to the file td.general/multipoles . This is required to, e.g., calculate optical absorption spectra of finite systems. The maximum value of lll can be set with the variable TDMultipoleLmax.


Outputs the orbital angular momentum of the system to td.general/angular , which can be used to calculate circular dichroism.


(Experimental) Outputs the expectation value of the spin, which can be used to calculate magnetic circular dichroism.


(Experimental) Outputs the projection of the time-dependent Kohn-Sham Slater determinant onto the ground state (or approximations to the excited states) to the file td.general/populations . Note that the calculation of populations is expensive in memory and computer time, so it should only be used if it is really needed. See TDExcitedStatesToProject.


If set (and if the atoms are allowed to move), outputs the coordinates, velocities, and forces of the atoms to the the file td.general/coordinates . On by default if MoveIons = yes .


When set, outputs the acceleration of the electronic dipole, calculated from the Ehrenfest theorem, in the file td.general/acceleration . This file can then be processed by the utility oct-harmonic-spectrum in order to obtain the harmonic spectrum.


If set, outputs the laser field to the file td.general/laser . On by default if TDExternalFields is set.


If set, octopus outputs the different components of the energy to the file td.general/energy . Will be zero except for every TDEnergyUpdateIter iterations.


(Experimental) If set, outputs the projections of the time-dependent Kohn-Sham wavefunctions onto the static (zero-time) wavefunctions to the file td.general/projections.XXX . Only use this option if you really need it, as it might be computationally expensive. See TDProjStateStart. The output interval of this quantity is controled by the variable TDOutputComputeInterval In case of states parallelization, all the ground-state states are stored by each task.


If set, outputs the local magnetic moments, integrated in sphere centered around each atom. The radius of the sphere can be set with LocalMagneticMomentsSphereRadius.


If set, outputs the vector gauge field corresponding to a spatially uniform (but time-dependent) external electrical potential. This is only useful in a time-dependent periodic run. On by default if GaugeVectorField is set.


If set, the ionic temperature at each step is printed. On by default if MoveIons = yes.


Write Fourier transform of the electron density to the file ftchd.X, where X depends on the kick (e.g. with sin-shaped perturbation X=sin). This is needed for calculating the dynamic structure factor. In the case that the kick mode is qbessel, the written quantity is integral over density, multiplied by spherical Bessel function times real spherical harmonic. On by default if TDMomentumTransfer is set.


When set, outputs the dipole velocity, calculated from the Ehrenfest theorem, in the file td.general/velocity . This file can then be processed by the utility oct-harmonic-spectrum in order to obtain the harmonic spectrum.


Write the KS eigenvalues.


Write the multiple-ionization channels using the KS orbital densities as proposed in C. Ullrich, Journal of Molecular Structure: THEOCHEM 501, 315 (2000).


Output the total current (average of the current density over the cell).


Bader and Hirshfeld partial charges. The output file is called td.general/partial_charges .


Project propagated Kohn-Sham states to the states at t=0 given in the directory restart_proj (see %RestartOptions). This is an alternative to the option td_occup, with a formating more suitable for k-points and works only in k- and/or state parallelization


Compute non-interacting Floquet bandstructure according to further options: TDFloquetFrequency, TDFloquetSample, TDFloquetDimension. This is done only once per td-run at t=0. works only in k- and/or state parallelization


Output the number of excited electrons, based on the projections of the time evolved wave-functions on the ground-state wave-functions. The output interval of this quantity is controled by the variable TDOutputComputeInterval


Writes geometries in a separate file.


Writes velocities in a separate file.


Writes forces in a separate file.


Output the total heat current (average of the heat current density over the cell).


Writes the total magnetization, where the total magnetization is calculated at the momentum defined by TDMomentumTransfer. This is used to extract the magnon frequency in case of a magnon kick.


Writes photons_q in a separate file.

Photoemission spectroscopy (PES) outputs


Outputs the photoelectron wavefunctions. The file name is pes_wfs- plus the orbital number.


Outputs the photolectron density. Output file is pes_dens- plus spin species if spin-polarized calculation is performed.


Outputs the time-dependent photoelectron spectrum.

Specialized outputs

DFT+U(+V) outputs

The following outputs are only available for DFT+U(+V) calculations.

Occupation matrices

Outputs the occupation matrices of LDA+U.

Input example:

In the case of a periodic system, the matrix elemets of the occupation matrices $n^{I,n,l,\sigma}_{mm'}$ are given by

$$ n^{I,n,l,\sigma}_{mm'} = \sum_{n}\sum_{\mathbf{k}}^{\mathrm{BZ}} w_\mathbf{k}f_{n\mathbf{k}}^\sigma \langle\psi_{n,\mathbf{k}}^{\sigma} |\hat{P}^{I,n,l}_{mm'}|\psi_{n,\mathbf{k}}^{\sigma} \rangle , $$

where $w_{\mathbf{k}}$ is the $\mathbf{k}$-point weight and $f_{n\mathbf{k}}^\sigma$ is the occupation of the Bloch state $|\psi_{n,\mathbf{k}}^{\sigma} \rangle$.
Here, $| \phi^{I,n,l}_{m}\rangle$ are the localized orbitals that form the basis used to describe electron localization and $\hat{P}^{I,n,l}_{mm'}$ is the projector associated with these orbitals, usually defined as $\hat{P}^{I,n,l}_{mm'}=| \phi^{I,n,l}_{m}\rangle \langle \phi^{I,n,l}_{m'}|$.


If one uses the ACBN0 functional, one can output the Hubbard U computed by Octopus. Outputs the value of the effectiveU for each atoms.


Outputs file containing structure and magnetization of the localized subspace on the atoms as a vector associated with each atom, which can be visualized. For the moment, it only works if a +U is added on one type of orbital per atom.


Outputs the localized orbitals that form the correlated subspace.


Outputs the Kanamori interaction parameters U, U', and J. These parameters are not determined self-consistently, but are taken from the occupation matrices and Coulomb integrals comming from a standard +U calculation.

ModelMB outputs


Triggers the ModelMB wavefunctions to be output for each state.


Triggers the ModelMB density matrix to be output for each state, and the particles specified by the DensitytoCalc block. Calculates, and outputs, the reduced density matrix. For the moment the trace is made over the second dimension, and the code is limited to 2D. The idea is to model N particles in 1D as an N-dimensional non-interacting problem, then to trace out N-1 coordinates.

External libraries


Output for a run with BerkeleyGW. See Output::BerkeleyGW for further specification. delta_perturbation: Outputs the “kick”, or time-delta perturbation applied to compute optical response in real time. external_td_potential: Outputs the (scalar) time-dependent potential.