Theoretical Basis of Quantum-Mechanical Modeling of Functional Nanostructures

Abstract

The paper presents an analytical review of theoretical methods for modeling functional nanostructures. The main evolutionary changes in the approaches of quantum-mechanical modeling are described. The foundations of the first-principal theory are considered, including the stationery and time-dependent Schrödinger equations, wave functions, the form of writing energy operators, and the principles of solving equations. The idea and specifics of describing the motion and interaction of nuclei and electrons in the framework of the theory of the electron density functional are presented. Common approximations and approaches in the methods of quantum mechanics are presented, including the Born–Oppenheimer approximation, the Hartree–Fock approximation, the Thomas–Fermi theory, the Hohenberg–Kohn theorems, and the Kohn–Sham formalism. Various options for describing the exchange–correlation energy in the theory of the electron density functional are considered, such as the local density approximation, generalized and meta-generalized gradient approximations, and hybridization of the generalized gradient method. The development of methods of quantum mechanics to quantum molecular dynamics or the dynamics of Car–Parrinello is shown. The basic idea of combining classical molecular modeling with calculations of the electronic structure, which is reflected in the potentials of the embedded atom, is described.

1. Introduction

Symmetry effects and the corresponding structure of objects are ubiquitous in the nanoworld. These effects can appear at the structure level of individual atoms or molecules, where a symmetric distribution of electron clouds, their stable energy state in a similar spatial arrangement, is observed. Symmetric phenomena are also required when considering large groups of atoms, which may form individual clusters or may arise within the structure of a material, e.g., within nanocomposites. A proper analysis of symmetric effects allows us to systematize the interacting mechanisms between atoms, molecules and nanoparticles, and lay the foundation for further research. An example is the introduction of the crystal lattice types classification [1,2], which in its time allowed us to considerably simplify the description of the solid materials structure and the analysis of their structure defects.

Developing modern science makes it possible to create new promising materials, design high-tech devices, and improve existing industrial and production facilities. As a result, science is currently in close integration with many branches of technology and human activity, and its applied significance is growing. Scientific research has now reached such a level that the individual atoms or molecules manipulating mechanisms [3,4,5] fall under observation, and sometimes the impact on materials reaches the level of single electrons [6,7], bosons [8,9], fermions [10,11], photons [12,13] or other elementary particles.

A particular promising area of research is nanocomposites, i.e., a multicomponent material with nanoscale dimensional, physical or structural characteristics. Because of improved and sometimes fundamentally new properties, nanocomposites are widely used in various industry fields. Finely dispersed composite materials are actively used as catalysts, adsorption and electrochemical agents [14,15,16]. The constituent elements of nanocomposites–nanoparticles also have similar reactive characteristics [17,18]. Nanoclusters or nanostructures are independent functional materials, but they require stabilization and isolated storage due to their increased agglomeration [19]. Medium-grained nanocomposites have effective optical, electrical and thermal properties [20,21,22]. Larger (about 100 nm) multicomponent samples are used for mechanical reinforcement [23], fracture resistance and simultaneous lightening of structures [24], and for obtaining new materials with improved magnetic characteristics [25,26]. The division into spheres of using nanocomposites, depending on the dispersion and internal structure should be understood as rather arbitrary. Often, there is a complex improvement in the properties of designed samples, in contrast to the components integrated in their composition.

Another popular area of implementing nano-objects is developing nanoelectromechanical systems [27,28]. In such systems, nanomaterials change their properties or create an electromagnetic response under the influence of an external force, electric field or luminous flux. With the computing devices evolution and the trend towards their miniaturization, nanoelectromechanical systems function as nanoscale recording devices: resonators [29], nanomotors [30], thermal sensors [31] and various other sensors. An example of the described devices is a nanomechanical resonator in the form of a multilayer gallium arsenide nanostructure proposed in [29]. The nano resonator is capable of converting the natural vibrations energy into radiation, thus observing the change in external forces. The developers [29] plan to create highly accurate sensors for medical and chemical analysis, based on such systems.

A decrease in the size of technical devices brings about a significant increase in their total surface area, and in the fraction of the friction force of the nanosystem relative general mechanical behavior. The increase in contributing the friction force causes appearing excessive heat release, which is not always useful, and often has no need. A similar problem also arises in nanoelectromechanical systems for electronic computing purposes [32]. The density of integrated circuit structural elements is currently reaching such a level that heat release limits their further miniaturization. This difficulty can be solved by scientific research in such an area of quantum electronics as spintronics. Spintronics devices and materials are based on the energy transfer using a spin current, rather than an electron current, which provides fundamentally new prospects for implementing these objects into microprocessors [33], memory cards [34] and sensors [35]. For example, Josephson contacts with a tunable critical current are extremely important for the creation of Josephson neurons and synapses. Such tunable nonlinear elements with memory can be created on the basis of hybrid superconductor–ferromagnetic (S/F) structures [36,37,38,39,40,41,42,43,44]. On the basis of F1-S-F2 structures, configurable linear elements with memory (controllable kinetic inductors) can also be manufactured [45,46]. Quantum-mechanical modeling of such artificial materials and functional nanostructures with specified properties remains one of the most urgent tasks of modern electronics and spintronics. To give an example: for efficient operation of multilayer structures in tunable Josephson and inductive elements it is extremely important to select materials and technological modes of their sputtering so as to ensure the required quality of boundaries and layers. Only a thorough first-principles simulation of the entire process of multilayer heterostructures formation gives chances to obtain devices with the specified characteristics on their basis. The required procedure is described in detail in Section 6 of this article.

Despite all the advances in modern nanotechnology, research at the molecular level remains an extremely challenging task. An undoubtedly important part is experiments and direct work with new materials. However, theoretical research methods play an important role in improving and optimizing the technology for manufacturing nanoscale objects. This work is an overview with elements of applying mathematical models to specific applied problems. The aim of the work is to demonstrate the relevance of using the family of quantum mechanics methods, their evolution and optimization to study various properties and structure features of nanomaterials, including symmetry effects, for example, in Josephson contacts in spintronics devices. The use of mathematical modeling methods for the study of specific technological processes and nanostructured materials is described in the authors’ previously published works [47,48,49,50,51,52,53,54].

2. First-Principal Modeling for Studying the Properties of Nanomaterials

Understanding, describing and modeling the structure and properties of nanomaterials is based on the quantum mechanics methods. These methods appeared some time ago, but they are still successfully used for research and analysis of nanotechnology objects [55,56]. The quantum mechanics mechanism describes the nature and electronic structure of matter in a physically reasonable way and is constantly being supplemented and improved, so that it has not lost its relevance. A review of the quantum mechanics theory, approximations and evolutionary approaches is given in [57,58,59].

The Schrödinger equation is the basis for describing the atom’s dynamics using quantum mechanics. In the general nonrelativistic case, the Schrödinger equation describes the motion of nuclei and electrons by the following equation

$$i\hslash \frac{\partial \Psi \left(R_1,\dots ,R_{N_n},r_1,\dots ,r_{N_e},t\right)}{\partial t}=\widehat{H}\Psi \left(R_1,\dots ,R_{N_n},r_1,\dots ,r_{N_e},t\right),$$

(1)

where i is unit imaginary number; $\hslash$ is Planck’s constant; $\Psi \left(R_1,\dots ,R_{N_n},r_1,\dots ,r_{N_e},t\right)$ is common wave function of nuclei and electrons; $N_n$ is number of nuclei in the system; $N_e$ is number of electrons; t is time; $R_i$ and $r_j$ are the radius-vectors of the nuclei and the electron at a certain point in time; $\widehat{H}$ is Hamilton operator.

The Hamilton operator includes the compound interaction operators, as well as the nuclei and electrons kinetic energy

$$\widehat{H}=-{\displaystyle \sum _i\frac{{\hslash }^2}{2M_i}{\nabla }_i^2}-{\displaystyle \sum _j\frac{{\hslash }^2}{2m}{\nabla }_j^2}-{\displaystyle \sum _{i,j}\frac{Z_i{e}^2}{\left|R_i-r_j\right|}}+\frac{1}{2}{\displaystyle \sum _{i,j,i\ne j}\frac{Z_i{Z}_j{e}^2}{\left|R_i-R_j\right|}}+\frac{1}{2}{\displaystyle \sum _{i,j,i\ne j}\frac{e^2}{\left|r_i-r_j\right|}},$$

(2)

where $M_i$ and $Z_i$ are masses and charges of nuclei; $m$ and $e$ is mass and charge of the electron, respectively. Hereinafter, the summation relative to nuclei is made up to the variable $N_n$, and the summation in the electron operators to the variable $N_e$.

To describe the atomic and electronic structure of a nanomaterial, it is necessary to solve an Equation (1). The Schrödinger equation solution involves a large amount of variables, limiting its exact solution. The separation of the variables into slow nuclei and fast electrons, called the adiabatic approximation or Born–Oppenheimer approximation, allows us to significantly simplify the equation [57]

$$\{\begin{array}{l}{M}_i{\ddot{R}}_i=-\frac{\partial E}{\partial R_i},\\ {\widehat{H}}_{\mathrm{BO}}\left(R_1,\dots ,R_{N_n},r_1,\dots ,r_{N_e},t\right)\hspace{0.17em}\Psi \left(r_1,\dots ,r_{N_e},t\right)=E\left(R_1,\dots ,R_{N_n},t\right)\hspace{0.17em}\Psi \left(r_1,\dots ,r_{N_e},t\right),\end{array}$$

(3)

where ${\widehat{H}}_{\mathrm{BO}}$ is Hamilton operator in the Born–Oppenheimer approximation.

The adiabatic approximation allows us to solve the Schrödinger equation gradually. At each time step, the motion of electrons in the stationary field of the nuclei is considered first. Then, based on the total electron distribution in the system, new nuclei positions are calculated. The primary interest in constructing quantum-mechanical fields is aimed at studying the contribution of the electronic component, so the stationary Schrödinger equation for the multielectron subsystem located in the field of stationary nuclei will be considered further, but taking into account the spin state of electrons

$${\widehat{H}}_{\mathrm{BO}}\Psi \left(r_1,s_1\dots ,r_{N_e},s_{N_e}\right)=E\Psi \left(r_1,s_1\dots ,r_{N_e},s_{N_e}\right),$$

(4)

$${\widehat{H}}_{\mathrm{BO}}=-{\displaystyle \sum _j\frac{{\hslash }^2}{2m}{\nabla }_j^2}-{\displaystyle \sum _{i,j}\frac{Z_i{e}^2}{\left|R_i-r_j\right|}}+\frac{1}{2}{\displaystyle \sum _{i,j,i\ne j}\frac{e^2}{\left|r_i-r_j\right|}},$$

(5)

where $E$ is energy of the considered system; $s_i$ is electron spin. The solution of Equation (4) makes it possible to calculate the ground state of the nanosystem corresponding to the lowest energy, as well as its excitation spectrum formed by solutions with higher energy.

The simplest assumption for solving an Equation (4) is the electrons’ noninteraction, which allows us to represent the Hamiltonian operator through the sum and total wave function through the product of the individual contributions

$$\widehat{H}={\displaystyle \sum _i\left(-\frac{{\hslash }^2}{2m}{\nabla }_i^2+v_{\mathrm{ion}}\left(r_i\right)\right)},\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{v}_{\mathrm{ion}}\left(r_i\right)=-{\displaystyle \sum _j\frac{Z_j{e}^2}{\left|R_j-r_i\right|}},\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\Psi \left(r_1,s_1\dots ,r_{N_e},s_{N_e}\right)={\displaystyle \prod _i\mathsf{\psi }\left(r_i,s_i\right)},$$

(6)

where $\hspace{0.17em}{v}_{\mathrm{ion}}\left(r_i\right)$ is the electron external field in the field of the nucleus; $\mathsf{\psi }\left(r_i,s_i\right)$ are one-electron wave functions.

Due to the electrons noninteraction hypothesis, the equations for particles can be separated

$$\left(-\frac{{\hslash }^2}{2m}{\nabla }_i^2+v_{\mathrm{ion}}\left(r_i\right)\right){\mathsf{\psi }}_{i\mathsf{\sigma }}\left(r\right)={\mathsf{\epsilon }}_{i\mathsf{\sigma }}{\mathsf{\psi }}_{i\mathsf{\sigma }}\left(r\right),$$

(7)

where ${\mathsf{\psi }}_{i\mathsf{\sigma }}\left(r\right)$ and ${\mathsf{\epsilon }}_{i\mathsf{\sigma }}$ is single-particle wave function and energy of the electron in the spin state $\mathsf{\sigma }$. The characteristics of the electron gas are found by elementary summation

$$E={\displaystyle \sum _{i,\mathsf{\sigma }}f\left({\mathsf{\epsilon }}_{i\mathsf{\sigma }}\right)\hspace{0.17em}{\mathsf{\epsilon }}_{i\mathsf{\sigma }},\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}n\left(r\right)}={\displaystyle \sum _{i,\mathsf{\sigma }}f\left({\mathsf{\epsilon }}_{i\mathsf{\sigma }}\right)\hspace{0.17em}{\mathsf{\psi }}_{i\mathsf{\sigma }}\left(r\right)\hspace{0.17em}{\mathsf{\psi }}_{i\mathsf{\sigma }}^{*}\left(r\right)},\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}f\left({\mathsf{\epsilon }}_{i\mathsf{\sigma }}\right)=\frac{1}{e^{\mathsf{\beta }\left({\mathsf{\epsilon }}_{i\mathsf{\sigma }}-\mathsf{\mu }\right)}+1},$$

(8)

where $E$ is average nanosystem energy; $n\left(r\right)$ is particles density; ${\mathsf{\psi }}_{i\mathsf{\sigma }}^{*}\left(r\right)$ is conjugate wave function of the electron; $f\left({\mathsf{\epsilon }}_{i\mathsf{\sigma }}\right)$ is Fermi function; $\mathsf{\mu }$ is chemical potential; $\mathsf{\beta }$ is reverse temperature.

The Hartree–Fock approximation [60,61] also simplifies the solution of the Schrödinger equation via single-particle wave functions, with the assumption that each electron moves in an averaged field. This field $U\left(r\right)$ is created by all electrons of the nanosystem and is gradually refined, so it is called self-consistent

$$\begin{gathered} \left(-\frac{{\hslash }^2}{2m}{\nabla }_i^2+U\left(r_i\right)\right)\hspace{0.17em}{\mathsf{\psi }}_{i\mathsf{\sigma }}\left(r\right)+{\displaystyle \sum _{j,{\mathsf{\sigma }}_i}{\displaystyle \int d{r}^{\prime }{\mathsf{\psi }}_{j{\mathsf{\sigma }}_j}^{*}\left(r^{\prime }\right)\frac{e^2}{\left|r-r^{\prime }\right|}}}{\mathsf{\psi }}_{j{\mathsf{\sigma }}_j}\left(r^{\prime }\right)\hspace{0.17em}{\mathsf{\psi }}_{i{\mathsf{\sigma }}_i}\left(r\right) \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}-{\displaystyle \sum _j{\displaystyle \int d{r}^{\prime }{\mathsf{\psi }}_{j\mathsf{\sigma }}^{*}\left(r^{\prime }\right)\hspace{0.17em}{\mathsf{\psi }}_{i\mathsf{\sigma }}\left(r^{\prime }\right)\frac{e^2}{\left|r-r^{\prime }\right|}\hspace{0.17em}}}{\mathsf{\psi }}_{j\mathsf{\sigma }}\left(r\right)={\mathsf{\epsilon }}_{i\mathsf{\sigma }}{\mathsf{\psi }}_{i\mathsf{\sigma }}\left(r\right). \end{gathered}$$

(9)

The first term in the left side of Equation (9) describes the one-electron contribution to the nanosystem energy, the second term is responsible for the Coulomb interaction and the third term corresponds to the exchange energy. The multielectron wave function is calculated in the form of the matrix determinant

$$\Psi \left(r_1,s_1\dots ,r_{N_e},s_{N_e}\right)=\frac{1}{\sqrt{N_e!}}\det \left[{\mathsf{\psi }}_j\left(r_i,s_i\right)\right].$$

(10)

The Hartree–Fock approximation allows one-particle wave functions to be used as variational quantities, simultaneously setting the value of the parameter ${\mathsf{\epsilon }}_{i\mathsf{\sigma }}$ equal to the one-electron excitations energy. Despite its importance for developing quantum-mechanical approaches, this method has not solved the problem of increasing the number of the system degrees of freedom with increasing numbers of atoms and electrons in it, which means that the question of rationalizing the numerical solution of the Schrödinger equation remains open.

3. The Family of Electron Density Functional Methods

Further evolution of the first-principal methods for modeling the properties and structure of nanomaterials made it possible to solve the problem of increasing complexity and dimensionality with increasing observation particles. An approach appeared where the nanosystem state energy for nuclei and electrons interacting in some external field began to be described in the form of a functional depending on the particle density. Since all properties of the system were exhaustively determined by its electron density, this method is called the electron density functional due to the principles of mathematical apparatus construction [62,63]. In this case, the electron density is a function of space but does not depend on the simultaneous arrangement of all nuclei and electrons, which significantly reduces the complexity of the equations to be solved. Hereinafter, unless otherwise indicated, mathematical calculations will be given in the atomic system of units for convenience.

The most primitive representative of the group of electron density functional methods is the Thomas–Fermi theory [64,65,66]. In the framework of this theory, the kinetic energy of interacting electrons is approximated by the corresponding expression for free particles

$$T\left[n\left(r\right)\right]=C_{\mathrm{F}}{\displaystyle \int n^{\raisebox{1ex}{$5$}\!\left/ \!\raisebox{-1ex}{$3$}\right.}\left(r\right)}\hspace{0.17em}dr=\frac{3}{10}{\left(3{\pi }^2\right)}^{\raisebox{1ex}{$2$}\!\left/ \!\raisebox{-1ex}{$3$}\right.}{\displaystyle \int n^{\raisebox{1ex}{$5$}\!\left/ \!\raisebox{-1ex}{$3$}\right.}\left(r\right)}\hspace{0.17em}dr,\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{\displaystyle \int n\left(r\right)}\hspace{0.17em}dr=N_e,$$

(11)

where $T\left[n\left(r\right)\right]$ is kinetic energy; $C_{\mathrm{F}}$ is numeric constant; $n\left(r\right)$ is electrons density; $N_e$ is total amount of electrons. The last expression in (11) is a normalization requirement.

When calculating the total energy of a nanosystem in the Thomas–Fermi model, the contribution of the exchange and correlation component is neglected, which gives the following expression

$$E_{\mathrm{TF}}\left[n\left(r\right)\right]=\frac{3}{10}{\left(3{\pi }^2\right)}^{\raisebox{1ex}{$2$}\!\left/ \!\raisebox{-1ex}{$3$}\right.}{\displaystyle \int n^{\raisebox{1ex}{$5$}\!\left/ \!\raisebox{-1ex}{$3$}\right.}\left(r\right)}\hspace{0.17em}dr+{\displaystyle \int n\left(r\right)\hspace{0.17em}{v}_{\mathrm{ion}}\left(r\right)\hspace{0.17em}dr}+\frac{1}{2}{\displaystyle \int {\displaystyle \int \frac{n\left(r\right)\hspace{0.17em}n\left(r^{\prime }\right)}{\left|r^{\prime }-r\right|}dr\hspace{0.17em}d{r}^{\prime }}},$$

(12)

where $v_{\mathrm{ion}}\left(r\right)$ is external potential field for a particular electron at a point in space. In expression (12), the second term describes the nuclei and electrons interaction and the last term is responsible for the Coulomb exchange of energy between the electrons. It should be noted that the Thomas–Fermi theory relies on the pseudo classical approximation, describes the kinetic energy of particles in a simplified way, and does not consider the exchange and correlation contributions in the interelectron interaction. Although this theory is not used to study the properties of real nanosystems, it is an important step in developing quantum mechanics methods.

The family of electron density functional methods was further developed in the works of Hohenberg–Kohn [67]. In these studies, it was proved that any system of interacting electron gas affected by an external potential can be clearly represented by the electron density, and there exists a universal functional with an extremum at the point of the electron density ground state to describe the external potential. The Hamiltonian for the system of interacting electrons in this case will have the form

$$\widehat{H}=-\frac{1}{2}{\displaystyle \sum _i{\nabla }_i^2+{\displaystyle \sum _i{v}_{\mathrm{ion}}\left(r_i\right)}}+\frac{1}{2}{\displaystyle \sum _{i,j,i\ne j}\frac{1}{\left|r_i-r_j\right|}},\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{v}_{\mathrm{ion}}\left(r_i\right)=-{\displaystyle \sum _j\frac{Z_j}{\left|r_i-R_j\right|}},$$

(13)

where $\hspace{0.17em}{v}_{\mathrm{ion}}\left(r_i\right)$ is the external potential created by the atomic nuclei; $Z_j$ are atomic nucleus charges.

All properties of a nanosystem are clearly defined by the electron density, which, in turn, is identified by the energy functional

$$E_{\mathrm{TF}}\left[n\left(r\right)\right]=F_{\mathrm{HK}}\left[n\left(r\right)\right]+{\displaystyle \int n\left(r\right)\hspace{0.17em}{v}_{\mathrm{ion}}\left(r\right)\hspace{0.17em}dr},\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{F}_{\mathrm{HK}}\left[n\left(r\right)\right]=T\left[n\left(r\right)\right]+V_{\mathrm{ee}}\left[n\left(r\right)\right],$$

(14)

where $n\left(r\right)$ is electron density; $F_{\mathrm{HK}}\left[n\left(r\right)\right]$ is system internal state functional; $T\left[n\left(r\right)\right]$ is kinetic energy functional; $V_{\mathrm{ee}}\left[n\left(r\right)\right]$ is interelectronic interaction functional.

The Hohenberg–Kohn theory can be generalized and extended to nonstationary nanosystems, systems with external magnetic and electric fields and spin-polarized configurations. In the case of spin-polarized systems, it is necessary to distinguish several groups of electrons with different densities. In this case, the external potential will be clearly identified as a function of these electron densities. For the variant of magnetic and electric fields, an additional vector field is introduced along with the scalar external potential.

The next development of the electron density theory was aimed at elaborating and investigating features of the system internal state functional $F_{\mathrm{HK}}\left[n\left(r\right)\right]$. The approach proposed by Kohn–Sham [68] implies a decomposition of the internal state functional into a number of components that have a physical justification

$$\begin{array}{c}{F}_{\mathrm{HK}}\left[n\left(r\right)\right]=T_{\mathrm{S}}\left[n\left(r\right)\right]+\frac{1}{2}{\displaystyle \int {\displaystyle \int \frac{n\left(r\right)\hspace{0.17em}n\left(r^{\prime }\right)}{\left|r^{\prime }-r\right|}dr\hspace{0.17em}d{r}^{\prime }}}+E_{\mathrm{XC}}\left[n\left(r\right)\right],\\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em} E\left[n\left(r\right)\right]=T_{\mathrm{S}}+V_{\mathrm{H}}+E_{\mathrm{XC}}+V_{\mathrm{ion}},\hfill \end{array}$$

(15)

where $T_{\mathrm{S}}$ is kinetic energy of non-interacting free electrons; $V_{\mathrm{H}}$ is Hartree or Coulomb interelectron interaction energy; $E_{\mathrm{XC}}$ is exchange–correlation contribution to energy; $V_{\mathrm{ion}}={\displaystyle \int n\left(r\right)\hspace{0.17em}{v}_{\mathrm{ion}}\left(r\right)\hspace{0.17em}dr}$ is the external potential of the nanosystem affecting the electrons. Hartree potential is used for the Coulomb repulsion energy of electrons

$$v_{\mathrm{H}}\left(r\right)={\displaystyle \int \frac{\hspace{0.17em}n\left(r^{\prime }\right)}{\left|r^{\prime }-r\right|}d{r}^{\prime }},\hspace{0.17em}\hspace{0.17em}{V}_{\mathrm{H}}={\displaystyle \int n\left(r\right)\hspace{0.17em}{v}_{\mathrm{H}}\left(r\right)\hspace{0.17em}dr}=\frac{1}{2}{\displaystyle \int {\displaystyle \int \frac{n\left(r\right)\hspace{0.17em}n\left(r^{\prime }\right)}{\left|r^{\prime }-r\right|}dr\hspace{0.17em}d{r}^{\prime }}}\hspace{0.17em}.$$

(16)

The Kohn–Sham formalism allows the electron density of the nanosystem ground state to be replaced by the density of free electrons and to search for a solution using single-particle orbitals ${\mathsf{\psi }}_{i\mathsf{\sigma }}\left(r\right)$. The index i runs through the one-electron states, the identifier $\mathsf{\sigma }$ is introduced to account for their spin. Then the nanosystem total electron density will be calculated according to the following expression

$$n\left(r\right)={{\displaystyle \sum _{i,\mathsf{\sigma }}\left|{\mathsf{\psi }}_{i\mathsf{\varpi }}\left(r\right)\right|}}^2.$$

(17)

The system of self-consistent Kohn–Sham equations is a set of one-particle Schrödinger equations

$$\left[-\frac{1}{2}{\nabla }_i^2+v_{\mathrm{KS}}\left(r\right)\right]{\mathsf{\psi }}_{i\mathsf{\sigma }}\left(r\right)={\mathsf{\epsilon }}_{i\mathsf{\sigma }}{\mathsf{\psi }}_{i\mathsf{\sigma }}\left(r\right),\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{v}_{\mathrm{KS}}\left(r\right)=v_{\mathrm{ion}}\left(r\right)+v_{\mathrm{H}}\left(r\right)+v_{\mathrm{XC}}\left(r\right),$$

(18)

where $V_{\mathrm{KS}}\left(r\right)$ is Kohn–Sham potential, also called effective potential; $v_{\mathrm{XC}}\left(r\right)$ is exchange–correlation potential. The general energy state of the nanosystem is determined by a set of wave functions ${\mathsf{\psi }}_{i\mathsf{\sigma }}\left(r\right)$ and Equation (18) have a variational character, i.e., they are solved by iterative methods. The energy functional has a single form for a given set of electrons. The difference between the Kohn–Sham theory and the Hartree method is the presence of exchange–correlation energy, which arises due to bound, not independent, motion of the electrons and the difference in energy due to their antisymmetry. There is no exact analytical description for interpreting the exchange–correlation energy and potential, so various approximations are used to estimate their value and represent the structure.

The simplest approximation for the nanosystem exchange–correlation energy is the Local-Density Approximation (LDA). This approximation is used for systems with slowly changing electron densities [69,70] and is based on the approximation of a homogeneous electron gas

$$E_{\mathrm{XC}}^{\mathrm{LDA}}\left[n\left(r\right)\right]={\displaystyle \int n\left(r\right)\hspace{0.17em}{\mathsf{\epsilon }}_{\mathrm{XC}}^{\mathrm{LDA}}\left[n\left(r\right)\right]}\hspace{0.17em}dr,$$

(19)

where ${\mathsf{\epsilon }}_{\mathrm{XC}}^{\mathrm{LDA}}\left[n\left(r\right)\right]$ is exchange–correlation energy of a single electron in a homogeneous field. There are several ways to parameterize the energy of the form (19). The approximation quite accurately describes the long-range Coulomb interaction, the kinetic component of the Hamiltonian, and the structural characteristics of real systems. The disadvantages of the local density approximation include the potential independent of the orbital type and the inclusion of the interaction of electrons with their own field, which in LDA and further modifications can be refined and eliminated.

Generalized Gradient Approximations (GGA), which take into account inhomogeneous electron density distribution, are more widely used in the theory of the electron density functional [71,72]. In GGA model, the exchange–correlation energy includes electron density gradient and is represented in the following form

$$E_{\mathrm{XC}}^{\mathrm{GGA}}\left[n\left(r\right)\right]={\displaystyle \int n\left(r\right)\hspace{0.17em}{\mathsf{\epsilon }}_{\mathrm{XC}}^{\mathrm{GGA}}\left[n\left(r\right),\nabla n\left(r\right)\right]}\hspace{0.17em}dr,$$

(20)

where $\nabla n\left(r\right)$ is electron density gradient.

If it is necessary to consider the spin state of the electrons, several groups of electrons with densities $n_{\uparrow }$ and $n_{\downarrow }$ are in for consideration. The expressions for the exchange–correlation energy (19) and (20) will be written accordingly

$$\begin{gathered} E_{\mathrm{XC}}^{\mathrm{LDA}}\left[n_{\uparrow },n_{\downarrow }\right]={\displaystyle \int n\left(r\right)\hspace{0.17em}{\mathsf{\epsilon }}_{\mathrm{XC}}^{\mathrm{LDA}}\left[n_{\uparrow },n_{\downarrow }\right]}\hspace{0.17em}dr, \\ {E}_{\mathrm{XC}}^{\mathrm{GGA}}\left[n_{\uparrow },n_{\downarrow }\right]={\displaystyle \int n\left(r\right)\hspace{0.17em}{\mathsf{\epsilon }}_{\mathrm{XC}}^{\mathrm{GGA}}\left[n_{\uparrow },n_{\downarrow },\nabla n_{\uparrow },\nabla n_{\downarrow }\right]}\hspace{0.17em}dr. \end{gathered}$$

(21)

Due to its construction, the GGA approximation only partially takes into account the inhomogeneity of the electron density distribution; therefore, the functionals of this method are also called semilocal. Improving the approach by introducing additional terms of the Taylor series energy expansion into consideration allows us to eliminate the drawback of the GGA approximation, which is transformed into the Meta-Generalized Gradient Approximations (MGGA) [73,74]. The exchange–correlation energy has a higher order of accuracy

$$E_{\mathrm{XC}}^{\mathrm{MGGA}}\left[n\left(r\right)\right]={\displaystyle \int n\left(r\right)\hspace{0.17em}{\mathsf{\epsilon }}_{\mathrm{XC}}^{\mathrm{MGGA}}\left[n\left(r\right),\nabla n\left(r\right),{\nabla }^{2}n\left(r\right),\tau \left(r\right),\mathsf{\mu }\left(r\right),\dots ,\gamma \left(r\right)\right]}\hspace{0.17em}dr,$$

(22)

where $\tau \left(r\right)$ is kinetic energy density of non-interacting particles; $\mathsf{\mu }\left(r\right),\dots ,\gamma \left(r\right)$ are local model factors. Functionals ${\mathsf{\epsilon }}_{\mathrm{XC}}^{\mathrm{GGA}}$ and ${\mathsf{\epsilon }}_{\mathrm{XC}}^{\mathrm{MGGA}}$ do not have clearly defined configuration and are chosen individually for different physical processes and problems. When constructing theoretical MGGA dependences, a sufficiently large number of experimental parameters are used. Finding and justifying these experimental quantities is a difficult task, which is a disadvantage of this approach and limits its application.

The next modification of the exchange–correlation energy in the electron density functional theory is the hybridization of the GGA method [75,76]. The general form of the expression for energy calculation is

$$E_{\mathrm{XC}}^{\mathrm{Hybrid}}\left[n\left(r\right)\right]=\mathsf{\alpha }{E}_{\mathrm{X}}^{\mathrm{EXX}}\left[n\left(r\right)\right]+\left(1-\mathsf{\alpha }\right)E_{\mathrm{X}}^{\mathrm{GGA}}\left[n\left(r\right)\right]+E_{\mathrm{C}}^{\mathrm{GGA}}\left[n\left(r\right)\right],$$

(23)

where $E_{\mathrm{X}}^{\mathrm{EXX}}\left[n\left(r\right)\right]$ is a precisely defined exchange energy potential, i.e., the Hartree–Fock-like potential; $E_{\mathrm{X}}^{\mathrm{GGA}}\left[n\left(r\right)\right]$ is exchange component in the form of a GGA functional; $\mathsf{\alpha }$ is proportionality coefficient regulating the share of exchange energies calculated by different approximations; $E_{\mathrm{C}}^{\mathrm{GGA}}\left[n\left(r\right)\right]$ is correlation energy functional.

Due to their high accuracy, hybrid functionals are popular for determining molecular configurations, atomization and ionization energies, and studies of harmonic and relaxation vibrational processes. Nevertheless, for all of the previously considered approaches in the theory of the electron density functional, the question of the choice of basis wave functions and the form of the potential remains open. For these, there are currently a large number of varieties (attached plane waves (APW), atomic sphere approximation (ASA), base construction based on Vannier functions (VF), base selection in the form of Muffin-tin orbitals (MTO) and linearized Muffin-tin orbitals (LMTO), Gaussian orbitals, Slater-type orbitals, and numerical orbital basis sets), but they will not be considered within this work.

In addition, the standard electron density functionals do not take into account the contribution of dispersion or Van der Waals forces. These forces are induced by the motion of charged particles and lead to electron density fluctuations. For a number of problems, the dispersion forces in a nanosystem can be neglected due to the smallness of their contribution, but Van der Waals forces have a significant effect on such parameters as the mechanical characteristics of ionic crystals, the cohesion energy and the surface properties of semiconductors.

Describing the energy of dispersion forces in quantum mechanics is variational [77]. Two alternative approaches are used to calculate them, which, in one form or another, bring about corrections to the interaction energies of particles. In the first case, the arising dipole fluctuations are considered as a result of the electrons’ motion, so corrections are introduced into the exchange–correlation part of the energy and the corresponding functionals

$$\begin{gathered} E_{\mathrm{XC}}\left[n\left(r\right)\right]=E_{\mathrm{X}}^{\mathrm{GGA}}\left[n\left(r\right)\right]+E_{\mathrm{C}}^{\mathrm{LDA}}\left[n\left(r\right)\right]+E_{\mathrm{C}}^{\mathrm{nl}}\left[n\left(r\right)\right]\hspace{0.17em}, \\ {E}_{\mathrm{C}}^{\mathrm{nl}}\left[n\left(r\right)\right]={\displaystyle \int {\displaystyle \int dr\hspace{0.17em}d{r}^{\prime }}n\left(r\right)\hspace{0.17em}\mathsf{\phi }\left(r,r^{\prime }\right)n\left(r^{\prime }\right)}\hspace{0.17em}, \end{gathered}$$

(24)

where $E_{\mathrm{X}}^{\mathrm{GGA}}\left[n\left(r\right)\right]$ is the energy exchange contribution obtained in the generalized gradient approximation is often written in the form of the Perdew–Burke–Ernzerhof functional $E_{\mathrm{X}}^{\mathrm{GGA}}\left[n\left(r\right)\right]=E_{\mathrm{X}}^{\mathrm{PBE}}\left[n\left(r\right)\right]$; $E_{\mathrm{C}}^{\mathrm{LDA}}\left[n\left(r\right)\right]$ is correlation energy fraction for the local density method; $E_{\mathrm{C}}^{\mathrm{nl}}\left[n\left(r\right)\right]$ is the long-range nonlinear energy correlations functional, which includes dispersion forces. The interaction function $\mathsf{\phi }\left(r,r^{\prime }\right)$ is called the core because it is the function that accounts for the differences and changes in the charge density gradient. The choice of the form of this function is important. Depending on the form $E_{\mathrm{C}}^{\mathrm{nl}}\left[n\left(r\right)\right]$, the scientific literature distinguishes between calculating the dispersion forces using the Dion method (vdW-DF) [78] and Random Phase Approximation (RPA) [79].

An alternative approach to describe the dispersion forces is their representation, not in terms of electron density, but through the atom’s interaction. The total energy of the nanosystem also contains an additional summand

$$E_{\mathrm{tot}}=E_{\mathrm{DFT}}+E_{\mathrm{disp}},\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{E}_{\mathrm{disp}}=-{\displaystyle \sum _{i,j}{\displaystyle \sum _{k=6,8,10,\dots }\frac{C_k}{{\left|R_i-R_j\right|}^k}}},$$

(25)

where $E_{\mathrm{tot}}$ is total energy; $E_{\mathrm{DFT}}$ is energy obtained by the electron density functional method within the Kohn–Sham formalism; $E_{\mathrm{disp}}$ is dispersion correction including interaction between instantaneous and induced dipoles, permanent and induced dipoles; $C_k$ are numerical coefficients; $R_i$ и $R_j$ are radius vectors of atomic nuclei. The summation in the second expression (25) by i and j is performed on all pairs of atoms in the nanosystem, the index k being responsible for the varieties and perturbations arising from Van der Waals interactions.

Due to the high degrees in the denominator of the dispersion correction, these forces decrease rapidly with distance. Therefore, quantum-mechanical calculations often make do with one or two terms in the summation (25). The first element of the sum, which is inversely proportional to distance to the sixth power, characterizes the dipole interaction. The description of dispersion forces through pair potentials using the damping function $f_{\mathrm{dmf}}$ is common in the scientific literature

$$\begin{gathered} E_{\mathrm{disp}}=-s_6{\displaystyle \sum _{i,j}\frac{C_{6,ij}}{{\left|R_i-R_j\right|}^6}}{f}_{\mathrm{dmf}}\left(\left|R_i-R_j\right|\right), \\ {f}_{\mathrm{dmf}}\left(\left|R_i-R_j\right|\right)={\left[1+\exp \left(-\mathsf{\alpha }\left(\frac{\left|R_i-R_j\right|}{s_R{R}_r}-1\right)\right)\right]}^{-1}, \end{gathered}$$

(26)

where $C_{6,ij}$ is pairwise dispersion coefficient, depending on the types of interacting atoms; $s_6$ is the scaling factor determined by the type of damping function.

The damping function is used to avoid singularities at small distances and to prevent double-counting of correlations. The form of the function is variational, depending on the method used to compute the dispersion forces. In (26) we present the expression formulated in the Grimme methods (DFT-D, DFT-D2, DFT-D3) [80], for which the parameters of the equations have the following valuations: $C_{6,ij}=\sqrt{C_{6,i}{C}_{6,j}}$, $R_r=R_{vdW,i}+R_{vdW,j}$ is the sum of the Van der Waals radii for the two atoms, $s_R=1,0$, $\mathsf{\alpha }=20.0$. An alternative form of damping function notation is used in the Tkachenko and Scheffler method (DFT+vdW) [81], bound cluster theory (DFT/CC) [82], Becke and Johnson approximation (XDM) [83]. There are also approaches for calculating the dispersion forces not only through pairwise interaction, but also taking into account the multi-particle contribution, e.g., the multi-electron dispersion method (MBD) [84].

4. Quantum Molecular Dynamics Method

A specific direction in quantum mechanics is developing theoretical methods for studying the properties of complex molecular systems called quantum molecular dynamics. The method is based on the Born–Oppenheimer approximation, which assumes the partitioning of variables into massive nuclei and fast electrons. The essence of this formalism has been previously described by relations (3)–(5). In the approach originally described, the solution of a partial differential equations system (3) was obtained through a decomposition into basis functions with a further diagonalization. The quantum molecular dynamics method proposed by Car–Parrinello [85] uses function minimization.

The idea of Car–Parrinello molecular dynamics (CPMD) is the combined calculation of electronic and atomic subsystems. The motion of nuclei is determined by the set of their coordinates $\left\{R_i|\hspace{0.17em}i=1,\dots ,N_n\right\}$, the electronic degrees of freedom are characterized by a set of quantum-mechanical wave functions $\left\{{\mathsf{\psi }}_j|\hspace{0.17em}j=1,\dots ,N_e\right\}$

$$M_i{\ddot{R}}_i=-\frac{\partial E}{\partial R_i}=F_i,\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\mathsf{\mu }{\ddot{\mathsf{\psi }}}_j\left(r,t\right)=-\widehat{H}{\mathsf{\psi }}_j\left(r,t\right)+{\displaystyle \sum _k{\mathsf{\Lambda }}_{jk}{\mathsf{\psi }}_k\left(r,t\right)},$$

(27)

where $F_i$ is the resultant force exerted on the atom; $E$ and $\widehat{H}$ is the energy functional and the Kohn–Sham Hamiltonian; $\mathsf{\mu }$ is fictitious electron mass; ${\mathsf{\Lambda }}_{jk}$ are indefinite Lagrange multipliers. The fictitious electron mass is chosen to organize the convergence of the algorithm as quickly as possible. The Lagrange multipliers are responsible for the additional forces that ensure the orthonormality of the wave functions.

In the Car–Parrinello theory, nanosystem total energy is a function of the total coefficients of the wave electron function expansion on some basis. When minimized by a given coefficient, there is a cooling and stabilization of the system state. The first part of expressions (27) is the classical Newton motion equations for a set of particles, which are solved iteratively until the equilibrium condition is reached. Consistency of the ion and electron subsystems is achieved by simultaneous minimization of the energy functional. It should be noted that this method does not describe the true nanosystem dynamics, but rather its fictitious evolution, leading to an equilibrium energy-stable state of a multielectron particle system.

The Car–Parrinello molecular dynamics is known now to reproduce the properties of semiconductor and dielectric materials quite well. For metallic systems near the band gap, there are many states with closely spaced eigenvalues, which brings about significant variations in the electron density even with a small change in the total energy [86]. In addition, the use of the Car–Parrinello approach is limited to a small number of nanosystem atoms due to the high complexity of solving the equations and significant computational costs.

There is another approach to quantum molecular dynamics, which is based on the Born-Oppenheimer approximation. This approach consists in separating the description of the nuclear and electron subsystems, hence it is called the Born-Oppenheimer molecular dynamics (BOMD) [87,88]. The motion of the nuclei is described by the classical equation of mechanics, but the energy and forces acting on them are calculated by solving the Schrödinger equation at each time step for the electronic wave functions. The essence of this approach results in the following expression

$$M_i{\ddot{R}}_i=-{\nabla }_i\left[\underset{{\mathsf{\psi }}_1,\dots ,{\mathsf{\psi }}_{N_e}}{\min }E\left(R_1,\dots ,R_{N_n},{\mathsf{\psi }}_1,,\dots ,{\mathsf{\psi }}_{N_e}\right)\right],$$

(28)

where $E$ is the total energy, which is minimized by the set of electron wave functions at fixed coordinates of the nuclei. Minimization is usually performed through the solution of the Kohn–Scham equations or by means of multidimensional optimization algorithms.

The BOMD method is highly accurate and can sufficiently describe, for example, the properties of dielectrics and metals. Nevertheless, the implementation of this approach, as well as the CPMD method, requires significant computational costs. Therefore, BOMD and CPMD are mainly used for nanosystems of small size and volume.

In order to reduce the complexity and increase the capacity of the system under study, methods of combining classical molecular modeling with electronic structure calculations have been actively developed. This idea has found application in the embedded atom potentials (Embedded Atom Model–EAM, Modified Embedded Atom Method–MEAM) [89,90]. The Embedded Atom Method derives from the theoretical positions of the electron density functional and describes the behavior of a nanosystem in the form of the following set of equations

$$M_i{\ddot{R}}_i=-\frac{\partial U_i}{\partial R_i}=F_i,\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{U}_i=F\left({\mathsf{\rho }}_i\right)+{\displaystyle \sum _{j,i\ne j}\mathsf{\varphi }\left(\left|R_i-R_j\right|\right)},\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{\mathsf{\rho }}_i={\displaystyle \sum _{j,i\ne j}\mathsf{\psi }\left(\left|R_i-R_j\right|\right)},$$

(29)

where $F\left({\mathsf{\rho }}_i\right)$ is immersion function characterizing the interaction of an individual atom with a medium with electron density ${\mathsf{\rho }}_i$; $\mathsf{\varphi }\left(\left|R_i-R_j\right|\right)$ is pair potential; $\mathsf{\psi }\left(\left|R_i-R_j\right|\right)$ is the electron density created by the j-th atom at the location of the i-th atom.

The use of immersion functions in the EAM and MEAM methods makes these potentials multi-particle, more accurately reflecting the properties of many nanomaterials, including metals. A limitation of this theory is the need to prepare potential parameters as there are a rather large number of them. Parameters are obtained mainly by experimental methods or by other theoretical approaches. Nevertheless, due to its speed and ability to analyze huge nanosystems (millions of atoms), the embedded atom method remains one of the popular methods for studying the properties of nanomaterials.

In conclusion of this section, we present a general scheme of methods for modeling nanostructures, shown in Figure 1, which allows us to evaluate the areas of different methods application in terms of mathematical complexity to the practical problems of nanosystems design and research.

5. Application of Quantum Mechanics Methods to Calculate the Optimal Structure of Molecules and Parameters of Force Potentials

The results presented in this section illustrate the application of quantum mechanics methods to optimize the geometry of molecular systems and find their interaction parameters. In applied problems, there are often special-purpose compositions where modeling parameters are not known. This problem is particularly relevant for metal salts found, for example, in condensation problems or technological processes. The above problem can be solved by using the quantum mechanics apparatus at the preliminary stage of modeling, which provides the search for interaction parameters and necessary constants for further investigation of mechanisms of agglomeration and structuring nanoparticles of the target composition.

We show the effectiveness of quantum-mechanical modeling of the functional nanostructure properties using the example of a relevant applied problem. We consider the effectiveness of preliminary preparation of nanoparticle composition modeling parameters used to fertilize plants with nanoparticles through the pores formed in the gas medium [49,91]. Fertilizing the soil is ineffective because most of it does not reach the root system, and hence does not get inside the crops. Therefore, developing new approaches of directed feeding to increase their yields is currently in demand. There is a fairly efficient and economically valid method of introducing minerals by burning a smoke composite pellet, which contains in its composition a complex of necessary nutrients. The expediency and economic validity of the approach of applying fertilizer through the gas medium is currently confirmed by agronomic studies, and the method is actively used in agriculture.

For the above-described problem of feeding plants from a gas medium, its molecular composition was determined. The study of the condensation processes of the mixture molecules was impossible without setting the parameters of the force potentials and molecules’ equilibrium states. Equilibrium optimal configurations were determined for all types of considered molecules. To calculate the equilibrium structure of the initial substance molecules, quantum mechanics was used, based on the Hartree–Fock method and the approximation of the linear combination of atomic orbitals for calculating molecular orbitals. The stationary Schrödinger equation in the Born–Oppenheimer formulation (1)–(3) is solved. The method was implemented in the HyperChem program. The 6-31G* basis was chosen as the basis, which sufficiently describes a wide class of chemical elements and compounds. This basis was used because it covers different types of electronic orbitals, and therefore is able to describe the energy states and behavior of metals and their salts, which is important for this applied problem.

The equilibrium geometry of the molecules was calculated by the Fletcher–Reeves conjugate gradient method [92,93]. As a result of the optimization and energy analysis of the nanosystem elements molecular structure, such data as standard bond lengths between atoms, values of equilibrium angles, intramolecular charge redistribution due to orbital rearrangement, and force-field interaction constants were obtained.

For simple and well-known molecules such as H₂O, N₂, O₂, and CO₂, the obtained geometrical data were used for comparison with the already known results in order to confirm the accuracy of the mathematical model construction and completeness of the chosen basis. The calculated values for comparison with the previously known confirmed data [94,95,96,97,98] are shown in

Table 1. Comparison of quantum-mechanical modeling results with previously known theoretical data.

Molecule	Optimization Parameters	Confirmed Data	Simulation Data	Error, %
H2O	H-O bond length, nm	0.0957	0.0947	1.05
	H-O-H angle size, deg	104.523	105.502	0.94
	H Atom charge, 10−19 C	+0.668	+0.695	4.04
	O Atom charge, 10−19 C	−1.336	−1.391	4.19
N2	N-N bond length, nm	0.1095	0.1078	1.55
O2	O-O bond length, nm	0.1210	0.1166	3.64
CO2	C-O bond length, nm	0.1163	0.1143	1.72
	O-C-O angle size, deg	180.00	180.00	0.00

. A graphical representation of molecules with optimized geometry for the example is demonstrated in Figure 2: (a) carbon dioxide molecule, (b) water molecule, and (c) potassium carbonate molecule.

The data comparison presented in Table 1 for nitrogen, oxygen, water, and carbon dioxide molecules shows a good correlation of the obtained results and a small error value. These types of molecules were chosen as test molecules because of their prevalence and availability of reliable information. The absolute value for the bond length did not differ from the previously known data by more than 0.005 nm. The relative value of the compression–expansion error of the interatomic bond was 3.64% in the maximum calculation and was obtained for an oxygen molecule. The effective charges in the molecules in question were determined without exceeding a deviation of 5%. The insignificant difference between the quantum-mechanical calculations and the information obtained by other authors testifies to the accuracy of the chosen orbital basis and the adequacy of the constructed mathematical model.

To summarize, we can conclude that the quantum mechanics exact methods can be and are used for the problems where detailed distributing electron clouds, accounting for their equilibrium state and full energy spectrum are important. Modeling data in such problems largely depend on the wave function basis, but at the moment there is a fairly wide set of suitable bases. Also, when using the exact methods of quantum mechanics, there is a strong correlation between the complexity of solving the problem and the number of atoms and electrons under consideration. Therefore, such modeling approaches are common for the study of small systems.

6. Application of the Methods of the Density Functional Theory to Analyze the Structure of Nanomaterials

The next applied problem in this work is the application of the density functional theory to the study of the nanomaterials structure. Structure analysis problems in the nanotechnology area are commonplace; often they need to be solved with adaptation to specific technological parameters and processes. The number of atoms considered in a nanosystem is determined individually in each problem, but structural defects often have an extended spatial character and cannot be investigated by direct methods of quantum mechanics. A way out of the situation can be, for example, the use of the molecular dynamics apparatus with force potentials that take into account the electron clouds nature and distribution in the system.

Below we will consider an example of modeling the relaxation of pure metal nanocrystals at constant normal temperature and pressure, and periodic boundary conditions. The force interaction in the problem was described using the modified embedded atom method. At the initial moment of time, the atoms of the considered metal were located in the nodes of the crystal lattice. The lattice structure and its dimensional properties depended on the type of the initial metal. For Cu and Au a cubic face-centered structure with crystal parameters of 3.615 Å and 4.0781 Å is typical. The appearance of the simulated crystal cube is demonstrated in Figure 3a. Computational experiments were aimed at studying the relaxation and self-organization processes of metal atoms in the solid state. Atoms moved near their crystal lattice nodes at rest and in the absence of external forces, but did not leave them.

The concept of short-range order helps to investigate the internal structure of a material in detail. Short-range order is a mutual ordering in the arrangement of particles in a material that is translated and repeated over insignificant distances commensurate with the size of the atom. Analysis and description of short-range order reveals patterns in structure based on the calculation of pairwise correlation and radial functions. The radial material function maps independent correlations between atoms in a nanosystem and characterizes the mutual arrangement of particles as a function of the central object radius. The radial function determines the probability that other atoms are placed in an elementary spherical layer. The principle of calculating the radial function is illustrated in Figure 3b.

For gold and copper, radial functions were calculated by modeling to compare with previously known theoretical or experimental data. Due to their similar crystal structure, the graphs of the obtained functions in Figure 4 are of the same type; the markers show a set of reference points of previously confirmed similar data by other authors [99,100]. Comparison of the radial functions for nanocrystals of the considered elements with the previously obtained data shows satisfactory agreement of the results. The radial functions have a harmonic character with a varying amplitude. Peaks and change of decreasing and increasing intervals on graphs coincide; points of extremums are close on the numerical values between them.

The obtained correspondence between the radial functions of the modeled system and the previously known reliable theoretical and experimental data allows us to speak of a properly chosen and adequately functioning potential. The structure of the nanocrystal in the modeling process is preserved, the atoms move vibrationally near the lattice nodes, as justified from a physical point of view and is completely accurate. The first and most significant maximum of radial functions corresponds to the number of nearest neighbors in the crystal structure.

Further analysis of the nanocrystal atom behavior was aimed at investigating the accuracy of the initial and current atom velocity distributions. The distribution of the nanosystem velocities was controlled by their autocorrelation function, for different moments of time presented in Figure 5.

The autocorrelation function shows the mutual dependence between the initial velocity of the atoms and the calculated velocity values at subsequent points in time. For the considered nanosystem, the instantaneous values in Figure 5 vary, but the autocorrelation function’s average value is zero. The limiting tendency to zero indicates the independence of the velocity distribution in relation to time, as it should be for correctly solved atom motion equations.

The results obtained in solving this applied problem show that it is possible to reliably describe the properties and structure of atomic systems of large size (millions of atoms and more) using potentials based on the theory of the electron density functional. These potentials are being actively improved and adapted for a more physically based reproduction of the wide range of element properties. In contrast to the direct quantum mechanics methods, the density functional theory enables the analysis of large-scale defects of material structure, such as pores, channels, cracks, and dislocation accumulations.

As an example to illustrate the modeling of a real material structure, we consider the structure of a multilayer Nb-Co spin valve nanosystem, the detailed formation modeling of which was performed by the molecular dynamics method using the embedded atom potential [101].

The images in Figure 6 characterize the qualitative view of the processes of spin valve formation from the niobium and cobalt layers and the structure of the layers. The structure of the layers formed by niobium atoms is close to crystalline. The groups of atoms are combined into nanocrystallites with vertical spatial orientation. In Figure 6, thin vertical solid lines highlight the boundaries of these nanocrystals. Cobalt nanofilms have an amorphous structure. The obtained results are in good correspondence with the experimental studies of the similar multilayer nanosystems structure [102].

7. Conclusions

The quantum mechanics methods undoubtedly occupy a key position in the study of the nanosystems’ behavior and their properties. Obtaining information on the crystal structure, their surface properties, energy characteristics, point effects, optical parametrons, activation energies and other important features is often based on first-principles calculations. High accuracy of quantum mechanics methods is closely connected with the complexity of used algorithms and approaches, which to some extent is a limitation of this theory.

Various assumptions and modifications are able to eliminate limitations of the quantum mechanics model, as well as reduce complexity and dimensionality of solved problems by reducing the number of investigated objects’ degrees of freedom. The basic approximations of the first-principles methods include the Born–Oppenheimer approximation, which assumes the separation of variables into massive nuclei and fast electrons, and the Hartree–Fock approximation, which transforms the dependence of the Schrödinger equations on the multi-electron wave function into a correlation relation of one-electron functions.

An important place in the quantum mechanics methods is occupied by the theory of the electron density functional, which allowed us to solve the problem of increasing complexity and dimensionality when the number of observation particles increases. In this approach, the nanosystem state energy for nuclei and electrons interacting in a certain external field is described in the form of a functional that depends on the particle density. Evolutionary changes and variants of describing the exchange–correlation energy, including Thomas–Fermi theory, Hohenberg–Kohn theorems, Kohn–Sham formalism, local density approximation, generalized and meta-generalized approximations and hybridization of the generalized gradient method demonstrate the relevance and importance of the electronic density functional theory, and show the development of this method.

There are approaches for introducing corrections to the quantum mechanics theory to account for dispersion or Van der Waals forces, which are essential in determining such parameters as mechanical characteristics of ionic crystals, cohesion energy and surface properties of semiconductors. Known methods for introducing dispersion corrections consider them as part of the exchange–correlation energy of the system and as an element of pair interaction of atoms. The development of quantum mechanics methods to quantum molecular dynamics or Car–Parrinello dynamics is shown. The basic idea of combining classical molecular modeling with calculations of the electronic structure, which is reflected in the embedded atom potentials, is described.

Despite the successful use of quantum mechanics methods to study the properties of real nanomaterials and nanosystems, the improvement of these approaches continues. At present, theoretical methods and practical implementations in the field of magnetic and superconducting materials, quasi-particle interactions, high-temperature and ultra-low-temperature behavior of nano cells are under close attention. The development of quantum mechanics approaches in these scientific fields will make it possible to obtain fundamentally new materials and discover promising properties of already known objects in the future.

To summarize, we can conclude that the applicability of a particular quantum mechanics method depends on the type of problem to be solved. Exact methods, including the Hartree–Fock approximation, are appropriate when the size of the nanosystem is small, for example, in the case of calculating the optimal configuration of molecules or obtaining energy constants for molecular dynamics. The theory of the electron density functional is indispensable when electronic, phonon, and spin effects are important, but the size of the nanosystem exceeds the limitations of the quantum mechanics exact methods, and additional approximations cannot be avoided. There are also a number of problems where the order and arrangement of atoms in a nanomaterial is under close consideration. The influence of edge effects and structuring mechanisms, defect migration and issues of their occurrence arising during the nanofilm deposition or nanocomposite deformation require considering a large number of atoms. Such problems have been successfully solved by embedded atom methods within the framework of molecular dynamics theory, since these approaches reproduce the behavior and describe the properties of metals and closely packed nanoobjects quite accurately.

Despite their rather long existence, the quantum mechanics methods are actively developed and successfully used to study the properties of atomistic systems, to analyze the mechanisms of their interaction, self-organization and structuring. Depending on the structural features, atomic or molecular symmetry, nanomaterials can significantly change their properties. For example, fullerenes and nanotubes can become semiconductors or have high-temperature superconductivity if their molecular symmetry is broken. Deviation from symmetrical contact between nanofilms can negatively affect their functional properties. Symmetry effects and cases of symmetry breaking make nanomaterials unique, which requires a comprehensive study using, among other things, mathematical modeling.