**MaxwellFunctions**
*Section*: Time-Dependent
*Type*: block

This block specifies the shape of a "spatial-dependent function", such as the
envelope needed when using the `MaxwellFunctions` block. Each line in the block
specifies one function. The first element of each line will be a string
that defines the name of the function. The second element specifies which type
of function we are using; in the following we provide an example for each of the
possible types:

*Options*:

**mxf_const_wave**:

`%MaxwellFunctions`

"function-name" | mxf_const_wave | kx | ky | kz | x0 | y0 | z0

%

The function is constant plane wave \( f(x,y,z) = a0 * cos( kx*(x-x0) + ky*(y-y0) + kz*(z-z0) ) \)

**mxf_const_phase**:

`%MaxwellFunctions`

"function-name" | mxf_const_phase | kx | ky | kz | x0 | y0 | z0

%

The function is a constant phase of \( f(x,y,z) = a0 * (kx*x0 + ky*y0 + kz*z0)

**mxf_gaussian_wave**:

`%MaxwellFunctions`

"function-name" | mxf_gaussian_wave | kx | ky | kz | x0 | y0 | z0 | width

%

The function is a Gaussian, \( f(x,y,z) = a0 * \exp( -( kx*(x-x0) + ky*(y-y0) + kz*(z-z0) )^2 / (2 width^2) ) \)

**mxf_cosinoidal_wave**:

`%MaxwellFunctions`

"function-name" | mxf_cosinoidal_wave | kx | ky | kz | x0 | y0 | z0 | width

%

\( f(t) = \cos( \frac{\pi}{2} \frac{kx*(x-x0)+ky*(y-y0)+kz*(z-z0)-2\width}{width} + pi ) \)

If \( | k_x*x+k_y*y+k_z*z - x_0 | > \xi_0 \), then \( f(x,y,z) = 0 \).

**mxf_logistic_wave**:

`%MaxwellFunctions`

"function-name" | mxf_logistic_wave | kx | ky | kz | x0 | y0 | z0 | growth | width

%

The function is a logistic function, \( f(x,y,z) = a0 * 1/(1+exp(growth*(kx*(x-x0)+ky*(y-y0)+kz*(kz*(z-z0))+width/2))) * 1/(1+exp(-growth*(kx*(x-x0)+ky*(y-y0)+kz*(kz*(z-z0))-width/2))) \)

**mxf_trapezoidal_wave**:

`%MaxwellFunctions`

"function-name" | mxf_trapezoidal_wave | kx | ky | kz | x0 | y0 | z0 | growth | width

%

The function is a logistic function, \( f(x,y,z) = a0 * ( ( 1-growth*(k*(r-r0)-width/2)*Theta(k*(r-r0)-width/2) ) * Theta(-(k*(r-r0)+width/2+1/growth)) + (-1+growth*(k*(r-r0)+width/2)*Theta(k*(r-r0)+width/2) ) * Theta(-(k*(r-r0)-width/2+1/growth)) \)

**mxf_from_expr**:

`%MaxwellFunctions`

"function-name" | mxf_from_expr | "expression"

%

The temporal shape of the field is given as an expression (e.g.,`cos(2.0*x-3*y+4*z)`. The letter*x*,*y*,*z*means spatial coordinates, obviously. The expression is used to construct the function*f*that defines the field.

**TDExternalFields**
*Section*: Time-Dependent
*Type*: block

The block `TDExternalFields` describes the type and shape of time-dependent
external perturbations that are applied to the system, in the form
\(f(x,y,z) \cos(\omega t + \phi (t)) g(t)\), where \(f(x,y,z)\) is defined by
by a field type and polarization or a scalar potential, as below; \(\omega\)
is defined by `omega`; \(g(t)\) is defined by
`envelope_function_name`; and \(\phi(t)\) is the (time-dependent) phase from `phase`.

These perturbations are only applied for time-dependent runs. If
you want the value of the perturbation at time zero to be
applied for time-independent runs, use `TimeZero = yes`.

Each line of the block describes an external field; this way you can actually have more
than one laser (*e.g.* a "pump" and a "probe").

There are two ways to specify \(f(x,y,z)\) but both use the same `omega | envelope_function_name [| phase]`
for the time-dependence.
The float `omega` will be the carrier frequency of the
pulse (in energy units). The envelope of the field is a time-dependent function whose definition
must be given in a `TDFunctions` block. `envelope_function_name` is a string (and therefore
it must be surrounded by quotation marks) that must match one of the function names
given in the first column of the `TDFunctions` block.
`phase` is optional and is taken to be zero if not provided, and is also a string specifying
a time-dependent function.

(A) type = `electric field, magnetic field, vector_potential`

For these cases, the syntax is:

`%TDExternalFields
type | nx | ny | nz | omega | envelope_function_name | phase
%`

The

(B) type =

scalar_potential | "spatial_expression" | omega | envelope_function_name | phase

%

The scalar potential is any expression of the spatial coordinates given by the string "spatial_expression", allowing a field beyond the dipole approximation.

A NOTE ON UNITS:

It is very common to describe the strength of a laser field by its intensity, rather than using the electric-field amplitude. In atomic units (or, more precisely, in any Gaussian system of units), the relationship between instantaneous electric field and intensity is: \( I(t) = \frac{c}{8\pi} E^2(t) \).

It is common to read intensities in W/cm\(^2\). The dimensions of intensities are [W]/(L\(^2\)T), where [W] are the dimensions of energy. The relevant conversion factors are:

Hartree / (\(a_0^2\) atomic_time) = \(6.4364086 \times 10^{15} \mathrm{W/cm}^2\)

eV / ( Å\(^2 (\hbar\)/eV) ) = \(2.4341348 \times 10^{12} \mathrm{W/cm}^2\)

If, in atomic units, we set the electric-field amplitude to \(E_0\), then the intensity is:

\( I_0 = 3.51 \times 10^{16} \mathrm{W/cm}^2 (E_0^2) \)

If, working with

\( I_0 = 1.327 \times 10^{13} (E_0^2) \mathrm{W/cm}^2 \)

**electric_field**: The external field is an electric field, the usual case when we want to describe a laser in the length gauge.**magnetic_field**: The external field is a (homogeneous) time-dependent magnetic field.**vector_potential**: The external field is a time-dependent homogeneous vector potential, which may describe a laser field in the velocity gauge.**scalar_potential**: The external field is an arbitrary scalar potential, which may describe an inhomogeneous electrical field.

**TDFreezeDFTUOccupations**
*Section*: Time-Dependent
*Type*: logical
*Default*: no

The occupation matrices than enters in the LDA+U potential
are not evolved during the time evolution.

**TDFreezeHXC**
*Section*: Time-Dependent
*Type*: logical
*Default*: no

The electrons are evolved as independent particles feeling the Hartree and
exchange-correlation potentials from the ground-state electronic configuration.

**TDFreezeOrbitals**
*Section*: Time-Dependent
*Type*: integer
*Default*: 0

(Experimental) You have the possibility of "freezing" a number of orbitals during a time-propagation.
The Hartree and exchange-correlation potential due to these orbitals (which
will be the lowest-energy ones) will be added during the propagation, but the orbitals
will not be propagated.
*Options*:

**sae**: Single-active-electron approximation. This option is only valid for time-dependent calculations (`CalculationMode = td`). Also, the nuclei should not move. The idea is that all orbitals except the last one are frozen. The orbitals are to be read from a previous ground-state calculation. The active orbital is then treated as independent (whether it contains one electron or two) -- although it will feel the Hartree and exchange-correlation potentials from the ground-state electronic configuration.

It is almost equivalent to setting`TDFreezeOrbitals = N-1`, where`N`is the number of orbitals, but not completely.

**TDFreezeU**
*Section*: Time-Dependent
*Type*: logical
*Default*: no

The effective U of LDA+U is not evolved during the time evolution.

**TDFunctions**
*Section*: Time-Dependent
*Type*: block

This block specifies the shape of a "time-dependent function", such as the
envelope needed when using the `TDExternalFields` block. Each line in the block
specifies one function. The first element of each line will be a string
that defines the name of the function. The second element specifies which type
of function we are using; in the following we provide an example for each of the
possible types:

*Options*:

**tdf_cw**:

`%TDFunctions`

"function-name" | tdf_cw | amplitude

%

The function is just a constant of value`amplitude`: \( f(t) \) = amplitude

**tdf_gaussian**:

`%TDFunctions`

"function-name" | tdf_gaussian | amplitude | tau0 | t0

%

The function is a Gaussian, \( f(t) = F_0 \exp( - (t-t_0)^2/(2\tau_0^2) ) \), where \(F_0\) = amplitude.

**tdf_cosinoidal**:

`%TDFunctions`

"function-name" | tdf_cosinoidal | amplitude | tau0 | t0

%

\( f(t) = F_0 \cos( \frac{\pi}{2} \frac{t-2\tau_0-t_0}{\tau0} ) \)

If \( | t - t_0 | > \tau_0 \), then \( f(t) = 0 \).

**tdf_trapezoidal**:

`%TDFunctions`

"function-name" | tdf_trapezoidal | amplitude | tau0 | t0 | tau1

%

This function is a trapezoidal centered around`t0`. The shape is determined by`tau0`and`tau1`. The function ramps linearly for`tau1`time units, stays constant for`tau0`time units, and then decays to zero linearly again for`tau1`time units.

**tdf_from_file**:

`%TDFunctions`

"function-name" | tdf_from_file | "filename"

%

The temporal shape of the function is contained in a file called`filename`. This file should contain three columns: first column is time, second and third column are the real part and the imaginary part of the temporal function*f*(*t*).

**tdf_from_expr**:

`%TDFunctions`

"function-name" | tdf_from_expr | "expression"

%

The temporal shape of the field is given as an expression (e.g.,`cos(2.0*t)`. The letter*t*means time, obviously. The expression is used to construct the function*f*that defines the field.

**TDGlobalForce**
*Section*: Time-Dependent
*Type*: string

If this variable is set, a global time-dependent force will be
applied to the ions in the x direction during a time-dependent
run. This variable defines the base name of the force, that
should be defined in the `TDFunctions` block. This force
does not affect the electrons directly.

**TDScissor**
*Section*: Time-Dependent
*Type*: float
*Default*: 0.0

(experimental) If set, a scissor operator will be applied in the
Hamiltonian, shifting the excitation energies by the amount
specified. By default, it is not applied.

**ABCapHeight**
*Section*: Time-Dependent::Absorbing Boundaries
*Type*: float
*Default*: -0.2 a.u.

When `AbsorbingBoundaries = cap`, this is the height of the imaginary potential.

**ABShape**
*Section*: Time-Dependent::Absorbing Boundaries
*Type*: block

Set the shape of the absorbing boundaries. Here you can set the inner
and outer bounds by setting the block as follows:

`%ABShape
inner | outer | "user-defined"
%`

The optional 3rd column is a user-defined expression for the absorbing boundaries. For example, \(r\) creates a spherical absorbing zone for coordinates with \({\tt inner} < r < {\tt outer}\), and \(z\) creates an absorbing plane. Note, values

**ABWidth**
*Section*: Time-Dependent::Absorbing Boundaries
*Type*: float

Specifies the boundary width. For a finer control over the absorbing boundary
shape use ABShape.

**AbsorbingBoundaries**
*Section*: Time-Dependent::Absorbing Boundaries
*Type*: flag
*Default*: not_absorbing

To improve the quality of the spectra by avoiding the formation of
standing density waves, one can make the boundaries of the simulation
box absorbing and use exterior complex scaling.
*Options*:

**not_absorbing**: Reflecting boundaries.**mask**: Absorbing boundaries with a mask function.**cap**: Absorbing boundaries with a complex absorbing potential.**exterior**: Exterior complex scaling (not yet implemented).

**MaxwellABPMLPower**
*Section*: Time-Dependent::Absorbing Boundaries
*Type*: float
*Default*: 3.5

Exponential of the polynomial profile for the non-physical conductivity of the PML.
Should be between 2 and 4

**MaxwellABPMLReflectionError**
*Section*: Time-Dependent::Absorbing Boundaries
*Type*: float
*Default*: 1.0e-16

Tolerated reflection error for the PML

**MaxwellABWidth**
*Section*: Time-Dependent::Absorbing Boundaries
*Type*: float

Width of the region used to apply the absorbing boundaries. The default value is twice
the derivative order.

**MediumElectricSigma**
*Section*: Time-Dependent::Absorbing Boundaries
*Type*: float
*Default*: 0.

Electric conductivity of the linear medium.

**MediumEpsilonFactor**
*Section*: Time-Dependent::Absorbing Boundaries
*Type*: float
*Default*: 1.0.

Linear medium electric susceptibility.

**MediumMagneticSigma**
*Section*: Time-Dependent::Absorbing Boundaries
*Type*: float
*Default*: 0.

Magnetic conductivity of the linear medium.

**MediumMuFactor**
*Section*: Time-Dependent::Absorbing Boundaries
*Type*: float
*Default*: 1.0

Linear medium magnetic susceptibility.

**MediumWidth**
*Section*: Time-Dependent::Absorbing Boundaries
*Type*: float
*Default*: 0.0 a.u.

Width of the boundary region with medium

**PESMask2PEnlargeFactor**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: float
*Default*: 1.0

Mask two points enlargement factor. Enlarges the mask box by adding two
points at the edges of the box in each direction (x,y,z) at a distance
L=Lb*`PESMask2PEnlargeFactor` where *Lb* is the box size.
This allows to run simulations with an additional void space at a price of
adding few points. The Fourier space associated with the new box is restricted
by the same factor.

Note: needs ` PESMaskPlaneWaveProjection = nfft_map or pnfft_map `.

**PESMaskEnlargeFactor**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: float
*Default*: 1

Mask box enlargement level. Enlarges the mask bounding box by a `PESMaskEnlargeFactor`.
This helps to avoid wavefunction wrapping at the boundaries.

**PESMaskFilterCutOff**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: float
*Default*: -1

In calculation with `PESMaskMode = fullmask_mode` and NFFT, spurious frequencies
may lead to numerical instability of the algorithm. This option gives the possibility
to filter out the unwanted components by setting an energy cut-off.
If `PESMaskFilterCutOff = -1` no filter is applied.

**PESMaskIncludePsiA**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: logical
*Default*: false

Add the contribution of \(\Psi_A\) in the mask region to the photo-electron spectrum.
Literally adds the Fourier components of:
\(\Theta(r-R1) \Psi_A(r)\)
with \(\Theta\) being the Heaviside step function.
With this option PES will contain all the contributions starting from the inner
radius \(R1\). Use this option to improve convergence with respect to the box size
and total simulation time.
Note: Carefully choose \(R1\) in order to avoid contributions from returning electrons.

**PESMaskMode**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: integer
*Default*: mask_mode

PES calculation mode.
*Options*:

**mask_mode**: Mask method.**fullmask_mode**: Full mask method. This includes a back action of the momentum-space states on the interaction region. This enables electrons to come back from the continuum.**passive_mode**: Passive analysis of the wf. Simply analyze the plane-wave components of the wavefunctions on the region*r*>*R1*. This mode employs a step masking function by default.

**PESMaskPlaneWaveProjection**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: integer
*Default*: fft_map

With the mask method, wavefunctions in the continuum are treated as plane waves.
This variable sets how to calculate the plane-wave projection in the buffer
region. We perform discrete Fourier transforms (DFT) in order to approximate
a continuous Fourier transform. The major drawback of this approach is the built-in
periodic boundary condition of DFT. Choosing an appropriate plane-wave projection
for a given simulation in addition to `PESMaskEnlargeFactor` and
`PESMask2PEnlargeFactor`will help to converge the results.

NOTE: depending on the value of `PESMaskMode` `PESMaskPlaneWaveProjection`,
may affect not only performance but also the time evolution of the density.
*Options*:

**fft_out**: FFT filtered in order to keep only outgoing waves. 1D only.**fft_map**: FFT transform.**nfft_map**: Non-equispaced FFT map.**pfft_map**: Use PFFT library.**pnfft_map**: Use PNFFT library.

**PESMaskShape**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: integer
*Default*: m_sin2

The mask function shape.
*Options*:

**m_sin2**: sin2 mask.

**PESMaskSize**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: block

Set the size of the mask function.
Here you can set the inner (R1) and outer (R2) radius by setting
the block as follows:

`%PESMaskSize
R1 | R2 | "user-defined"
%`

The optional 3rd column is a user-defined expression for the mask function. For example,

**PESMaskSpectEnergyMax**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: float
*Default*: maxval(mask%Lk)\(^2/2\)

The maximum energy for the PES spectrum.

**PESMaskSpectEnergyStep**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: float

The PES spectrum energy step.

**PESMaskStartTime**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: float
*Default*: -1.0

The time photoelectrons start to be recorded. In pump-probe simulations, this allows
getting rid of an unwanted ionization signal coming from the pump.
NOTE: This will enforce the mask boundary conditions for all times.

**PES_Flux_ARPES_grid**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: logical

Use a curvilinear momentum space grid that compensates the transformation
used to obtain ARPES. With this choice ARPES data is laid out on a Cartesian
regular grid.
By default true when `PES_Flux_Shape = pln`.

**PES_Flux_AvoidAB**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: logical
*Default*: yes

For `PES_Flux_Shape = cub`, checks whether surface points are inside the
absorbing zone and discards them if set to yes (default).

**PES_Flux_BZones**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: block

The block `PES_Flux_BZones` specifies the number of Brillouin
zones along each periodic dimensions. This value will be used
to generate the mesh in reciprocal space.

%PES_Flux_EnergyGrid

NBZx | NBZy

%

The option can also be specified with the same value for each direction
as `PES_Flux_BZones = NBZ`.
By default `PES_Flux_BZones = 2`; this means that we will have,
on top of the 1st zone, two additional zones at positive and negative
reciprocal lattice translation for each periodic dimension.

**PES_Flux_DeltaK**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: float
*Default*: 0.002

Spacing of the k-mesh in |k| (equidistant).

**PES_Flux_EnergyGrid**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: block

The block `PES_Flux_EnergyGrid` specifies the energy grid to
be used. Only for ` PES_Flux_Shape = pln`.

%PES_Flux_EnergyGrid

Emin | Emax | DeltaE

%

**PES_Flux_Gpoint_Upsample**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: integer

Increase resolution in momentum by adding Gpoints in between adjacent
Kpoints. For each additional Gpoint G an entire Kpoint grid centered at
G is added to the momentum grid.

**PES_Flux_Kmax**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: float
*Default*: 1.0

The maximum value of |k|.

**PES_Flux_Lmax**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: integer
*Default*: 1

Maximum order of the spherical harmonic to be integrated on an equidistant spherical
grid (to be changed to Gauss-Legendre quadrature).

**PES_Flux_Lsize**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: block

For `PES_Flux_Shape = cub` sets the dimensions along each direction. The syntax is:

`%PES_Flux_Lsize
xsize | ysize | zsize
%
`

where

**PES_Flux_Offset**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: block

Shifts the surface for `PES_Flux_Shape = cub`. The syntax is:

`%PES_Flux_Offset
xshift | yshift | zshift
%
`

**PES_Flux_Radius**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: float

The radius of the sphere, if `PES_Flux_Shape == sph`.

**PES_Flux_Shape**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: integer

The shape of the surface.
*Options*:

**cub**: Uses a parallelepiped as surface. All surface points are grid points. Choose the location of the surface with variable`PES_Flux_Lsize`(default for 1D and 2D).**sph**: Constructs a sphere with radius`PES_Flux_Radius`. Points on the sphere are interpolated by trilinear interpolation (default for 3D).**pln**: This option is for periodic systems. Constructs planes perpendicular to the non-periodic dimension symmetrically placed at positive and negative values of`PES_Flux_Lsize`.

**PES_Flux_StepsPhiK**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: integer
*Default*: 90

Number of steps in \(\phi\) (\(0 \le \phi \le 2 \pi\)) for the spherical grid in k.

**PES_Flux_StepsPhiR**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: integer

Number of steps in \(\phi\) (\(0 \le \phi \le 2 \pi\)) for the spherical surface.

**PES_Flux_StepsThetaK**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: integer
*Default*: 45

Number of steps in \(\theta\) (\(0 \le \theta \le \pi\)) for the spherical grid in k.

**PES_Flux_StepsThetaR**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: integer

Number of steps in \(\theta\) (\(0 \le \theta \le \pi\)) for the spherical surface.

**PES_Flux_UseMemory**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: logical

Use memory to tabulate Volkov plane-wave components on the surface.
This option speeds up calculations by precomputing plane-wave phases
on the suface.
By default true when `PES_Flux_Shape = cub`.

**PES_spm_DeltaOmega**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: float

The spacing in frequency domain for the photoelectron spectrum (if `PES_spm_OmegaMax > 0`).
The default is `PES_spm_OmegaMax/500`.

**PES_spm_OmegaMax**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: float
*Default*: 0.0

If `PES_spm_OmegaMax > 0`, the photoelectron spectrum is directly calculated during
time-propagation, evaluated by the PES_spm method. `PES_spm_OmegaMax` is then the maximum frequency
(approximate kinetic energy) and `PES_spm_DeltaOmega` the spacing in frequency domain of the spectrum.

**PES_spm_Radius**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: float

The radius of the sphere for the spherical grid (if no `PES_spm_points`
are given).

**PES_spm_StepsPhiR**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: integer
*Default*: 90

Number of steps in \(\phi\) (\(0 \le \phi \le 2 \pi\)) for the spherical grid (if no
`PES_spm_points` are given).

**PES_spm_StepsThetaR**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: integer
*Default*: 45

Number of steps in \(\theta\) (\(0 \le \theta \le \pi\)) for the spherical grid (if no
`PES_spm_points` are given).

**PES_spm_points**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: block

List of points at which to calculate the photoelectron spectrum by the sample point
method. If no points are given, a spherical grid is generated automatically.
The exact syntax is:

`%PES_spm_points
x1 | y1 | z1
%
`

**PES_spm_recipe**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: integer
*Default*: phase

The type for calculating the photoelectron spectrum in the sample point method.
*Options*:

**raw**: Calculate the photoelectron spectrum according to A. Pohl, P.-G. Reinhard, and E. Suraud,*Phys. Rev. Lett.***84**, 5090 (2000).**phase**: Calculate the photoelectron spectrum by including the Volkov phase (approximately), see P. M. Dinh, P. Romaniello, P.-G. Reinhard, and E. Suraud,*Phys. Rev. A.***87**, 032514 (2013).

**PhotoElectronSpectrum**
*Section*: Time-Dependent::PhotoElectronSpectrum
*Type*: integer
*Default*: none

This variable controls the method used for the calculation of
the photoelectron spectrum. You can specify more than one value
by giving them as a sum, for example:
`PhotoElectronSpectrum = pes_spm + pes_mask`
*Options*:

**none**: The photoelectron spectrum is not calculated. This is the default.**pes_spm**: Store the wavefunctions at specific points in order to calculate the photoelectron spectrum at a point far in the box as proposed in A. Pohl, P.-G. Reinhard, and E. Suraud,*Phys. Rev. Lett.***84**, 5090 (2000).**pes_mask**: Calculate the photo-electron spectrum using the mask method. U. De Giovannini, D. Varsano, M. A. L. Marques, H. Appel, E. K. U. Gross, and A. Rubio,*Phys. Rev. A***85**, 062515 (2012).**pes_flux**: Calculate the photo-electron spectrum using the t-surff technique,*i.e.*, spectra are computed from the electron flux through a surface close to the absorbing boundaries of the box. (Experimental.) L. Tao and A. Scrinzi,*New Journal of Physics***14**, 013021 (2012).

**ArnoldiOrthogonalization**
*Section*: Time-Dependent::Propagation
*Type*: integer

The orthogonalization method used for the Arnoldi procedure.
Only for TDExponentialMethod = lanczos.
*Options*:

**cgs**: Classical Gram-Schmidt (CGS) orthogonalization. The algorithm is defined in Giraud et al., Computers and Mathematics with Applications 50, 1069 (2005).**drcgs**: Classical Gram-Schmidt orthogonalization with double-step reorthogonalization. The algorithm is taken from Giraud et al., Computers and Mathematics with Applications 50, 1069 (2005). According to this reference, this is much more precise than CGS or MGS algorithms.

**IonsConstantVelocity**
*Section*: Time-Dependent::Propagation
*Type*: logical
*Default*: no

(Experimental) If this variable is set to yes, the ions will
move with a constant velocity given by the initial
conditions. They will not be affected by any forces.

**IonsTimeDependentDisplacements**
*Section*: Time-Dependent::Propagation
*Type*: block

(Experimental) This variable allows you to specify a
time-dependent function describing the displacement of the ions
from their equilibrium position: \(r(t) = r_0 + \Delta
r(t)\). Specify the displacements dx(t), dy(t), dz(t) as
follows, for some or all of the atoms:

`%IonsTimeDependentDisplacements
atom_index | "dx(t)" | "dy(t)" | "dz(t)"
%`

The displacement functions are time-dependent functions and should match one of the function names given in the first column of the

**MaxwellAbsorbingBoundaries**
*Section*: Time-Dependent::Propagation
*Type*: block

Type of absorbing boundaries used for Maxwell propagation in each direction.

Example:

`%MaxwellAbsorbingBoundaries
cpml | cpml | cpml
%`

**not_absorbing**: No absorbing boundaries.**mask**: A mask equal to the wavefunctions mask is applied to the Maxwell states at the boundaries**maxwell_mask**: A different mask than the wavefunctions mask is applied on Maxwell states**cpml**: Perfectly matched layer absorbing boundary**mask_zero**: Absorbing boundary region is set to zero

**MaxwellBoundaryConditions**
*Section*: Time-Dependent::Propagation
*Type*: block

Defines boundary conditions for the electromagnetic field propagation.

Example:

`%MaxwellBoundaryConditions
zero | mirror_pec | consant
%`

**zero**: Boundaries are set to zero.**constant**: Boundaries are set to a constant.**mirror_pec**: Perfect electric conductor.**mirror_pmc**: Perfect magnetic conductor.**plane_waves**: Boundaries feed in plane waves.**periodic**: Periodic boundary conditions (not yet implemented).**medium**: Boundaries as linear medium (not yet implemented).

**MaxwellMediumBox**
*Section*: Time-Dependent::Propagation
*Type*: block

Defines parameters for the box.

Example:

`%MaxwellMediumBox
center_x | center_y | center_z | x_length | y_length | z_length | epsilon_factor | mu_factor | sigma_e | sigma_m | edged/smooth
%`

Position of center (three components) and length (three components), followed by permittivity factor, electric conductivity and magnetic conductivity, and finally type of edge.

**edged**: Medium box edges are steep.**smooth**: Medium box edged and softened.

**MaxwellTDETRSApprox**
*Section*: Time-Dependent::Propagation
*Type*: integer
*Default*: no

Whether to perform aproximations to the ETRS propagator.
*Options*:

**no**: No approximations.**const_steps**: Use constant current density.

**MaxwellTDIntervalSteps**
*Section*: Time-Dependent::Propagation
*Type*: integer

This variable determines how many intervall steps the Maxwell field propagation
does until it reaches the matter time step. In case that MaxwellTDIntervalSteps is
equal to one, the Maxwell time step is equal to the matter one. The default value is 1.

**MaxwellTDOperatorMethod**
*Section*: Time-Dependent::Propagation
*Type*: integer
*Default*: op_fd

The Maxwell Operator e.g. the curl operation can be obtained by
two different methods, the finid-difference or the fast fourier
transform.
*Options*:

**op_fd**: Maxwell operator calculated by finite differnce method**op_fft**: Maxwell operator calculated by fast fourier transform

**MaxwellTDSCFThreshold**
*Section*: Time-Dependent::Propagation
*Type*: float
*Default*: 1.0e-6

Since the Maxwell-KS propagator is non-linear, each propagation step
should be performed self-consistently. In practice, for most
purposes this is not necessary, except perhaps in the first
iterations. This variable holds the number of propagation steps
for which the propagation is done self-consistently.

This variable controls the accuracy threshold for the self consistency.

**MoveIons**
*Section*: Time-Dependent::Propagation
*Type*: logical

This variable controls whether atoms are moved during a time
propagation run. The default is yes when the ion velocity is
set explicitly or implicitly, otherwise is no.

**RecalculateGSDuringEvolution**
*Section*: Time-Dependent::Propagation
*Type*: logical
*Default*: no

In order to calculate some information about the system along the
evolution (e.g. projection onto the ground-state KS determinant,
projection of the TDKS spin-orbitals onto the ground-state KS
spin-orbitals), the ground-state KS orbitals are needed. If the
ionic potential changes -- that is, the ions move -- one may want
to recalculate the ground state. You may do this by setting this
variable.

The recalculation is not done every time step, but only every
`RestartWriteInterval` time steps.

**TDDynamics**
*Section*: Time-Dependent::Propagation
*Type*: integer
*Default*: ehrenfest

Type of dynamics to follow during a time propagation.
For BO, you must set `MoveIons = yes`.
*Options*:

**ehrenfest**: Ehrenfest dynamics.**bo**: Born-Oppenheimer (Experimental).

**TDEnergyUpdateIter**
*Section*: Time-Dependent::Propagation
*Type*: integer

This variable controls after how many iterations Octopus
updates the total energy during a time-propagation run. For
iterations where the energy is not updated, the last calculated
value is reported. If you set this variable to 1, the energy
will be calculated in each step.

**TDExpOrder**
*Section*: Time-Dependent::Propagation
*Type*: integer
*Default*: 4

For `TDExponentialMethod` = `standard` or `chebyshev`,
the order to which the exponential is expanded. For the Lanczos approximation,
it is the Lanczos-subspace dimension.

**TDExponentialMethod**
*Section*: Time-Dependent::Propagation
*Type*: integer
*Default*: taylor

Method used to numerically calculate the exponential of the Hamiltonian,
a core part of the full algorithm used to approximate the evolution
operator, specified through the variable `TDPropagator`.
In the case of using the Magnus method, described below, the action of the exponential
of the Magnus operator is also calculated through the algorithm specified
by this variable.
*Options*:

**lanczos**: Allows for larger time-steps. However, the larger the time-step, the longer the computational time per time-step. In certain cases, if the time-step is too large, the code will emit a warning whenever it considers that the evolution may not be properly proceeding -- the Lanczos process did not converge. The method consists in a Krylov subspace approximation of the action of the exponential (see M. Hochbruck and C. Lubich,*SIAM J. Numer. Anal.***34**, 1911 (1997) for details). Two more variables control the performance of the method: the maximum dimension of this subspace (controlled by variable`TDExpOrder`), and the stopping criterion (controlled by variable`TDLanczosTol`). The smaller the stopping criterion, the more precisely the exponential is calculated, but also the larger the dimension of the Arnoldi subspace. If the maximum dimension allowed by`TDExpOrder`is not enough to meet the criterion, the above-mentioned warning is emitted.**taylor**: This method amounts to a straightforward application of the definition of the exponential of an operator, in terms of its Taylor expansion.

\(\exp_{\rm STD} (-i\delta t H) = \sum_{i=0}^{k} {(-i\delta t)^i\over{i!}} H^i.\)

The order*k*is determined by variable`TDExpOrder`. Some numerical considerations from Jeff Giansiracusa and George F. Bertsch suggest the 4th order as especially suitable and stable.**chebyshev**: In principle, the Chebyshev expansion of the exponential represents it more accurately than the canonical or standard expansion. As in the latter case,`TDExpOrder`determines the order of the expansion.

There exists a closed analytic form for the coefficients of the exponential in terms of Chebyshev polynomials:

\(\exp_{\rm CHEB} \left( -i\delta t H \right) = \sum_{k=0}^{\infty} (2-\delta_{k0})(-i)^{k}J_k(\delta t) T_k(H),\)

where \(J_k\) are the Bessel functions of the first kind, and H has to be previously scaled to \([-1,1]\). See H. Tal-Ezer and R. Kosloff,*J. Chem. Phys.***81**, 3967 (1984); R. Kosloff,*Annu. Rev. Phys. Chem.***45**, 145 (1994); C. W. Clenshaw,*MTAC***9**, 118 (1955).

**TDIonicTimeScale**
*Section*: Time-Dependent::Propagation
*Type*: float
*Default*: 1.0

This variable defines the factor between the timescale of ionic
and electronic movement. It allows reasonably fast
Born-Oppenheimer molecular-dynamics simulations based on
Ehrenfest dynamics. The value of this variable is equivalent to
the role of \(\mu\) in Car-Parrinello. Increasing it
linearly accelerates the time step of the ion
dynamics, but also increases the deviation of the system from the
Born-Oppenheimer surface. The default is 1, which means that both
timescales are the same. Note that a value different than 1
implies that the electrons will not follow physical behaviour.

According to our tests, values around 10 are reasonable, but it
will depend on your system, mainly on the width of the gap.

Important: The electronic time step will be the value of
`TDTimeStep` divided by this variable, so if you have determined an
optimal electronic time step (that we can call *dte*), it is
recommended that you define your time step as:

`TDTimeStep` = *dte* * `TDIonicTimeScale`

so you will always use the optimal electronic time step
(more details).

**TDLanczosTol**
*Section*: Time-Dependent::Propagation
*Type*: float
*Default*: 1e-5

An internal tolerance variable for the Lanczos method. The smaller, the more
precisely the exponential is calculated, and also the bigger the dimension
of the Krylov subspace needed to perform the algorithm. One should carefully
make sure that this value is not too big, or else the evolution will be
wrong.

**TDMaxSteps**
*Section*: Time-Dependent::Propagation
*Type*: integer
*Default*: 1500

Number of time-propagation steps that will be performed. You
cannot use this variable together with `TDPropagationTime`.

**TDMaxwellKSRelaxationSteps**
*Section*: Time-Dependent::Propagation
*Type*: integer
*Default*: 0

Time steps in which the coupled Maxwell-matter system relax the Maxwell states evolves under
free dynamics conditions. After these many steps, the external fields and currents are
switched on. The full requested simulation effectively states after this value.

**TDMaxwellTDRelaxationSteps**
*Section*: Time-Dependent::Propagation
*Type*: integer
*Default*: 0

Time steps needed to relax the Maxwell states in the presence of a matter system, to avoid
spurious relaxation effects. After these steps, the Maxwell-matter coupling can be switched on.
of the relaxation dynamics.

**TDPropagationTime**
*Section*: Time-Dependent::Propagation
*Type*: float

The length of the time propagation. You cannot set this variable
at the same time as `TDMaxSteps`. By default this variable will
not be used.

The units for this variable are \(\hbar\)/Hartree (or \(\hbar\)/eV if you
selected `ev_angstrom` as input units). The approximate conversions to
femtoseconds are 1 fs = 41.34 \(\hbar\)/Hartree = 1.52 \(\hbar\)/eV.

**TDPropagator**
*Section*: Time-Dependent::Propagation
*Type*: integer
*Default*: etrs

This variable determines which algorithm will be used to approximate
the evolution operator \(U(t+\delta t, t)\). That is, given
\(\psi(\tau)\) and \(H(\tau)\) for \(\tau \le t\),
calculate \(t+\delta t\). Note that in general the Hamiltonian
is not known at times in the interior of the interval \([t,t+\delta t]\).
This is due to the self-consistent nature of the time-dependent Kohn-Sham problem:
the Hamiltonian at a given time \(\tau\) is built from the
"solution" wavefunctions at that time.

Some methods, however, do require the knowledge of the Hamiltonian at some
point of the interval \([t,t+\delta t]\). This problem is solved by making
use of extrapolation: given a number \(l\) of time steps previous to time
\(t\), this information is used to build the Hamiltonian at arbitrary times
within \([t,t+\delta t]\). To be fully precise, one should then proceed
*self-consistently*: the obtained Hamiltonian at time \(t+\delta t\)
may then be used to interpolate the Hamiltonian, and repeat the evolution
algorithm with this new information. Whenever iterating the procedure does
not change the solution wavefunctions, the cycle is stopped. In practice,
in `Octopus` we perform a second-order extrapolation without a
self-consistency check, except for the first two iterations, where obviously
the extrapolation is not reliable.

The proliferation of methods is certainly excessive. The reason for it is that
the propagation algorithm is currently a topic of active development. We
hope that in the future the optimal schemes are clearly identified. In the
mean time, if you do not feel like testing, use the default choices and
make sure the time step is small enough.
*Options*:

**qoct_tddft_propagator**: WARNING: EXPERIMENTAL**caetrs**: (experimental) Corrected Approximated Enforced Time-Reversal Symmetry (AETRS), this is the previous propagator but including a correction step to the exponential.**runge_kutta4**: WARNING: EXPERIMENTAL. Implicit Gauss-Legendre 4th order Runge-Kutta.**runge_kutta2**: WARNING: EXPERIMENTAL. Implicit 2nd order Runge-Kutta (trapezoidal rule). Similar, but not identical, to Crank-Nicolson method.**expl_runge_kutta4**: WARNING: EXPERIMENTAL. Explicit RK4 method.**cfmagnus4**: WARNING EXPERIMENTAL**etrs**: The idea is to make use of time-reversal symmetry from the beginning:

\( \exp \left(-i\delta t H_{n} / 2 \right)\psi_n = \exp \left(i\delta t H_{n+1} / 2 \right)\psi_{n+1}, \)

and then invert to obtain:

\( \psi_{n+1} = \exp \left(-i\delta t H_{n+1} / 2 \right) \exp \left(-i\delta t H_{n} / 2 \right)\psi_{n}. \)

But we need to know \(H_{n+1}\), which can only be known exactly through the solution \(\psi_{n+1}\). What we do is to estimate it by performing a single exponential: \(\psi^{*}_{n+1}=\exp \left( -i\delta t H_{n} \right) \psi_n\), and then \(H_{n+1} = H[\psi^{*}_{n+1}]\). Thus no extrapolation is performed in this case.**aetrs**: Approximated Enforced Time-Reversal Symmetry (AETRS). A modification of previous method to make it faster. It is based on extrapolation of the time-dependent potentials. It is faster by about 40%. The only difference is the procedure to estimate \(H_{n+1}\): in this case it is extrapolated via a second-order polynomial by making use of the Hamiltonian at time \(t-2\delta t\), \(t-\delta t\) and \(t\).**exp_mid**: Exponential Midpoint Rule (EM). This is maybe the simplest method, but it is very well grounded theoretically: it is unitary (if the exponential is performed correctly) and preserves time-reversal symmetry (if the self-consistency problem is dealt with correctly). It is defined as: \( U_{\rm EM}(t+\delta t, t) = \exp \left( -i\delta t H_{t+\delta t/2}\right)\,. \)**crank_nicolson**: Classical Crank-Nicolson propagator. \( (1 + i\delta t H_{n+1/2} / 2) \psi_{n+1} = (1 - i\delta t H_{n+1/2} / 2) \psi_{n} \)**crank_nicolson_sparskit**: Classical Crank-Nicolson propagator. Requires the SPARSKIT library. \( (1 + i\delta t H_{n+1/2} / 2) \psi_{n+1} = (1 - i\delta t H_{n+1/2} / 2) \psi_{n} \)**magnus**: Magnus Expansion (M4). This is the most sophisticated approach. It is a fourth-order scheme (a feature which it shares with the ST scheme; the other schemes are in principle second-order). It is tailored for making use of very large time steps, or equivalently, dealing with problem with very high-frequency time-dependence. It is still in a experimental state; we are not yet sure of when it is advantageous.

**TDSCFThreshold**
*Section*: Time-Dependent::Propagation
*Type*: float
*Default*: 1.0e-6

Since the KS propagator is non-linear, each propagation step
should be performed self-consistently. In practice, for most
purposes this is not necessary, except perhaps in the first
iterations. This variable holds the number of propagation steps
for which the propagation is done self-consistently.

The self consistency has to be measured against some accuracy
threshold. This variable controls the value of that threshold.

**TDStepsWithSelfConsistency**
*Section*: Time-Dependent::Propagation
*Type*: integer
*Default*: 0

Since the KS propagator is non-linear, each propagation step
should be performed self-consistently. In practice, for most
purposes this is not necessary, except perhaps in the first
iterations. This variable holds the number of propagation steps
for which the propagation is done self-consistently.

The special value `all_steps` forces self-consistency to
be imposed on all propagation steps. A value of 0 means that
self-consistency will not be imposed. The default is 0.
*Options*:

**all_steps**: Self-consistency is imposed for all propagation steps.

**TDSystemPropagator**
*Section*: Time-Dependent::Propagation
*Type*: integer
*Default*: verlet

A variable to set the propagator in the multisystem framework.
This is a temporary solution, and should be replaced by the
TDPropagator variable.
*Options*:

**verlet**: (Experimental) Verlet propagator.**beeman**: (Experimental) Beeman propagator without predictor-corrector.**beeman_scf**: (Experimental) Beeman propagator with predictor-corrector scheme.**exp_mid**: (Experimental) Exponential midpoint propagator without predictor-corrector.**exp_mid_scf**: (Experimental) Exponential midpoint propagator with predictor-corrector scheme.

**TDTimeStep**
*Section*: Time-Dependent::Propagation
*Type*: float

The time-step for the time propagation. For most propagators you
want to use the largest value that is possible without the
evolution becoming unstable.

The default value is the maximum value that we have found
empirically that is stable for the spacing \(h\):
\(dt = 0.0426 - 0.207 h + 0.808 h^2\)
(from parabolic fit to Fig. 4 of http://dx.doi.org/10.1021/ct800518j,
probably valid for 3D systems only).
However, you might need to adjust this value.

**TemperatureFunction**
*Section*: Time-Dependent::Propagation
*Type*: integer
*Default*: "temperature"

If a thermostat is used, this variable indicates the name of the
function in the `TDFunctions` block that will be used to control the
temperature. The values of the temperature are given in
degrees Kelvin.

**Thermostat**
*Section*: Time-Dependent::Propagation
*Type*: integer
*Default*: none

This variable selects the type of thermostat applied to
control the ionic temperature.
*Options*:

**none**: No thermostat is applied. This is the default.**velocity_scaling**: Velocities are scaled to control the temperature.**nose_hoover**: Nose-Hoover thermostat.

**ThermostatMass**
*Section*: Time-Dependent::Propagation
*Type*: float
*Default*: 1.0

This variable sets the fictitious mass for the Nose-Hoover
thermostat.

**TDDeltaKickTime**
*Section*: Time-Dependent::Response
*Type*: float
*Default*: 0.0

The delta-perturbation that can be applied by making use of the `TDDeltaStrength` variable,
can be applied at the time given by this variable. Usually, this time is zero, since one wants
to apply the delta pertubation or "kick" at a system at equilibrium, and no other time-dependent
external potential is used. However, one may want to apply a kick on top of a laser field,
for example.

**TDDeltaStrength**
*Section*: Time-Dependent::Response
*Type*: float
*Default*: 0

When no laser is applied, a delta (in time) perturbation with
strength `TDDeltaStrength` can be applied. This is used to
calculate, *e.g.*, the linear optical spectra. If the ions are
allowed to move, the kick will affect them also.
The electric field is \(-(\hbar k / e) \delta(t)\) for a dipole with
zero wavevector, where *k* = `TDDeltaStrength`, which causes
the wavefunctions instantaneously to acquire a phase \(e^{ikx}\).
The unit is inverse length.

**TDDeltaStrengthMode**
*Section*: Time-Dependent::Response
*Type*: integer
*Default*: kick_density

When calculating the density response via real-time propagation,
one needs to perform an initial kick on the KS system, at
time zero. Depending on what kind of response property one wants to obtain,
this kick may be done in several modes. For use to calculate triplet excitations,
see MJT Oliveira, A Castro, MAL Marques, and A Rubio, *J. Nanoscience and Nanotechnology* **8**, 3392 (2008).
*Options*:

**kick_density**: The total density of the system is perturbed. This mode is appropriate for electric dipole response, as for optical absorption.**kick_spin**: The individual spin densities are perturbed oppositely. Note that this mode is only possible if the run is done in spin-polarized mode, or with spinors. This mode is appropriate for the paramagnetic dipole response, which can couple to triplet excitations.**kick_spin_and_density**: A combination of the two above. Note that this mode is only possible if the run is done in spin-polarized mode, or with spinors. This mode is intended for use with symmetries to obtain both of the responses at once, at described in the reference above.**kick_magnon**: Rotates the magnetization. Only works for spinors. Can be used in a supercell or my making use of the generalized Bloch theorem. In the later case (see`SpiralBoundaryConditions`) spin-orbit coupling cannot be used.

**TDDeltaUserDefined**
*Section*: Time-Dependent::Response
*Type*: string

By default, the kick function will be a dipole. This will change if (1) the variable
`TDDeltaUserDefined` is present in the inp file, or (2) if the block `TDKickFunction`
is present in the `inp` file. If both are present in the `inp` file, the `TDKickFunction`
block will be ignored. The value of `TDDeltaUserDefined` should be a string describing
the function that is going to be used as delta perturbation.

**TDKickFunction**
*Section*: Time-Dependent::Response
*Type*: block

If the block `TDKickFunction` is present in the input file, and the variable
`TDDeltaUserDefined` is not present in the input file, the kick function to
be applied at time zero of the time-propagation will not be a "dipole" function
(*i.e.* \(\phi \rightarrow e^{ikx} \phi\), but a general multipole in the form \(r^l Y_{lm}(r)\).

Each line has three columns: integers *l* and *m* that defines the
multipole, and a weight. Any number of lines may be given, and the kick will be the sum of those
multipoles with the given weights.

This feature allows calculation of quadrupole, octupole, etc., response functions.

**TDMomentumTransfer**
*Section*: Time-Dependent::Response
*Type*: block

Momentum-transfer vector for the calculation of the dynamic structure factor.
When this variable is set, a non-dipole field is applied, and an output file
`ftchd` is created (it contains the Fourier transform of the charge density
at each time). The type of the applied external field can be set by
an optional last number. Possible options are `qexp` (default), `qcos`,
`qsin`, or `qcos+qsin`. In the formulae below,
\(\vec{q}\) is the momentum-transfer vector.
*Options*:

**qexp**: External field is \(e^{i \vec{q} \cdot \vec{r}}\).**qcos**: External field is \(\cos \left( i \vec{q} \cdot \vec{r} \right)\).**qsin**: External field is \(\sin \left( i \vec{q} \cdot \vec{r} \right)\).**qbessel**: External field is \(j_l \left( \vec{q} \cdot \vec{r} \right) Y_{lm} \left(\vec{r} \right)\). In this case, the block has to include two extra values (*l*and*m*).

**TDMultipleMomentumTransfer**
*Section*: Time-Dependent::Response
*Type*: block

For magnon kicks only.
A simple way to specify momentum-transfer vectors for the calculation of
the magnetization dynamics. This variable should be used for a supercell.
For each reciprocal lattice vectors, the code will kick the original magnetization
using all the multiples of it.
The syntax reads:

`%TDMultipleMomentumTransfer
N_x | N_y | N_z
%`

and will include the (2N_x+1)*(2N_y+1)*(2N_z+1) multiples vectors of the reciprocal lattice vectors of the current cell.

**TDEasyAxis**
*Section*: Time-Dependent::Response::Dipole
*Type*: block

For magnon kicks only.
This variable defines the direction of the easy axis of the crystal.
The magnetization is kicked in the plane transverse to this vector

**TDPolarization**
*Section*: Time-Dependent::Response::Dipole
*Type*: block

The (real) polarization of the delta electric field. Normally
one needs three perpendicular polarization directions to calculate a
spectrum (unless symmetry is used).
The format of the block is:

`%TDPolarization
pol1x | pol1y | pol1z
pol2x | pol2y | pol2z
pol3x | pol3y | pol3z
%`

The default value for

Note that the directions do not necessarily need to be perpendicular when symmetries are used.

WARNING: If you want to obtain the cross-section tensor, the

**TDPolarizationDirection**
*Section*: Time-Dependent::Response::Dipole
*Type*: integer

When a delta potential is included in a time-dependent run, this
variable defines in which direction the field will be applied
by selecting one of the lines of `TDPolarization`. In a
typical run (without using symmetry), the `TDPolarization` block
would contain the three Cartesian unit vectors (the default
value), and one would make 3 runs varying
`TDPolarization` from 1 to 3.
If one is using symmetry, `TDPolarization` should run only from 1
to `TDPolarizationEquivAxes`.

**TDPolarizationEquivAxes**
*Section*: Time-Dependent::Response::Dipole
*Type*: integer
*Default*: 0

Defines how many of the `TDPolarization` axes are equivalent. This information is stored in a file and then
used by `oct-propagation_spectrum` to rebuild the full polarizability tensor from just the
first `TDPolarizationEquivAxes` directions. This variable is also used by `CalculationMode = vdw`.

**TDPolarizationWprime**
*Section*: Time-Dependent::Response::Dipole
*Type*: block

This block is needed only when
`TDPolarizationEquivAxes` is set to 3. In such a case,
the three directions (*pol1*, *pol2*, and *pol3*) defined in
the `TDPolarization` block should be related by symmetry
operations. If *A* is the symmetry operation that takes you
from *pol1* to *pol2*, then `TDPolarizationWprime`
should be set to the direction defined by *A*\(^{-1}\)*pol3*.
For more information see MJT Oliveira
*et al.*, *J. Nanoscience and Nanotechnology* **8**,
3392 (2008).

**MaxwellTDOutput**
*Section*: Time-Dependent::TD Output
*Type*: flag
*Default*: maxwell_energy

Defines what should be output during the time-dependent
Maxwell simulation. Many of the options can increase the computational
cost of the simulation, so only use the ones that you need. In
most cases the default value is enough, as it is adapted to the
details of the TD run.
*Options*:

**b_field_surface**: Output of the B field sliced along the planes x=0, y=0, z=0 for each field component**maxwell_energy**: Output of the electromagnetic field energy into the folder`td.general/maxwell`.**maxwell_fields**: Output of the electromagnetic field at the origin of the simulation box into the folder`td.general/fields`**mean_poynting**: Output of the mean Poynting vector**e_field_surface**: Output of the E field sliced along the planes x=0, y=0, z=0 for each field component

**TDExcitedStatesToProject**
*Section*: Time-Dependent::TD Output
*Type*: block

**[WARNING: This is a *very* experimental feature]**
To be used with `TDOutput = populations`.
The population of the excited states
(as defined by *t*, and |Phi_I> is the excited state of interest) can be approximated -- it is not clear
how well -- by substituting for those real many-body states the time-dependent Kohn-Sham
determinant and some modification of the Kohn-Sham ground-state determinant (*e.g.*,
a simple HOMO-LUMO substitution, or the Casida ansatz for excited states in linear-response
theory. If you set `TDOutput` to contain `populations`, you may ask for these approximated
populations for a number of excited states, which will be described in the files specified
in this block: each line should be the name of a file that contains one excited state.

This file structure is the one written by the casida run mode, in the files in the directory `*_excitations`.
The file describes the "promotions" from occupied
to unoccupied levels that change the initial Slater determinant
structure specified in ground_state. These promotions are a set
of electron-hole pairs. Each line in the file, after an optional header, has four
columns:

*i a \(\sigma\) weight*

where *i* should be an occupied state, *a* an unoccupied one, and \(\sigma\)
the spin of the corresponding orbital. This pair is then associated with a
creation-annihilation pair \(a^{\dagger}_{a,\sigma} a_{i,\sigma}\), so that the
excited state will be a linear combination in the form:

\(\left|{\rm ExcitedState}\right> =
\sum weight(i,a,\sigma) a^{\dagger}_{a,\sigma} a_{i,\sigma} \left|{\rm GroundState}\right>\)

where *weight* is the number in the fourth column.
These weights should be normalized to one; otherwise the routine
will normalize them, and write a warning.

**TDFloquetDimension**
*Section*: Time-Dependent::TD Output
*Type*: integer
*Default*: -1

Order of Floquet Hamiltonian. If negative number is given, downfolding is performed.

**TDFloquetFrequency**
*Section*: Time-Dependent::TD Output
*Type*: float
*Default*: 0

Frequency for the Floquet analysis, this should be the carrier frequency or integer multiples of it.
Other options will work, but likely be nonsense.

**TDFloquetSample**
*Section*: Time-Dependent::TD Output
*Type*: integer
*Default*: 20

Number of points on which one Floquet cycle is sampled in the time-integral of the Floquet analysis.

**TDMultipoleLmax**
*Section*: Time-Dependent::TD Output
*Type*: integer
*Default*: 1

Maximum electric multipole of the density output to the file `td.general/multipoles`
during a time-dependent simulation. Must be non-negative.

**TDOutput**
*Section*: Time-Dependent::TD Output
*Type*: flag
*Default*: multipoles + energy (+ others depending on other options)

Defines what should be output during the time-dependent
simulation. Many of the options can increase the computational
cost of the simulation, so only use the ones that you need. In
most cases the default value is enough, as it is adapted to the
details of the TD run. If the ions are allowed to be moved, additionally
the geometry and the temperature are output. If a laser is
included it will output by default.

Note: the output files generated by this option are updated
every `RestartWriteInterval` steps.
*Options*:

**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 \(l\) 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.**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**spin**: (Experimental) Outputs the expectation value of the spin, which can be used to calculate magnetic circular dichroism.**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.**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`.

**TDOutputComputeInterval**
*Section*: Time-Dependent::TD Output
*Type*: integer
*Default*: 50

The TD output requested are computed
when the iteration number is a multiple of the `TDOutputComputeInterval` variable.
Must be >= 0. If it is 0, then no output is written.
Implemented only for projections and number of excited electrons for the moment.

**TDOutputDFTU**
*Section*: Time-Dependent::TD Output
*Type*: flag
*Default*: none

Defines what should be output during the time-dependent
simulation, related to LDA+U.

Note: the output files generated by this option are updated
every `RestartWriteInterval` steps.
*Options*:

**effective_u**: Writes the effective U for each orbital set as a function of time.

**TDProjStateStart**
*Section*: Time-Dependent::TD Output
*Type*: integer
*Default*: 1

To be used with `TDOutput = td_occup`. Not available if `TDOutput = populations`.
Only output projections to states above `TDProjStateStart`. Usually one is only interested
in particle-hole projections around the HOMO, so there is no need to calculate (and store)
the projections of all TD states onto all static states. This sets a lower limit. The upper limit
is set by the number of states in the propagation and the number of unoccupied states
available.