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

Geometry

  %Output
    geometry
  %
  

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.

Forces

  %Output
    forces
  %
  

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.

Potential

  %Output
    potential
  %
  

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

• v0 for the local part of the ionic potential
• vc for the classical potential (if it exists)
• vh for the Hartree potential
• vks for the local part of the Kohn-Sham potential,
• vxc for the exchange-correlation potentials.

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

  %Output
    density  
  %
  

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

  %Output
    kinetic_energy_density
  %
  

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

  %Output
    current
  %
  

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

  %Output
    energy_density
  %
  

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.

Wavefunctions

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.

  %Output
    wfs
  %
  

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.

  %Output
    wfs_sqmod
  %
  

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

  %Output
    wfs_fourier
  %
  

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); ].

  %Output
    elf
  %
  

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} ,$$

where

$$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.

  %Output
    elf_basins
  %
  

Bader population analysis

  %Output
    bader
  %
  

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

  %Output
    el_pressure
  %
  

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.

  %Output
    matrix_elements
  %

  %OutputMatrixElements
    momentum
    ks_multipoles
  %
  

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

• momentum
• ang_momentum
• one_body
• two_body
• two_body_exc_k
• ks_multipoles
• dipole Please refer to the variable description of OutputMatrixElements for more details.

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

  %Output
    tpa
  %
  

Density of states

  %Output
    dos
  %
  

Dipole-moment, dipole-moment density and Polarizability density

Octopus can compute ground-state and excited dipoles. The ground-state dipole moment is computed by default at the end of a ground-state calculation, see SCFCalculateDipole. The values are written in the static/info file.

The excited-state dipoles can be obtained using the Casida/Petersilka equation, see the corresponding tutorial for more details.

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.

  %Output
    pol_density
  %
  

Grid-point coordinates

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

  %Output
    mesh_r
  %
  

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 .

  %Output
    xc_density
  %
  

Heat current density

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

  %Output
    heat_current
  %
  

Electron-photon correlation function

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

  %Output
    photon_correlator
  %
  

J_flow

todo: document J_flow option!

Exchange correlation torque

  %Output
    xc_torque
  %
  

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.

current_kpt

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

density_kpt

Outputs the electronic density resolved in momentum space.

eigenval_kpt

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.

multipoles

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.

angular

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

spin

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

populations

(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.

geometry

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 .

dipole_acceleration

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.

laser

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

energy

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

td_occup

(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.

local_mag_moments

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

gauge_field

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.

temperature

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

ftchd

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.

dipole_velocity

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.

eigenvalues

Write the KS eigenvalues.

ionization_channels

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

total_current

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

partial_charges

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

td_kpoint_occup

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

td_floquet

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

n_excited_el

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

coordinates_sep

Writes geometries in a separate file.

velocities_sep

Writes velocities in a separate file.

forces_sep

Writes forces in a separate file.

total_heat_current

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

total_magnetization

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.

photons_q

Writes photons_q in a separate file.

Photoemission spectroscopy (PES) outputs

PES_wfs

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

PES_density

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

PES

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.

  %Output
    occ_matrices
  %
  

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'}|$.

effectiveU

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

magnetization

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.

local_orbitals

Outputs the localized orbitals that form the correlated subspace.

kanamoriU

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

mmb_wfs

Triggers the ModelMB wavefunctions to be output for each state.

mmb_den

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

BerkeleyGW

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.