Bond Dissociation Energy

Calculate the bond dissociation energy (BDE) of ligands attached to the surface of the core. The calculation consists of five distinct steps:

1. Dissociate all combinations of \({n}\) ligands (\(Y\), see optional.qd.dissociate.lig_count) a nd an atom from the core (\(X\), see optional.qd.dissociate.core_atom) within a radius \(r\) from aforementioned core atom (see optional.qd.dissociate.lig_core_dist and optional.qd.dissociate.core_core_dist). The dissociated compound has the general structure of \(XY_{n}\).

2. Optimize the geometry of \(XY_{n}\) at the first level of theory (\(1\)). Default: ADF MOPAC [1, 2, 3].

3. Calculate the “electronic” contribution to the BDE (\(\Delta E\)) at the first level of theory (\(1\)): ADF MOPAC [1, 2, 3]. This step consists of single point calculations of the complete quantum dot, \(XY_{n}\) and all \(XY_{n}\)-dissociated quantum dots.

4. Calculate the thermalchemical contribution to the BDE (\(\Delta \Delta G\)) at the second level of theory (\(2\)). Default: ADF UFF [4, 5]. This step consists of geometry optimizations and frequency analyses of the same compounds used for step 3.

  1. \(\Delta G_{tot} = \Delta E_{1} + \Delta \Delta G_{2} = \Delta E_{1} + (\Delta G_{2} - \Delta E_{2})\).

Default Settings

optional:
    qd:
        dissociate:
            core_atom: Cd
            core_index: null
            lig_count: 2
            core_core_dist: 5.0  # Ångström
            lig_core_dist: 5.0  # Ångström
            lig_core_pairs: 1
            topology: {}

            keep_files: True
            job1: AMSJob
            s1: True
            job2: AMSJob
            s2: True

Arguments

optional.qd.dissociate
optional:
    qd:
        dissociate:
            core_atom: Cd
            core_index: null
            lig_count: 2
            lig_pairs: 1
            core_core_dist: null  # Ångström
            lig_core_dist: 5.0  # Ångström
            topology:
                7: vertice
                8: edge
                10: face

optional.qd.dissociate.core_atom
Parameter

The atomic number or atomic symbol of the core atoms (\(X\)) which are to be dissociated. The core atoms are dissociated in combination with \(n\) ligands (\(Y\), see dissociate.lig_count). Yields a compound with the general formula \(XY_{n}\).

Atomic indices can also be manually specified with dissociate.core_index

If one is interested in dissociating ligands in combination with a molecular species (e.g. \(X = {NR_4}^+\)) the atomic number (or symbol) can be substituted for a SMILES string represting a poly-atomic ion (e.g. tetramethyl ammonium: C[N+](C)(C)C).

If a SMILES string is provided it must satisfy the following 2 requirements:

  1. The SMILES string must contain a single charged atom; unpredictable behaviour can occur otherwise.

  2. The provided structure (including its bonds) must be present in the core.

Warning

This argument has no value be default and thus must be provided by the user.

optional.qd.dissociate.lig_count
Parameter

The number of ligands, \(n\), which is to be dissociated in combination with a single core atom (\(X\), see dissociate.core_atom).

Yields a compound with the general formula \(XY_{n}\).

Warning

This argument has no value be default and thus must be provided by the user.

optional.qd.dissociate.core_index
Parameter
  • Type - int or list [int], optional

  • Default valueNone

Alternative to dissociate.lig_core_dist and dissociate.core_atom. Manually specify the indices of all to-be dissociated atoms in the core. Core atoms will be dissociated in combination with the \(n\) closest ligands.

Note

Atom numbering follows the PLAMS [1, 2] convention of starting from 1 rather than 0.

Note

The yaml format uses null rather than None as in Python.

optional.qd.dissociate.core_core_dist
Parameter
  • Type - float or int, optional

  • Default valueNone

The maximum to be considered distance (Ångström) between atoms in dissociate.core_atom. Used for determining the topology of the core atom

(see dissociate.topology) and whether it is exposed to the surface of the core or not. It is recommended to use a radius which encapsulates a single (complete) shell of neighbours.

If not specified (or equal to 0.0) CAT will attempt to guess a suitable value based on the cores’ radial distribution function.

optional.qd.dissociate.lig_core_dist
Parameter
  • Type - float or int, optional

  • Default valueNone

Dissociate all combinations of a single core atom (see dissociate.core_atom) and the \(n\) closests ligands within a user-specified radius.

Serves as an alternative to dissociate.lig_core_dist, which removes a set number of combinations rather than everything withing a certain radius.

The number of ligands dissociated in combination with a single core atom is controlled by dissociate.lig_count.

_images/BDE_XY2.png

optional.qd.dissociate.lig_pairs
Parameter
  • Type - int, optional

  • Default valueNone

Dissociate a user-specified number of combinations of a single core atom (see dissociate.core_atom) and the \(n\) closests ligands.

Serves as an alternative to dissociate.lig_core_dist, removing a preset number of (closest) pairs rather than all combinations within a certain radius.

The number of ligands dissociated in combination with a single core atom is controlled by dissociate.lig_count.

optional.qd.dissociate.topology
Parameter
  • Type - dict

  • Default value{}

A dictionary which translates the number neighbouring core atoms (see dissociate.core_atom and dissociate.core_core_dist) into a topology. Keys represent the number of neighbours, values represent the matching topology.

Example

Given a dissociate.core_core_dist of 5.0 Ångström, the following options can be interpreted as following:

optional:
    qd:
        dissociate:
            7: vertice
            8: edge
            10: face

Core atoms with 7 other neighbouring core atoms (within a radius of 5.0 Ångström) are marked as "vertice", the ones with 8 neighbours are marked as "edge" and the ones with 10 neighbours as "face".

optional.qd.dissociate.qd_opt
Parameter
  • Type - bool

  • Default valueFalse

Whether to optimize the quantum dot and \(XY_{n}\) -dissociated quantum dot.


Arguments - Job Customization

optional.qd.dissociate
optional:
    qd:
        dissociate:
            keep_files: True
            job1: AMSJob
            s1: True
            job2: AMSJob
            s2: True

optional.qd.dissociate.keep_files
Parameter
  • Type - bool

  • Default valueTrue

Whether to keep or delete all BDE files after all calculations are finished.

optional.qd.dissociate.xyn_pre_opt
Parameter
  • Type - bool

  • Default valueTrue

Pre-optimize the \(XY_{n}\) fragment with UFF.

Note

Requires AMS.

optional.qd.dissociate.job1
Parameter

A type object of a Job subclass, used for calculating the “electronic” component (\(\Delta E_{1}\)) of the bond dissociation energy. Involves single point calculations.

Alternatively, an alias can be provided for a specific job type (see Type Aliases).

Setting it to True will default to AMSJob, while False is equivalent to optional.qd.dissociate = False.

optional.qd.dissociate.s1
Parameter
s1:
    input:
        mopac:
            model: PM7
        ams:
            system:
                charge: 0

The job settings used for calculating the “electronic” component (\(\Delta E_{1}\)) of the bond dissociation energy.

Alternatively, a path can be provided to .json or .yaml file containing the job settings.

Setting it to True will default to the ["MOPAC"] block in CAT/data/templates/qd.yaml, while False is equivalent to optional.qd.dissociate = False.

optional.qd.dissociate.job2
Parameter

A type object of a Job subclass, used for calculating the thermal component (\(\Delta \Delta G_{2}\)) of the bond dissociation energy. Involves a geometry reoptimizations and frequency analyses.

Alternatively, an alias can be provided for a specific job type (see Type Aliases).

Setting it to True will default to AMSJob, while False will skip the thermochemical analysis completely.

optional.qd.dissociate.s2
Parameter
s2:
    input:
        uff:
            library: uff
        ams:
            system:
                charge: 0
                bondorders:
                    _1: null

The job settings used for calculating the thermal component (\(\Delta \Delta G_{2}\)) of the bond dissociation energy.

Alternatively, a path can be provided to .json or .yaml file containing the job settings.

Setting it to True will default to the the MOPAC block in CAT/data/templates/qd.yaml, while False will skip the thermochemical analysis completely.

Index

dissociate_ligand(mol, lig_count[, ...])

Remove \(XY_{n}\) from mol with the help of the MolDissociater class.

MolDissociater(mol, core_idx, ligand_count)

The MolDissociater class; serves as an API for dissociate_ligand().

MolDissociater.remove_bulk([max_vec_len])

Remove out atoms specified in MolDissociater.core_idx which are present in the bulk.

MolDissociater.assign_topology()

Assign a topology to all core atoms in MolDissociater.core_idx.

MolDissociater.get_pairs_closest(lig_idx[, ...])

Create and return the indices of each core atom and the \(n\) closest ligands.

MolDissociater.get_pairs_distance(lig_idx[, ...])

Create and return the indices of each core atom and all ligand pairs with max_dist.

MolDissociater.combinations(cor_lig_pairs[, ...])

Create a list with all to-be removed atom combinations.

MolDissociater.__call__(combinations)

Start the dissociation process.

API

nanoCAT.bde.dissociate_xyn.dissociate_ligand(mol, lig_count, lig_core_pairs=1, lig_core_dist=None, core_atom=None, core_index=None, core_smiles=None, core_core_dist=None, topology=None, **kwargs)[source]

Remove \(XY_{n}\) from mol with the help of the MolDissociater class.

The dissociation process consists of 5 general steps:

Examples

>>> from typing import Iterator

>>> import numpy as np
>>> from scm.plams import Molecule

# Define parameters
>>> mol = Molecule(...)
>>> core_idx = [1, 2, 3, 4, 5]
>>> lig_idx = [10, 20, 30, 40]
>>> lig_count = 2

# Start the workflow
>>> dissociate = MolDissociater(mol, core_idx, lig_count)
>>> dissociate.assign_topology()
>>> pairs: np.ndarray = dissociate.get_pairs_closest(lig_idx)
>>> combinations: Iterator[tuple] = dissociate.get_combinations(pairs)

# Create the final iterator
>>> mol_iterator: Iterator[Molecule] = dissociate(cor_lig_combinations)
Parameters
  • mol (plams.Molecule) – A molecule.

  • lig_count (int) – The number of to-be dissociated ligands per core atom/molecule.

  • lig_core_pairs (int, optional) – The number of to-be dissociated core/ligand pairs per core atom. Core/ligand pairs are picked based on whichever ligands are closest to each core atom. This option is irrelevant if a distance based criterium is used (see lig_dist).

  • lig_core_dist (float, optional) – Instead of dissociating a given number of core/ligand pairs (see lig_pairs) dissociate all pairs within a given distance from a core atom.

  • core_index (int or Iterable [int]) – An index or set of indices with all to-be dissociated core atoms. See core_atom to define core_idx based on a common atomic symbol/number.

  • core_atom (int or str, optional) – An atomic number or symbol used for automatically defining core_idx. Core atoms within the bulk (rather than on the surface) are ignored.

  • core_smiles (str, optional) – A SMILES string representing molecule containing core_idx. Provide a value here if one wants to disociate an entire molecules from the core and not just atoms.

  • core_core_dist (float, optional) – A value representing the mean distance between the core atoms in core_idx. If None, guess this value based on the radial distribution function of mol (this is generally recomended).

  • topology (Mapping [int, str], optional) – A mapping neighbouring of atom counts to a user specified topology descriptor (e.g. "edge", "vertice" or "face").

  • **kwargs (Any) – For catching excess keyword arguments.

Returns

A generator yielding new molecules with \(XY_{n}\) removed.

Return type

Generator [plams.Molecule]

Raises

TypeError – Raised if core_atom and core_idx are both None or lig_core_pairs and lig_core_dist are both None.

class nanoCAT.bde.dissociate_xyn.MolDissociater(mol, core_idx, ligand_count, max_dist=None, topology=None)[source]

The MolDissociater class; serves as an API for dissociate_ligand().

Parameters
mol

A PLAMS molecule consisting of cores and ligands.

Type

plams.Molecule

core_idx

An iterable with (1-based) atomic indices of all core atoms valid for dissociation.

Type

int or Iterable [int]

ligand_count

The number of ligands to-be dissociation with a single atom from MolDissociater.core_idx.

Type

int

max_dist

The maximum distance between core atoms for them to-be considered neighbours. If None, this value will be guessed based on the radial distribution function of MolDissociater.mol.

Type

float, optional

topology

A mapping of neighbouring atom counts to a user-specified topology descriptor.

Type

dict [int, str], optional

MolDissociater.remove_bulk(max_vec_len=0.5)[source]

Remove out atoms specified in MolDissociater.core_idx which are present in the bulk.

The function searches for all neighbouring core atoms within a radius MolDissociater.max_dist. Vectors are then constructed from the core atom to the mean positioon of its neighbours. Vector lengths close to 0 thus indicate that the core atom is surrounded in a (nearly) spherical pattern, i.e. it’s located in the bulk of the material and not on the surface.

Performs in inplace update of MolDissociater.core_idx.

Parameters

max_vec_len (float) – The maximum length of an atom vector to-be considered part of the bulk. Atoms producing smaller values are removed from MolDissociater.core_idx. Units are in Angstroem.

MolDissociater.assign_topology()[source]

Assign a topology to all core atoms in MolDissociater.core_idx.

The topology descriptor is based on:

If no topology description is available for a particular neighbouring atom count, then a generic f"{i}_neighbours" descriptor is used (where i is the neighbouring atom count).

Performs an inplace update of all Atom.properties.topology values.

MolDissociater.get_pairs_closest(lig_idx, n_pairs=1)[source]

Create and return the indices of each core atom and the \(n\) closest ligands.

Parameters
  • lig_idx (int or Iterable [int]) – The (1-based) indices of all ligand anchor atoms.

  • n_pairs (int) – The number of to-be returned pairs per core atom. If n_pairs > 1 than each successive set of to-be dissociated ligands is determined by the norm of the \(n\) distances.

Returns

A 2D array with the indices of all valid ligand/core pairs.

Return type

2D numpy.ndarray [int]

MolDissociater.get_pairs_distance(lig_idx, max_dist=5.0)[source]

Create and return the indices of each core atom and all ligand pairs with max_dist.

Parameters
  • lig_idx (int or Iterable [int]) – The (1-based) indices of all ligand anchor atoms.

  • max_dist (float) – The radius (Angstroem) used as cutoff.

Returns

A 2D array with the indices of all valid ligand/core pairs.

Return type

2D numpy.ndarray [int]

MolDissociater.combinations(cor_lig_pairs, lig_mapping=None, core_mapping=None)[source]

Create a list with all to-be removed atom combinations.

Parameters
  • cor_lig_pairs (numpy.ndarray) – An array with the indices of all core/ligand pairs.

  • lig_mapping (Mapping, optional) – A mapping for translating (1-based) atomic indices in cor_lig_pairs[:, 0] to lists of (1-based) atomic indices. Used for mapping ligand anchor atoms to the rest of the to-be dissociated ligands.

  • core_mapping (Mapping, optional) – A mapping for translating (1-based) atomic indices in cor_lig_pairs[:, 1:] to lists of (1-based) atomic indices. Used for mapping core atoms to the to-be dissociated sub structures.

Returns

A set of 2-tuples. The first element of each tuple is a frozenset with the (1-based) indices of all to-be removed core atoms. The second element contains a frozenset with the (1-based) indices of all to-be removed ligand atoms.

Return type

set [tuple]

MolDissociater.__call__(combinations)[source]

Start the dissociation process.