GrainBoundaryScattering

class GrainBoundaryScattering(configuration, filename, object_id, grain_boundary_generators, kpoint_density=None, optimization_region_length=None, optimize_geometry_parameters=None, energies=None, transmission_kpoint_density=None, self_energy_calculator=None, infinitesimal=None, passivate_electrode_surfaces=None, device_algorithm_parameters=None, contour_parameters=None, relaxation_calculator=None, log_filename_prefix=None, number_of_processes_per_task=None)

Class for performing a grain boundary scattering study.

Parameters:
  • configuration (DeviceConfiguration) – The device configuration with attached calculator for which to perform the study. The configuration must include at least one metallic region acting as the gate.

  • filename (str) – The full or relative filename path the Study object should be saved to. See nlsave().

  • object_id (str) – The name of the study that the Study object should be saved to within the file. This needs to be a unique name in this file. See nlsave().

  • grain_boundary_generators (list or Table of GrainBoundaryGenerator) – List or Table of GrainBoundaryGenerator to generate the device grain boundary structures.

  • kpoint_density (PhysicalQuantity of type length.) – The k-point density to be used for the self-consistent DFT-NEGF calculation.
    Default: 4.0 * Ang

  • optimization_region_length (PhysicalQuantity of type length) – The length along the transport direction of the region to be optimized within the central region of the device. Note that the length of the optimization region is automatically capped to keep the given electrode extensions fixed on either side.
    Default: 10.0 * Ang

  • optimize_geometry_parameters (OptimizeGeometryParameters) – The parameters to use for optimizing the geometry. Note that constraints, pre_step_hook and post_step_hook must be left unset. Also note that max_stress, target_stress, restart_strategy and enable_optimization_stop_file do not have any effect since stress is not optimized, and the restart mechanism within OptimizeGeometry() is disabled.
    Default: OptimizeGeometryParameters()

  • energies (list of PhysicalQuantity of type energy) – A list of the energies for which the TransmissionSpectrum at each voltage should be calculated.
    Default: Energy range that covers the bias window, plus \(30 k_B T\).

  • transmission_kpoint_density (PhysicalQuantity of type length.) – The k-point density to be used for the transmission spectrum calculations.
    Default: 7.0 * Ang

  • self_energy_calculator (DirectSelfEnergy | RecursionSelfEnergy | SparseRecursionSelfEnergy | KrylovSelfEnergy) – The SelfEnergyCalculator to be used for the TransmissionSpectrum calculation at each voltage.
    Default: RecursionSelfEnergy(storage_strategy=NoStorage())

  • infinitesimal (PhysicalQuantity of type energy) – Small positive energy, used to move the TransmissionSpectrum calculation at each voltage away from the real axis. This is only relevant for recursion-style self-energy calculators.
    Default: 1.0e-6 * eV

  • passivate_electrode_surfaces (bool) – Whether the left and right surfaces of the central region of the device along the transport direction should be passivated when extracting the central region for optimization.
    Default: False.

  • device_algorithm_parameters (DeviceAlgorithmParameters) – The DeviceAlgorithmParameters used for the device simulation
    Default: The default of the calculator.

  • contour_parameters (ContourParameters) – The parameters used for the complex contour integration for the device calculation.
    Default: The default of the calculator.

  • relaxation_calculator (Calculator) – The calculator used for relaxing the atomic coordinates of the bulk configuration and of the grain boundary configurations. The calculator must contain a basis set for the elements in the configuration. Note that the calculator is taken as a reference corresponding to the bulk unit cell given in bulk_configuration; the density_mesh_cutoff and k_point_sampling parameters of the NumericalAccuracyParameters of the calculator will be scaled consistently with the supercell size. The k_point_sampling will be scaled to maintain the k-point density approximately equal to the reference.

  • log_filename_prefix (str | LogToStdOut) – Filename prefix for the logging output of the calculations, each to be stored in a separate file. If LogToStdOut, all logging will instead be sent to standard output.
    Default: 'grain_boundary_scattering_'

  • number_of_processes_per_task (int) – The number of processes that will be used to execute each task. If this value is greater than or equal to the total number of available processes, each single task will be executed collaboratively over all processes. Otherwise, a delegator-worker scheme is used; in this case, one process will be set aside as the delegator, and the remaining ones will be grouped into workers and execute tasks concurrently.
    Default: All available processes execute each task collaboratively.

addGrainBoundaryGenerator(grain_boundary_generator=None)

Add a grain boundary generator to the list of grain boundary generators.

Parameters:

grain_boundary_generator (GrainBoundaryGenerator) – The grain boundary generator to generate the device grain boundary structures.

calculateMayadasShatzkesResistivityVsGrainBoundarySize(bulk_resistivity, bulk_mean_free_path, grain_boundary_lengths, grain_boundary_weights=None, temperature=PhysicalQuantity(300.0, K), spin=<class 'NL.ComputerScienceUtilities.NLFlag.Spin.Sum'>, fermi_shift=PhysicalQuantity(0.0, eV))

Calculate the resistivity vs. grain boundary size based on the Mayadas and Shatzkes model. Returns a list of resistivities corresponding to the valus of the GB lengths.

Parameters:
  • bulk_resistivity (PhysicalQuantity of type Meter / Siemens) – The bulk resistivity of the material without grain boundaries.

  • bulk_mean_free_path (PhysicalQuantity of type length) – The bulk resistivity of the material without grain boundaries.

  • grain_boundary_lengths (PhysicalQuantity of type length) – List of grain boundary lengths. The GB lengths are the average characteristic dimension of the grain boundaries.

  • grain_boundary_weights – The weight factors for the individual grain boundary specific resistances.
    Default: list of equal weights, i.e. numpy.ones(number_of_grain_boundaries)

  • temperature (PhysicalQuantity of type Kelvin) – The electrode temperatures.

  • spin (Spin.Up | Spin.Down | Spin.Sum) – The spin component to query result for.
    Default: Spin.Sum

  • fermi_shift (PhysicalQuantity of type energy) – Manuel shift of the Fermi level shift.
    Default: 0.0 * eV

Returns:

The resistivity vs. grain boundary size based on the Mayadas and Shatzkes model.

Return type:

PhysicalQuantity of type Meter/Siemens

calculatedGrainBoundaries()
Returns:

The list GrainBoundaryGenerator descriptions for which all the calculation results are present in the cache.

Return type:

list of GrainBoundaryGenerator

dependentStudies()
Returns:

The list of dependent studies.

Return type:

list of Study

electrodeTotalEnergies(contribution=None)

Retrieve the total energies for the specified electrode.

Parameters:

contribution (Left | Right) – The electrode for which to get the energy.

Returns:

The total energies. If not available, returns list of None.

Return type:

list of TotalEnergy | None

electrodeTransmissionSpectra(contribution=None)

Retrieve the transmission spectra for the specified electrode.

Parameters:

contribution (Left | Right) – The electrode for which to get the transmission.

Returns:

The transmission spectra. If not available, returns list of None.

Return type:

list of TransmissionSpectrum | None

energies()
Returns:

The list of energies used for the TransmissionSpectrum calculation.

Return type:

list of PhysicalQuantity of type energy

filename()
Returns:

The filename where the study object is stored.

Return type:

str

formationEnergies()
Returns:

A list of GB formation energies for all the grain boundaries. If a grain boundary calculation is not finished, a None will be returned in the list for that grain boundary. The GB formation energy calculated as (E_device / N_device - E_elec / N_elec), where E_device is the total energy of the device configuration, E_elec is the energy of the electrodes, and N_device, N_elec are the number of atoms in the device and electrodes.

Return type:

PhysicalQuantity of type energy | None

formationEnergiesBulk()
Returns:

A list of GB formation energies for all the grain boundaries calculated from the relaxed bulk grain boundary configuration. If a grain boundary calculation is not finished, a None will be returned in the list for that grain boundary. The GB formation energy calculated as (E_device / N_device - E_elec / N_elec), where E_device is the total energy of the device configuration, E_elec is the energy of the electrodes, and N_device, N_elec are the number of atoms in the device and electrodes.

Return type:

PhysicalQuantity of type energy | None

grainBoundaryGenerators()
Returns:

All the grain boundary generators given as input.

Return type:

GrainBoundaryGenerator

infinitesimal()
Returns:

The small positive energy used to move the TransmissionSpectrum calculation away from the real axis for the case of a recursion-style self-energy calculator.

Return type:

PhysicalQuantity of type energy

kpointDensity()
Returns:

The kpoint density given as input.

Return type:

PhysicalQuantity of type length.

logFilenamePrefix()
Returns:

The filename prefix for the logging output of the study.

Return type:

str | LogToStdOut

nlprint(stream=None)

Print a string containing an ASCII table useful for plotting the Study object.

Parameters:

stream (python stream) – The stream the table should be written to.
Default: NLPrintLogger()

numberOfProcessesPerTask()
Returns:

The number of processes to be used to execute each task. If None, all available processes should execute each task collaboratively.

Return type:

int | None | ProcessesPerNode

numberOfProcessesPerTaskResolved()
Returns:

The number of processes to be used to execute each task. Default values are resolved based on the current execution settings.

Return type:

int

objectId()
Returns:

The name of the study object in the file.

Return type:

str

optimizationRegionLength()
Returns:

The length along the transport direction of the optimization region. Note that the optimization region might be capped in order to keep at least one repeat of the minimal electrode fixed on either side.

Return type:

PhysicalQuantity of type length

optimizeGeometryParameters()
Returns:

The parameters used for optimizing the geometry.

Return type:

OptimizeGeometryParameters

passivateElectrodeSurfaces()
Returns:

Whether the left and right surfaces of the central region of the device along the transport direction are passivated when extracting the central region for optimization.

Return type:

bool

reflectionCoefficients(temperature=PhysicalQuantity(300.0, K), spin=<class 'NL.ComputerScienceUtilities.NLFlag.Spin.Sum'>, fermi_shift=PhysicalQuantity(0.0, eV))

Retrieve the GB reflection coefficient for all the grain boundary generators.

Parameters:
  • temperature (PhysicalQuantity of type Kelvin) – The electrode temperatures.

  • spin (Spin.Up | Spin.Down | Spin.Sum) – The spin component to query result for.
    Default: Spin.Sum

  • fermi_shift (PhysicalQuantity of type energy) – Manuel shift of the Fermi level shift.
    Default: 0.0 * eV

Returns:

The GB reflection coefficients calculated as r = 1 - R_elec / R_device, where R_device, R_elec = (R_left + R_right)/2 are the resistances of the device, and the average of left and right electrodes, respectively.

Return type:

list of float | None

relaxedGrainBoundaryDevice(grain_boundary_generator)

Retrieve the relaxed device configuration.

Parameters:

grain_boundary_generator (GrainBoundaryGenerator) – The grain boundary generator for which to get the device configuration.

Returns:

The relaxed device configuration. If not available, returns None.

Return type:

DeviceConfiguration | None

resistances(temperature=PhysicalQuantity(300.0, K), spin=<class 'NL.ComputerScienceUtilities.NLFlag.Spin.Sum'>, fermi_shift=PhysicalQuantity(0.0, eV))

Retrieve the GB scattering resistance for all the grain boundary generators.

Parameters:
  • temperature (PhysicalQuantity of type Kelvin) – The electrode temperatures.

  • spin (Spin.Up | Spin.Down | Spin.Sum) – The spin component to query result for.
    Default: Spin.Sum

  • fermi_shift (PhysicalQuantity of type energy) – Manuel shift of the Fermi level shift.
    Default: 0.0 * eV

Returns:

The GB scattering resistance calculated as R_device - 0.5*(R_left + R_right), where R_device, R_left, R_right are the resistances of the device, and left and right electrodes, respectively.

Return type:

list PhysicalQuantity of type 1/Siemens | None

saveToFileAfterUpdate()
Returns:

Whether the study is automatically saved after it is updated.

Return type:

bool

selfEnergyCalculator()
Returns:

The SelfEnergyCalculator used for the TransmissionSpectrum calculation.

Return type:

DirectSelfEnergy | RecursionSelfEnergy | SparseRecursionSelfEnergy | KrylovSelfEnergy

specificResistivities(temperature=PhysicalQuantity(300.0, K), spin=<class 'NL.ComputerScienceUtilities.NLFlag.Spin.Sum'>, fermi_shift=PhysicalQuantity(0.0, eV))

Retrieve the GB area specific resistivity for all the grain boundary generators.

Parameters:
  • temperature (PhysicalQuantity of type Kelvin) – The electrode temperatures.

  • spin (Spin.Up | Spin.Down | Spin.Sum) – The spin component to query result for.
    Default: Spin.Sum

  • fermi_shift (PhysicalQuantity of type energy) – Manuel shift of the Fermi level shift.
    Default: 0.0 * eV

Returns:

A list with the specific resistivity for each grain boundary, calculated as (R_device - 0.5*(R_left + R_right)) * A, where A is the cross sectional area, and R_device, R_left, R_right are the resistances of the device, and left and right electrodes, respectively.

Return type:

list PhysicalQuantity of type Meter^2/Siemens | None

totalEnergies()

Retrieve the total energies for the device configurations.

Returns:

The total energies. If not available, returns list of None.

Return type:

list of TotalEnergy | None

totalEnergiesBulk()

Retrieve the total energies for the bulk GB configurations.

Returns:

The total energies. If not available, returns list of None.

Return type:

list of TotalEnergy | None

totalEnergyInitialBulk()

Retrieve the total energies for the initial bulk configuration.

Returns:

The total energy If not available, returns list of None.

Return type:

list of TotalEnergy | None

transmissionKpointDensity()
Returns:

The kpoint density for the transmissson calculations given as input.

Return type:

PhysicalQuantity of type length.

transmissionSpectra()

Retrieve the transmission spectra

Returns:

The transmission spectra. If not available, returns list of None.

Return type:

list of TransmissionSpectrum | None

uniqueString()

Return a unique string representing the state of the object.

unrelaxedGrainBoundaryDevice(grain_boundary_generator)

Retrieve the unrelaxed device configuration.

Parameters:

grain_boundary_generator (GrainBoundaryGenerator) – The grain boundary generator for which to get the device configuration.

Returns:

The unrelaxed device configuration. If not available, returns None.

Return type:

DeviceConfiguration | None

update()

Run the calculations for the study object.

updatedGrainBoundaryDevice(grain_boundary_generator)

Retrieve the updated, and optimized device configuration

Parameters:

grain_boundary_generator (GrainBoundaryGenerator) – The grain boundary generator for which to get the device configuration.

Returns:

The optimized device configuration. If not available, returns None.

Return type:

DeviceConfiguration | None

weightedReflectionCoefficient(grain_boundary_weights=None, temperature=PhysicalQuantity(300.0, K), spin=<class 'NL.ComputerScienceUtilities.NLFlag.Spin.Sum'>, fermi_shift=PhysicalQuantity(0.0, eV))

Retrieve the weighted reflection coefficient.

Parameters:
  • grain_boundary_weights – The weight factors for the individual grain boundary specific resistances.
    Default: list of equal weights, i.e. numpy.ones(number_of_grain_boundaries)

  • temperature (PhysicalQuantity of type Kelvin) – The electrode temperatures.

  • spin (Spin.Up | Spin.Down | Spin.Sum) – The spin component to query result for.
    Default: Spin.Sum

  • fermi_shift (PhysicalQuantity of type energy.) – Manuel shift of the Fermi level shift.
    Default: 0.0 * eV

Returns:

The weighted GB specific resistance, calculated as sum(r_i * w_i) / sum(w_i), where r_i is the specific resistances and w_i are the weight factors.

Return type:

PhysicalQuantity of type Meter^2/Siemens

weightedSpecificResistivity(grain_boundary_weights=None, temperature=PhysicalQuantity(300.0, K), spin=<class 'NL.ComputerScienceUtilities.NLFlag.Spin.Sum'>, fermi_shift=PhysicalQuantity(0.0, eV))

Retrieve the weighted specific resistivity.

Parameters:
  • grain_boundary_weights – The weight factors for the individual grain boundary specific resistivities.
    Default: list of equal weights, i.e. numpy.ones(number_of_grain_boundaries)

  • temperature (PhysicalQuantity of type Kelvin) – The electrode temperatures.

  • spin (Spin.Up | Spin.Down | Spin.Sum) – The spin component to query result for.
    Default: Spin.Sum

  • fermi_shift (PhysicalQuantity of type energy) – Manuel shift of the Fermi level shift.
    Default: 0.0 * eV

Returns:

The weighted GB specific resistivity, calculated as sum(r_i * w_i) / sum(w_i), where r_i is the specific resistivities and w_i are the weight factors.

Return type:

PhysicalQuantity of type Meter^2/Siemens

Usage Example

This example script shows how to perform a grain boundary scattering calculation of Copper with two different grain boundaries considered. The grain boundaries are specified using the GrainBoundaryGenerator objects.

from QuantumATK import *
from SMW import *

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

# Set up lattice
lattice = FaceCenteredCubic(3.61496*Angstrom)

# Define elements
elements = [Copper]

# Define coordinates
fractional_coordinates = [[ 0.,  0.,  0.]]

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

# -------------------------------------------------------------
# Calculator
# -------------------------------------------------------------
k_point_sampling = MonkhorstPackGrid(
    na=11,
    nb=11,
    nc=11,
    )
numerical_accuracy_parameters = NumericalAccuracyParameters(
    density_mesh_cutoff=85.0*Hartree,
    k_point_sampling=k_point_sampling,
    )

calculator = LCAOCalculator(
    numerical_accuracy_parameters=numerical_accuracy_parameters,
    )

bulk_configuration.setCalculator(calculator)

# -------------------------------------------------------------
# Grain Boundary Scattering
# -------------------------------------------------------------
# Grain boundary generators.
grain_boundary_generators = []

# Grain Boundary Generator 0
grain_boundary_generator = GrainBoundaryGenerator(
    rotation_axis=[1, 1, 1],
    rotation_angle=60.00*Degrees,
    sigma=3,
    boundary_plane=[-1, -1, -1],
    minimum_electrode_length=5.0*Angstrom,
    buffer_layer=15.0*Angstrom,
    overlap_distance=0.5*Angstrom,
)
grain_boundary_generators.append(grain_boundary_generator)

# Grain Boundary Generator 1
grain_boundary_generator = GrainBoundaryGenerator(
    rotation_axis=[0, 0, 1],
    rotation_angle=36.87*Degrees,
    sigma=5,
    boundary_plane=[-2, -1, 0],
    minimum_electrode_length=5.0*Angstrom,
    buffer_layer=15.0*Angstrom,
    overlap_distance=0.5*Angstrom,
)
grain_boundary_generators.append(grain_boundary_generator)

# Optimization parameters.
optimize_geometry_parameters = OptimizeGeometryParameters(
     max_forces=0.05*eV/Angstrom,
     max_steps=200,
     max_step_length=0.2*Angstrom,
     optimizer_method=LBFGS(),
)

# Calculate scattering.
grain_boundary_scattering = GrainBoundaryScattering(
    configuration=bulk_configuration,
    filename='Copper-Grain-Boundary-Scattering.hdf5',
    object_id='grainboundaryscattering',
    grain_boundary_generators=grain_boundary_generators,
    kpoint_density=4.0*Angstrom,
    optimization_region_length=10.0*Angstrom,
    optimize_geometry_parameters=optimize_geometry_parameters,
    energies=numpy.linspace(-0.1, 0.1, 21)*eV,
    transmission_kpoint_density=12.0*Angstrom,
    self_energy_calculator=RecursionSelfEnergy(),
    infinitesimal=1e-06*eV,
    log_filename_prefix='grainboundaryscattering_',
    number_of_processes_per_task=1,
)
grain_boundary_scattering.update()

grain_boundary_scattering.py

Notes

Note

Study objects behave differently from analysis objects. See the Study object overview for more details.

The GrainBoundaryScattering object allows the user to study electonic transport and scattering properties of metal grain boundaries in an automated work flow. The purpose is to bundle a number of individual tasks into a single work flow. The grain boundaries are specified using GrainBoundaryGenerator objects, which can generate a device configuration from a bulk configuration of a cubic lattice type. As seen in the example above, the GrainBoundaryScattering is thus constructed with a bulk configuration and a list of GrainBoundaryGenerator objects.

The work flow is the following:

  • Initially the bulk configuration is optimized in order to minimize the forces and stress.

For each grain boundary structure the following tasks will be performed:

  • Construct the device configuration using the GB generator and the relaxed bulk configuration.

  • Convert the device configuration to a bulk configuration as in OptimizeDeviceConfiguration.

  • Perform structural relaxation of the bulk configuration as in OptimizeDeviceConfiguration. Only the central atoms around the grain boundary will be relaxed. The length of the region with structural relaxations is determined by the input parameters optimization_region_length.

  • Convert the relaxed bulk configuration to a device configuration.

  • Perform a self-consistent DFT-NEGF calculation for the relaxed device configuration.

  • Total energy of relaxed DeviceConfiguration.

  • Total energy of electrode configuration.

  • TransmissionSpectrum calculation for the device configuration, \(T_{dev}(E)\)

  • TransmissionSpectrum calculation for the pristine (electrode) structure (‘Sharvin resistance’), \(T_{elec}(E)\)

The grain boundary resistance is calculated as [1]

\[R_{GB}=R_{dev}-R_{elec} = \frac{1}{G_{dev}} - \frac{1}{G_{elec}}\]

where

\[G_{dev} = \int T_{dev}(E) \left(-\frac{\partial f(E)}{\partial E}\right)dE\]

is the device conductance obtained from the device transmission function, \(T_{dev}(E)\), and the energy derivative of the Fermi-Dirac distribution function. Similarly the electrode conductance and resistance can be obtained from the electrode transmission function. The electrode transmission represents the transmission through a pure metal without any defects or grain boundaries.

The specific resistivity is determined as

\[\gamma = R_{GB}\cdot A\]

where A is the cross-sectional area.

The grain boundary reflection coefficient is defined as

\[r = 1 - \frac{G_{dev}}{G_{elec}}\]

A weighted reflection coefficient can be obtained as

\[\langle r \rangle = \frac{\sum_i w_i r_i}{\sum w_i}\]

where the sum runs over all the considered grain boundaries. \(r_i\) and \(w_i\) are the reflection coefficient and relative weight of the i’th grain boundary. The weighted resistivity is defined by a similar equation. The weights, or distribution, of different defects depends on the experimental conditions (metal type, growth conditions etc.) and will be user determined.

The conductivity as function of grain boundary size can finally be obtained from the Mayadas and Shatzkes model as

\[\sigma = \sigma_0 \cdot 3 \left(\frac{1}{3} - \frac{\alpha}{2} + \alpha^2 - \alpha^3 \ln(1 + \frac{1}{\alpha} ) \right)\]
\[\alpha =\frac{l_0}{L_g}\frac{\langle r \rangle}{1 - \langle r \rangle}\]

where \(l_0\) is the bulk mean free path, \(L_g\) the average grain size, and \(\sigma_0\) is the bulk conductivity in the absence of grain boundaries. The bulk conductivity and mean free path are user inputs, but can be calculated from the Mobility object.

References