ElectricFieldConstraint

class ElectricFieldConstraint(electric_field=None, born_effective_charge_parameters=None, charge_sum_rule_method=None, polarization_parameters=None, update_strategy=None, optimize_for_field_direction=None, enable_stress_correction=None)

Add a response to an applied electic field to energy, forces and stress.

Parameters:
  • electric_field (PhysicalQuantity of type Volt / length) – A vector representing the screened electric field inside the material.
    Default: [0.0, 0.0, 0.0] * Volt / Meter

  • born_effective_charge_parameters (BornEffectiveChargeParameters) – A container with the parameters to be passed to BornEffectiveCharge when evaluating the electric field correction.
    Default: BornEffectiveChargeParameters()

  • charge_sum_rule_method (Disabled | DistributeEvenly | DistributeRelative.) – Keyword specifying if a charge sum-rule should be applied to the Born effective charges. Can be Disabled for no sum-rule application. Alternatively the charge error is either distributed evenly (DistributeEvenly) or relative (DistributeRelative) to the Born effective charges on the atoms.
    Default: Disabled

  • polarization_parameters (PolarizationParameters) – A container with the parameters to be passed to Polarization when evaluating the electric field enthalpy.
    Default: PolarizationParameters()

  • update_strategy (FixElectricFieldCorrection() | UpdateElectricFieldCorrection) – The strategy used to update Polarization and Born Effective Charge during the optimization or molecular dynamics. FixElectricFieldCorrection() will calculate Polarization and BornEffectiveCharge upfront based on the initial configuration and keep the value fixed. UpdateElectricFieldCorrection will calculate Polarization and BornEffectiveCharge at each step, unless the change in configuration coordinates is below a given threshold.
    Default: FixElectricFieldCorrection()

  • optimize_for_field_direction (bool) – For an electric field along a Cartesian direction and orthorhombic cell, it is possible to speed up the calculation by evaluating only Polarization components along the direction of the field. If True, this optimization is enforced when possible.
    Default: True

  • optimize_for_field_direction – For an electric field along a Cartesian direction and orthorhombic cell, it is possible to speed up the calculation by evaluating only Polarization components along the direction of the field. If True, this optimization is enforced when possible.
    Default: True

  • enable_stress_correction (bool) – Whether also the stress tensor should include the electric field correction.
    Default: True

bornEffectiveChargeParameters()
Returns:

The container with the BornEffectiveCharge parameters.

Return type:

BornEffectiveChargeParameters

calculatedBornEffectiveCharge()

Retrieve the BornEffectiveCharge analysis available on this object, if any.

Returns:

The born effective charge and the configuration for which it was evaluated. In the case of NEB, a list is returned.

Return type:

None | (AtomicConfiguration, BornEffectiveCharge) | list of (AtomicConfiguration, BornEffectiveCharge)

calculatedPolarization()

Retrieve the last Polarization calculated to evaluate the correction, if any.

Returns:

The polarization and the configuration for which it was evaluated. In the case of NEB, a list is returned.

Return type:

None | (AtomicConfiguration, Polarization) | list of AtomicConfiguration, list of Polarization

chargeSumRuleMethod()
Returns:

Whether a charge sum-rule should be applied to the Born effective charges.

Type:

Disabled | DistributeEvenly | DistributeRelative.

electricField()
Returns:

A vector representing the screened electric field inside the material.
Default: [0.0, 0.0, 0.0] * Volt / Meter

Return type:

PhysicalQuantity of type Volt / length

enableStressCorrection()
Returns:

Whether the stress tensor should include the electric field correction.

Type:

bool

frozenDegreesOfFreedom(local_atoms=None)
Parameters:

local_atoms (list of int | None) – The group of atoms from which the frozen degrees of freedom should be calculated, e.g. a thermalized group of atoms.
Default: All atoms.

Returns:

The number of degrees of freedom that are frozen by this constraint object.

Return type:

int

optimizeForFieldDirection()
Returns:

Whether the optimization with respect to Polarization components is enabled or not.

Return type:

bool

polarizationParameters()
Returns:

The container with the Polarization parameters.

Return type:

Polarization

uniqueString()

Return a unique string representing the state of the object.

updateStrategy()
Returns:

The strategy used to update Polarization and Born Effective Charges during the optimization or molecular dynamics.

Return type:

FixPolarizationToInitialValue | UpdatedPolarizationDuringTrajectory

Usage Examples

Optimize the geometry of a tetragonal BaTiO3 with an applied electric field opposite to the initial polarization:

# -*- coding: utf-8 -*-
setVerbosity(MinimalLog, ELECTRIC_FIELD_CONSTRAINT=True)

# -------------------------------------------------------------
# Bulk Configuration
# -------------------------------------------------------------

# Set up lattice
lattice = SimpleTetragonal(3.97816645571*Angstrom, 4.07906435409*Angstrom)

# Define elements
elements = [Barium, Titanium, Oxygen, Oxygen, Oxygen]

# Define coordinates
fractional_coordinates = [[ 0.00000000008 ,  0.000000000035, -0.010189450414],
                          [ 0.500000000154,  0.500000000008,  0.465308176297],
                          [ 0.500000000067,  0.50000000004 ,  1.011196206321],
                          [ 0.500000000058,  0.000000000049,  0.498397341501],
                          [ 0.000000000038,  0.500000000045,  0.498397341501]]

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

# -------------------------------------------------------------
# Calculator
# -------------------------------------------------------------
#----------------------------------------
# Basis Set
#----------------------------------------
basis_set = [
    LDABasis.Oxygen_SingleZetaPolarized,
    LDABasis.Titanium_SingleZetaPolarized,
    LDABasis.Barium_SingleZetaPolarized,
    ]

#----------------------------------------
# Exchange-Correlation
#----------------------------------------
exchange_correlation = LDA.PZ

k_point_sampling = MonkhorstPackGrid(
    na=2,
    nb=2,
    nc=2,
    )
numerical_accuracy_parameters = NumericalAccuracyParameters(
    density_mesh_cutoff=80.0*Hartree,
    k_point_sampling=k_point_sampling
    )

iteration_control_parameters = IterationControlParameters(
    tolerance=1e-02,
    )

calculator = LCAOCalculator(
    basis_set=basis_set,
    exchange_correlation=exchange_correlation,
    numerical_accuracy_parameters=numerical_accuracy_parameters,
    iteration_control_parameters=iteration_control_parameters,
    )

bulk_configuration.setCalculator(calculator)
nlprint(bulk_configuration)
bulk_configuration.update()
nlsave('bulk_configuration.hdf5', bulk_configuration)

# -------------------------------------------------------------
# Optimize Geometry
# -------------------------------------------------------------

# Use default parameters for the Born effective charges and
# Polarization.
born_effective_charge_parameters = BornEffectiveChargeParameters()
polarization_parameters = PolarizationParameters()

# Setup an electric field constraint object.
efield = ElectricFieldConstraint(
    [0.0, 0.0, 0.1] * Volt / Ang,
    born_effective_charge_parameters=born_effective_charge_parameters,
    polarization_parameters=polarization_parameters,
    update_strategy=UpdateElectricFieldCorrection(max_displacement=0.*Angstrom))

constraints = [efield]

bulk_configuration = OptimizeGeometry(
    bulk_configuration,
    max_forces=0.05*eV/Ang,
    max_steps=200,
    max_step_length=0.2*Ang,
    constraints=constraints,
    trajectory_filename='optimize_bulk_configuration.hdf5',
    trajectory_interval=1*Second,
    restart_strategy=RestartFromTrajectory(),
    optimize_cell=False,
    enable_optimization_stop_file=True,
    )
nlsave('optimize_bulk_configuration.hdf5', bulk_configuration)
nlprint(bulk_configuration)

efield_optimization.py

In the example above the parameter update_strategy=UpdateElectricFieldCorrection(max_displacement=0.*Angstrom) implies that electric field derived force and stress terms (see section below) are recalculated at each optimization step. The polarization and Born effective charges along the optimization trajectory are printed in the output log, as long as the log region ELECTRIC_FIELD_CONSTRAINT is active (refer to setVerbosity documentation).

Note

The order of the constraint specified in input to OptimizeGeometry, OptimizeNudgedElasticBand and MolecularDynamics matters. Constraint objects modify energies, forces or stress in the order they are specified. Constraints which fix energy, forces or stress to a given value, like FixAtomConstraints, should always be specified after ElectricFieldConstraint.

ElectricFieldConstraint can be used in the same way in OptimizeNudgedElasticBand and MolecularDynamics.

Notes

ElectricFieldConstraint is used to model the effect of an external electric field in DFT calculations with periodic boundary conditions.

Within the DFT formalism, the electric enthalpy of an insulator under a uniform electric field is expressed as [1]:

\[E^{\varepsilon}[\mathbf{R}, \varepsilon] = E^{0}[\mathbf{R}] - \mathbf{P}[\mathbf{R}] \cdot \mathbf{\varepsilon}\]

where \(\mathbf{\varepsilon}\) represents the electric field, \(\mathbf{R}\) the atomic coordinates, \(E^{0}\) the Kohn-Sham energy functional and \(\mathbf{P}\) the macroscopic polarization [2].

Only the dominant first-order response to the field is taken into account, i.e. the field-induced polarization of the electronic wave-function, which introduces an explicit dependency of polarization and Kohn-Sham energy on the electric field, is neglected.

The force on an atom \(i\) is:

\[F_{i}^{\varepsilon}[\mathbf{R}] = F_{i}^{0} + Z_{i}^{*}[\mathbf{R}] \cdot \mathbf{\varepsilon}\]

where the second term of the equation is a field-induced force dependent on the Born Effective Charge tensor \(Z_{i}^{*}\).

Similarly, the field-induced stress is:

\[\sigma^{\varepsilon}_{\alpha \beta} = -\sum_{i} F_{i, \alpha} R_{i, \beta}\]

The field-induced correction to energy, forces and stress can then be determined by calculating a Polarization and BornEffectiveCharge analysis at a given atomic positions.

The ElectricFieldConstraint evaluates such corrections and can be used in all methods which accept constraint objects, i.e. OptimizeGeometry, OptimizeNudgedElasticBand and MolecularDynamics.

Note

ElectricFieldConstraint adds a large overhead to calculation time, mainly due to the evaluation of BornEffectiveCharge. The Born effective charges and the polarization must in principle be recalculated for any change in atomic coordinates. However, ElectricFieldConstraint allows to limit the computational burden by performing such calcualtions more sparsely via the update_strategy parameter. The validity of such approach will depend on the system and the quantities under evaluation, and should be verified case by case. For reasonable static geometries, a single evaluation of the Born effective charges may be sufficient.