Modules of methods of calculation

Tight-Binding method

Semiempirical Tight Binding method was realized.

class kvazar.core.tb.TB(struct, parm)[source]

Class that realizes Tight-Binding method.

Constructor of class TB.

Args:

struct (Structure): nanostructure.

lib (module): dynamic library, it realize filling Hamiltonian matrix, calculating repulsive energie, forces acting on atoms.

parms (dict): parameters of method.

Kwargs:

functions (list): list of functions that will be loaded from lib library.
compute()[source]

Calculating full energy of nanostructure.

Returns:

energy_total(float): full energy of nanostructure
compute_ebond()[source]

Calculating bond energy using formula

$E_{val} = \sum\limits_{\gamma}{n_\gamma \epsilon_\gamma}$

where \epsilon_{\gamma} is a value on diagonal of diagonalised hamiltonian

Returns:

energy_bond(float): bond energy of atoms
compute_erep()[source]

Repulsive energy between two atoms

$E_{core}(r) = E_{core}(r_0) \left(\frac{r_0}{r}\right)^{m_a}exp\left(-m_b\left(\frac{r}{r_c}\right)^{m_c} + m_b\left( \frac{r_0}{r_c}\right)^{m_c}\right)$

Total repulsive energy

$E_{rep} = \sum\limits_{i < j}{E_{core}(r_{ij})}$

Returns:
energy_rep(float): repulsive energy
ensure_forces()[source]

Function changes method for calculating forces, acting on atoms of structure, to ensure, that hamiltonian matrix is not zero before calculation start and eigenvectors can be calculated.

load_force_function(lib)[source]

Loading of function that calculate forces and creating wrap-methods from dynamic library.

Args:

functions (list): list of functions, that will be loaded from lib

lib (module): dynamic library

load_tb_functions(functions, lib)[source]

Loading of functions from dynamic library and creating wrap-methods.

Args:

functions (list): list of functions, that will be loaded from lib

lib (module): dynamic library

preprocessing(struct)[source]

Preprocessing for TB method, that:

  1. Checks existing of parameters for different atom types, involved in structure.
  2. Creates Hamiltonian matrix.
  3. Calculates number of orbitals of each type and shift indices for orbitals in hamiltonian.
  4. Calculates number of occupied energy levels
Args:
struct (Structure): nanostructure
  1. Updates data of tight bonds in structure
  2. Zeros hamiltonian matrix
Args:

tb_links(numpy.array): 2D arrays of bonds.

num_tb_links(numpy.array): 1D array, it contains numbers of bonds for each atom

Molecular mechanics

Molecular mechanics method was realized.

class kvazar.core.mm.MM(struct, lib, ff_functions, ff_params)[source]

Wrap class for computational modules of molecular mechanical method.

Constructor of class:

Args:

struct(structure): nanostructure

lib(module): dynamic library

ff_functions(list): list of names of force field functions

ff_params(list): list of names of parameters of force field functions

compute()[source]

Returns energy of structure

compute_i(i)[source]

Returns summ of energies, calculated by functions from list ‘functions_i’

compute_ij(i, j)[source]

Returns summ of energies, calculated by functions from list ‘functions_ij’

loadForceFunction(lib)[source]

Adds function ‘forces’ to module

Args:

lib(module): dynamic library
loadMMFunctions(ff_functions, lib)[source]

Args:

ff_functions(list): list of names of atom types

lib(module): dynamic library for molecular mechanics computation

preprocessing(struct)[source]

It checks matching of atom types and atom names of structure with the same parameters of force field. If there are no atom types in structure function rebuilds it from atom names. If there are no atom names, it prints, that it’s not possible to check out correctness of atom types.

Reactive empirical bond order potential

Method REBO was realised according to article: Donald W. Brenner, Olga A. Shenderova, Judith A. Harrison, Steven J. Stuart, Boris Ni, and Susan B. Sinnott, “A second-generation reactive empirical bond order (REBO) potential energy expression for hydrocarbons” // J. Phys.: Condens. Matter 14(2002), P. 783–802. Description for method was taken from this article.

Chemical binding energy E_b can be written as a sum over nearest neighbours in the form:

(1)E_b = \sum_{i}\sum_{j(>i)}[V^R(r_{ij}) - b_{ij}V^A(r_{ij})]

The equation (1) is used for the total potential energy. Following the earlier hydrocarbon bonding expression, the empirical bond order function used here is written as a sum of terms:

(2)\overline{b}_{ij} = \frac{1}{2}[b_{ij}^{\sigma-\pi} + b_{ji}^{\sigma-\pi}] + b_{ij}^\pi

Values for the function b_{ij}^{\sigma-\pi} and b_{ji}^{\sigma-\pi} depend on the local coordination and bond angles for atoms i and j respectively. The function b_{ij}^\pi is further written as a sum of two terms:

(3)b_{ij}^\pi = \Pi_{ij}^{RC} + b_{ij}^{DH}

The first generation hydrocarbon expression used Morse-type terms for the pair interactions in equation (1). It was determined that this form is too restrictive to simultaneously fit equilibrium distances, energies, and force constants for carbon-carbon bonds. This form has the further disadvantage that both terms go to finite values as the distance between atoms decreases, limiting the possibility of modeling process involving energetic atomic collisions. In the second-generation potential next forms are used for attractive ad repulsive potentials:

(4)V^{R}(r) = f^{c}(r)(1+Q/r)Ae^{-\alpha r}

and

(5)V^{A}(r) = f^{c}(r)\sum_{n = 1,3}B_{n}e^{-\beta_{n} r}

Analytic forms for each term in equation (2) are fitted to both the discrete values of the bond orders above, and to other properties of solid-state and molecular carbon as described below. The first term in equation (2):

(6)b_{ij}^{\sigma-\pi} = \left[ 1 + \sum_{k(\ne i,j)}f^{c}_{ik}(r_{ik})G(cos(\theta_{ijk}))e^{\lambda_{ijk}} + P_{ij}(N^{C}_i , N^{H}_i)\right]^{-\frac{1}{2}}

where f^{c}(r) is a function to ensure that interations include only nearest neighbours. P represents a bicubic spline and quantities N^{C}_i and N^{H}_i represent the number of carbon and hydrogen atoms, respectively, that are neighbours of atom i. It’s defined by next equations:

(7)N^{C}_i = \sum^{carbon\;atoms}_{k (\ne i,j)}f^{c}_{ik}(r_{ik})

and

(8)N^{H}_i = \sum^{hydrogen\;atoms}_{l (\ne i,j)}f^{c}_{il}(r_{il})

The term \pi^{RC}_{ij} is taken as tricubic spline F:

(9)\pi^{RC}_{ij} = F_{ij}\left(N_{i}^{t}, N_{j}^{t}, N_{ij}^{conj}\right)

N_{ij}^{conj} is given by next expression:

(10)N_{ij}^{conj} = 1 + \left[\sum_{k(\ne i,j)}^{carbon}f_{ik}^c (r_{ik})F(X_{ik}) \right] ^ 2 + \left[\sum_{l(\ne i,j)}^{carbon}f_{il}^c (r_{il})F(X_{il}) \right] ^ 2

where

(11)$F(\chi_{ik}) = \begin{cases}
                                   1 & \chi_{ik} < 2 \\
        [1+cos(2\pi(\chi_{ik}-2))]/2 & 2 < \chi_{ik} < 3 \\
                                   0 & \chi_{ik} > 3
\end{cases}$

and

(12)\chi_{ik} = N_{k}^t - f_{ik}^c (r_{ik})

The term b_{ij}^{DH} in equation (3) :

(13)b_{ij}^{DH} = T_{ij}(N_{i}^t ,N_{j}^t ,N_{ij}^{conj})\left[\sum_{k(\ne i,j)}\sum_{l(\ne i,j)}(1 - cos^{2}(\theta_{ijkl}))f_{ik}^{c}(r_{ik})f_{jl}^{c}(r_{il})\right]

where

(14)cos(\theta_{ijkl}) = e_{ijk}e_{ijl}

The function T_{ij}(N_{i}^t ,N_{j}^t ,N_{ij}^{conj}) is a tricubic spline, e_{ijk} and e_{ijl} are unit vectors in the direction of the cross products R_{ji} \times R_{ik} and R_{ij} \times R_{jl}, respecrively, where R are vectors connecting subscripted atoms.

The value of f^{c}(r) is defined by a switching function of the form:

(15)$f_{ij}^{c}(r) = \begin{cases}
1 & r < D_{ij}^{min} \\
\left[ 1 + cos(\pi((r - D_{ij}^{min})/(D_{ij}^{max} - D_{ij}^{min}))\right]/2 & D_{ij}^{min} < r < D_{ij}^{max} \\
0 & r > D_{ij}^{max}
\end{cases}$

Adaptive intermolecular reactive empirical bond order potential

AIREBO method was realized according to article: Steven J. Stuart, Alan B. Tutein and Judith A. Harrison. “A reactive potential for hydrocarbons with intermolecular interactions” // Journ. of Chem. Phys., V 112,14 (2000), p. 6472-6486.

AIREBO potential can be represented by a sum over pairwise interactions, including covalent bonding REBO interactions, LJ terms, and torsion interactions:

(16)E = \frac{1}{2}\sum_{i}\sum_{j\neq{i}}\left(E_{ij}^{REBO} + E_{ij}^{LJ} +\sum_{k\neq{i,j}}\sum_{l\neq{i,j,k}}E_{kijl}^{tors} \right)

Van-der-Vaals interaction is described through the using of the Lennard-Jones potential:

(17)E_{ij}^{LJ} = S(t_r (r_{ij}))S(t_b (b_{ij}^* ))C_{ij}V_{ij}^{LJ}(r_{ij}) + [1 - S(t_r(r_{ij}))]C_{ij}V_{ij}^{LJ}(r_{ij})

V_{ij}^{LJ} is the traditional LJ term:

(18)V_{ij}^{LJ} = 4\epsilon_{ij}\left( \left(\frac{\sigma_{ij}}{r_{ij}}\right)^{12} - \left(\frac{\sigma_{ij}}{r_{ij}}\right)^6 \right)

It was modified by several sets of switching functions. S(t) can be represented in the next form:

(19)S(t) = \Theta(-t)+\Theta(t)\Theta(1-t)[1 - t^2(3-2t)]

Magnitude of LJ interactions depends on bonding environment. Gradual exclusion of Lennard-Jones interactions with changings of r_{ij} is controlled by scalling function t_{b}:

(20)t_b(b_{ij}) = \frac{b_{ij} - b_{ij}^{min}}{b_{ij}^{max} - b_{ij}^{min}}

If atoms i and j are not connected by two or fewer intermediate atoms, LJ interactions between them are controlled by next switching function:

(21)C_{ij} = 1 - max\{w_{ij}(r_{ij}),w_{ik}(r_{ik})w_{kj}(r_{kj}), \forall k, \\
w_{ik}(r_{ik})w_{kl}(r_{kl})w_{lj}(r_{lj}), \forall k,l \}

where:

w_{ij}(r_{ij}) = S^\prime(t_c(r_{ij}))

S^\prime(t_c(r_{ij})) = \Theta(-t)+\Theta(t)\Theta(1-t)\frac{1}{2}[1+cos(\pi t)]

The torsional part of equation (16) for the dihedral angle determined by atoms i, j, k, l has the next form:

(22)E^{tors}_{kijl} = w_{ki}(r_{ki})w_{ij}(r_{ij})w_{jl}(r_{jl})V^{tors}(\omega_{kijl})

where V^{tors}(\omega_{kijl}) is the torsional potential:

(23)V^{tors}(\omega_{kijl}) = \frac{256}{405}\epsilon_{kijl}cos^{10}(\omega_{kijl}/2)-\frac{1}{10}\epsilon_{kijl}

class kvazar.core.airebo.AIREBO(struct, ff_params, comm=None)[source]

This class is the successor of class MM

Constructor of class AIREBO:

Args:

struct(structure): nanostructure

ff_params(list): list of names of potentials functions

init_cutoff(r_in, r_out, cutoff_neighbours)[source]

Args:

r_in(float): internal cutoff radius

r_out(float): external cutoff radius

cutoff_neighbours(list): list of neighbours get into the field with r_out radius

preprocessing(struct)[source]

It corrects atom type and atom boxnames, group them and find out. Masses are taken from parameters.

update_cutoff(cutoff_neighbours)[source]

Updates list of neighbours in cutoff area

Args:

cutoff_neighbours(list): list of neighbours

Coarse-grained modeling

MARTINI force field for coarse-grained modeling was realized.

class kvazar.core.martini.MARTINI(struct, ff_params)[source]

This class is the successor of class MM

Constructor of class MARTINI:

Args:

struct(structure): nanostructure

ff_params(list): list of names of force field parameters(functions)

init_cutoff(r_in, r_out, cutoff_neighbours)[source]

Args:

r_in(float): internal cutoff radius

r_out(float): external cutoff radius

cutoff_neighbours(list): list of neighbours get into the field with r_out radius

preprocessing(struct)[source]

It corrects atom type and atom boxnames, group them and find out, if bead is the part of protein, take mass of bead from parameters.

update_cutoff(cutoff_neighbours)[source]

Updates list of neighbours in cutoff area

Args:

cutoff_neighbours(list): list of neighbours