Descriptor

The descriptor module contains methods for describing the local atomic environment; that is, feature fectors that can be fed to machine-learning modules.

Gaussian

class amp.descriptor.gaussian.FingerprintCalculator(neighborlist, Gs, cutoff, fortran)[source]

For integration with .utilities.Data

Parameters:
  • neighborlist (list of str) – List of neighbors.
  • Gs (dict) –

    Dictionary of symbols and lists of dictionaries for making symmetry functions. Either auto-genetrated, or given in the following form, for example:

    >>> Gs = {"O": [{"type":"G2", "element":"O", "eta":10.},
    ...             {"type":"G4", "elements":["O", "Au"],
    ...              "eta":5., "gamma":1., "zeta":1.0}],
    ...       "Au": [{"type":"G2", "element":"O", "eta":2.},
    ...              {"type":"G4", "elements":["O", "Au"],
    ...               "eta":2., "gamma":1., "zeta":5.0}]}
    
  • cutoff (float) – Radius above which neighbor interactions are ignored.
  • fortran (bool) – If True, will use fortran modules, if False, will not.
calculate(image, key)[source]

Makes a list of fingerprints, one per atom, for the fed image.

Parameters:
  • image (object) – ASE atoms object.
  • key (str) – key of the image after being hashed.
get_fingerprint(index, symbol, neighborsymbols, neighborpositions)[source]

Returns the fingerprint of symmetry function values for atom specified by its index and symbol.

neighborsymbols and neighborpositions are lists of neighbors’ symbols and Cartesian positions, respectively.

Parameters:
  • index (int) – Index of the center atom.
  • symbol (str) – Symbol of the center atom.
  • neighborsymbols (list of str) – List of neighbors’ symbols.
  • neighborpositions (list of list of float) – List of Cartesian atomic positions.
Returns:

symbol, fingerprint – fingerprints for atom specified by its index and symbol.

Return type:

list of float

class amp.descriptor.gaussian.FingerprintPrimeCalculator(neighborlist, Gs, cutoff, fortran)[source]

For integration with .utilities.Data

Parameters:
  • neighborlist (list of str) – List of neighbors.
  • Gs (dict) –

    Dictionary of symbols and lists of dictionaries for making symmetry functions. Either auto-genetrated, or given in the following form, for example:

    >>> Gs = {"O": [{"type":"G2", "element":"O", "eta":10.},
    ...             {"type":"G4", "elements":["O", "Au"],
    ...              "eta":5., "gamma":1., "zeta":1.0}],
    ...       "Au": [{"type":"G2", "element":"O", "eta":2.},
    ...              {"type":"G4", "elements":["O", "Au"],
    ...               "eta":2., "gamma":1., "zeta":5.0}]}
    
  • cutoff (float) – Radius above which neighbor interactions are ignored.
  • fortran (bool) – If True, will use fortran modules, if False, will not.
calculate(image, key)[source]

Makes a list of fingerprint derivatives, one per atom, for the fed image.

Parameters:
  • image (object) – ASE atoms object.
  • key (str) – key of the image after being hashed.
get_fingerprintprime(index, symbol, neighborindices, neighborsymbols, neighborpositions, m, l)[source]

Returns the value of the derivative of G for atom with index and symbol with respect to coordinate x_{l} of atom index m.

neighborindices, neighborsymbols and neighborpositions are lists of neighbors’ indices, symbols and Cartesian positions, respectively.

Parameters:
  • index (int) – Index of the center atom.
  • symbol (str) – Symbol of the center atom.
  • neighborindices (list of int) – List of neighbors’ indices.
  • neighborsymbols (list of str) – List of neighbors’ symbols.
  • neighborpositions (list of list of float) – List of Cartesian atomic positions.
  • m (int) – Index of the pair atom.
  • l (int) – Direction of the derivative; is an integer from 0 to 2.
Returns:

fingerprintprime – The value of the derivative of the fingerprints for atom with index and symbol with respect to coordinate x_{l} of atom index m.

Return type:

list of float

class amp.descriptor.gaussian.Gaussian(cutoff=<Cosine cutoff with Rc=6.500 from amp.descriptor.cutoffs>, Gs=None, dblabel=None, elements=None, version=None, fortran=True, mode='atom-centered')[source]

Bases: object

Class that calculates Gaussian fingerprints (i.e., Behler-style).

Parameters:
  • cutoff (object or float) – Cutoff function, typically from amp.descriptor.cutoffs. Can be also fed as a float representing the radius above which neighbor interactions are ignored; in this case a cosine cutoff function will be employed. Default is a 6.5-Angstrom cosine cutoff.
  • Gs (dict) –

    Dictionary of symbols and lists of dictionaries for making symmetry functions. Either auto-genetrated, or given in the following form, for example:

    >>> Gs = {"O": [{"type":"G2", "element":"O", "eta":10.},
    ...             {"type":"G4", "elements":["O", "Au"],
    ...              "eta":5., "gamma":1., "zeta":1.0}],
    ...       "Au": [{"type":"G2", "element":"O", "eta":2.},
    ...              {"type":"G4", "elements":["O", "Au"],
    ...               "eta":2., "gamma":1., "zeta":5.0}]}
    
  • dblabel (str) – Optional separate prefix/location for database files, including fingerprints, fingerprint derivatives, and neighborlists. This file location can be shared between calculator instances to avoid re-calculating redundant information. If not supplied, just uses the value from label.
  • elements (list) – List of allowed elements present in the system. If not provided, will be found automatically.
  • version (str) – Version of fingerprints.
  • fortran (bool) – If True, will use fortran modules, if False, will not.
  • mode (str) – Can be either ‘atom-centered’ or ‘image-centered’.
Raises:

RuntimeError

calculate_fingerprints(images, parallel=None, log=None, calculate_derivatives=False)[source]

Calculates the fingerpints of the images, for the ones not already done.

Parameters:
  • images (list or str) – List of ASE atoms objects with positions, symbols, energies, and forces in ASE format. This is the training set of data. This can also be the path to an ASE trajectory (.traj) or database (.db) file. Energies can be obtained from any reference, e.g. DFT calculations.
  • parallel (dict) – Configuration for parallelization. Should be in same form as in amp.Amp.
  • log (Logger object) – Write function at which to log data. Note this must be a callable function.
  • calculate_derivatives (bool) – Decides whether or not fingerprintprimes should also be calculated.
tostring()[source]

Returns an evaluatable representation of the calculator that can be used to restart the calculator.

amp.descriptor.gaussian.Kronecker(i, j)[source]

Kronecker delta function.

Parameters:
  • i (int) – First index of Kronecker delta.
  • j (int) – Second index of Kronecker delta.
Returns:

The value of the Kronecker delta.

Return type:

int

class amp.descriptor.gaussian.NeighborlistCalculator(cutoff)[source]

For integration with .utilities.Data

For each image fed to calculate, a list of neighbors with offset distances is returned.

Parameters:cutoff (float) – Radius above which neighbor interactions are ignored.
calculate(image, key)[source]

For integration with .utilities.Data

For each image fed to calculate, a list of neighbors with offset distances is returned.

Parameters:
  • image (object) – ASE atoms object.
  • key (str) – key of the image after being hashed.
amp.descriptor.gaussian.calculate_G2(neighborsymbols, neighborpositions, G_element, eta, cutoff, Ri, fortran)[source]

Calculate G2 symmetry function.

Ideally this will not be used but will be a template for how to build the fortran version (and serves as a slow backup if the fortran one goes uncompiled). See Eq. 13a of the supplementary information of Khorshidi, Peterson, CPC(2016).

Parameters:
  • neighborsymbols (list of str) – List of symbols of all neighbor atoms.
  • neighborpositions (list of list of float) – List of Cartesian atomic positions.
  • G_element (dict) – Symmetry functions of the center atom.
  • eta (float) – Parameter of Gaussian symmetry functions.
  • cutoff (float) – Radius above which neighbor interactions are ignored.
  • Ri (int) – Index of the center atom.
  • fortran (bool) – If True, will use the fortran subroutines, else will not.
Returns:

ridge – G2 fingerprint.

Return type:

float

amp.descriptor.gaussian.calculate_G2_prime(neighborindices, neighborsymbols, neighborpositions, G_element, eta, cutoff, i, Ri, m, l, fortran)[source]

Calculates coordinate derivative of G2 symmetry function for atom at index i and position Ri with respect to coordinate x_{l} of atom index m.

See Eq. 13b of the supplementary information of Khorshidi, Peterson, CPC(2016).

Parameters:
  • neighborindices (list of int) – List of int of neighboring atoms.
  • neighborsymbols (list of str) – List of symbols of neighboring atoms.
  • neighborpositions (list of list of float) – List of Cartesian atomic positions of neighboring atoms.
  • G_element (dict) – Symmetry functions of the center atom.
  • eta (float) – Parameter of Behler symmetry functions.
  • cutoff (object or float) – Cutoff function, typically from amp.descriptor.cutoffs. Can be also fed as a float representing the radius above which neighbor interactions are ignored; in this case a cosine cutoff function will be employed. Default is a 6.5-Angstrom cosine cutoff.
  • i (int) – Index of the center atom.
  • Ri (float) – Position of the center atom.
  • m (int) – Index of the atom force is acting on.
  • l (int) – Direction of force.
  • fortran (bool) – If True, will use the fortran subroutines, else will not.
Returns:

ridge – Coordinate derivative of G2 symmetry function for atom at index a and position Ri with respect to coordinate x_{l} of atom index m.

Return type:

float

amp.descriptor.gaussian.calculate_G4(neighborsymbols, neighborpositions, G_elements, gamma, zeta, eta, cutoff, Ri, fortran)[source]

Calculate G4 symmetry function.

Ideally this will not be used but will be a template for how to build the fortran version (and serves as a slow backup if the fortran one goes uncompiled). See Eq. 13c of the supplementary information of Khorshidi, Peterson, CPC(2016).

Parameters:
  • neighborsymbols (list of str) – List of symbols of neighboring atoms.
  • neighborpositions (list of list of float) – List of Cartesian atomic positions of neighboring atoms.
  • G_elements (dict) – Symmetry functions of the center atom.
  • gamma (float) – Parameter of Gaussian symmetry functions.
  • zeta (float) – Parameter of Gaussian symmetry functions.
  • eta (float) – Parameter of Gaussian symmetry functions.
  • cutoff (object or float) – Cutoff function, typically from amp.descriptor.cutoffs. Can be also fed as a float representing the radius above which neighbor interactions are ignored; in this case a cosine cutoff function will be employed. Default is a 6.5-Angstrom cosine cutoff.
  • Ri (int) – Index of the center atom.
  • fortran (bool) – If True, will use the fortran subroutines, else will not.
Returns:

ridge – G4 fingerprint.

Return type:

float

amp.descriptor.gaussian.calculate_G4_prime(neighborindices, neighborsymbols, neighborpositions, G_elements, gamma, zeta, eta, cutoff, i, Ri, m, l, fortran)[source]

Calculates coordinate derivative of G4 symmetry function for atom at index i and position Ri with respect to coordinate x_{l} of atom index m.

See Eq. 13d of the supplementary information of Khorshidi, Peterson, CPC(2016).

Parameters:
  • neighborindices (list of int) – List of int of neighboring atoms.
  • neighborsymbols (list of str) – List of symbols of neighboring atoms.
  • neighborpositions (list of list of float) – List of Cartesian atomic positions of neighboring atoms.
  • G_elements (dict) – Symmetry functions of the center atom.
  • gamma (float) – Parameter of Behler symmetry functions.
  • zeta (float) – Parameter of Behler symmetry functions.
  • eta (float) – Parameter of Behler symmetry functions.
  • cutoff (object or float) – Cutoff function, typically from amp.descriptor.cutoffs. Can be also fed as a float representing the radius above which neighbor interactions are ignored; in this case a cosine cutoff function will be employed. Default is a 6.5-Angstrom cosine cutoff.
  • i (int) – Index of the center atom.
  • Ri (float) – Position of the center atom.
  • m (int) – Index of the atom force is acting on.
  • l (int) – Direction of force.
  • fortran (bool) – If True, will use the fortran subroutines, else will not.
Returns:

ridge – Coordinate derivative of G4 symmetry function for atom at index i and position Ri with respect to coordinate x_{l} of atom index m.

Return type:

float

amp.descriptor.gaussian.dCos_theta_ijk_dR_ml(i, j, k, Ri, Rj, Rk, m, l)[source]

Returns the derivative of Cos(theta_{ijk}) with respect to x_{l} of atomic index m.

See Eq. 14f of the supplementary information of Khorshidi, Peterson, CPC(2016).

Parameters:
  • i (int) – Index of the center atom.
  • j (int) – Index of the first atom.
  • k (int) – Index of the second atom.
  • Ri (float) – Position of the center atom.
  • Rj (float) – Position of the first atom.
  • Rk (float) – Position of the second atom.
  • m (int) – Index of the atom force is acting on.
  • l (int) – Direction of force.
Returns:

dCos_theta_ijk_dR_ml – Derivative of Cos(theta_{ijk}) with respect to x_{l} of atomic index m.

Return type:

float

amp.descriptor.gaussian.dRij_dRml(i, j, Ri, Rj, m, l)[source]

Returns the derivative of the norm of position vector R_{ij} with respect to coordinate x_{l} of atomic index m.

See Eq. 14c of the supplementary information of Khorshidi, Peterson, CPC(2016).

Parameters:
  • i (int) – Index of the first atom.
  • j (int) – Index of the second atom.
  • Ri (float) – Position of the first atom.
  • Rj (float) – Position of the second atom.
  • m (int) – Index of the atom force is acting on.
  • l (int) – Direction of force.
Returns:

dRij_dRml – The derivative of the noRi of position vector R_{ij} with respect to x_{l} of atomic index m.

Return type:

list of float

amp.descriptor.gaussian.dRij_dRml_vector(i, j, m, l)[source]

Returns the derivative of the position vector R_{ij} with respect to x_{l} of itomic index m.

See Eq. 14d of the supplementary information of Khorshidi, Peterson, CPC(2016).

Parameters:
  • i (int) – Index of the first atom.
  • j (int) – Index of the second atom.
  • m (int) – Index of the atom force is acting on.
  • l (int) – Direction of force.
Returns:

The derivative of the position vector R_{ij} with respect to x_{l} of atomic index m.

Return type:

list of float

amp.descriptor.gaussian.make_symmetry_functions(elements)[source]

Makes symmetry functions as in Nano Letters function by Artrith.

Elements is a list of the elements, as in: [“C”, “O”, “H”, “Cu”]. G[0] = {“type”:”G2”, “element”: “O”, “eta”: 0.0009} G[40] = {“type”:”G4”, “elements”: [“O”, “Au”], “eta”: 0.0001, “gamma”: 1.0, “zeta”: 1.0}

If G (a list) is fed in, this will add to it and return an expanded version. If not, it will create a new one.

Parameters:elements (list of str) – List of symbols of all atoms.
Returns:G – Symmetry functions if not given by the user.
Return type:dict of lists

Zernike

class amp.descriptor.zernike.FingerprintCalculator(neighborlist, Gs, nmax, cutoff, fortran)[source]

For integration with .utilities.Data

calculate(image, key)[source]

Makes a list of fingerprints, one per atom, for the fed image.

Parameters:
  • image (object) – ASE atoms object.
  • key (str) – Key of the image after being hashed.
get_fingerprint(index, symbol, n_symbols, Rs)[source]

Returns the fingerprint of symmetry function values for atom specified by its index and symbol.

n_symbols and Rs are lists of neighbors’ symbols and Cartesian positions, respectively.

Parameters:
  • index (int) – Index of the center atom.
  • symbol (str) – Symbol of the center atom.
  • n_symbols (list of str) – List of neighbors’ symbols.
  • Rs (list of list of float) – List of Cartesian atomic positions of neighbors.
Returns:

symbols, fingerprints – Fingerprints for atom specified by its index and symbol.

Return type:

list of float

class amp.descriptor.zernike.FingerprintPrimeCalculator(neighborlist, Gs, nmax, cutoff, fortran)[source]

For integration with .utilities.Data

calculate(image, key)[source]

Makes a list of fingerprint derivatives, one per atom, for the fed image.

Parameters:
  • image (object) – ASE atoms object.
  • key (str) – Key of the image after being hashed.
get_fingerprintprime(index, symbol, n_indices, n_symbols, Rs, p, q)[source]

Returns the value of the derivative of G for atom with index and symbol with respect to coordinate x_{i} of atom index m. n_indices, n_symbols and Rs are lists of neighbors’ indices, symbols and Cartesian positions, respectively.

Parameters:
  • index (int) – Index of the center atom.
  • symbol (str) – Symbol of the center atom.
  • n_indices (list of int) – List of neighbors’ indices.
  • n_symbols (list of str) – List of neighbors’ symbols.
  • Rs (list of list of float) – List of Cartesian atomic positions.
  • p (int) – Index of the pair atom.
  • q (int) – Direction of the derivative; is an integer from 0 to 2.
Returns:

fingerprint_prime – The value of the derivative of the fingerprints for atom with index and symbol with respect to coordinate x_{i} of atom index m.

Return type:

list of float

amp.descriptor.zernike.Kronecker(i, j)[source]

Kronecker delta function.

i : int
First index of Kronecker delta.
j : int
Second index of Kronecker delta.
Returns:Kronecker delta
Return type:int
class amp.descriptor.zernike.NeighborlistCalculator(cutoff)[source]

For integration with .utilities.Data

For each image fed to calculate, a list of neighbors with offset distances is returned.

Parameters:cutoff (object or float) – Cutoff function, typically from amp.descriptor.cutoffs. Can be also fed as a float representing the radius above which neighbor interactions are ignored; in this case a cosine cutoff function will be employed. Default is a 6.5-Angstrom cosine cutoff.
calculate(image, key)[source]

For integration with .utilities.Data For each image fed to calculate, a list of neighbors with offset distances is returned.

Parameters:
  • image (object) – ASE atoms object.
  • key (str) – Key of the image after being hashed.
class amp.descriptor.zernike.Zernike(cutoff=<Cosine cutoff with Rc=6.500 from amp.descriptor.cutoffs>, Gs=None, nmax=5, dblabel=None, elements=None, version='2016.02', mode='atom-centered', fortran=True)[source]

Bases: object

Class that calculates Zernike fingerprints.

Parameters:
  • cutoff (object or float) – Cutoff function, typically from amp.descriptor.cutoffs. Can be also fed as a float representing the radius above which neighbor interactions are ignored; in this case a cosine cutoff function will be employed. Default is a 6.5-Angstrom cosine cutoff.
  • Gs (dict) –

    Dictionary of symbols and dictionaries for making symmetry functions. Either auto-genetrated, or given in the following form, for example:

    >>> Gs = {"Au": {"Au": 3., "O": 2.}, "O": {"Au": 5., "O": 10.}}
    
  • nmax (integer or dict) – Maximum degree of Zernike polynomials that will be included in the fingerprint vector. Can be different values for different species fed as a dictionary with chemical elements as keys.
  • dblabel (str) – Optional separate prefix/location for database files, including fingerprints, fingerprint derivatives, and neighborlists. This file location can be shared between calculator instances to avoid re-calculating redundant information. If not supplied, just uses the value from label.
  • elements (list) – List of allowed elements present in the system. If not provided, will be found automatically.
  • version (str) – Version of fingerprints.
  • mode (str) – Can be either ‘atom-centered’ or ‘image-centered’.
  • fortran (bool) – If True, will use fortran modules, if False, will not.
Raises:

RuntimeError, TypeError

calculate_fingerprints(images, parallel=None, log=None, calculate_derivatives=False)[source]

Calculates the fingerpints of the images, for the ones not already done.

Parameters:
  • images (list or str) – List of ASE atoms objects with positions, symbols, energies, and forces in ASE format. This is the training set of data. This can also be the path to an ASE trajectory (.traj) or database (.db) file. Energies can be obtained from any reference, e.g. DFT calculations.
  • parallel (dict) – Configuration for parallelization. Should be in same form as in amp.Amp.
  • log (Logger object) – Write function at which to log data. Note this must be a callable function.
  • calculate_derivatives (bool) – Decides whether or not fingerprintprimes should also be calculated.
tostring()[source]

Returns an evaluatable representation of the calculator that can be used to restart the calculator.

amp.descriptor.zernike.binomial(n, k, factorial)[source]

Returns C(n,k) = n!/(k!(n-k)!).

amp.descriptor.zernike.calculate_R(n, l, rho, factorial)[source]

Calculates R_{n}^{l}(rho) according to the last equation of wikipedia.

amp.descriptor.zernike.calculate_Z(n, l, m, x, y, z, factorial)[source]

Calculates Z_{nl}^{m}(x, y, z) according to the unnumbered equation afer Eq. (11) of “3D Zernike Descriptors for Content Based Shape Retrieval”, Computer-Aided Design 36 (2004) 1047-1062.

amp.descriptor.zernike.calculate_Z_prime(n, l, m, x, y, z, p, factorial)[source]

Calculates dZ_{nl}^{m}(x, y, z)/dR_{p} according to the unnumbered equation afer Eq. (11) of “3D Zernike Descriptors for Content Based Shape Retrieval”, Computer-Aided Design 36 (2004) 1047-1062.

amp.descriptor.zernike.calculate_q(nu, k, l, factorial)[source]

Calculates q_{kl}^{nu} according to the unnumbered equation afer Eq. (7) of “3D Zernike Descriptors for Content Based Shape Retrieval”, Computer-Aided Design 36 (2004) 1047-1062.

amp.descriptor.zernike.der_position(m, n, Rm, Rn, l, i)[source]

Returns the derivative of the norm of position vector R_{mn} with respect to x_{i} of atomic index l.

Parameters:
  • m (int) – Index of the first atom.
  • n (int) – Index of the second atom.
  • Rm (float) – Position of the first atom.
  • Rn (float) – Position of the second atom.
  • l (int) – Index of the atom force is acting on.
  • i (int) – Direction of force.
Returns:

der_position – The derivative of the norm of position vector R_{mn} with respect to x_{i} of atomic index l.

Return type:

list of float

amp.descriptor.zernike.generate_coefficients(elements)[source]

Automatically generates coefficients if not given by the user.

Parameters:elements (list of str) – List of symbols of all atoms.
Returns:G
Return type:dict of dicts

Bispectrum

class amp.descriptor.bispectrum.Bispectrum(cutoff=<Cosine cutoff with Rc=6.500 from amp.descriptor.cutoffs>, Gs=None, jmax=5, dblabel=None, elements=None, version='2016.02', mode='atom-centered')[source]

Bases: object

Class that calculates spherical harmonic bispectrum fingerprints.

Parameters:
  • cutoff (object or float) – Cutoff function, typically from amp.descriptor.cutoffs. Can be also fed as a float representing the radius above which neighbor interactions are ignored; in this case a cosine cutoff function will be employed. Default is a 6.5-Angstrom cosine cutoff.
  • Gs (dict) –

    Dictionary of symbols and dictionaries for making fingerprints. Either auto-genetrated, or given in the following form, for example:

    >>> Gs = {"Au": {"Au": 3., "O": 2.}, "O": {"Au": 5., "O": 10.}}
    
  • jmax (integer or half-integer or dict) – Maximum degree of spherical harmonics that will be included in the fingerprint vector. Can be also fed as a dictionary with chemical species as keys.
  • dblabel (str) – Optional separate prefix/location for database files, including fingerprints, fingerprint derivatives, and neighborlists. This file location can be shared between calculator instances to avoid re-calculating redundant information. If not supplied, just uses the value from label.
  • elements (list) – List of allowed elements present in the system. If not provided, will be found automatically.
  • version (str) – Version of fingerprints.
  • Raises
  • ------- – RuntimeError, TypeError
calculate_fingerprints(images, parallel=None, log=None, calculate_derivatives=False)[source]

Calculates the fingerpints of the images, for the ones not already done.

Parameters:
  • images (list or str) – List of ASE atoms objects with positions, symbols, energies, and forces in ASE format. This is the training set of data. This can also be the path to an ASE trajectory (.traj) or database (.db) file. Energies can be obtained from any reference, e.g. DFT calculations.
  • parallel (dict) – Configuration for parallelization. Should be in same form as in amp.Amp.
  • log (Logger object) – Write function at which to log data. Note this must be a callable function.
  • calculate_derivatives (bool) – Decides whether or not fingerprintprimes should also be calculated.
tostring()[source]

Returns an evaluatable representation of the calculator that can be used to restart the calculator.

amp.descriptor.bispectrum.CG(a, alpha, b, beta, c, gamma, factorial)[source]

Clebsch-Gordan coefficient C_{a alpha b beta}^{c gamma} is calculated acoording to the expression given in Varshalovich Eq. (3), Section 8.2, Page 238.

class amp.descriptor.bispectrum.FingerprintCalculator(neighborlist, Gs, jmax, cutoff)[source]

For integration with .utilities.Data

calculate(image, key)[source]

Makes a list of fingerprints, one per atom, for the fed image.

Parameters:
  • image (object) – ASE atoms object.
  • key (str) – key of the image after being hashed.
get_fingerprint(index, symbol, n_symbols, Rs)[source]

Returns the fingerprint of symmetry function values for atom specified by its index and symbol.

n_symbols and Rs are lists of neighbors’ symbols and Cartesian positions, respectively.

Parameters:
  • index (int) – Index of the center atom.
  • symbol (str) – Symbol of the center atom.
  • n_symbols (list of str) – List of neighbors’ symbols.
  • Rs (list of list of float) – List of Cartesian atomic positions of neighbors.
Returns:

symbols, fingerprints – fingerprints for atom specified by its index and symbol.

Return type:

list of float

class amp.descriptor.bispectrum.NeighborlistCalculator(cutoff)[source]

For integration with .utilities.Data

For each image fed to calculate, a list of neighbors with offset distances is returned.

calculate(image, key)[source]
amp.descriptor.bispectrum.U(j, m, mp, omega, theta, phi, factorial)[source]

Calculates rotation matrix U_{MM’}^{J} in terms of rotation angle omega as well as rotation axis angles theta and phi, according to Varshalovich, Eq. (3), Section 4.5, Page 81. j, m, mp, and mpp here are J, M, M’, and M’’ in Eq. (3).

amp.descriptor.bispectrum.WignerD(j, m, mp, alpha, beta, gamma, factorial)[source]

Returns the Wigner-D matrix. alpha, beta, and gamma are the Euler angles.

amp.descriptor.bispectrum.binomial(n, k, factorial)[source]

Returns C(n,k) = n!/(k!(n-k)!).

amp.descriptor.bispectrum.calculate_B(j1, j2, j, G_element, cutoff, cutofffn, factorial, n_symbols, rs, psis, thetas, phis)[source]

Calculates bi-spectrum B_{j1, j2, j} according to Eq. (5) of “Gaussian Approximation Potentials: The Accuracy of Quantum Mechanics, without the Electrons”, Phys. Rev. Lett. 104, 136403.

amp.descriptor.bispectrum.calculate_c(j, mp, m, G_element, cutoff, cutofffn, factorial, n_symbols, rs, psis, thetas, phis)[source]

Calculates c^{j}_{m’m} according to Eq. (4) of “Gaussian Approximation Potentials: The Accuracy of Quantum Mechanics, without the Electrons”, Phys. Rev. Lett. 104, 136403

amp.descriptor.bispectrum.generate_coefficients(elements)[source]

Automatically generates coefficients if not given by the user.

Parameters:elements (list of str) – List of symbols of all atoms.
Returns:G
Return type:dict of dicts
amp.descriptor.bispectrum.m_values(j)[source]

Returns a list of m values for a given j.

Cutoff functions

This script contains different cutoff function forms.

Note all cutoff functions need to have a “todict” method to support saving/loading as an Amp object.

All cutoff functions also need to have an Rc attribute which is the maximum distance at which properties are calculated; this will be used in calculating neighborlists.

class amp.descriptor.cutoffs.Cosine(Rc)[source]

Bases: object

Cosine functional form suggested by Behler.

Parameters:Rc (float) – Radius above which neighbor interactions are ignored.
__call__(Rij)[source]
Parameters:Rij (float) – Distance between pair atoms.
Returns:The vaule of the cutoff function.
Return type:float
prime(Rij)[source]

Derivative of the Cosine cutoff function.

Parameters:Rij (float) – Distance between pair atoms.
Returns:The vaule of derivative of the cutoff function.
Return type:float
todict()[source]
class amp.descriptor.cutoffs.Polynomial(Rc, gamma=4)[source]

Bases: object

Polynomial functional form suggested by Khorshidi and Peterson.

Parameters:
  • gamma (float) – The power of polynomial.
  • Rc (float) – Radius above which neighbor interactions are ignored.
__call__(Rij, gamma)[source]
Parameters:Rij (float) – Distance between pair atoms.
Returns:value – The vaule of the cutoff function.
Return type:float
prime(Rij, gamma)[source]

Derivative of the Cosine cutoff function.

Parameters:
  • Rc (float) – Radius above which neighbor interactions are ignored.
  • Rij (float) – Distance between pair atoms.
Returns:

The vaule of derivative of the cutoff function.

Return type:

float

todict()[source]
amp.descriptor.cutoffs.dict2cutoff(dct)[source]

This function converts a dictionary (which was created with the to_dict method of one of the cutoff classes) into an instantiated version of the class. Modeled after ASE’s dict2constraint function.