OnsagerCalc

OnsagerCalc:

The OnsagerCalc module defines the Intersitial class (for computation of interstitial-mediated diffusion), and VacancyMediated class (for computation of vacancy-mediated diffusion).

Onsager calculator module: Interstitialcy mechanism and Vacancy-mediated mechanism

Class to create an Onsager “calculator”, which brings two functionalities: 1. determines what input is needed to compute the Onsager (mobility, or L) tensors 2. constructs the function that calculates those tensors, given the input values.

This class is designed to be combined with code that can, e.g., automatically run some sort of atomistic-scale (DFT, classical potential) calculation of site energies, and energy barriers, and then in concert with scripts to convert such data into rates and probabilities, this will allow for efficient evaluation of transport coefficients.

This implementation will be for vacancy-mediated solute diffusion assumes the dilute limit. The mathematics is based on a Green function solution for the vacancy diffusion. The computation of the GF is included in the GFcalc module.

Now with HDF5 write / read capability for VacancyMediated module

class OnsagerCalc.Interstitial(crys, chem, sitelist, jumpnetwork)

A class to compute interstitial diffusivity; uses structure of crystal to do most of the heavy lifting in terms of symmetry.

Takes in a crystal that contains the interstitial as one of the chemical elements, to be specified by chem, the sitelist (list of symmetry equivalent sites), and jumpnetwork. Both of the latter can be computed automatically from crys methods, but as they are lists, can also be editted or constructed by hand.

__init__(crys, chem, sitelist, jumpnetwork)

Initialization; takes an underlying crystal, a choice of atomic chemistry, a corresponding Wyckoff site list and jump network.

Parameters

• crys – Crystal object

• chem – integer, index into the basis of crys, corresponding to the chemical element that hops

• sitelist – list of lists of indices, site indices where the atom may hop; grouped by symmetry equivalency

• jumpnetwork – list of lists of tuples: ( (i, j), dx ) symmetry unique transitions; each list is all of the possible transitions from site i to site j with jump vector dx; includes i->j and j->i

__str__()

Human readable version of diffuser

__weakref__

list of weak references to the object (if defined)

diffusivity(pre, betaene, preT, betaeneT, CalcDeriv=False)

Computes the diffusivity for our element given prefactors and energies/kB T. Also returns the negative derivative of diffusivity with respect to beta (used to compute the activation barrier tensor) if CalcDeriv = True The input list order corresponds to the sitelist and jumpnetwork

Parameters

• pre – list of prefactors for unique sites

• betaene – list of site energies divided by kB T

• preT – list of prefactors for transition states

• betaeneT – list of transition state energies divided by kB T

Return D[3,3]

diffusivity as a 3x3 tensor

Return DE[3,3]

diffusivity times activation barrier (if CalcDeriv == True)

elastodiffusion(pre, betaene, dipole, preT, betaeneT, dipoleT)

Computes the elastodiffusion tensor for our element given prefactors, energies/kB T, and elastic dipoles/kB T The input list order corresponds to the sitelist and jumpnetwork

Parameters

• pre – list of prefactors for unique sites

• betaene – list of site energies divided by kB T

• dipole – list of elastic dipoles divided by kB T

• preT – list of prefactors for transition states

• betaeneT – list of transition state energies divided by kB T

• dipoleT – list of elastic dipoles divided by kB T

Return D[3,3]

diffusivity as 3x3 tensor

Return dD[3,3,3,3]

elastodiffusion tensor as 3x3x3x3 tensor

generateJumpGroupOps()

Generates a list of group operations that transform the first jump in the jump network into all of the other members; one group operation for each.

Return siteGroupOps

list of list of group ops that mirrors the structure of jumpnetwork.

generateJumpSymmTensorBasis()

Generates a list of symmetric tensor bases for the first representative transition in our jump network

Return TensorSet

list of list of symmetric tensors

generateSiteGroupOps()

Generates a list of group operations that transform the first site in each site list into all of the other members; one group operation for each.

Return siteGroupOps

list of list of group ops that mirrors the structure of site list

generateSiteSymmTensorBasis()

Generates a list of symmetric tensor bases for the first representative site in our site list.

Return TensorSet

list of symmetric tensors

generatetags()

Create tags for unique interstitial states, and transition states.

Return tags

dictionary of tags; each is a list-of-lists

Return tagdict

dictionary that maps tag into the index of the corresponding list.

Return tagdicttype

dictionary that maps tag into the key for the corresponding list.

jumpDipoles(dipoles)

Returns a list of the elastic dipole for each transition, given the dipoles for the representatives. (“populating” the full set of dipoles)

Parameters

dipoles – list of dipoles for the first representative transition

Return dipolelist

list of lists of dipole for each jump[site][3][3]

static jumpnetworkYAML(jumpnetwork, dim=3)

Dumps a “sample” YAML formatted version of the jumpnetwork with data to be entered

losstensors(pre, betaene, dipole, preT, betaeneT)

Computes the internal friction loss tensors for our element given prefactors, energies/kB T, and elastic dipoles/kB T The input list order corresponds to the sitelist and jumpnetwork

Parameters

• pre – list of prefactors for unique sites

• betaene – list of site energies divided by kB T

• dipole – list of elastic dipoles divided by kB T

• preT – list of prefactors for transition states

• betaeneT – list of transition state energies divided by kB T

Return lambdaL

list of tuples of (eigenmode, L-tensor) where L-tensor is a 3x3x3x3 loss tensor L-tensor needs to be multiplied by kB T to have proper units of energy.

makesupercells(super_n)

Take in a supercell matrix, then generate all of the supercells needed to compute site energies and transitions (corresponding to the representatives).

Parameters

super_n – 3x3 integer matrix to define our supercell

Return superdict

dictionary of states, transitions, transmapping, and indices that correspond to dictionaries with tags.

• superdict[‘states’][i] = supercell of site;

• superdict[‘transitions’][n] = (supercell initial, supercell final);

• superdict[‘transmapping’][n] = ((site tag, groupop, mapping), (site tag, groupop, mapping))

• superdict[‘indices’][tag] = index of tag, where tag is either a state or transition tag.

ratelist(pre, betaene, preT, betaeneT)

Returns a list of lists of rates, matched to jumpnetwork

siteDipoles(dipoles)

Returns a list of the elastic dipole on each site, given the dipoles for the representatives. (“populating” the full set of dipoles)

Parameters

dipoles – list of dipoles for the first representative site

Return dipolelist

array of dipole for each site [site][3][3]

static sitelistYAML(sitelist, dim=3)

Dumps a “sample” YAML formatted version of the sitelist with data to be entered

siteprob(pre, betaene)

Returns our site probabilities, normalized, as a vector

symmratelist(pre, betaene, preT, betaeneT)

Returns a list of lists of symmetrized rates, matched to jumpnetwork

class OnsagerCalc.VacancyMediated(crys, chem, sitelist, jumpnetwork, Nthermo=0, NGFmax=4)

A class to compute vacancy-mediated solute transport coefficients, specifically L_vv (vacancy diffusion), L_ss (solute), and L_sv (off-diagonal). As part of that, it determines what quantities are needed as inputs in order to perform this calculation.

Based on crystal class. Also now includes its own GF calculator and cacheing, and storage in HDF5 format.

Requires a crystal, chemical identity of vacancy, list of symmetry-equivalent sites for that chemistry, and a jumpnetwork for the vacancy. The thermodynamic range (number of “shells” – see crystalStars.StarSet for precise definition).

GFcalculator(NGFmax=0)

Return the GF calculator; create a new one if NGFmax is being changed

Lij(bFV, bFS, bFSV, bFT0, bFT1, bFT2, large_om2=100000000.0)

Calculates the transport coefficients: L0vv, Lss, Lsv, L1vv from the scaled free energies. The Green function entries are calculated from the omega0 info. As this is the most time-consuming part of the calculation, we cache these values with a dictionary and hash function.

Parameters

• bFV[NWyckoff] – beta*eneV - ln(preV) (relative to minimum value)

• bFS[NWyckoff] – beta*eneS - ln(preS) (relative to minimum value)

• bFSV[Nthermo] – beta*eneSV - ln(preSV) (excess)

• bFT0[Nomega0] – beta*eneT0 - ln(preT0) (relative to minimum value of bFV)

• bFT1[Nomega1] – beta*eneT1 - ln(preT1) (relative to minimum value of bFV + bFS)

• bFT2[Nomega2] – beta*eneT2 - ln(preT2) (relative to minimum value of bFV + bFS)

• large_om2 – threshold for changing treatment of omega2 contributions (default: 10^8)

Return Lvv[3, 3]

vacancy-vacancy; needs to be multiplied by cv/kBT

Return Lss[3, 3]

solute-solute; needs to be multiplied by cv*cs/kBT

Return Lsv[3, 3]

solute-vacancy; needs to be multiplied by cv*cs/kBT

Return Lvv1[3, 3]

vacancy-vacancy correction due to solute; needs to be multiplied by cv*cs/kBT

__init__(crys, chem, sitelist, jumpnetwork, Nthermo=0, NGFmax=4)

Create our diffusion calculator for a given crystal structure, chemical identity, jumpnetwork (for the vacancy) and thermodynamic shell.

Parameters

• crys – Crystal object

• chem – index identifying the diffusing species

• sitelist – list, grouped into Wyckoff common positions, of unique sites

• jumpnetwork – list of unique transitions as lists of ((i,j), dx)

• Nthermo – range of thermodynamic interaction (in successive jumpnetworks)

• NGFmax – parameter controlling k-point density of GF calculator; 4 seems reasonably accurate

__str__()

Human readable version of diffuser

__weakref__

list of weak references to the object (if defined)

addhdf5(HDF5group)

Adds an HDF5 representation of object into an HDF5group (needs to already exist).

Example: if f is an open HDF5, then VacancyMediated.addhdf5(f.create_group(‘Diffuser’)) will (1) create the group named ‘Diffuser’, and then (2) put the VacancyMediated representation in that group.

Parameters

HDF5group – HDF5 group

clearcache()

Clear out the GF cache values

generate(Nthermo)

Generate the necessary stars, vector-stars, and jump networks based on the thermodynamic range.

Parameters

Nthermo – range of thermodynamic interactions, in terms of “shells”, which is multiple summations of jumpvect

generatematrices()

Generates all the matrices and “helper” pieces, based on our jump networks. This has been separated out in case the user wants to, e.g., prune / modify the networks after they’ve been created with generate(), then generatematrices() can be rerun.

generatetags()

Create tags for vacancy states, solute states, solute-vacancy complexes; omega0, omega1, and omega2 transition states.

Return tags

dictionary of tags; each is a list-of-lists

Return tagdict

dictionary that maps tag into the index of the corresponding list.

Return tagdicttype

dictionary that maps tag into the key for the corresponding list.

interactlist()

Return a list of solute-vacancy configurations for interactions. The points correspond to a vector between a solute atom and a vacancy. Defined by Stars.

Return statelist

list of PairStates for the solute-vacancy interactions

classmethod loadhdf5(HDF5group)

Creates a new VacancyMediated diffuser from an HDF5 group.

Parameters

HDFgroup – HDF5 group

Return VacancyMediated

new VacancyMediated diffuser object from HDF5

makeLIMBpreene(preS, eneS, preSV, eneSV, preT0, eneT0, **ignoredextraarguments)

Generates corresponding energies / prefactors for corresponding to LIMB (Linearized interpolation of migration barrier approximation). Returns a dictionary. (we ignore extra arguments so that a dictionary including additional entries can be passed)

Parameters

• preS[NWyckoff] – prefactor for solute formation

• eneS[NWyckoff] – solute formation energy

• preSV[Nthermo] – prefactor for solute-vacancy interaction

• eneSV[Nthermo] – solute-vacancy binding energy

• preT0[Nomeg0] – prefactor for vacancy jump transitions (follows jumpnetwork)

• eneT0[Nomega0] – transition energy for vacancy jumps

Return preT1[Nomega1]

prefactor for omega1-style transitions (follows om1_jn)

Return eneT1[Nomega1]

transition energy/kBT for omega1-style jumps

Return preT2[Nomega2]

prefactor for omega2-style transitions (follows om2_jn)

Return eneT2[Nomega2]

transition energy/kBT for omega2-style jumps

makesupercells(super_n)

Take in a supercell matrix, then generate all of the supercells needed to compute site energies and transitions (corresponding to the representatives).

Note: the states are lone vacancy, lone solute, solute-vacancy complexes in our thermodynamic range. Note that there will be escape states are endpoints of some omega1 jumps. They are not relaxed, and have no pre-existing tag. They will only be output as a single endpoint of an NEB run; there may be symmetry equivalent duplicates, as we construct these supercells on an as needed basis.

We’ve got a few classes of warnings (from most egregious to least) that can issued if the supercell is too small; the analysis will continue despite any warnings:

1. Thermodynamic shell states map to different states in supercell

2. Thermodynamic shell states are not unique in supercell (multiplicity)

3. Kinetic shell states map to different states in supercell

4. Kinetic shell states are not unique in supercell (multiplicity)

The lowest level can still be run reliably but runs the risk of errors in escape transition barriers. Extreme caution should be used if any of the other warnings are raised.

Parameters

super_n – 3x3 integer matrix to define our supercell

Return superdict

dictionary of states, transitions, transmapping, indices that correspond to dictionaries with tags; the final tag reference is the basesupercell for calculations without defects.

• superdict[‘states’][i] = supercell of state;

• superdict[‘transitions’][n] = (supercell initial, supercell final);

• superdict[‘transmapping’][n] = ((site tag, groupop, mapping), (site tag, groupop, mapping))

• superdict[‘indices’][tag] = (type, index) of tag, where tag is either a state or transition tag.

• superdict[‘reference’] = supercell reference, no defects

maketracerpreene(preT0, eneT0, **ignoredextraarguments)

Generates corresponding energies / prefactors for an isotopic tracer. Returns a dictionary. (we ignore extra arguments so that a dictionary including additional entries can be passed)

Parameters

• preT0[Nomeg0] – prefactor for vacancy jump transitions (follows jumpnetwork)

• eneT0[Nomega0] – transition energy state for vacancy jumps

Return preS[NWyckoff]

prefactor for solute formation

Return eneS[NWyckoff]

solute formation energy

Return preSV[Nthermo]

prefactor for solute-vacancy interaction

Return eneSV[Nthermo]

solute-vacancy binding energy

Return preT1[Nomega1]

prefactor for omega1-style transitions (follows om1_jn)

Return eneT1[Nomega1]

transition energy for omega1-style jumps

Return preT2[Nomega2]

prefactor for omega2-style transitions (follows om2_jn)

Return eneT2[Nomega2]

transition energy for omega2-style jumps

omegalist(fivefreqindex=1)

Return a list of pairs of endpoints for a vacancy jump, corresponding to omega1 or omega2 Solute at the origin, vacancy hopping between two sites. Defined by om1_jumpnetwork

Parameters

fivefreqindex – 1 or 2, corresponding to omega1 or omega2

Return omegalist

list of tuples of PairStates

Return omegajumptype

index of corresponding omega0 jumptype

static preene2betafree(kT, preV, eneV, preS, eneS, preSV, eneSV, preT0, eneT0, preT1, eneT1, preT2, eneT2, **ignoredextraarguments)

Read in a series of prefactors (\(\exp(S/k_\text{B})\)) and energies, and return \(\beta F\) for energies and transition state energies. Used to provide scaled values to Lij(). Can specify all of the entries using a dictionary; e.g., preene2betafree(kT, **data_dict) and then send that output as input to Lij: Lij(*preene2betafree(kT, **data_dict)) (we ignore extra arguments so that a dictionary including additional entries can be passed)

Parameters

• kT – temperature times Boltzmann’s constant kB

• preV – prefactor for vacancy formation (prod of inverse vibrational frequencies)

• eneV – vacancy formation energy

• preS – prefactor for solute formation (prod of inverse vibrational frequencies)

• eneS – solute formation energy

• preSV – excess prefactor for solute-vacancy binding

• eneSV – solute-vacancy binding energy

• preT0 – prefactor for vacancy transition state

• eneT0 – energy for vacancy transition state (relative to eneV)

• preT1 – prefactor for vacancy swing transition state

• eneT1 – energy for vacancy swing transition state (relative to eneV + eneS + eneSV)

• preT2 – prefactor for vacancy exchange transition state

• eneT2 – energy for vacancy exchange transition state (relative to eneV + eneS + eneSV)

Return bFV

beta*eneV - ln(preV) (relative to minimum value)

Return bFS

beta*eneS - ln(preS) (relative to minimum value)

Return bFSV

beta*eneSV - ln(preSV) (excess)

Return bFT0

beta*eneT0 - ln(preT0) (relative to minimum value of bFV)

Return bFT1

beta*eneT1 - ln(preT1) (relative to minimum value of bFV + bFS)

Return bFT2

beta*eneT2 - ln(preT2) (relative to minimum value of bFV + bFS)

tags2preene(usertagdict, VERBOSE=False)

Generates energies and prefactors based on a dictionary of tags.

Parameters

• usertagdict – dictionary where the keys are tags, and the values are tuples: (pre, ene)

• VERBOSE – (optional) if True, also return a dictionary of missing tags, duplicate tags, and bad tags

Return thermodict

dictionary of ene’s and pre’s corresponding to usertagdict

Return missingdict

dictionary with keys corresponding to tag types, and the values are lists of lists of symmetry equivalent tags that are missing

Return duplicatelist

list of lists of tags in usertagdict that are (symmetry) duplicates

Return badtaglist

list of all tags in usertagdict that aren’t found in our dictionary

OnsagerCalc.arrays2vTKdict(vTKarray, valarray, vTKsplits)

Takes two arrays of vTK keys and values, and the splits to separate vTKarray back into vTK and returns a dictionary indexed by the vTK.

Parameters

• vTKarray – array of vTK entries

• valarray – array of values

• vTKsplits – split placement for vTK entries

Return vTKdict

dictionary, indexed by vTK objects, whose entries are arrays

OnsagerCalc.vTKdict2arrays(vTKdict)

Takes a dictionary indexed by vTK objects, returns two arrays of vTK keys and values, and the splits to separate vTKarray back into vTK

Parameters

vTKdict – dictionary, indexed by vTK objects, whose entries are arrays

Return vTKarray

array of vTK entries

Return valarray

array of values

Return vTKsplits

split placement for vTK entries

class OnsagerCalc.vacancyThermoKinetics

Class to store (in a hashable manner) the thermodynamics and kinetics for the vacancy

Parameters

• pre – prefactors for sites

• betaene – energy for sites / kBT

• preT – prefactors for transition states

• betaeneT – transition state energy for sites / kBT

__eq__(other)

Return self==value.

__hash__()

Return hash(self).

__ne__(other)

Return self!=value.

__repr__()

Return a nicely formatted representation string

static vacancyThermoKinetics_constructor(loader, node)

Construct a GroupOp from YAML

static vacancyThermoKinetics_representer(dumper, data)

Output a PairState