Elsevier Science Home
Computer Physics Communications Program Library
Full text online from Science Direct
Programs in Physics & Physical Chemistry
CPC Home

[Licence| Download | New Version Template] acjd_v1_0.gz(28 Kbytes)
Manuscript Title: Equation of motion method for the electronic structure of disordered transition metal oxides.
Authors: M.T. Michalewicz, H.B. Shore, N. Tit, J.W. Halley
Program title: Eq_of_Motion
Catalogue identifier: ACJD_v1_0
Distribution format: gz
Journal reference: Comput. Phys. Commun. 71(1992)222
Programming language: Fortran.
Computer: CRAY-2.
Operating system: UNICOS 5.1, UNICOS 6.1.5.
Word size: 64
Keywords: Solid state physics, Band structure, Electronic structure, Disordered systems, Transition metal oxides, Insulators, Semiconductors, Tight binding Hamiltonian.
Classification: 7.3.

Nature of problem:
The equation of motion method, as implemented in our program, consists of the following sequence of steps. A solid is described by a tight binding Hamiltonian. The time evolution of a system is determined by the Schrodinger equation for the amplitude of the Green's function F. It is formally solved by polynomial expansion of the exponent. The time evolved amplitude is then Fourier transformed to energy domain. The negative imaginary part of this quantity, divided by Pi gives the density of states. Depending on the initial conditions it is possible to calculate the local density of states at any selected site(s) j (Fi(t=0) = delta i,j) or the total density of states (random F(t=0)). The sample is periodic in x and y directions and bound by surfaces normal to the z axis. It can be 'grown' in three dimensions, so 2D, 3D and the surface density of states (DOS) can be obtained. The defects (oxygen vacancies) are described by a screened Coulomb potential and a positive soft core. The tables of neighbours, direction cosines and distances between the pairs of ions and the Hamiltonian matrix including the orientational corrections obtained from the Slater-Koster prescript- ion are saved. The expansion of the exponent is an efficient computational strategy for the following reason. For a given energy resolution, Delta E in the final spectral function (here usually the density of states), and a bandwidth B, one only needs the value of the corresponding function of time at N=B/Delta E points separated by intervals h/N Delta E. (In the application in question, for example, N=256 time points are often adequate.) On the other hand, if one uses a low order expansion of the exponent, a much shorter time step will be required in the integration of the equations of motion. Thus one can either use a low order expansion and many time steps, using only a few of the resulting time values in the resulting Fourier transform or, alternatively, a high order expansion which takes one directly from one point used in the Fourier transform to the next. The computational cost of the two approachs would be same if the number of terms required in an expansion of the exponent increased linearly with the time interval over which the expansion is interpolating. However this is not the case. The number of required terms in the exponential expansion increases more slowly than linearly. As a result it is more cost efficient to use the second method. (An idea of the advantage of the second method may be easily gleaned by expanding e**(i Pi) to a needed accuracy and comparing the number of terms required to the number of terms required to expand e**((i Pi)/2) twice. That difference (12 terms versus 16 terms) is of order 30 per cent at the accuracies to which we are working.) An additional speed-up of 20 per cent was achieved by using Chebyshev fit to an expotential function.

The program is linear IordMx both in central computer memory and in cpu time, hence the only limitation is the available computer time and the memory. Note, however that the program only vectorizes on computers with hard wired gather scatter operations as discussed in section 7 of the long writeup. The size of the system can be from 1 rutile unit cell up to about 10**3 rutile unit cells. However the periodic boundary conditions require that the number of unit cells in x and y directions be even numbers. This condition can be relaxed for very small (non- periodic) systems.

Running time:
The approximate formula for the cpu time on Cray Y-MP is Nions x 0.1 seconds. The system consisting of approximately 4 x 10**3 ions requires approximately 3 x 10**2 seconds on the Cray-2 or Cray Y-MP.