In July 1998, we organized a CECAM workshop on Local Orbital Methods for Large Scale Atomistic Simulations. At that time, there was renewed interest in the methods based in localized orbitals, from empirical TB to accurate ab-initio calculations. Part of the motivation of this interest was the realization that localization concepts provide a road towards low complexity electronic structure algorithms (by that meaning a low scaling of the computational workload involved). The workshop was devoted to a wide range of aspects, including TB parametrizations, applications of TB formulations, approximate DFT methods, fully selfconsistent DFT calculations, quantum chemistry methods, linear scaling algorithms and non-standard localized basis sets.
However, one of the key questions that arose was the issue of the optimal basis sets for exploiting the localization properties in ab-initio electronic structure methods. Since that date, much work has been done in different lines, and this is the subject of the present meeting.
As a matter of introduction to this workshop, I will review some of the ideas that were discussed in the previous one, and advances made in the field since then.
During the past about 20 years enormous progress has been made in the development of efficient methods of computations of molecules, clusters and solids within Density-Functional Theory (DFT). There is a lot of experience on the performance of exchange-correlation energy (EXC) functionals, though from a 'first principle' point of view the status is still unsatisfactory. Although only approximate model density functionals are known, even for their application work has to be done to improve the efficiency and numerical accuracy of the solution of the Kohn-Sham (KS) equations. Nowadays, local orbital basis sets for the KS-orbitals are widely used in DFT methods. Atomic orbitals (AO's) and their modifications, though they have been used since the early days of the quantum mechanical treatment of molecules and solids, are by far not obsolete in the age of high performance computing of quantum mechanical many particle systems. Their advantages are clearly their physical meaning, allowing the smallest basis sets of all basis set based methods and their general suitability and flexibility, making them easily applicable for systems containing any type of atom of the periodic table. These advantages are contrasted by the cumbersome numerical handling and by the practical problems in basis set extensions for high energy excited states. For the highest accuracy demands the AO's are represented numerically on a radial mesh [1]. Numerical handling of the AO's can be much improved by the use of linear combinations of analytical functions as Gauss type functions (GTO) or Slater type functions (STO). The advantage of the STO's is their correct asymptotic behaviour on cost of a higher computational effort for the integral calculations, compared to GTO's. But usually the basis sets of STO expansions can be considerbly smaller than GTO expansions at the same level of accuracy.
We often use STO expansions of AO's (
)
in the form [2]:
are the spherical harmonics,
the 
are chosen as a geometric
series with 
and
(Z - nuclear charge
of the atom), 
.
Such representation gives AO's with a quality comparable to numerical AO's.
The integrals can be calculated to a large extent analytically - see e.g. [3].
The representation of the 
in molecular and
solid state calculations allows small basis sets:
):
are used.
For close packed solids 
(
- Wigner-Seitz radius) and for molecules
and less dense packed solids
(
- covalent radius of the atom)
have been found as optimal choices. In principle 
can be determined by variation to minimize the final total energy.
This concept of confined AO's can be applied for 'first principle' solid state as well as molecular calculations, but also for semiempirical DF based schemes. Order-N type realizations have been realized also, as well as relativistic or semirelativistic formulations for heavy atoms.
Plane wave basis sets have a number of very useful features, one of these is their transferabilty across a range of bonding enviroments and another the fact that a single parameter, the cutoff energy, determines the convergence of the basis set in a straightforward and stable fashion. One of the great challenges for localised orbital basis sets is to achieve compact transferable basis sets which can also be scaled to give greater accuracy. Ideally we would like a suite of basis sets that take us from minimal basis sets, with reasonable accuracy through to fully converged basis sets. To do this we need to find some form of "optimal" basis set for a given accuracy which gives a similar accuracy when used in a wide range of bonding enviroments. We will present work [1] on producing optimal localised atomic type basis sets using the concept of spillage [2]. We have done this for C, Si and Al, which allows us to investigate these issues for systems that range from wide band gap semiconductor through conventional semiconductor through to nearly free electron metal. To test the transferability of the orbitals we looked at both a number of different structures within these materials and also at the description of SiC using the orbitals generated for bulk Si and C. We will also discuss the implementation of these localised orbital methods for transition metal systems.
The true test of any methodology is when it is employed to perform real calculations. We will look at how well the method performs when used to for calculations on two different systems, the problems that were encountered and the implications for the orbital optimisation. The two systems that we have investigated are the bonding of C60 to Si (100) [3] and the adhesion of Co on Si (100) [4]. Both these systems are large enough that the advantages of the localised basis sets over the plane wave basis sets both in cpu time and memory are considerable.
A requirement for linear-scaling calculations is the use of
localized basis sets. Bessel functions, gaussians [1-3],
and real space grids [4] have been
proposed, as well as strictly localized
numerical pseudo-atomic orbitals (NAOs) [5].
The latter are the solutions of the one-dimensional radial
Kohn-Sham pseudoatomic Hamiltonian under some confining
constraint that makes them zero or negligible beyond some radii.
These radii control the quality and the efficiency of the calculations.
It is therefore important to find the radial shapes of the
orbitals that optimize both the quality and the
convergence cutoff-radius. Different confinement options proposed in
the literature are analyzed, including confining potentials as
,
[5-7] or direct confinement of the atomic wavefunctions [8,9].
They are compared with a new proposal
that overcomes previous drawbacks [10].
Variational optimizations have been performed for
the different schemes in different systems.
The results are compared with fully converged plane-wave
calculations on a significant variety of
systems, including covalent, ionic and metallic.
Optimal double-
plus
polarization basis sets show very good agremeent with the plane-wave
structural parameters and cohesion energies.
The transferability of the obtained basis sets is tested in several
cases and it is found to be satisfactory.
A balanced way to define globally the confinement conditions for
all the orbitals in a system is also presented, based on applying a
fictitious pressure on the orbitals and then minimizing
the corresponding enthalpy.
The Amsterdam Density Functional (ADF) program [1] is based on Slater type orbitals and uses numerical integrations throughout. This has its consequences for the used algorithms, such as a density fit for the Coulomb potential, and subsequent evaluation of the Coulomb potential in the grid.
By straightforward distribution of the grid and atom pairs over processors an efficient parallelized version of ADF became available around 1997 [2]. A later project was concerned with removing the computational bottlenecks in the most time-consuming parts of the code: the calculation of the Kohn-Sham potential in the grid, the Fock matrix build, density fits, etc. This led to a near-linear scaling implementation [3], even though bottlenecks such as matrix diagonalization are still present. The description of the used techniques and applications to DNA can be found in C Fonseca Guerra's thesis [4].
A third step involved the combination of the near linear scaling code with the parallelization. The problem here is the disturbed load balance, due to the fact that the computational weight of a block of grid points depends on its position in the molecule. To counteract this, the points were first localized using a space-filling curve, and then statically distributed according to an estimate for their weight. Using time-dependent DFT calculations as an example, it is demonstrated that good scaling with the number of CPUs and system size is obtained both for the SCF and for property calculations [5]. The talk summarizes these previous developments.