PlaneWaveCalculator

class PlaneWaveCalculator(wave_function_cutoff=None, basis_set=None, exchange_correlation=None, numerical_accuracy_parameters=None, iteration_control_parameters=None, poisson_solver=None, algorithm_parameters=None, charge=None, fixed_spin_moment=None, correction_extension=None, checkpoint_handler=None, store_wave_functions=None, processes_per_kpoint=None, parallel_parameters=None)

Class for representing calculations using the ATK-PW plane wave method for BulkConfigurations.

Parameters:

  • wave_function_cutoff (PhysicalQuantity of type energy) – The energy cutoff of the plane wave basis functions to include the basis set. This value determines the size of the basis set. See Wave function cutoff.
    Default: 0.25 * density_cutoff from numerical_accuracy_parameters

  • basis_set (BasisSet | PAWDataSet) – An object describing the basis set to use. Note that this parameter is only used to determine the pseudopotentials to use.
    Default: BasisGGAPseudoDojo.Medium

  • exchange_correlation (ExchangeCorrelation) – The exchange correlation to use for this calculation.
    Default: GGA.PBE

  • numerical_accuracy_parameters (NumericalAccuracyParameters) –

    The NumericalAccuracyParameters used for the calculation.
    Default:

    NumericalAccuracyParameters(
    density_cutoff=1.0e-6,
    density_mesh_cutoff=80.0*Hartree,
    interaction_max_range=20.0*Angstrom,
    k_point_sampling=MonkhorstPackGrid(1, 1, 1),
    number_of_reciprocal_points=1024,
    occupation_method=FermiDirac(1000.0*Kelvin),
    radial_step_size=0.001*Bohr,
    reciprocal_energy_cutoff=1250.0*Hartree
    bands_per_electron=1.2)

  • iteration_control_parameters (IterationControlParameters) –

    The IterationControlParameters used for the calculation. For non-self-consistent, set this to NonSelfconsistent.
    Default:

    IterationControlParameters(
    tolerance=1.0e-5,
    max_steps=100,
    algorithm=PulayMixer(),
    damping_factor=0.1,
    number_of_history_steps=20,
    start_mixing_after_step=0,
    mixing_variable=HamiltonianVariable,
    linear_dependence_threshold=0.0,
    preconditioner=Preconditioner.Off)

  • poisson_solver (FastFourierSolver) – The Poisson solver to calculate the electrostatic potential.
    Default: A default FastFourierSolver object

  • algorithm_parameters (AlgorithmParameters) –

    The AlgorithmParameters used for calculating the density matrix.
    Default:

    AlgorithmParameters(
    density_matrix_method=GeneralizedDavidsonSolver(),
    store_grids=True,
    store_basis_on_grid=Automatic,
    store_energy_density_matrix=Automatic,
    scf_restart_step_length=0.1*Angstrom
    use_symmetries=Automatic)

  • charge (float) – The charge of the system, a charge of -1 corresponds to one extra electron.
    Default: 0.0

  • fixed_spin_moment (float | False) – Total spin moment (per unit cell) to use, defined as \(\Delta N = N_{\uparrow} - N_{\downarrow}\), where \(N_{\uparrow}\) and \(N_{\downarrow}\) are the number of electrons in the Up and Down spin channels, respectively. When specified the spin moment will be fixed at the given value by introducing separate Fermi levels for the Up and Down spin channels. Can only be specified when doing a calculation with polarized spin. If set to False the spin moment will not be fixed - a single Fermi level is used.
    Default: False

  • correction_extension (GrimmeDFTD2 | GrimmeDFTD3) – The correction extension to used, when calculating energy, forces and stress.
    Default: None.

  • checkpoint_handler (CheckpointHandler) – The checkpoint handler used for specifying the save-file and the time interval.
    Default: No checkpoint handling.

  • store_wave_functions (bool) –

    Flag that decides whether or not the wave functions should be saved along with the potential.


    Default: True For Hybrid calculations

    False All other cases.

  • processes_per_kpoint (Automatic | int) – The number of processes to use per kpoint. When set to Automatic the number will be determined automatically from the total number of k-points and processes such as to keep the number as small as possible and at the same time minimize the number of idle processes. One may set this number manually in order to reduce the memory requirements for each process.
    Default: Automatic

  • parallel_parameters (ParallelParameters) –

    The parameters used to control parallelization options.
    Default:

    ParallelParameters(
    processes_per_neb_image=None,
    processes_per_individual=None,
    processes_per_bias_point=None,
    processes_per_saddle_search=1)

algorithmParameters()

Query method for the algorithm parameters.

Returns:

The algorithm parameters object.

Return type:

AlgorithmParameters

basisSet()

Get the basis set.

Returns:

The basis set set on the calculator.

Return type:

BasisSet

charge()

Get the charge of the system.

Returns:

The charge in units of the electron charge.

Return type:

float

checkpointHandler()

Return the checkpoint handler.

correctionExtension()
Returns:

The correction extension or None if not set.

Return type:

GrimmeDFTD2 | GrimmeDFTD3 | None

exchangeCorrelation()

Get the exchange-correlation.

Returns:

The exchange correlation set on the calculator.

Return type:

ExchangeCorrelation

fixedSpinMoment()

Get the fixed spin moment.

Returns:

The fixed spin moment or False if the spin moment is not held fixed.

Return type:

float | False

isConverged()
Returns:

True when the call to “update()” resulted in a converged SCF loop.

Return type:

bool

iterationControlParameters()

Query method for the IterationControlParameters.

Returns:

The iteration control parameters set on the calculator.

Return type:

IterationControlParameters

metatext()
Returns:

The metatext of the object or None if no metatext is present.

Return type:

str | None

numericalAccuracyParameters()

Query method for the NumericalAccuracyParameters.

Returns:

The numerical accuracy parameters set on the calculator.

Return type:

NumericalAccuracyParameters

parallelParameters()
Returns:

The parallel parameters object.

Return type:

ParallelParameters

poissonSolver()
Returns:

The Poisson solver set on the calculator.

Return type:

DirectSolver | MultigridSolver | FastFourierSolver | FastFourier2DSolver

processesPerKpoint()
Returns:

The number of processes to use per kpoint.

Return type:

int

setBasisSet(basis_set)

Set the basis set.

Parameters:

basis_set (BasisSet) – The basis set to check and set.

setCheckpointHandler(checkpoint_handler)

Set the checkpoint handler. :param checkpoint_handler: The CheckpointHandler to set on the calculator. :type checkpoint_handler: CheckpointHandler

setCorrectionExtension(correction_extension)

Set the iteration correction extension.

Parameters:

correction_extension (GrimmeDFTD2 | GrimmeDFTD3 | None) – The correction extension to use.

setExchangeCorrelation(exchange_correlation)

Set the exchange-correlation.

Parameters:

exchange_correlation (ExchangeCorrelation) – The exchange correlation to set.

setIterationControlParameters(iteration_control_parameters)

Set the iteration control parameters.

Parameters:

iteration_control_parameters (IterationControlParameters) – The iteration control parameters to set.

setMetatext(metatext)

Set a given metatext string on the object.

Parameters:

metatext (str | None) – The metatext string that should be set. A value of “None” can be given to remove the current metatext.

setNumericalAccuracyParameters(numerical_accuracy_parameters)

Set the numerical accuracy parameters.

Parameters:

numerical_accuracy_parameters (NumericalAccuracyParameters) – The numerical accuracy parameters to set.

setParallelParameters(parallel_parameters)

Method for setting the parallel parameters.

Parameters:

parallel_parameters (ParallelParameters | None) – The parallel parameters to set. If None a default version of the parameters is used.

setPoissonSolver(poisson_solver)

Set the Poisson solver.

Parameters:

poisson_solver (DirectSolver | MultigridSolver | FastFourierSolver | FastFourier2DSolver) – The Poisson solver to set on the calculator.

storeWaveFunctions()
Returns:

The flag that decides whether or not the wave functions should be stored on disk.

Return type:

bool

uniqueString()

Return a unique string representing the state of the object.

upgrade(configuration)

Private method for updating the calculator from the configuration, if it is possible. @private

versionUsed()
Returns:

The version of ATK used to update the calculator.

Return type:

str

waveFunctionCutoff()

Get the plane wave cutoff energy.

Returns:

The plane wave cutoff energy.

Return type:

PhysicalQuantity of type energy

Notes

Wave function cutoff

A calculation with the PlaneWaveCalculator allows to represent the Bloch states in a systematically improvable discrete Fourier expansion:

\[\psi_{\mathbf{k},n}(\mathbf{r}) = \sum_\mathbf{G}^{|\mathbf{k} + \mathbf{G}| < Q_{\rm max}} C_{\mathbf{k} + \mathbf{G},n}~e^{i\left(\mathbf{k} + \mathbf{G}\right) \cdot \mathbf{r}},\]

where \(\mathbf{G}\) is a reciprocal lattice vector. The number of plane waves to include in the expansion is determined by the reciprocal lattice vector cutoff length \(Q_{\rm max}=\frac{1}{2} E_{\rm cut}\) where \(E_{\rm cut}\) is the kinetic energy cutoff of the wave functions set by the keyword argument wave_function_cutoff. By increasing the wave function cutoff one increases the number of plane wave basis functions used in the calculation. This leads to a better approximation of the exact wave function, but at an increased computational cost.

Density mesh cutoff

Since the density is defined as the square of the wave function amplitudes:

\[\rho(\mathbf{r}) = \sum_{\mathbf{k},n} f(\epsilon_{\mathbf{k}, n}) |\psi_{\mathbf{k}, n}(\mathbf{r})|^2\]

we, in principle, need twice the wave function resolution for the density. By default, we set the density_mesh_cutoff on NumericalAccuracyParameters to 4 times the wave_function_cutoff. In some cases this may be unnecessary, and can be overwritten by the user by specifying it explicitly.

Bands per electron

Since the plane wave basis is typically very large, a direct diagonalization of the Hamiltonian will be very inefficient and very time consuming. Instead iterative eigensolvers (GeneralizedDavidsonSolver and PPCGSolver) are used to solve for the lowest eigenstates that are needed to describe the ground state. The number of bands to solve for is set by the bands_per_electron keyword in NumericalAccuracyParameters. One has to be sure to solve for all bands that are fully or partly occupied through the SCF cycle, but keeping the number as low as possible ensures the fastest calculation.

By default bands_per_electron is set to 1.2 which corresponds to adding upwards of 20% more bands than there are bands below the Fermi level (for a non-metal). This is a safe setting for semi-conductors. For metals one might need to increase the number of bands, especially when using a large smearing. For insulators fewer bands will likely suffice. For closed shell systems, one can set it to 1.0.

Initialization methods

The iterative eigensolvers generally need a starting guess for the wave functions and the density used in constructing the initial Hamiltonian. The initialization method is set on the eigensolver supplied to the AlgorithmParameters and can be one of the following two:

Usage Examples

Si FCC calculation with default parameters

The following script is an example of how to set up a PlaneWaveCalculator for a Si FCC system with all defaults.

# Set up a Bravais lattice.
lattice = FaceCenteredCubic(5.4306*Angstrom)

# Define elements.
elements = [Silicon, Silicon]

# Define coordinates.
fractional_coordinates = [[0.0, 0.0, 0.0],
                     [0.25, 0.25, 0.25]]

# Set up configuration.
bulk_configuration = BulkConfiguration(
    bravais_lattice=lattice,
    elements=elements,
    fractional_coordinates=fractional_coordinates)

# Set up calculator.
calculator = PlaneWaveCalculator()
bulk_configuration.setCalculator(calculator)

# Run the scf loop.
bulk_configuration.update()

Accurate Si FCC calculation

If we want to perform a more accurate calculation, with a custom density_mesh_cutoff, an example script could be:

# Define a non-default cutoff.
wave_function_cutoff = 80 * Hartree

# Make some custom numerical accuracy parameters.
numerical_accuracy_parameters = NumericalAccuracyParameters(
    # Large k-point sampling.
    k_point_sampling=MonkhorstPackGrid(11, 11, 11),
    # Take the density mesh cutoff as 2.5 times the wave function cutoff instead of 4.
    density_mesh_cutoff=2.5 * wave_function_cutoff)

# Converge to a higher than default tolerance.
iteration_control_parameters = IterationControlParameters(tolerance=1.0e-7)

# Set up calculator.
calculator = PlaneWaveCalculator(
    # Take a higher number of bands (default for Si-fcc is 4).
    bands_above_fermi_level=10,
    wave_function_cutoff=wave_function_cutoff,
    numerical_accuracy_parameters=numerical_accuracy_parameters,
    iteration_control_parameters=iteration_control_parameters)