New Orbital Free Simulation Method Based on the Density Functional Theory

A practical way to simulate multi-atomic systems without using of wave functions (orbitals) is proposed. Kinetic functionals for each type of atoms are constructed and then are used for complex systems. On examples of clusters containing Al, Si, C, and O it is shown that this method can describe structures and energies of multi-atomic systems not worse than the Kohn-Sham method but faster. Besides, it is demonstrated that the orbital-free version of the density functional theory may be used for finding equilibrium configurations of multi-atomic systems with covalent bonding. The equilibrium interatomic distances, interbonding angles and binding energies for Si3 and C3 clusters are found in good accordance with known data.


Introduction
The orbital-free (OF) approach is an alternative to the Kohn-Sham (KS) [1] method to simulate multi-atomic systems in the framework of the density functional theory (DFT) [2]. The OF approach operates with the electron density only (without wave functions) and being developed enough can be applied for simulation of very large systems: up to millions atoms [3]. Several groups [3][4][5][6][7][8][9][10][11] are working in this area with different success, and the calculation of the kinetic energy is noted as a main problem. In the previous papers [12,13] it was suggested that there is no universal way to describe the kinetic energy of different atoms and compounds. In the present work it is described how it is possible to extend this approach to systems with more large systems.

A General Description of the of Approach
As it is known the DFT claims that the energy E of the ground state of any quantum system can be found by minimization of the some functional depending only on the electronic density of this system ρ(r): (1) where V(r) is an external potential, | | There are some realistic approximations for exchangecorrelation potential µ ex− c ( ρ) there; the potential Hartree ( ) φ r may be calculated using Fourier transformations or Poisson equations; the external potential V( ) r usually consists of atomic potentials. The only real problem is the kinetic potential µ kin .
Pseudopotential Approach for Dimers In practice, the DFT calculations are simpler if one uses pseudopotentials instead of full electron potentials. Therefore let us rewrite the above equations in the pseudopotential approach, and, for simplicity, let us limit yourself by twoatomic systems and s-and p-components of pseudopotentials.
We present the total density 12 ρ as a sum of partial ε is the electron energy per electron for the two-atomic system with the total number of electrons N 12 . In the other words we have to find the density ρ 12 that satisfies the system of two equations:  (6) where где (4-s), (4-p) for a two-atomic system as it follows below: Obviously, for two isolated atoms it possible to write equations similar to (7) and (8): (13) From which it is followed: (17) Let us introduce notations: Now equations (5-s) and (5-p) can be transformed to simple forms: 0 12 s 12-12 12 12 (23) It is possible to solve these equations using some probe functions for 12

12-s
The probe functions must lead to the equilibrium interatomic distances and binding energy for dimers.
The binding energy for a dimer (per one atom) would be calculated as follows: Dimers with identical atoms Al, Si, and C were taken as test elements. The FHI98pp [14] package was used as a generator of pseudo-potentials and equilibrium partial electron densities. Exchange and correlation potentials were calculated in the local density approach [15,16]. Studied atoms were located in a cubic cell of the L size (L=30 a.u.; 1 a.u. = 0.529 Å). The cell was divided on 150×150×150 elementary sub-cells for the integration with the step ∆L of 0.406 a.u. Results of calculations were compared with published data.  Namely, for Al they are: Calculated values of interatomic distances and binding energies for the Al 2 , Si 2 , and C 2 dimers are collected in Table  1 in comparison with known published data. It is clear that agreement is rather good. Notation: "KS" means the calculations used the FHI96md package [17] based on the Kohn-Sham method.

Dimers with different atoms
As different atoms have different functions for kinetic energy, some procedure to calculate the total kinetic functions kin 12 ∆µ in the space of the atomic system has to be developed. Near each atom it has to be approximately equal to its atomic function, but it has to be equal to mixture of the specific atomic functions between atoms. It seems that the simple way to construct this total function kin 12 ∆µ is to summarize the specific atomic functions with some weights: The weights have to be determined through gauss functions those are fitted to atomic densities: An example of fitting of densities is demonstrated, Figure  1 the fitting of atomic densities for silicon and oxygen. Value A and B for Si, Al, C and O are collected in Table 2.   Table 3. It is clear from Table 3  Interbonding angles To describe the angle depending of interatomic bonding it is useful to analyze the reasons of this depending in the standard quantum-mechanical approach, which uses wave functions and electronic states. As it is specified in the work [25] the angle peculiarities of the cluster Si 3 are defined by the Yang-Teller effect caused by existence of the energy gap between occupied and empty states. In other words, the difference of structures of semiconductor and metal small clusters is connected with the difference of their bond wave functions: namely, covalent atoms have localized functions orientated between nearest atoms, while metallic atoms have dispersed functions without orientation in the space.
In the orbital-free case wave functions are absent, electronic states are absent too, and, therefore, one cannot speak about any energy gap. In the OF approach there is only the electronic density which defines all energy and structure of the polyatomic system. However the main quantummechanical rules remain fair and in this case. Besides the Schrödinger's (or Kohn-Sham) equations out of which wave functions and electron states are brought, there is the Paulie's principle specifying that in one quantum state there can be only two electrons (without taking into account a spin). In the OF case this principle may be paraphrased by the following way: a covalent bond is formed by two electrons, the common wave function of which is localized in the space between two nearest atoms. It is obvious that the quantity of the electrons which are responsible for this bond doesn't change as the distance between atoms changes (if, of course, the bond isn't broken at all and the electronic structure isn't reconstructed completely). In case of metals the conduction states are close each other and electrons can easily "flow" from one state to another during the changing of the atomic geometry.
Talking about covalent bonding in the language of the electronic density one obliges to base on the results of wave functions calculations. The main result of such calculations is that an atomic system with saturated covalent bonds (having two electrons in each bond) has extremely clear energy gap between occupied and empty states.
Let us consider some multi-atomic system which density has local maximums ρ ij between the nearest neighbors i,j. Let a part of these maximums overcome the value of ρ 0 in some etalon system (for example, in dimer or trimer). Then the excess of energy appears: C is an unknown coefficient. The density between the nearest atoms may be changed by two ways: 1) changing the interatomic distance; 2) changing the interbonding angle. The first way is caused by quantum forces calculated from interatomic interactions described above. The second way may be characterized by some additional, so called Paulie forces.
The electronic density of the trimer ρ trim may be found as follows: where R A , R B , and R C are coordinates of points in which the A, B, and C atoms with densities ρ at are situated. The binding energy is   Figure 2 presents dependence of the energy of the Si 3 cluster on the angle between bonds. It demonstrates a minimum at the angle of 78º. For C 3 a minimum was found at 180º (a linear chain), and for Al 3 (where covalent bonds absent) all angles were found of 60º (an equilateral triangle).
Equilibrium values of interatomic distances, angles α, and binding energies for the trimers Al 3 , Si 3 , and C 3 are collected in Table 4 in comparison with known data and results of Kohn-Sham calculations. Parameter C of 25000 was used for silicon and carbon clusters; for aluminum C was equal zero. To experience the OF approach for ability to describe correctly three-dimensional systems we investigated fouratomic clusters Si4, Al4 and С4 with structures of straight lines, rhombuses, and trigonal pyramids. It is known that the Si4 cluster is built in the form of a rhombus [18,19]; and four carbon atoms form a straight linear chain [23,24]. As for aluminum, there is no consensus [21,22,31]. Some authors claim that the most favorable configuration is rhombus, others favor pyramid. What gives our method?
Our Calculations show that for the four-atomic cluster the most favorable configuration is the linear chain with the interatomic distance of 1.27 Å and the binding energy of 6.3 eV per atom. The Kohn-Sham results are 1.30 Å and 6.76 eV.
It had been found that the rhombus configuration is favorable for the Si4 cluster with the interatomic distance of 2.35 Å and the binding energy of 2.7 eV (the gain is 0.7 eV in comparison with a pyramid). The Kohn-Sham results are 2.32 Å and 4.1 eV (the gain is 0.6 eV).
In the Al case our calculations demonstrate benefit of the pyramid; however the difference with the rhombic configuration is only 0.2 eV per atom. The Kohn-Sham method yields the opposite results: 0.13 eV in favor of a rhombus.

Conclusions
The basic principles of modeling of atomic interactions within of an orbital-free version of the density functional theory were formulated; and modeling of some clusters was carried out. The possibility to simulate interactions of atoms of non-identical types in the framework of the orbital-free version of the density functional theory it is shown. A rather simple technique was used for this purpose, namely: first, the atomic kinetic functions were found for homo-atomic dimers Si 2 , Al 2 , C 2 and for the SiO dimer; second, some atomic weights were proposed using gaussians associated with atomic densities; third, kinetic functions for hetero-atomic dimers were constructed. Equilibrium interatomic distances and dissociation energies for the SiC, SiAl, AlC, SiO and CO dimers were found in good comparison with other data.
Restriction principle for the interatomic density (following from Paulie's principle) allows us to describe angular dependences of the interatomic bonding in polyatomic clusters. In particular, it had been shown that for the Al 3 cluster the equilateral triangle is favorable; the Si 3 trimer is characterized by the isosceles triangle with angles of 80 and 50 degrees, and three atoms of carbon built the linear chain. Calculated equilibrium interatomic distances and the values of binding energy are well compared with known data.
As calculation of kinetic energy and interbonding angles are key points in modeling of polyatomic systems in the orbital-free approach, it is possible to consider that this work opens a direct way to design an effective method of modeling of complicated nanosystems and supermolecules with a high number of atoms.