Sicut omnes homines naturaliter scire desiderant veritatem, ita naturale desiderium inest hominibus fugiendi errores, et eos cum facultas adfuerit confutandi[].
Molecular modeling of transition metal complexes (TMC), reproducing characteristic features of their stereochemistry and electronic structure, is in high demand in relation with studies and development of various processes of complex formation with an accent on ion extraction, ion exchange, isotope separation, neutralization of nuclear waste, and also when studying structure and reactivity of metal-containing enzymes. Solving these technological problems requires modeling methods allowing massive simulations of potential energy surfaces (PES) of TMCs in a wide range of molecular geometries including (in the case of, say, complexation processes) internuclear separations corresponding to dissociation of coordination bonds between metal ions and ligands' donor atoms. The tools generally available for performing a required modeling range from fully empirical molecular mechanics (MM) to quantum mechanical (QM) or quantum chemical (QC) methods of different degree of refinement and sophistication. From a bird view perspective all mentioned approaches seem to be rather successful in a technical sense since normally it is possible to find a suitable and inexpensive method for modeling at least certain classes of molecules. However, when it goes about TMCs the methods of each family may turn to be inconsistent with the problem at hand. The sources of these seemingly unpredictable and nonsystematic failures of all mentioned modeling methods will be one of the main topics of the present paper. In the next Section we briefly review the existing methodologies as applied to TMCs and provide explanations of their limitations in this context. Further Sections present a formal point of view on constructing a QM (likely ab initio, semiempirical, and DFT) description of TMCs and apply it to analysis of the corresponding difficulties. The rest of the paper is devoted to description of the effective Hamiltonian of crystal field (EHCF) method and of its semiempirical implementation satisfying which can be used for TMCs' modeling and to its application to analysis of spin-active complexes of iron (II). Finally discussion and conclusions are given.
The elementary empirical tool for the molecular modeling of polyatomic systems is the method of molecular mechanics (MM) BurkertAllinger1982,Dashevskii-eng. It explicitly employs intuitively transparent features of molecular electronic structure like localization of chemical bonds and groups. The basic assumption of the MM is the possibility to directly parameterize molecular PES in the form of a sum of contributions (force fields) relevant to bonds, their interactions, and to interactions of non-bonded atoms:
| (0) |
In the literature Rappe1992,Rappe1993,Hay1993,Hay1998,Comba1995,Comba2003,Erras-Hanauer2003 various MM constructions are considered as effective methods for modeling PES of TMCs. It is noteworthy that in the case of metal complexes in general and of TMCs in particular the very basic characteristics of electronic structure comprising the basis of MM may be questioned. In fact when it goes about the metal ion in a complex it is not possible to single out transferable two-center bonds involving the metal. Also the number of bonds formed by a metal atom (the coordination number) may be variable and namely these variabilities may be the main topic requiring the modeling as indicated in the Introduction. Also the great variety of accessible coordination polyhedra makes it difficult to set the preferential valence angles. In review [] the extensive summary of results of calculations on coordination compounds of a wide variety of metals by the MM methods (as of 1993) is given. During the following decade, numerous subsequent works quoted in reviews [,,] were performed, in which PES of special classes of metal complexes in that or another manner is parameterized by some MM-like force fields. As it can be seen from the recent review [] the situation did not change too much since then. The conceptual problems mentioned above manifest themselves in extremely cumbersome and awkward appearence of the set of force fields in case of metal atoms as compared to traditional 'organic' force field systems. For example, it becomes necessary to introduce a double set of optimal valence angles for octahedral (or plane squared) complexes to assure these important molecular shapes are reproduced in the calculation as are the relative energies of the cis- and trans-isomers Hay1993,Comba1995. The number of other bonding parameters also rapidly grows, and it is difficult either to assign any clear physical sense to all these, or to restrict reasonable interval of parameter values and thus to separate probable ones from improbable.
An alternative to the valence force field approach based on the concept of preferable valence angles is to reload the responsibility for the description of the shapes of coordination polyhedra in TMCs to the non-bonding interactions. Such approach exists in the literature in two versions. The first is represented by the Kepert model, termed also as one of the 'points on a sphere' (POS) [,] in which the terms responsible for the metal-ligand bond stretching energy are taken into account by harmonic terms as previously. Everything that concerns the dependence of the bending energy of the valence angles at the metal atom, is replaced by the terms representing an effective interaction (van-der-Waals-like) either between the donor atoms, or between effective repulsion centers placed somewhere on the metal-donor atom bond Pletnev1997,Razumov2001. The second version is called the electrostatic model [,]. It completely neglects any specific bonding interactions in the nearest coordination sphere and substitutes the Coulomb interaction between effective charges residing on atoms of the complex for the overall interaction energy. The repulsive part of the metal-ligand van-der-Waals potential acts to prevent collapse of the system. The basic weakness of this approach is, certainly, lack of the reasonable method allowing independent estimates for effective charges.
Despite the considerable progress achieved in MM modeling of TMCs (for example, MM models of the complexes of cyclic polyamines with metals like Cu(II), Co(III), Ni(II) are reported Comba1999,Lindoy2000,Niketic2001,Zimmer2003), many questions including those of practical importance, remain unanswered. The first one is the problem of consistent modeling of metal complexes with variable number of the ligands. The need for such description arises in context of molecular-dynamic studies (see, for example, []) of metal ions soluted in complexation solvents containing chelating ligands (crown-ethers, cyclic polyamines etc.). In such systems one may expect formation of numerous complexes with different number of ligands or degree of coordination (the chelate number), which should be considered at one level of accuracy to keep uniform energy scale. Obviously, the harmonic approximation for stretching energy of metal-donor atom bond usually employed in MM, as, for example, in [], can not describe such effects. A direct replacement of the harmonic potential by another one, with more suitable asymptotic behaviour (for example, by the Morse potential), does not solve the problem, since it neglects many other factors, that apparently matter (different mutual influence effects just to give an example).
Another important point specific namely for TMCs is the presence of the partially filled d-shell on the metal ion which produces a whole set of electronic states of the complex of different total spin and spatial symmetry in a narrow energy range close to the ground state energy. Geometry dependence of these energies may be rather confusing which results in existence of the areas in the nuclear coordinate space where the PESs belonging to different electronic terms, closely approach each other and even intersect, leading to experimentally observed spin transitions Guetlich1981,Koenig1985,Koenig1991,TopCurrChemSpinTran or Jahn-Teller distortions []. Thus, the very problem of including the transition metals in the MM context implies certain contradiction: in the presence of several close in energy (or even crossing) electronic terms there is no object for the MM modeling in a strict sense, since there is no uniform (and single) PES of the complex. This specificity of the electronic structure of TMCs can be clearly observed in the results on blue copper proteins with aproximately trigonal-bipyrmidal coordination of the coper ion as reviewed in []. The Cu2+ cation is known to be a Jahn-Teller ion due to the spatial degeneracy of its respective 2Eg and 2T2g ground state terms in the octahedral and tetrahedral environments. The latter Jahn-Teller instability is inherited also by the trigonal bipyramidal environment where the ground state is 2E due to the electron count in the d-shell of the Cu2+ cation. Clearly the spatial degeneracy of the ground state is the limiting case of the closeness of electronic terms on the energy scale. This degeneracy is lifted when the molecular geometry deviates from the symmetrical arrangement and this is the content of the Jahn-Teller theorem (see for details []; an original and what is important - a concise proof is given in Pupyshev2005). Technically the Jahn-Teller instability manifests itself in the presence of multiple minima on the PES, having a close total energy. It must be understood, however, that these minima are a result of the sufficiently quantum behaviour of the d-shell of the Cu2+ cation which as it has been noticed previously in a certain sense prevent the usage of the classical MM picture. Indeed, as it is mentioned in DT1999eng,DRPT2002, the physical pre-condition of successful use of MM theories for common organic molecules is that their electronic excited states are well separated from the respective ground states on the energy scale. Only one quantum state of their electronic system is experimentally observed in `organics' at ambient conditions and the MM (a sort of classical) description becomes valid. By contrast, the behavior of the metal valence d-shell is sufficiently quantum: several electronic states may appear in a narrow energy range close to its ground state and this quantum feature requires a special care, not reducible to a simplistic adjustment of the form and parameters of no matter how sophisticated force fields.
A plausible way out of this situation has been proposed by R. Deeth ( Deeth2001 and references therein). In order to handle quantum beavior of the d-shell the ligand filed stabilisation energy (LFSE) term is added to the MM energy expression Eq. (1). The LFSE is written as a sum of the orbital energies of the d-orbitals in its turn calculated in the angular overlap approximation (see below) whose parameters are taken to be linearly dependent on the internuclear separation betwen the metal and donor atoms. Applying such a model solves many complications inherent to the MM of TMCs since the LFSE is a pure quantum contribution to the energy. For exmple the Jahn-Teller in Cu2+ compounds must be perfectly covered within such a setting. On the other hand the LFSE is by construction a sum of one-electron energy contributions whereas the energy of the d-shell is very much dependent of the two-electron contributions to the energy particularly when it goes about relative energies of the stets of different total spins and spatial symmetries. Bringing the latter into the MM context requires much more evolved and refined theory which will be explained below.
Turning in this context to a main topic of our interest, namely to modeling of the spin active TMCs we notice that the above considerations apply to them in a large extent. The change of the spin state of a complex is possible if at least two different electronic states (differing by the value of the total spin) have their respective minima at quite similar geometries of the complex at hand so that their respective total energies become equal at some intermediate geometry. As in the case of the Jahn-Teller Cu2+ cation in that of the spin-active ions (e.g. d6 Fe2+) the unique PES of the complex does not exist and at least two of them (the low-spin - LS - for S = 0 and the high-spin - HS - for S = 2) must be considered. Previously the MM force fields using different parameter sets for different spin states of the central atom were in use [], but due to no predictive force they are considered to be obsolete by now. However, the basic principles of their construction do not differ from those which explicitly use different parameter sets for say axial and equatorial ligands in the Cu2+ complexes [] since the latter are as well designed to imitate by means of a classical potential sufficiently quantum characteristics of the TMC's electronic structure. On this way one can expect pretty different sets of parameters say for four-coordinate complexes of the Ni2+ ion which must be tetrahedral in their triplet states and square planar in the singlet states. In this respect the recent paper [] is very remarkable. The authors try to construct the MM potential capable to describe transfromation between the square pyramidal and two trigonal bipyramidal forms of the pentacoordinate [Ni(acac)2py] complex (acac stands for acetylacetone, py - for pyridine ligands). To do so these authors propose to employ a specially designed force field dependent on L-Ni-L¢ angle posessing two minima at 90° and 120° separated by a barrier of the hight larger than 5 eV (500 kJ/mole). This clearly indicates some problems which can be clearly revealed by a simple analysys: The trigonal bipyramidal forms of the complex are obviously [] triplet (two d-levels degenerate in the trigonal field filled by two electrons) whereas the square pyramidal from may well be singlet. This spin switch has to take place somewhere along the rearragnement reaction coordinate but is by no means reflected in the MM picture. There remains a question whether the aproach emplying the LFSE is capable to describe such a low-symmetric and potentially correlation dpendent situation.
We see, that in the case of TMC any description of PES by the MM methods may be as well rather successful and rather poor. The borderline between potentially successful and unsuccessful cases looks out rather peculiar from the point of view of standard chemical nomenclature. Why Ni(II) is sometimes successful, and sometimes not, Cu(II), Fe(II), and Co(II) are very difficult, whereas Co(III) brings no special problem? A general conclusion is that a more detailed description of the electronic structure of TMCs than one implicitly put in the base of the entire MM picture is necessary. It must take into account all the important features of the former. Physically it is rather clearly formulated in terms of experimental accessibility of low-lying excited states, whose spectrum is responsible for the observed i.e. quantum, behavior of the TMCs. Quite naturally a quantum description is given by methods of quantum chemistry (QC). The latter further subdivides into ab initio, DFT, and semiempirical domains. Although the semiempirical methods are nowadays frequently treated as obsolete, taking into account the number of atoms in TMCs and thus incurred computational costs (see [] and below) which may become prohibitive despite considerable progress of the computational hardware they still deserve attention as a pragmatic tool for the massive PES simulations. On the other hand the DFT based methods which are almost unianimously considered to be the method of choice for TMCs Ziegler1991 will be shown to suffer of basically the same structural deficiency as do the semiempirical ones. We shall analyse briefly these extended classes of QC methods used for evaluation of PESs and of other properties of TMCs from a common point of view allowing to further consider specific difficulties pertinent to each of these classes of methods when applied to numerical modeling of TMCs. The main problem referenced in the literature in relation to QC description of TMCs is that of electron correlation []. It is a general belief that the correlations are possibly reproduced only by high-quality ab initio methods, but also can be adequately modeled by the DFT based methods. On the other hand the common opinion is that the semi-empirical methods are not suitable for modeling correlation effects at all. We shall show that two latter opinions are largely an exaggeration and that the potential of the DFT when applied to TMCs is considerably oversold whereas broadly understood semi-empirical methods by contrast still may be useful.
The relevant formal treatment starts from the notion that all the quantities related to electronic structure of molecules can be calculated with use of only one- and two-electron density matrices for the relevant electronic state of the system under study []. Taking the energy for the sake of definiteness we get:
| (0) |
| (0) |
The density matrices are by definition partial integrals of the corresponding trial wave functions YCQGS(xw | x1,x2,¼,xN) obtained for the given composition C and nuclear configuration Q so that they have the specified total spin S and spatial symmetry G:
| (0) |
The expressions Eqs. (2) ,(4) are completely general. To address the aspects important for the TMCs' modeling we notice that the statement that the motion of electrons is correlated can be given an exact sense only with use of the two-electron density matrix Eq. (4). Generally, it looks like [] (with subscripts and variables' notations xw omitted for brevity):
| (0) |
| (0) |
In the following Sections we analyse the previously listed classes of QC methods of electronic structure modeling in terms of the density matrices.
The modeling by ab initio QC methods bases on complete description of electronic structure for which it is necessary to consider a set of one-electron states (basis functions), number of electrons in the system and nuclear charges. All consequent modeling is the computer work which involves calculating the matrix components h(1), h(2) of the electronic Hamiltonian, for the set of selected basis functions (whose parameters are above denoted as w).
In ab initio methods the HFR approximation is used for build-up of initial estimate for r(1) and r(2) which have to be further improved by methods of configurational interaction in the complete active space (CAS) [], or by Mø ller-Plesset perturbation theory (MPn) of order n, or by the coupled clusters' Purvis1982,Chiles1981 methods. In fact, any reasonable result within the ab initio QC requires at least minimal involvement of electron correlation. All the technical tricks invented to go beyond the HFR calculation scheme in terms of different forms of the trial wave function or various perturbative procedures represent in fact attempts to estimate somehow the second term of Eq. (5) - the cumulant c of the two-particle density matrix.
In application of non-empirical methods to TMCs there exist specific difficulties caused by the correlation strength. This can be formulated as essential deviation of r(2) from the HFR approximate form Eq.( 6) which makes it necessary to take it into account at the initial stage of calculation. Meanwhile, the listed (systematic) methods of taking the correlation into account are based on the assumption that the correlations appear as a smaller corrections to the mainly HFR approximate wave function at least when it goes about the ground state. This is an assumption leading to the whole variety of the single-reference (SR) perturbative and coupled cluster methods, where by the single refence state to be improved is assumed to be a single Slater determinant. The actual physics of TMC's is sometimes much more complex. Even obtaining of the approximate solutions of the electronic problem within the HFR approximation although they are required only as starting points (reference states) for further improvements in case of TMCs may represent a serious problem. It is known that for TMCs the HFR methods in many cases yield the electronic structure breaking the Aufbauprinzip, according to which MOs are filled by electrons beginning from the lowest energy levels. However, any variational function of the HFR approximation giving the minimum of energy with respect to relevant variational parameters must satisfy this requirement. Another problem well known to practical workers in the field is the slow convergency of the HFR iterations or nonrare cases of being trapped into oscillatory regime. These problems are numerical manifestations of electron correlation. In this situation the HFR solution even if it is obtained may lay too far from the correct ground state of the TMC. The latter can not be derived from this approximation by those homeopathic medication which is provided by the perturbation or coupled cluster theories. The problem is that for example the formally 'excited' configuration may have the same energy as the 'ground state' one thus preventing proper treatment by either MP or CC methods. From the general point of view the situation is comepletely clear - one has to use configuration interaction (CI - multireference - MR) or CAS methods. Pragmatically, however, there remains the question: what amount and which configurations have to be included. In any case the poor initial approximation requires for curing a large number of configurations.
The ab initio calculations on TMCs date back to late sixties when the first examples of such calculations in the HFR approximation as applied to simple NiF64- ion had been published Gladney1969,Richardson1969,Moskowitz1970. In those early times the optimistic belief [] was that `It is mostly computational limitations which have in the past more or less prevented a wide application of the ab initio techniques to the chemistry of transition metal compounds ... with technical developments which may be forecast for the next few years, this type of calculations will probably become much more common'. It, however, happened that within ten years a collection of papers edited by one of the previous authors came out where the description of TMCs has been recognized as a 'challenge' Veillard1986. The above analysis shows the reason. Nevertheless, within two subsequent decades the hardware improved significantly so that the TMCs of modest size became available for direct more or less complete numerical ab inito study. Examples of such approach are numerous in the literature. Their range is extremely wide: from studies of structure and properties of 'helide' molecules HeM2+, where M is the doubly charged cation of the first transition row metal [] by various ab initio methods (including those with relativistic corrections). Of course, this is not the topic of our main interest. On the other hand the wide usage of the ab initio methods as a molecular modeling tool for TMCs is still prevented by enormous computational costs. In ab initio HFR MO LCAO methods used as zero-approximation calculations of correlation corrections required to make the result somehow acceptable are so complex that the dependence of time and other necessary computing resources on the size of the molecular system (N as number of AOs) scales up as N5¸N7. Therefore, at larger sizes of systems under study calculations of TMCs electronic structure become very expensive. This prompted an approach which hardly can be called methodologically sound, but which is widely represented in the literature: considering at a (currently acquired) ab initio level only a smaller part of the molecule of interest ignoring the rest of the system as it is demonstrated in the representative collection of reviews []. It is clear that the ab initio methods do not provide any tool for adjusting the presumably ``exact'' result obtained for a nonexisting model to the needs of analysis of an experimental situation. Nevertheless, when cautiously applied this model approach can be useful. In most studies on TMCs of say biological or industrial interest only a model small compound of it is actually considered. The following provides representative examples. In paper [] various models of the active site of metal-containing enzyme glutathiontransferase having structures with five and six coordinated ions Mn2+ and Fe2+ are explored. In both cases it is assumed, that the respective metal ion is in its HS state. It essentially facilitates calculation since such a ground state can be reasonably modeled by a single-configuration (HFR) wave function. Structural studies of the Jahn-Teller effect in TMCs by the ab initio CC methods are performed in [,]. An attractive testbed for testing various QC methods as applied to TMCs is provided by metal-porphyrins being rather interesting from various practical points of view but simulataneously polyatomic enough to raise the efficiency issues and also well studied experimentally. The early attempts to apply ab initio QC methods are reviewed in Ref. DedieuRohmerVeillard1982. With use of models of various extent of realism (including those with exact number of atoms and electrons) it was shown that HFR MO LCAO turns out to be good in the extremal cases of Co(II) and Mn(II) whose ground state spin were reproduced correctly to be 1/2 and 5/2, respectively, and fails in the practically most important case of Fe(II) prophirine which was known to be of intermediate total spin S = 1 in its ground state. Despite 20 years of development this result quite well established experimentally many times escaped from ab initio workers. Even the most recent results [,] do not allow to make a definitive conclusion on the capacity of the ab initio methods.
Unfortunately the ab initio workers are not cautious enough as is exemplified by Ref. []. In it a wide set of experimental data on catalytic reactions taking place in the presence of Pd complexes with substituted phosphine ligands PR3 is modeled by ab initio methods applied to models where the whole variety of the ligands is represented by the unsubstituted phosphine PH3 molecule. The main problem with such an approach is that the sensitivity of the processes under study to the number and nature of organic substituents at the phosphorous atom is well known in the literature.
The above mentioned computational costs lead to a necessity to find a fragile compromise between the requirements of precision and feasibility of a calculation. That of course raises interest to applying the hybrid QM/MM methods to TMCs. The more traditional version of this approach consists in taking rather large portion of the TMC (including the metal atom) to the QM subsystem and in treating the distant groups on a lower level of the theory. Approaches of that type are quite frequently applied to the TMCs and the recent review of it is given in []. The general problem of this approach is the nonsystematic character of the treatment of the intersubsystem border (junction) accepted in most standard packages. The detailed discussion of these problems is given in our recent review TT2003. More specific argument can be borrowed from the metal porphyrin problem as well. In the literature there are reliable experimental data concerning the influence of the periferal substitutents on the electronic structure of the central ion in the case of tetraphenylporhyrinates of iron (III) additionally substituted in the phenyl rings (see [] and references therein). This type of effect cannot be attributed to any steric hinrdance or whatever of that sort and must be ascribed to the influence of the substituent upon the electronic structure of the transition metal atom. In the standard setting employing the QM/MM technique is precluded for such a problem and a total QC calculation is required. The latter is, however, very difficult since although chemically tetraphenyl porphyrins do not seem to be very different from the nonsubstituted ones the total number of basis functions is approximately as twice as large for the substituted species thus increasing the required computational resources may by factors 25 to 27.
In the previous Section we briefly described the problems arising when the ab initio QC methods are applied to the modeling of TMCs. These problems may be considered largely as technical ones: if the computer power is sufficient the required solution of the many electron problem can be obtained by brute force even if the initial guess for the wave function is poor. Pragmatically, however, the resource requirements may become prohibitively high for using the ab initio QC techniques as a tool for massive PES modeling. In this situation the semi-empirical methods can again come into play as 40 years since the pioneer works Clack1972,Clack1974,Clack1974a where the CNDO and INDO parameterizations by Pople and Beveridge [] were extended to transition metal compounds. Now there is an extensive sector of semiempirical methods differing by expedients of parametrizations of the HFR approximation in the valence basis. In many of them the parametrization at least is formally extended to the transition metal atoms, for example, in methods ZINDO/1, SAM1, PM3(tm), PM3* etc. BaconZerner1979,Boehm1981a,Boehm1981b,Dewar1993-SAM1,Cundari1998,Cundari1999,Mohr2003, although, principles of parametrization may differ as stipulated by the need to reproduce different experimental characteristics. The attempts to construct an acceptable parameterization for TMCs are all undertaken within the framework of the HFR MO LCAO paradigm. It is easy to understand that the nature of failures which accompany this direction of research as long as it exists lays precisely in the inadequate treatment of the cumulant of the two-electron density matrix by the HFR MO LCAO.
The procedure of developing a semi-empirical parameterization can be generally formalized in terms of Eq. (2) as follows. A set of experimental energies E(CQGS) corresponding to different chemical compositions C, molecular geometries Q, and electronic states with specific values of S and G is given. In the case when a response to an external field is to be reproduced the latter can be included into the coordinate set Q. Developing a parameterization means to find certain (sub)set of parameters w which minimizes the norm of the deviation vector dEw with the components E(CQGS)-E(CQGS | w) numbered by the tuples CQGS:
| (0) |
| (0) |
| (0) |
|
| (0) |
| (0) |
We do not intend to further elaborate on characteristics of the semiempirical methods. It it enough to say, that all of them which are restricted to the HFR approximation suffer from the shortcoming described above and, hence, one has not to have too much hope to reach a consistent description of TMCs within the HFR framework. The recent semiempirical attempt to parameterize the TMCs in the PM3(tm) method [] is very instructive in this respect. The calculations carried out in Ref. Lindoy1996 show that the method is not capable to reproduce even very simple characteristics in a series of TMCs having similar structure, though other authors [,] state that in some cases reasonable estimates of geometrical characteristics may be received, nevertheless. This situation can be understood by thorough consideration of the sets of objects chosen for analysis in different works. In Lindoy1996 authors study a uniform set of about 30 Ni2+ complexes with the ligands bound by nitrogen donor atoms. The analysis of this series performed there clearly shows that the PM3(tm) method fails for these Ni2+ complexes for the now understandable reason. However, in Cundari1999,Basiuk2001 the authors try to explore a comparable number of complexes but much more dispersed over the range of classes, which includes compounds of the first and second transition row atom, HS and LS ones, those having ``ionic'' and ``covalent'' bonds etc. For that reason in the test sets [,] the problematic classes are represented by a couple of examples each, which look out as completely isolated exceptions. This can serve as an example of how trying to test the method on a wide and apparently ''random'' selection of objects may lead to a smeared picture due to absence of clear criteria of any adequate classification of the chosen set. That said above does not mean that a semiempirical parameterization based on the HFR MO LCAO scheme and valid for a certain narrow class of compounds or even for a specific purpose cannot be built. It is done for example in [] for iron(II) porphyrins. But in a more general case there is no way to arrive to any definite conclusion [] about the validity of a semi-empirical parameterization in the HFR context. On the other hand we have to mention that the semiempirical method ZINDO/1 [] which allows for some true correlation by taking into account the configuration interaction may be considered as a prospective setting for further parameterization, provided the HFR solution required by this method as a zero approximation can be obtained. This will be discussed in a more detail below.
Methods of density functional theory (DFT) originate from the Xa
method originally proposed by Slater [] on the base of
statistical description of atomic electron structure within the Thomas-Fermi
theory []. From the point of view of Eq. (
3), fundamental idea of the DFT based methods consist first
of all in approximate treatment of the electron-electron interaction energy
which is represented as:
| ||||||||||
| (0) |
| ||||||||||||||||||||||||||||
Further pragmatic moves are described in details in numerous books and reviews of which we cite the most consize and recent Ref. []. Two furhter hypotheses are an important complement to the above cited theorems. One is the locality hypothesis another is the Kohn-Sham orbital trick. The locality has been seriously questionned by Nesbet in recent papers [,], however, it remains the only practically inplemented solution for the DFT. The simplest one is the local density approximation (LDA):
| (0) |
During last decades the DFT based methods have received a wide circulation in calculations on TMCs' electronic structure Ziegler1991,Noodleman1997,Chermette1998,Anthon2002,Siegbahn2003. It is, first of all, due to widespread use of extended basis sets, allowing to improve the quality of the calculated electronic density, and, second, due to development of successful (so called - hybrid) parameterizations for the exchange-correlation functionals (vide infra for discussion). It is generally believed, that the DFT-based methods give in case of TMCs more reliable results, than the HFR non-empirical methods and that their accuracy is comparable to that which can be achieved after taking into account perturbation theory corrections to the HFR at the MP2 or some limited CI level [,,].
As in all axiomatic theories relying upon existence theorems partcular
attention has to be paid to the consequences of these theorems which
sometimes can be rather peculiar. Although the general validity of the
Hohenberg-Kohn theorems cannot be questioned the example Eq. (
8) obviously presents two different wave functions
with different energies yielding the same one-electron density
matrix and thus of course the density itself. For that reason the
qualitative effects of electron correlation which are crucially important
for correct TMC modeling and for which the term c in Eq. (
5) takes care, can not be reproduced by the DFT based methods at
all since these do not contain the necessary elements of the theory for it
(although MP2 and even limited CI do). Whatever attempt to do that is going
to have a restricted character due to restricted treatment of the cumulant.
In that respect the situation is analogous to that in the HFR-based
semi-empirical methods. In terms of the ``data-fit'' model Eq. (
7) the DFT methods can be understood as ones with the fitting
model of the form:
|
In a more general setting the recipe [] can be considered as
an implementation of another suggestion by Gunnarsson and Lundqvist
Gunnarsson1976 and von Barth [] known also at a pretty early
stage of the development of the DFT technique of employing different
functionals to describe different spin or symmetry states. In other words
the simplified model for the data fit Eq.(11) changes to:
|
In the remarkable series of papers [,,]
the authors attempted to reproduce the relative energies of the LS and HS
terms in a series of pseudo-octahedral Fe(II) complexes with ligands bound
to the metal atom through sulfur and nitrogen donor atoms. It was achieved
only by the cost of adjustment of the weight of the Hartree-Fock exchange
energy in the hybrid B3LYP functional leading to development of the B3LYP* functional (see below). These authors conclude, that the hybrid
functionals (B3LYP) in the DFT based methods favor mainly the HS states. As
one can see namely the Fock exact exchange is responsible for the Hund's
rule conformity. Indeed the HFR estimates of the LS-HS splitting given in
all the cited papers amount several eV with the wrong sign (HS state below
the LS state). The reason was transparent yet in the year 1993: unbalanced
account of correlations and exchange in the DFT schemes. The Hartree-Fock
exchange strongly stabilizes the HS state of the open d-shell even if the
single-determinant wave function is used, whereas the correlations which can
potentially stabilize multi-determinant LS states are absent. The LSDA
estimate of exchange gap taken together with some contribution of
correlation by contrast strongly underestimates the latter and by this
favors the LS states to become lower on the energy scale. This discrepancy
has been tried to remove by developing hybrid functionals which reduces to
elimination of the disbalance between different contributions to the
exchange-correlation terms by taking different estimates of the latter and
ascribing them different weights fit to reproduce some data (ourdays taken
largely from a numerical experiment performed on an ab initio level,
but for sure whatever experimental data could be used). For the specific
example of the B3LYP functional - the most popular one:
|
In order to reach an agreement with similar data on iron(II) complexes with ligands containing nitrogen donor atoms [] the value of c3 = 0.12 had had to be introduced. Similar measures are necessary even for very simple metals cations like Cu2+. In [] it was found that the most commonly used DFT functionals give a too covalent ground state of D4h [CuCl4]2-. A novel hybrid functional with 38% HF exchange (c3 = 0.38) can give the good agreement between the calculated and experimental ligand field and ligand-to-metal charge transfer excited state energies. Incidentally, the EHCF theory (see below) gives as much as accurate results on the d-d transitions as all the standard DFT procedures whereas the much better agreement in [] is achieved by the cost of unusually high weight of the HF exchange, which basically has no other substantiation. The entire situation thus looks out to be rather comic: after 40 years of claims of existence of the unique density functional and after 30 years of similar claims of an extraordinary power of the DFT-based theories in the realm of TMCs it turns out that namely for TMCs no single functional could be found so far. Whatever implementation reasonable from a practical point of view requires specific nonuniversal functionals dependent on spin, symmetry and chemical nature not only of the metal, but also of the donor atoms in the ligands.
We see from the above discussion that staying within the DFT it is not possible to describe the multiplet structure of the d-shell (incidentally all the success stories reported so far are limited to the p-shells Nagy1998,Nagy1998a whereas in the d-shells only average energies of several multiplet states can be reproduced []). In this context it seems to be necessary to analyse the attempts to achieve it which are available in the literature (see [,]). These attempts, however, have a nature absolutely different from the DFT itself so they will be described in an appropriate place below.
Currently the time dependent DFT methods are becoming popular among the workers in the area of molecular modeling of TMCs. A comprehensive review of this area is recently given by renown workers in this field []. From this review one can clearly see [] that the equations used for the density evolution in time are formally equivalent to those known in the time dependent Hartree-Fock (TDHF) theory Dirac1929,Frenkel1934,Thouless1961 or in its equivalent - the random phase approximation (RPA) both well known for more than three quarters of a century (more recent references can be found in Rowe1968,McWeeny1992,BlaizotRipka). This allows to use the analysis performed for one of these equivalent theories to understand the features of others.
According to analysis of Ref. [] the excitation energies evaluated by TDDFT correspond to taking into account interactions between configurations obtained from the original single-determinant ground state by single electron excitations (CIS). This is obviously equivalent to the so called Tamm-Dankoff approximation in the energy domain []. For the latter it is known that in sertain situations when the HOMO-LUMO gap becomes small as compared to the Coulomb interaction matrix elements or in other words in a vicinity of the stability loss by the corresponding HFR solution the excitation energies thus obtained may become negative thus indicating some serious problems. The reason is quite transparent: the electron correlation (interaction of the configurations) is taken into account in an unbalanced manner; it is accounted in the singly excited manyfold but is completely neglected for the ground state. If the bare gap (orbital energy difference) is not too small (othervise the problem becomes evident) the unbalanced correlation account manifests itself in that the excitation energies estimated in the Tamm-Dankoff approximation are somewhat lower than necessary. That is precisely which is reported in [] on example of Pd complexes (and for other systems in Ref. []): the TDDFT excitation energies are systematically lower than the experimental ones. In this context it becomes clear that the TDDFT may be quite useful for obtaining the excitation energies in those cases when the ground state is well separated from the lower excited states and can be reasonably represented by a single determinant wave function may be for somehow renormalized quasiparticles interacting according to some effective law, but shall definitely fail when such a (basically the Fermi-liquid) picture is not valid. In a way this is what the other prominent authors in the field of the TDDFT recognize as inherent drawbacks of this approach Furche2004,Burke2005. According to these authors nothing can be done if the one-electron states used have wrong one-electron energies or if the ground state is not a single determinant built of the Kohn-Sham orbitals. But these are the situations we must be ready for when addressing the TMCs.
This brief analysis allows to conclude that the fact that ``the superiority of the TDDFT method ¼ has not been unequivocally established ... in particular for d® d transitions'' [] is not an unfortunate accident but a logical consequence of deeply rooted deficiencies inherent to the underlying single-determinant nature of the TDDFT method and the announced proof of superiority will hardly whenever take place.
The above review of the methods of molecular modeling (both QC - including DFT - and MM), given above, has shown that none of them is completely suitable for molecular modeling of TMCs. The MM methods do not allow to consider multiple PES corresponding to several energetically close electronic states of TMCs. Ab initio QC methods appear to be too demanding to computational resources when employed to model chemically interesting TMCs with bulky organic ligands; HFR-based semiempirical methods and even the DFT-based methods suffer from the same deficiencies as MM methods, since within their respective frameworks it is not possible to reproduce relative energies of electronic states of different spin multiplicity without serious ad hoc assumptions.
We shall note, that the difficulties arise precisely when modeling is to be applied to molecules involving transition metal atoms mainly of the second half of the first transition row. Moreover, even among the TMCs formed by these atoms the problems are not uniformly distributed: the normal chemical nomenclature does not provide here an adequate classification. When it goes about metal carbonyls or about metals of the second or even third transition row, the DFT methods seem to be able to do the job quite decently. However, turning to compounds of the first transition row metals with open d-shells raises many problems. Intuitively distinction in behaviour of two types of the metal compounds is clear to any chemist. In a row of isoelectronic species Ni(CO)4, Co(CO)4-, Ni(CN)42-, Fe(CO)42- they readily recognize the ``not a family member'', but probably fail to give a reason. In terms of the traditional theory of chemical bonding the classification is rather vaguely formulated in terms of covalent, polar covalent, ionic, metallic, coordination, donor-acceptor and other types of chemical bonds. Clearly enough such a classification is not relevant in the case of interest. Remarkable alternative systematic of types of chemical bonds is given in []. We reproduce it in Table partially abridging.
| Bond type | Electronic structure | Compounds example | Typical properties |
| Valence | MOs are localized between | CH4 | Distinct character of |
| pairs of atoms and occupied | NH4+ | bond energy, dipole moments, | |
| by two paired electrons | Diamond | frequencies, polarizabilities etc. | |
| C2H4 | |||
| Orbital | MOs are delocalized in one | Benzene | There are no distinct |
| or two dimensions | Graphite | characteristics; conductivity, | |
| cycles aromaticity | |||
| Coordination | MOs are delocalized | CuCl42- | There are no distinct |
| in space - three-dimensional | CoCl2 (crystal) | characteristics; variable | |
| coordination number and magnetic moment, | |||
| strong mutual influence of ligands |
Although this classification also is not particularly satisfying it can be used as a starting point for further discussion. The classification of Table 1 relies upon the MO LCAO, i.e. ultimately the HFR, picture of molecular electronic structure. As it is discussed above, the HFR is not very much reliable when it goes about TMCs. Nevertheless we can observe that the lack of regularity both in bond lengths, and in oscillatory frequencies of bonds in complexes is associated, according to Bersuker1986eng, with a three-dimensional delocalization of one-electronic states involved. As an example though, the metal-ligand bonds in TMCs are given and the optical spectra and magnetic moment distinguishing them from all other compounds are given as specific characteristics. However, the non-characteristicity of the bond lengths and valence angles leading to flexibility of shapes of coordination polyhedra and the coordination numbers themselves are equally common for the complexes of nontransition metals (for example, alkaline or alkali earths). This shows a necessity to turn to somewhat more formal description of molecular electronic structure than it can be provided by the traditional theory of chemical bonding, but still more qualitative than it appears from the numerical experiments arranged in the framework of no-matter-how-precise QC methods.
The qualitative description of the electronic structure can be given in terms of even older concept of ``chromophores''. According to the IUPAC definition the chromophore is an atom or group of atoms in the molecule that gives color to the molecule. This definition unites two aspects - one related to the system's response to an external perturbation: the spectrum. By this the concept of chromophore is related to experimental behavior of molecular systems. Another aspect relates to the structure understood as a localization of the excited states controlling the tentative response to that perturbation. The examples of chromophores are well known from the textbooks. On the part of the modern theory of electronic structure the concept of chromophore formalizes in the McWeeny's theory of electron groups. (The analogy between chromophore concept and McWeeny's theory for the special case of TMCs has been early noticed also in a remarkable work []). Within this theory the zero approximation to the system's electronic wave function is taken as an antisymmetrized product of rather local group multipliers referring to relatively isolated elements of molecular electronic structure. These elements - electron groups - are physically identified as two-electron two-center bonds, conjugated p-systems etc. Of course, these groups are not totally isolated and ascribing excitations to only one of them is an idealization. Nevertheless, the effective Hamiltonian technique is available to reduce manifestations of the intergroup interactions to renormalizations of the effective group Hamiltonians which allows to interpret the response of the system to any external perturbation in terms of excitations localized in the groups.
Further analysis is based on the idea that the characteristic experimental behavior of different classes of compounds and the suitability of those or other models used to describe this behavior is ultimately related to the extent to which the chromophores or electron groups physically present in the molecular system are reflected in these models. It is easy to notice, that the MM methods work well in case of molecules with local bonds designated in Table 1 as valence bonds; the QC methods apply both to the valence bonded systems, and for the systems with delocalized bonds (referred as ``orbital bonds'' in Table 1). The TMCs of interest, however, not covered either by MM or by standard QC techniques can be physically characterized as those bearing the d-shell chromophore. The magnetic and optical properties characteristic for TMCs are related to d- or f-states of metal ions. The basic features in the electronic structure of TMCs of interest, distinguishing these compounds from others are the following:
These properties of the d-shell chromophore (group) prove the necessity of the localized description of d-electrons of transition metal atom in TMCs with explicit account for effects of electron correlations in it. Incidentally, during the time of QC development (more than three quarters of century) there was a period when two directions based on two different approximate descriptions of electronic structure of molecular systems coexisted. This reproduced division of chemistry itself to organic and inorganic and took into account specificity of the molecules related to these classical fields. The organic QC was then limited by the Hückel method, the elementary version of the HFR MO LCAO method. The description of inorganic compounds - mainly TMCs,- within the QC of that time was based on the crystal field theory (CFT) [,]. The latter allowed qualitatively correct description of electronic structure, magnetism and optical absorption spectra of TMCs by explicitly addressing the d-shell chromophore. Let us consider the CFT in more detail.
Basics of the CFT were introduced in the classical work by Bethe [] devoted to the description of splitting of atomic terms in crystal environments of various symmetry. The splitting pattern itself is established by considering the change of symmetry properties of atomic wave functions while spatial symmetry goes down from the spherical (in the case of a free atom) to that of a point group of the crystal environment. The energies of the d-d-excitations in this model are obtained by diagonalizing the matrix of the Hamiltonian constructed in the basis of nd-electronic wave functions (nd is the number of d-electrons). Matrix elements of the Hamiltonian are expressed through the parameters describing the crystal field and those of the Coulomb repulsion of d-electrons, that is Slater-Condon parameters Fk, k = 0,2,4, or the Racah parameters A, B, and C. In the simplest version of the CFT these quantities are considered as empirical parameters and determined by fitting the calculated excitation energies to the experimental ones. This approach devoids any predictive force (except for the splitting pattern itself) due to presence of empirical parameters in the theory, which are specific for each compound. The CFT gives a description to the characteristic properties of TMCs at the phenomenological level. All important features of their electronic structure are fixed by this theory and the perpetual problem remains obtaining consistent estimates of its parameters (strength of the crystal field). All further development of the CFT was concentrated on attempts to obtain independent estimate of its parameters []. Within the standard CFT this problem, however, has no solution due to oversimplified picture of the transition metal ion environment (ligands). Indeed the CFT theory uses the ionic model of the environment and calculates the splitting of the initial term of the free metal ion as if it were a pure electrostatic effect. The symmetry of the environment is correctly reproduced even in this simplistic model, whereas all the chemical specifics of this environment gets lost. For this reason it is not surprising that the heaviest strike upon the CFT from the (semi)quantitative side was given by TMC spectroscopy yet in 30-ties. Spectroscopic experiments allowed to range the strengths of the crystal fields exerted by different ligands to the so called spectrochemical series [,,]:
| ||||||||||||||
The problem of estimating crystal field parameters can be solved by considering the CFT/LFT as a special case of the effective Hamiltonian theory for one group of electrons of the whole N-electronic system in the presence of other groups of electrons. The standard CFT ignores all electrons outside the d-shell and takes into account only the symmetry of the external field and the electron-electron interaction inside the d-shell. The sequential deduction of the effective Hamiltonian for the d-shell, carried out in the work [] is based on representation of the wave function of TMC as antisymmetrized product of group functions of d-electrons and other (valence) electrons of a complex. This allows to express the CFT's (LFT's or AOM's) parameters through characteristics of electronic structure of the environment of the metal ion. Further we shall characterize the effective Hamiltonian of crystal field (EHCF) method and its numerical results.
The TMCs' electronic wave function formalizing the CFT ionic model is one with a fixed number of electrons in the d-shell. In the EHCF method it is used as a zero approximation. The interactions responsible for electron transfers between the d-shell and the ligands are treated as perturbations. Following the standards semiempirical setting we restrict the AO basis for all atoms of the TMC by the valence orbitals. All the AOs of the TMC are then separated into two subsets from which one (the d-system) contains 3d-orbitals of the transition metal atom, and another (the ``ligand subsystem'', or the l-system) contains the 4s- and 4p-orbitals of the transition metal atom and the valence AOs of all ligand atoms. We shall try to cover within the present theory only the complexes, where excitation energies in the l-system are by far larger than the excitation energies in the d-shell of the metal atom. This singles out a subset: the Werner-type TMCs which from the point of view of chemical nomenclature can be characterized as ones with the closed electronic shell ligands, such as F-, Cl-, Br-, I-, saturated organic molecules with donor atoms etc.
Formaly the theory evolves in a following way. The low-energy d-d-spectrum of the TMC can be obtained if the Hamiltonian is rewritten in the form:
| (0) |
| (0) |
| (0) |
| (0) |
| (0) |
| (0) |
| (0) |
Within the EHCF method [] the single Slater determinant Fl has to be obtained from semiempirical HFR procedure. Solving the HFR
problem for the l-system yields the one-electron density matrix Pab, orbital energies ei, and the MO-LCAO coefficients cia. These quantities completely define the electronic structure
of the l-system and are used to calculate the effective Hamiltonian Eq. (18) by Eqs. (20)-(22), where QL = åa Î LPaa-ZL is the effective charge
of the ligand atom L; ZL is the core charge of the ligand atom L; VmnL is the matrix element of the potential energy operator
describing the interaction between a d-electron and a unit charge placed
on the ligand atom L; ni is the occupation number of the i-th l-MO (ni = 0 or 1); DEdi (DEid) is the energy
necessary to transfer an electron from the d-shell (from the i-th l-MO) to the i-th l-MO (to the d-shell):
| (0) |
|
This comprises the EHCF picture for the electronic structure of the Werner-type TMCs.
In the context of the EHCF construct described in the previous Section the problem of semiempirical modeling of TMCs' electronic structure is seen in a perspective somewhat different from that of the standard HFR MO LCAO-based setting. The EHCF provides a framework which implicitly contains the crucial element of the theory: the block of the two-electron density matrix cumulant related to the d-shell. Instead of hardly systematizeable attempts to extend a parameterization to the transition metals it is now possible to check in a systematic way the value of different parameterization schemes already developed in the ``organic'' context for the purpose of estimating the quantities necessary to calculate the crystal field according to prescriptions Eqs. (20) - (22) of the EHCF theory. Solving the eigenvalue problem with the effective Hamiltonian for the d-subsystem (Hdeff) with the matrix elements which are estimated with use of any ``organic'' semiempirical scheme with the CI wave function constructed in the basis of the d-system, one obtains the complete description of the many-electron states of the d-shell of the metal ion in the complex. In such a formulation the EHCF method was parameterized for calculations of various complexes of metals of the first transition row, with mono- and polyatomic ligands. In papers STM1992,STM1996a,STM1996b,SDTM2002 the parameters for the compounds with donor atoms C, N, O, F, Cl and for doubly and triply charged ions of V, Cr, Mn, Fe, Co, Ni and Cu are fitted. These parameters do not depend on details of chemical structure of the ligands, rather they are characteristic for each pair metal-donor atom. The dependence of the excerted effective field on details of geometry and chemical composition of the ligands is to be reproduced in a frame of a standard HFR-based semiempirical procedure. The further evaluations [,] have shown applicability of the fitted system of parameters for calculations of electronic structure and spectra of numerous complexes of divalent cations with use of merely the CNDO parameterization for the l-system. In [,] the EHCF method is also extended for calculations of the ligands by the INDO and MINDO/3 parameterizations. In all calculations the experimental multiplicity (spin) and spatial symmetry of the corresponding ground states was reproduced correctly. The summit of this aproach had been reached in Ref. [] by calculations on the cis-[Fe(NCS)2(bipy)2]2+ complex. Its molecular geometry is known both for the HS and LS isomers of the said compound. The calculation for the both reproduces the respective ground state spins and the spectra of low lying d-d-excitations in a remarkable agreement with experimental data. Another good example is the treatment of metal porphyrins with use of the EHCF method. As already said above, for the decades the ab initio methods fail to reproduce the experimental ground state of Fe(II) porphyrine. It is really a complex case since it is an intermediate spin (S = 1- i.e. neither HS nor LS) and spatially degenerate state (3E). However, applying even very sophisticated methods (including CASPT2 which is considered to be a method of choice for TMCs in the ab initio area) has not yet led to the desired success. According to [] the HS forms are ground states and the hope to get a correct result is rather meagre since the gap amounts up to 1 eV in favor of the HS state (although the interpretation given in Ref. [] is quite different). Meanwhile the EHCF method in its simplest setting (CNDO type of parameterization employed for the l-system) yields the experimental ground state 3E without any further adjustment of parameters.
The success of the EHCF method in reproducing the crystal field from
geometry data and ligand electronic structure as described by semiempirical
QC procedure poses a question on possible relation between the EHCF method
and the successful parameterization scheme for the LFT, the already
mentioned AOM. As it is shown below, a local version of the EHCF method
EHCF(L) derived and tested in our papers [,] represents
an effective tool allowing to estimate the AOM parameters with a good
precision. The derivation consists of two unitary transformations. The first
one is from the basis of canonical MOs (CMOs) of the l-system used in Eq. (22) to the basis of localized one-electron states representing
characteristic features of the ligand electroic structure - like presence
of lone pairs on the donor atoms. These are obtained by the maxY4
localization procedure []. This leads to the approximate
formula for the covalent contribution Eq. (22) to the effective
crystal field:
| (0) |
| (0) |
The second transformation is that of the d-orbitals from the global (laboratory) coordinate frame (GCF) to the diatomic coordinate frame (DCF) related to the ligand L , defined so that its z-axis is the straight line connecting the metal atom with the ligand donor atom, so that the resonance integrals bmL in Eq. (24) can be expressed through the tL vector of the resonance integrals between the metal d-AO's and the L-th LMO in the DCF:
| (0) |
| (0) |
| (0) |
The EHCF methodology allowed to perform systematic calculations of the crystal field for various ligand environments. The results of these calculations are in fair agreement with the experimental data, particularly with respect to the spin multiplicity of the ground states of the complexes. In the respective simple versions the EHCF/X methods treat the electronic structure of the ligands within a semiempirical approximation X. These methods are not, however, well suitable to conduct the systematic studies on PESs of TMCs. Further application of the EHCF methodology would be to develop a method for the calculation of PESs of TMCs. To do so we notice that the CNDO or INDO parameterizations for the ligands are probably enough accurate when it goes about the charge distribution in the ligands and the orbital energies at fixed experimental geometries, although, they do not suit for geometry optimizations (or more generally for searching PESs) of TMCs. Nevertheless, the EHCF method can be adapted for the PES search in a more general framework of the hybrid QM/MM methodology (standard reference here is []; for recent review see []). This finally allows to ``incorporate'' quantum and correlated behavior of TMC into the ``classical'' methodology of MM and to provide necessary flexibility for quantum/classical interface (see below).
This is done as follows. According to [] the total
electronic energy of the n-th state of a system with the wave function Eq.
(17) is
| (0) |
| (0) |
Appropriate test objects for this approach is provided by the spin isomers of TMCs already addressed in the context of the attempts to apply the DFT-based methods to them. As during spin transition (ST) variation of the Fe-N bond lengths makes up more than 10% of the bond length itself in the LS complex, the harmonic approximation does not suffice for the MM part of the energy. Thus, for Fe-N bond stretching potentials the Morse potential was used. By Eq. (30) terms of the singlet, triplet and quintet lowest states of the considered complex are constructed. Parameters of MM-potentials at which metal atom influence (angle bending and bond stretching) are fitted so that positions of minima on terms of singlet and quintet have as much as possible coincided with experimental distances Fe-N in HS- and LS-structures. The calculation has been carried out in our work []. The general scheme Eq. (30) of energy evaluation using the EHCF method for Edeff(n) requires an HFR semiempirical calculation of the l-system for each geometry of a complex. To clear this, the local version of the EHCF method which allows to calculate the crystal field at each geometry without repeating HFR calculations can be employed.
According to [,,,,] the covalent term Eq. (22) gives the main contribution (up to 90%) to the splitting of d-electron levels. Remaining 10-20% of the splitting comes from the Coulomb interaction with effective charges residing on the ligand atoms. The problem is how to calculate the covalent contribution to the splitting without recalculation the one-electron states of the l-system at each geometry. In Section 0.1 we reviewed the EHCF(L) theory which allows to estimate the crystal field in terms of local electronic structure parameters (ESP) of the ligands. By this method it can be done for arbitrary geometry of the complex, which is prerequisite for developing a hybrid QM/MM method.
The proposed approach in certain respects is resemblant to the general QM/MM techniques which are invented with the general purpose to treat different parts of polyatomic systems at different levels of theory. The general setting of this theory is discussed in detail in []. The main difference between the setting of the standard QM/MM technique and the present one is that the majority of authors working in the area of QM/MM see as a desirable feature a possibility to extend the subsystem to be treated on a quantum level as much as possible. This is seen as a medicine against the uncontrolable errors introduced by uncautious cutting the entire electronic system in parts treated by the QM and MM techniques respectively. The hybrid EHCF/MM technique uses somewhat opposite approach: it tries not to extend but to reduce the QM subsystem as much as possible and to treat the intersubsystem frontier in such a way that the interactions between the quantically and classically treated parts are sequentially taken into account. Since physically the true quantum effects - the low-energy excited states in TMCs, - are located in the d-shell we restrict the true quantum description to these latter. This is related to the very understanding of the notion ``quantum'' relevant to the present problem which we have already mentioned at the beginning: in organic chemistry one normally deals with the ground state only which on the energy scale is well separated from the lowest excited state. This is the physical reason why the classical (MM) description is possible for organics. The TMCs differ from that picture obviously due to low-energy excitations in the d-shell accessible in experiment, thus it must be treated on a quantum level.
The technical problem was to develop an adequate form of the intersubsystem junction precisely for the case when the quantum subsystem is represented by the d-shell. The source of the problem here is that as it is shown by clear advantage of the LFT taking into account the ligands' electronic structure over the CFT the former must be somehow economically reproduced in the otherwise MM calculation.
Since in the EHCF(L) the effective crystal field is given in terms of the l-system Green's function, the natural way to go further with this technique
is to apply the perturbation theory to obtain estimates of the l-system
Green's function entering Eqs. (24) and/or (27). It was
assumed and reasoned in [] that the bare Green's function for
the l-system has a block-diagonal form:
| (0) |
| (0) |
The Coulomb interaction between the ligands themselves and between each of them and the metal ion when turned on does not break the block diagonal structure of the bare Green's function G00l. Then the approximate Green's function for the l-system conserves the form eq. ( 31) but with the poles now corresponding to the orbital energies of the ligand molecules in the Coulomb field induced by the central ion and by other ligands (L ¢ ¹ L ) rather than to those of the free ligands.
The simplest picture of the effect of the central ion on the surrounding
ligands reduces to that of the Coulomb field affecting the positions of the
poles of the Green's function (orbital energies) of the free ligand. The
form of the CMO's of each ligand is left unchanged which corresponds to the
rigid ligands' MO's (RLMO) picture Ref. []. According to
Abrikosov1965, the effect of the Coulomb field upon the orbital energies
is represented by:
| (0) |
| |||||||||||||
| (0) |
| (0) |
This model can be improved by taking into account polarization effects in the ligand sphere. For this end the metal ion is considered as a point charge equal to its oxidation degree or formal charge, which is the sparkle model deAndrade1994.
Within models of the sparkle family the effect of the external Coulomb field does not reduce to the renormalization of the orbital energies as it is within the RLMO model (see above). By contrast, the electron distribution also changes when the ligand molecules are put into the field. We model this by classical polarizability. Accordingly the difference between effective charge on atom A in the complex (polarized) and that in the free ligand (non-polarized) is:
| (0) |
| (0) |
Though procedures of that sort are admitted in modern MM schemes directed to the systems with significant charge redistribution [] we consider such a procedure to be too resource consuming and restrict ourselves by several lower orders with respect to P in the expansion. Then the term Pdh0 corresponds to the first order perturbation by the Coulomb field induced by the metal ion and bare (non-polarized) ligand charges. The second order term corresponds to the perturbation due to the Coulomb field induced by the mutually polarized (upto the first order) charges:
| (0) |
The details on calculating mutual polarizabilities relevant to the EHCF/MM context can be found in Ref. []. The charges thus obtained are used for calculation of the Sii(f) term according to Eq. ( 34). This model can be termed as PS model (PS stands for Perturbative Sparkle). Specifically, PSn approximation level of the PS model stands for the charge corrections employing the series Eq. (38) up to the n-th order, while PS itself stands for the exact expression with the inverse matrix in the second row of the same equation. Then, Eqs. (37)-(39) comprise the perturbative form of the Sparkle model of the l-system's electronic structure (the PS model).
The proposed procedure improves the junction between the EHCF(L) method playing role of the QM procedure and the MM part, as shown below, where details of the calculations performed within this approximation are given.
The local version of EHCF method was implemented and used for the analysis
of the molecular geometries of complexes of iron (II) in works
DRPT2002,DPT2003,DT2004. The satisfactory agreement in the description of
complexes geometry with different total spins is achieved when the effect of
electrostatic field of the metal ion on the ligands is taken into account
through the electrostatic polarization of the ligands. Satisfactory
estimates of parameters of the crystal field for series of complexes of iron
(II) and cobalt (II) (both LS and HS ground states) are achieved. Totally 35
six-coordinated iron complexes with mono- and polydentate ligands,
containing both aliphatic, and aromatic donor nitrogen atoms (mixed
complexes with different types of donor nitrogen atoms and different spin
isomers of one complex are included to this number) and ten cobalt complexes
also with different types of donor nitrogen atoms and coordination numbers
ranging from four to six have been considered. Deviations of calculated bond
lengths Fe-N and Co-N from the experimental values are well enough
described by the normal distribution. Parameters of that distributions were
the following: the mean value (average deviation over the dataset, m = -0.037
Å and s=0.054 Å in the case of Fe(II) complexes, and
m=0.017 Å and s=0.044 Å in the case of Co(II) complexes.
The above values are quite acceptable for the
entire set of data but it turned out that they mask an inherent bias of
the proposed approach. In the iron(II) complexes the Fe-N bond lengths
of the HS complexes are systematically
underestimated whereas those in the LS come out
slightly overestimated. In fact the parameters of the fit of the empirical
distrbution function of deviations restricted to the LS complexes are
m=0.011Å and s=0.034Å and those restricted
to the HS complexes are m = -0.023Å and s=0.054Å\.
The reason seemed to be in the inherent ``stiffness'' of the Morse potential.
In order to avoid this, we tested another MM bond stretching potential for the
metal-ligand bonds in Fe(II) complexes:
| (0) |
| ||||||||||||||||||||||||||||||||||||||
The Fig. 1 presents the curves of the equivalent Morse and NR potentials. One can see that the NR potential increases much slower in the asymptotic region than the Morse one. The parameters of the NR potential equivalent to the Morse potential fitted in our paper [] are given in Table . The value of the A parameter can be identified with the interaction of some Coulomb charges. Extracting these effective values in Fe(II) complexes with nitrogen-containing complexes we get QFe=1.757 e and QN = -0.293 e; the latter is close to the CNDO charges on the donor nitrogen atoms obtained in the EHCF calculations for hexaammine Fe(II) complex [].
Then the NR potential with the parameters of Table is tested on the set of Fe(II) HS and LS molecules described in the paper [].
| Parameter | Bond | |
| Fe-NA | Fe-N3 | |
| A, kcal/mol·Å | -189.3 | -161.3 |
| B, kcal/mol·Å5 | -1084.4 | -1940.9 |
| C, kcal/mol·Å9 | 10817.4 | 19803.2 |
| D0 | 110.0 | 102.9 |
| a | 1.49 | 1.59 |
| r0 | 1.88 | 1.96 |
The corresponding empirical distribution functions for the Fe-N distances' deviations for the overall data set and for the HS and LS subsets separately, are given in Figs. 2-4 We get the following characteristics of the normal distribution of the deviations for the EHCF/MM calculations with NR potential on the overall data set: m = -0.031 Å, s=0.052Å. The both parameters are smaller than those obtained for the Morse potential thus indicating some improvement both in terms of systematic errors and the scattering of data. It is clearly seen from the EDF plots however that the systematic error remains; on the other hand, it is also seen that the mentioned difference in descriptions of the LS and HS complexes is now removed so the structures of both types of complexes are now reproduced accordingly with almost the same systematic error.
As it is known, the minima of the HS and LS terms lie on the right and on the left positions relative to the crossection point of the pure electronic quintet and singlet terms of the d-shell in Fe(II) complexes. Systematic error may be due to some shift of that point from its true position. We obtained a negative value of the mean error which is an indication that the crossection point is shifted towards shorter bond lengths. Thus, moving the crossection point to larger Fe-N distances, it is possible to remove the systematic error. To get rid of the systematic error its possible source has to be identified. As discussed in our paper [], the position of the crossection point in the EHCF/MM method depends on the the Racah parameters B and C, which we have previously accepted as those for the free ion. It can be achieved by slightly reducing of the Racah parameters as compared to the free ion values, to the values B=850 cm-1 and C=3400 cm-1. The set of MM parameters should be also changed in this case. We have done it first for the NR potential. The parameters for the Morse potential are also calibrated independently for the scaled Racah parameters. For the new values see Table .
| Parameter | Bond | |
| Fe-NA | Fe-N3 | |
| A, kcal/mol·Å | -149.0 | -136.0 |
| B, kcal/mol·Å5 | -1385.0 | -1660.0 |
| C, kcal/mol·Å9 | 13650.0 | 18000.0 |
| D0 | 73.0 | 65.0 |
| a | 1.69 | 1.73 |
| r0 | 1.94 | 1.97 |
Thus, a set of semiempirical methods based on EHCF approach allows with good precision to calculate geometrical characteristics (structure) and spectral transitions (splitting, electronic and Mössbauer spectra) of Fe(II) and Co(II) complexes, which is hardly accessible by existent QC methods or can be done only by enormous computational cost.
In the present paper we tried to demonstrate that the problems faced by most empirical and by (actual and so called) ab initio techniques when applied to modeling TMCs have deep roots in the specific features of the electronic structure of the latter and in approximations which tacitly drop the necessary elements of the theory required to reproduce these features of the former. Of course, the EHCF approach whose success story is described here in details is not completely isolated from other methods. In general picture, the various CAS techniques must be mentioned in relation to it, first of all. The characteristic feature uniting these two otherwise very different approaches is the selection of a small subset of one electron states followed by performing adequately complete correlation calculation restriced to this smaller subset. The general problem with such aproaches is that ususally it is taken for granted that the HRF MO LCAO is a good source for obtaining the states to be used in the correlated calculation. Two pitfalls can be expected and actually occur on this route. The first is that in the TMCs the HFR MO LCAOs can be difficult to obtain or those obtained are of a poor quality. The second is that even if the MO LCAOs are obtained correctly, they provide too much delocalized picture of electron distribution. In terms first proposed by J.-P. Malrieu and then extensively used by P. Fulde it is equivalent to saying that in the HFR solution for the TMC the number of electrons in the d-shell too much fluctuates around may be correct average (integer) value. In both cases the limited CI (CAS) techniques are applied to improve a very poor zero approximation. Taking only five MOs of appropriate symmetry to model the d-shell may be too naïve since the number of states to be included in the CI formation to reduce the excessive fluctuations must be much larger. Going to the one-electron states obtained from the canonical MO LCAOs by some localization technique, may be useful, but numerically expensive. The EHCF here advantageously uses the fact that the exact wave function of the TMCs most probably corresponds to very high localization of electrons in the d-shell which enables taking their delocalization into account as a perturbation. Among other approaches based on a similar vision of the situation in TMCs ones of Refs. Nieuwpoort1988,Seijo1999 must be mentioned.
When it comes to analysis of similar approaches stemming from the DFT the numerous attempts to cope with the multiplet states must be mentioned Goerling1993,Nagy1998. In these papers an attempt is made to construct symmetry dependent functionals capable to distinguish different multiplet states in a general direction proposed by [,]. It tuns out however, that the result [] is demonstrated for the lower multiplets of the C atom which are all Roothaan terms. It is not clear that this methodology is not going to work when applied to the d-shell multiplets which may be either non-Roothaan ones or even nontrivially correlated multiple terms.
Another group of appoaches can be qualified as an attempt of using the DFT in order to evaluate the parameters of the CFT/LFT theory. In this respect the papers [,] must be mentioned. The latter in a sence follows the same line as the old semiempirical implementation Zerner1996 where the MOs for the TMC molecule are first obtained by an approximate SCF-like procedure and then a CI is done in some restricted subspace of the latter. In some sence this approach is similar to the EHCF model too with the general difference that the one-electron states used to construct the complete CFT/LFT manifold are taken ''as is'' from the KS calculation. In this case one can expect some difficulties while selecting the MOs into the set of those to be used in constructing the CI (it is not obvious whether simple energy/symmetry criteria allow to select the necessary manifold of the KS orbitals to reproduce the states in the d-shell; and what shall be done when the symmetry is low?). Also the degree of delocalization of the KS orbitals may interfere in evaluation of the CFT/LFT parameters from the results of the DFT calculation. It looks like that it is precisely what happened in [] where the values of the Racah parametrs turned out to be strongly underestimated as compared to the values known to fit the experiment within the CFT/LFT model by this indicating the excess of delocalization of the KS orbitals as compared to that necessary to reproduce the experimental data.
Generally one can notice that almost whatever review on computational chemistry of TMCs starts from a sort of ''triple denial'' of the old CFT/LFT approaches as being pertinent to something which was happening ''once upon a time''. Our point of view on the CFT/LFT picture is absolutely different. It more or less corresponds to that given in the brilliant introduction to the paper []. The clearcut conclusion to be derived from there is that the CFT/LFT picture keeps track of very physical picture of the low-energy spectrum of the TMCs. Whatever discrepance between the results obtained by no matter how refined QC methods and those appearing from the CFT/LFT must be considered as failures of the QC rather than ``age effects'' of the CFT/LFT. It is the purpose for a QC study to reproduce results obtained within the CFT/LFT paradigm and it is not easily reacheable and in many cases has not been reached yet. This idea was the leading one in our studies on TMCs from very beginning and its adequate formal representation in terms of the group functions and the Löwdin partition technique provided a crucial step forward which allowed the numerical implementation of the EHCF method []. It immediately solved the problem of constructing semi-empirical description of the TMCs which otherwise remained unaccessible for 30 years. The cost of this was rejecting the HFR from of the wave function of the TMC which in the present context cannot be considered as a big loss. Further development of this approach and realizing its deeper relation to the general QM/MM setting helped in evolving the corresponding EHCF/MM hybrid scheme. The latter is in relation with those proposed by Deeth [] and Berne []. Both involve the d-shell energy as an additional contribution to that of the MM scheme and use the AOM model with interpolated parameters to estimate the latter. In the case of the approach [] there are two main problems. First is that the AOM parameters involved are assumed to depend only on the interatomic separation between the metal and donor atoms. This is obviously an oversimplification since from the formulae Eq. (27) it is clear that the lone pair orientation is of crucial importance. This is taken into account in the EHCF/MM method Second important flaw is the absence of any correlation in describing the d-shell in the model []. This precludes correct desciption of the switch between different spin states of the open d-shell, although in some situations different spin states can be described uniformly.
In the present paper we tried to demonstrate the feasibility of a semiempirical description of electronic structure and properties of the Werner TMCs on a series of examples. The main feature of the proposed approach was the careful following to the structural aspects of the theory in order to preclude the loss of its elements responsible for description of qualitative physical behavior of the objects under study, in oiur case of TMCs. If it is done the subsequent parameterization becomes sensible and succesful solutions of two long lasting problems: semi-empirical parameterization of transition metals complexes and of extending the MM description to these atoms can be suggested.
This work is supported by the RFBR grants Nos 04-03-32146, 04-03-32206. The authors are thankful to Profs. Paul Ziesche, Lothar Fritsche, Francesc Illas, and Marc Casida for sending (p)reprints of their works, to Profs. I.G. Kaplan, E. Ludeña, J.-P. Julien and to Dr. V.I. Pupyshev for valuable discussions. Organizers of the 9-th European Symposium on Quantum Systems in Chemistry and Physics (QSCP-9) in Les Houches are acknowledged for a kind support extended to A.L.T.