document: there has been recent progress in understanding physics of strongly correlated electronic systems and their electronic structure near a localization  delocalization transition through the development of dynamical mean 
   field theory ( dmft ) @xcite . merging this computationally tractable many  body technique with realistic local  density 
   approximation ( lda )  @xcite based electronic structure calculations of strongly correlated solids is promising due to its simplicity and correctness in both band and atomic limits . at present , much effort is being made in this direction including the developments of a lda+dmft method anisimovkotliar , lda++ approach @xcite , combined gw and dmft theory @xcite , spectral density functional theory @xcite as well as applications to various systems such as la@xmath0sr@xmath1tio@xmath2 latio3 , v@xmath3o@xmath2 @xcite , fe and ni @xcite , ce @xcite , pu @xcite , transition metal oxides@xcite , and many others . for a review , see ref . . 
   such _ ab initio _ dmft based self 
   consistent electronic structure algorithms should be able to explore all space of parameters where neither dopings nor even degeneracy itself is kept fixed as different states may appear close to the fermi level during iterations towards self  consistency . 
   this is crucial if one would like to calculate properties of realistic solid state system where bandwidth and the strength of the interaction is not known at the beginning . 
   it is very different from the ideology of model hamiltonians where the input set of parameters defines the regime of correlations , and the corresponding many 
   body techniques may be applied afterwards . 
   realistic dmft  simulations of material properties require fast scans of the entire parameter space to determine the interaction for a given doping , degeneracy and bandwidth via the solution of the general multiorbital anderson impurity model ( aim ) @xcite . 
   unfortunately , present approaches based on either non  crossing approximation ( nca )  or iterative perturbation theory ( ipt ) are unable to provide the solution to that problem due to a limited number of regimes where these methods can be applied @xcite . 
   the quantum monte carlo ( qmc ) technique dmft , jarrell is very accurate and can cope with multiorbital situation but not with multiplet interactions . 
   also its applicability so far has been limited either to a small number of orbitals or to unphysically large temperatures due to its computational cost . 
   recently some progress has been achieved using impurity solvers that improve upon the nca approximation rotors , jeschke , haule:2001 , but it has not been possible to retrieve fermi liquid behavior at very low temperatures with these methods in the orbitally degenerate case .    as universal impurity solvers have not yet being designed in the past we need to explore other possibilities , and this paper proposes interpolative approach for the self  energy in general multiorbital situation . 
   we stress that this is not an attempt to develop an alternative method for solving the impurity problem , but follow  up of the ideology of lda theory where approximations were designed by analytical fits @xcite to the quantum monte carlo simulations for homogeneous electron gas @xcite . 
   numerically very expensive qmc calculations for the impurity model display smooth self 
   energies at imaginary frequencies for a wide range of interactions and dopings , and it is therefore tempting to design such an interpolation . we also keep in mind that for many applications a high precision in reproducing the self  energies may not be required . 
   one of such applications is , for example , the calculation of the total energy ce , nature , science , nioprl which , as well known from lda based experience , may not be so sensitive to the details of the one  electron spectra . as a result , we expect that even crude evaluations of the self 
   energy shapes on imaginary axis may be sufficient for solving many realistic total energy problems , some of which have appeared already @xcite . 
   another point is a computational efficiency and numerical stability . bringing full self  consistent loops with respect to charge densities nature and 
   other spectral functions require many iterations towards the convergency which may not need too accurate frequency resolutions at every step . 
   however , the procedure which solves the impurity model should smoothly connect various regions of the parameter space . 
   this is a crucial point if one would like to have a numerically stable algorithm and our new interpolational approach ideally solves this problem .    in the calculations of properties such as the low  energy spectroscopy and 
   especially transport more delicate distribution of spectral weight is taken place at low energies , and the imaginary part of the analytically continued self  energy needs to be computed with a greater precision . 
   here we expect that our obtained spectral functions should be used with care . 
   also , in a few well distinct regimes , such , e.g. , as very near the mott transition , the behavior maybe much more complicated and more difficult to interpolate . 
   for the cases mentioned above  extensions of the interpolative methods should be implemented and its beyond the scope of the present work . 
   we can achieve a fast interpolative algorithm for the self 
   energy by utilizing a rational representation . 
   the coefficients in this interpolation can be found by forcing the self  energy to obey several limits and constrains . 
   for example , if infinite frequency ( hartree  fock ) limit , positions of the hubbard bands , low  frequency mass renormalization @xmath4 , mean number of particles @xmath5 as well as the value of the self  energy at zero frequency @xmath6 are known from independent calculation , the set of interpolating coefficients is well defined . in this work , 
   we explore the slave  boson mean  field ( sbmf ) approach @xcite and the hubbard i approximation @xcite to determine the functional dependence of these coefficients upon doping , degeneracy and the strength of the interaction @xmath7 . 
   we verify all trends produced by this interpolative procedure in the regimes of weak , intermediate and strong interactions and at various dopings conditions . 
   these trends are compared with known analytical limits as well as against calculations using the quantum monte carlo method . 
   also , compared with qmc are self  energies and spectral functions on both imaginary and real axes for selective values of dopings . 
   they indicate that the sbmf approach can predict such parameters of interpolation as @xmath8and @xmath4 with a good accuracy while the hubbard i method fails in a number of regimes . 
   however , the functional form of the atomic green function which appears within hubbard i can be used to determine positions of atomic satellites which helps to impose additional constraints on our procedure .    giving an extraordinary computational speed of this approach we generally find a very satisfactory accuracy in comparisons with the numerically more accurate qmc calculations . 
   if an increased accuracy is desired our method can be naturally extended by imposing more constraints and by implementing more refined impurity solvers other than the ones explored in this work . 
   the paper is organized as follows . in section 
   ii we discuss rational interpolation for the self  energy and list the constraints . in section iii 
   we discuss methods for solving anderson impurity model such as the slave  boson mean field and the hubbard i approximations which can be used to find these constraints . 
   a brief survey of the qmc method used to benchmark our algorithm is also given . 
   we present numerical comparisons of sbmf and hubbard i techniques against the qmc simulations for such quantities as quasiparticle residue and multiple occupancies . in section 
   iv we report the results of the interpolative method and compare the obtained spectral functions with the qmc . 
   in section v we discuss possible improvements of the algorithm . 
   section vi is the conclusion . 
   to be specific , we concentrate on the anderson impurity hamiltonian @xmath9 ,   \label{ham}\end{aligned}\]]describing the interaction of the impurity level @xmath10 with bands of conduction electrons @xmath11 via hybridization @xmath12 . 
   @xmath7 is the coulomb repulsion between different orbitals in the @xmath13band . 
   inspired by the success of the iterative perturbation theory @xcite , in order to solve the anderson impurity model in general multiorbital case , we use a rational interpolative formula for the self  energy . 
   this can be encoded into a form@xmath14}{\prod\limits_{m=1}^{m}[\omega -p_{m}^{(\sigma ) } ] } , \label{rat}\]]the coefficients @xmath15 , @xmath16 , or , alternatively , the poles @xmath17 , zeroes @xmath18 and @xmath19 in this equation are to be determined . the form ( [ rat ] ) 
   can be also viewed as a continuous fraction expansion but continuous fraction representation will not be necessary for the description of the method . 
   our basic assumption is that only a well distinct set of poles in the rational representation ( [ rat ] ) is necessary to reproduce an overall frequency dependence of the self  energy . 
   extensive experience gained from solving hubbard and periodic anderson model within dmft at various ratios of the on  site coulomb interaction @xmath7 to the bandwidth @xmath20 shows the appearance of lower and upper hubbard bands as well as renormalized quasiparticle peak in the spectrum of one  electron excitations @xcite . 
   it is clear that the hubbard bands are damped atomic excitations and to the lowest order approximation appear as the positions of the poles of the atomic green function . in the @xmath21 symmetry case which is described by the hamiltonian ( [ ham ] ) 
   , these energies are numerated by the number of electrons occupying impurity level , i.e. @xmath22 and the atomic green function takes a simple functional form @xmath23where @xmath24 are the probabilities to find an atom in configuration with @xmath25 electrons while combinatorial factor @xmath26 arrives due to equivalence of all states with @xmath25 electrons in @xmath21 . 
   we can represent the atomic green function ( [ leh ] ) using the rational representation ( [ rat ] ) , i.e. @xmath27}{\prod\limits_{n=1}^{n}[\omega 
   -p_{n}^{(g ) } ] } ,   \label{gzp}\]]where @xmath28 are all @xmath29 atomic poles , while @xmath30 denote @xmath31 zeroes with @xmath29 being the total degeneracy . 
   the centers of the hubbard bands are thus located at the atomic excitations @xmath32 . using standard expression for the atomic green function 
   @xmath33 $ ] , we arrive to a desired representation for the atomic self  energy    @xmath34}{\prod\limits_{n=1}^{n-1}[\omega -z_{n}^{(g ) } ] } . 
   \label{sat}\ ] ]    using this functional form for finite @xmath35 modifies the positions of poles and zeroes via recalculating probabilities @xmath24 which is equivalent to the famous hubbard i approximation ( discussed in more detail in the next section ) . 
   we now concentrate on the description of the quasiparticle peak which is present in metallic state of the system . for this 
   an extra pole and zero have to be added in eq ( [ sat ] ) . to see this , 
   let us consider the hubbard model for the @xmath21 case where the local green function can be written by the following hilbert transform @xmath36.$ ] if self  energy lifetime effects are ignored , the local spectral function becomes @xmath37 $ ] where @xmath38 is the non  interacting density of states . 
   the peaks of the spectral functions thus appear as zeroes in eq . 
   ( [ sat ] ) and in order to add the quasiparticle peak , one needs to add one extra zero ( denoted hereafter as @xmath39 to the numerator in eq . 
   ( [ sat ] ) . to make the self  energy 
   finite at @xmath40 one has to also add one more pole ( denoted hereafter as @xmath41 ) which should appear in the denominator or eq . 
   ( [ sat ] ) . 
   furthermore , frequently the hartree fock value for the self  energy can be computed separately and it is desirable to have a parameter in the functional form ( [ sat ] ) which will allow us to fix @xmath19 . an obvious candidate to be changed is that self  energy pole in ( [ sat ] ) which is closest to @xmath42 equal zero . 
   let us denote this parameter as @xmath43 and rewrite the denominator of ( [ sat ] ) as @xmath44 $ ] where the product is now extended over all zeroes of the atomic green functions except the one closest to zero and two extra poles @xmath41 and @xmath43 can control the width of the quasiparticle peak and @xmath45 thus , we arrive to the functional form for the self  energy    @xmath46}{(\omega -p_{1}^{(\sigma ) } ) ( \omega -p_{2}^{(\sigma ) } ) \prod\limits_{n=1}^{n-2}[\omega -z_{n}^{(g)}]}. \label{szp}\ ] ]    we are now ready to list all constrains of our interpolative scheme . to fix the unknown coefficients @xmath47 @xmath48 @xmath49 @xmath50 in eq . 
   ( [ szp ] ) and to write down the linear set of equations for the coefficients @xmath15 , @xmath16 in eq . 
   ( [ rat ] ) . 
   we use the following set of conditions . 
   _ a ) _ _ hartree fock value _ @xmath45 in the limit @xmath40 the self  energy takes its hartree  fock form @xmath51mean level occupancy _ _  _ _ @xmath52is defined as a sum over all matsubara frequencies for the green s function , i.e.@xmath53where @xmath54defines the impurity green function and @xmath35 is the hybridization function .    _ 
   b ) _ _ zero  frequency value _ 
   @xmath55 the so called friedel sum rule establishes the relation between the total density and the real part of the self  energy at zero frequency@xmath56    _ c ) quasiparticle mass renormalization value _ @xmath57 the slope of the self  energy at zero frequency controls the quasiparticle residue , @xmath4 using the following relationship @xmath58formally , constraints ( _ b _ ) and ( _ c _ ) hold for zero temperature only but we expect no significant deviations in many regions of parameters as long as we stay at low enough temperatures 
   . the behavior may be more elaborated in the vicinity of mott transition @xcite .    _ d ) _ _ positions of hubbard bands . 
   _    as we discussed , in order that the self  energy obeys the atomic limit and places the centers of the hubbard bands at the positions of the atomic excitations , we demand that    @xmath59    this condition fixes almost all self  energy zeroes @xmath18 in eq.([rat ] ) to the poles @xmath60 however , it alone does not ensure that the weight is correctly distributed among the hubbard bands and that the very distant hubbard bands disappear . for this to occur 
   , distant poles of green function have to be canceled out by nearby zeroes of the green function . 
   it is clear that each pole @xmath28 far from the fermi level has to be accompanied by a nearby zero @xmath30 in order the weight of the pole be small . 
   thus , the self  energy has poles at the position of green s function zeroes which can be encoded into the constrain @xmath61^{-1}=0 . 
   \label{zer}\]]we want to keep this property of the self 
   energy for finite @xmath62 and thus demand that self 
   energy diverges ( when lifetime effects are kept , it only reaches a local maximum ) at the zeros of the functional form ( [ leh ] ) of @xmath63 note that the relationship ( [ zer ] ) holds ( approximately ) for frequency @xmath42 larger than the renormalized bandwidth @xmath64 therefore the information about one @xmath65 which lies close to @xmath66 is omitted and replaced by the information about @xmath67 @xmath6 and @xmath4 as it is done by separating @xmath41 and @xmath43 in the denominator of eq . 
   ( szp ) . 
   we can now write down  a set of linear equations for all unknown coefficients in the expression ( [ rat ] ) . there is total @xmath68 of parameters @xmath15 and @xmath69 where we can always set @xmath70 the conditions _ 
   a)_,_b),c ) _ give @xmath71according to condition _ 
   d ) _ we can use all @xmath29 poles @xmath28 and @xmath72 zeroes @xmath73 the zero @xmath30 closest to @xmath66 is dropped out . 
   this brings additional @xmath74 equations for the coefficients and makes @xmath75 as the degree of the rational interpolation which are written below    @xmath76^{m}-(p_{n}^{(g)}+\mu -\epsilon _ { f})\sum_{m=0}^{n}b_{m}[p_{n}^{(g)}]^{m } & = & 0\text { for } n=1 ... n , \label{co5 } \\ 
   \sum_{m=0}^{n}b_{m}[z_{n}^{(g)}]^{m } & = & 0\text { for } n=1 ... n-2 .\label{co6}\end{aligned}\ ] ]    note that while @xmath77 may be rather large , the actual number of poles contributing to the self  energy behavior is indeed very small . 
   we can directly see this from eq . 
   ( [ sat ] ) which uses all @xmath29 poles @xmath28 fulfilling eq . 
   ( [ pol ] ) and uses @xmath72 zeroes @xmath30 directly related to @xmath72 poles @xmath78 clearly , when the spectral weight of the atomic excitation becomes small , the corresponding @xmath79 becomes close to @xmath30 and the cancellation occurs . 
   therefore in realistic situations when only the upper and lower hubbard bands have significant spectral weight along with the quasiparticle peak , the actual degree of the polynomial expansion is either two or three . 
   however , it is advantageous numerically and cheap computationally to keep all poles and zeroes in eq . 
   ( [ szp ] ) because the formula automatically distributes spectral weight over all existing hubbard bands .    in the limit when @xmath80 the self 
   energy automatically translates to the non  interacting one . 
   the atomic poles get close to each other but , most importantly , their spectral weight goes rapidly to zero as it gets accumulated within the quasiparticle band .    in the mott insulating regime , 
   the conditions _ b ) _ and _ c ) _ drop out while all poles @xmath28 and zeroes @xmath30 can be used to determine the interpolation . 
   however , in this regime it does not matter whether one of @xmath30 closest to @xmath66 is dropped out or kept , since we can always replace this information by information about @xmath81 therefore the mott transition can be studied without changing the constraints . 
   we thus see that in the insulating case the self  energy correctly reproduces the well  known result of the hubbard i  method where the green function is computed after eq . 
   ( [ imp ] ) with atomic self  energy . if the lifetime effects are computed , the parameters @xmath28 and @xmath30 become complex and the hubbard bands will acquire an additional bandwidth . 
   this effect is evident from the simulations using various perturbative or qmc impurity solvers and can be naturally incorporated into the interpolative formulas ( [ rat ] ) or ( [ szp ] ) . 
   however , in practical implementation below we will omit it for illustrative purposes . 
   let us now discuss the quality of interpolation from the perspective of the high  frequency behavior for the self  energy . 
   the latter can be viewed moments as expansion in terms of the moments @xmath82 , i.e. , @xmath83 . 
   most important for us is to look at highest moments which are given by the hartree fock value , eq . 
   ( [ inf ] ) involving single occupancy matrix @xmath84 , as well as the first moment    @xmath85u^{2 } ,   \label{ms1}\ ] ]    containing a double occupancy matrix @xmath86 we see that the interpolation in part relies on the accuracy in computing multiple occupancies which are the functionals of both atomic excitations and the hybridization function . in this regard , using exact atomic green function to find poles @xmath28 and zeroes @xmath30 as part of the constrained procedure may not be as accurate since it would assume the use of _ _  atomic _ _ multiple occupancies which do _ not _ carry information about @xmath87on the other hand , we can also use only a functional form of the atomic green function where the multiple occupancies are computed in a more rigorous manner . in the next section 
   we will show how this can be implemented using the sbmf multiple occupancies which will be found to be in better agreement with the quantum monte carlo data . 
   note that the moments @xmath82 themselves can be used in establishing the constraints for interpolation coefficients . 
   this would involve independent evaluations of @xmath88 , etc . as well as various integrals involving hybridization function @xmath35 . 
   however , we may run into ill  defined numerical problem since high  frequency information will be used to extract the low  frequency behavior . 
   therefore , it is more advantageous numerically to use some poles and zeros of @xmath89 as given by condition _ 
   d ) _ above . 
   we thus see that the interpolational scheme is defined completely once a prescription for obtaining parameters such as @xmath90 @xmath91 as well as poles @xmath28 , and zeroes @xmath30 is given . for this purpose 
   we will test two commonly used methods : sbmf method due to gutzwiller @xcite as described by kotliar and ruckenstein ruck and the well  known hubbard i approximation @xcite . 
   we compare these results against more accurate but computationally demanding quantum monte carlo calculations and establish the procedure to extract all necessary parameters . 
   note that once the constraints such as @xmath4 are computed from a given approximate method , some of the quantities such as the total number of particles , @xmath84 , and the value of the self  energy at zero frequency , @xmath92 can be computed fully self  consistently . 
   they can be compared with their non  self  consistent values . 
   if the approximate scheme already provides a good approximation for @xmath84 and satisfies the friedel sum rule , the self  consistency can be avoided hence accelerating the calculation . 
   indeed we found that inclusion of the self 
   consistency improves the results only marginally except when we are in a close vicinity to the mott transition but here we do not expect that our simple interpolative algorithm is very accurate . 
   we now give the description of the approximate methods for solving the impurity model and then present the comparisons of our interpolative procedure with the qmc calculations . 
   the quantum monte carlo method is a powerful and manifestly not perturbative approach in either interaction @xmath7 or the bandwidth @xmath20 . in the qmc method one 
   introduces a hubbard  stratonovich field and averages over it using the monte carlo sampling . 
   this is a controlled approximation using different expansion parameter , the size of the mesh for the imaginary time discretization . 
   unfortunately it is computationally very expensive as the number of time slices and the number of hubbard  stratonovich fields increases . 
   also the method works best at imaginary axis while analytical continuation is less accurate and has to be done with a great care . 
   extensive description of this method can be found in ref . . 
   we will use this method to benchmark our calculations using approximate algorithms described later in this section .      a fast approach 
   to solve a general impurity problem is the slave  boson method  @xcite . at the mean field level 
   , it gives the results similar to the famous gutzwiller approximation @xcite . 
   however , it is improvable by performing fluctuations around the saddle point . 
   this approach is accurate as it has been shown recently to give the exact critical value of @xmath7 in the large degeneracy limit  at half  filling florens . 
   the main idea is to rewrite atomic states consisting of @xmath25 electrons @xmath93 , @xmath94 with help of a set of slave 
   bosons @xmath95 . in the following , 
   we assume @xmath21 symmetric case , i.e. , equivalence between different states @xmath96 for fixed @xmath25 . formulae corresponding to a more general crystal  field case are given in appendix b. the creation operator of a physical electron is expressed via slave particles in the standard manner   @xcite . in order to recover the correct non  interacting limit at the mean  field level , the bose fields 
   @xmath97 can be considered as classical values found from minimizing a lagrangian @xmath98 corresponding to the hamiltonian ( [ ham ] ) . 
   two lagrange multipliers @xmath99 and @xmath100 should be introduced in this way , which correspond to the following two constrains : @xmath101    the numbers @xmath102 are similar to the probabilities @xmath24 discussed in connection to the formula for the atomic green function ( leh ) . 
   we thus see the physical meaning of the first constrain which is the sum of probabilities to find atom in any state is equal to one , and the second constrain gives the mean number of electrons coinciding with that found from @xmath103 . a combinatorial factor 
   @xmath104 arrives due to assumed equivalence of all states with @xmath25 electrons . 
   minimization of @xmath98 with respect to @xmath105 leads us to the following set of equations to determine the quantities @xmath105 : @xmath106\psi _ { n}+nbt\sum_{i\omega } \delta ( i\omega ) g_{g}(i\omega ) [ lr\psi _ { n-1}+\psi _ { n}bl^{2}]+\]]@xmath107=0 ,   \label{guz}\]]where @xmath108 , determines the mass renormalization , and the coefficients @xmath109 , @xmath110 are normalization constants as in refs . . 
   @xmath111 is the total energy of the atom with @xmath25 electrons in @xmath21 approximation . 
   ( [ guz ] ) , along with the constrains ( [ one ] ) , ( [ two ] ) constitute a set of non  linear equations which have to be solved iteratively . in practice , we consider eq . 
   ( [ guz ] ) as an eigenvalue problem with @xmath100 being the eigenvalue and @xmath105 being the eigenvectors of the matrix . 
   the physical root corresponds to the lowest eigenvalue of @xmath100 which gives a set of @xmath105 determining the mass renormalization @xmath112 since the matrix to be diagonalized depends non 
   linearly on @xmath105 via the parameters @xmath113 and @xmath114 and also on @xmath99 , the solution of the whole problem assumes the self  consistency : i ) we build an initial approximation to @xmath105 ( for example the hartree 
   fock solution ) and fix some @xmath99 , ii ) we solve eigenvalue problem and find new normalized @xmath105 , iii ) we mix new @xmath105 with the old ones using the broyden method  @xcite and build new @xmath115 and @xmath114 . 
   steps ii ) and iii ) are repeated until the self  consistency with respect to @xmath105 is reached . during the iterations we also vary @xmath99 to obey the constrains . 
   the described procedure provides a stable computational algorithm for solving aim and gives us an access to the low  frequency green s function and the self  energy of the problem via knowledge of the slope of @xmath116 and the value @xmath117 at zero frequency .    , 
   ( b ) dependence of the spectral weight @xmath118 on concentration , and ( c ) density  density correlation function , @xmath119 versus filling , @xmath5 , in the two  band hubbard model in @xmath120 and @xmath121,title="fig:",scaledwidth=50.0% ] +    the described slave  boson method gives the following expression for the self  energy : @xmath122the impurity green function @xmath123 in this limit is given by the expression @xmath124    as an illustration , we now give the solution of eq . 
   ( [ guz ] ) for non  degenerate case ( @xmath125 and at the particle  hole symmetry point with @xmath126 . consider a dynamical mean  field theory for the hubbard model which reduces the problem to solving the impurity model subject to the self  consistency condition with respect to @xmath35 . 
   starting with the semicircular density of states ( dos ) , the self  consistency condition is given by eq . 
   ( [ guz])@xmath127 we obtain the following simplifications : @xmath128 , @xmath129 , @xmath130 and @xmath131^{-1}.$ ] the sum @xmath132 appeared in eq . 
   ( [ guz ] ) scales as @xmath133 with the constant @xmath134 being the characteristic of a particular density of states and approximately equal to 0.2 in the semicircular dos case . 
   the self  consistent solution of eq . 
   ( [ guz ] ) is therefore possible and simply gives @xmath135 . 
   the mott transition occurs when no sites with double occupancies can be found , i.e. when @xmath136 a critical value of @xmath137 . for @xmath138 , 
   this gives @xmath139 and reproduces the result for @xmath140 known from the qmc calculation within a few percent accuracy . 
   as degeneracy increases , critical @xmath7 is shifted towards higher values @xcite . from numerical calculations we obtained the following values of the critical interactions in the half 
   filled case @xmath141 for @xmath142 ( @xmath143level)@xmath144 @xmath145 for @xmath146 ( @xmath147level)@xmath144 and @xmath148 for @xmath149 ( @xmath13level ) . 
   density  density correlation function @xmath119 for local states with @xmath25 electrons is proportional to the number of pairs formed by @xmath25 particles @xmath150 . since the probability for @xmath25 electrons to be occupied is given by : @xmath151 , the physical density 
   density correlator can be deduced from : @xmath152 . 
   similarly , the triple occupancy can be calculated from @xmath153 . 
   let us now check the accuracy of this method by comparing its results with the qmc data . 
   we consider the two  band hubbard model in @xmath154 orbitally degenerate case . 
   hybridization @xmath155 satisfies the dmft  self  consistency condition of the hubbard model on a bethe lattice @xmath156the coulomb interaction is chosen to be @xmath157 which is sufficiently large to open the mott gap at integer fillings . 
   all calculations are done for the temperature @xmath158 . 
   we first compare the average number of electrons vs. chemical potential determined from the slave bosons which is plotted in fig . 
   [ figsbmf](a ) . 
   this quantity is sensitive to the low 
   frequency part of the green function which should be described well by the present method . 
   we see that it reproduces the qmc data with a very high accuracy and only differs by 20 per cent very near the jump of @xmath159 at @xmath160 . 
   the quasiparticle residue @xmath4 versus filling @xmath5 is plotted in fig . [ 
   figsbmf](b ) . the slave  boson method gives the fermi liquid and provides estimates for the quasiparticle residue with the overall discrepancy of the order of 30% . 
   in fact , we have performed several additional calculations for other degeneracies ( @xmath161 and @xmath162 ) and for various parameters regimes . 
   the trend to overestimate mass renormalization can be seen in many cases . 
   it disappears only when @xmath7 approaches zero . 
   we need to point out , however , that i ) the extractions of zero frequency self 
   energy slopes from the high  temperature qmc is by itself numerically not well grounded procedure , as information for the self  energy is known at the matsubara points only , which is then extrapolated to @xmath66 , ii ) other methods for solving impurity model , such as nca or ipt display similar discrepancies and iii ) recent findings @xcite suggest that at least at half  filling quasiparticle residues deduced from slave bosons becomes exact when @xmath163 . 
   most importantly for our interpolative method is that the entire functional dependence of @xmath4 vs. filling , interaction and degeneracy is correctly captured . 
   its overall accuracy is acceptable as it is evident from our comparisons of the spectral functions presented in the next section and well within the main goal of our work to provide fast scans of the entire parameter space necessary for simulating real materials . 
   this is important as , for example , for general @xmath13electron material , the qmc method is prohibitively time consuming , but we expect from the sbmf method the results for mass renormalization not worse than 50% for such delicate regime as the vicinity of the mott transition . 
   [ figsbmf](c ) shows the density 
   density correlation function @xmath119 as a function of average occupation @xmath5 . 
   the discrepancy is most pronounced for fillings @xmath164 [ see the inset of fig . 
   [ figsbmf](c ) ] where the absolute values of @xmath119 are rather small . 
   although our slave  boson technique captures only the quasiparticle peak , it gives the correlation function in reasonable agreement with the qmc for dopings not too close to the mott transition . 
   now we turn to the hubbard  i approximation  @xcite which is closely related to the moments expansion method  @xcite . 
   consider many  body atomic states @xmath165 which in @xmath21 are all degenerate with index @xmath166 numerating these states for a given number of electrons @xmath167 the impurity green function is defined as the average@xmath168and becomes diagonal with all equal elements in @xmath169it is convenient to introduce the hubbard operators @xmath170and represent the one  electron creation and destruction operators as follows@xmath171the impurity green function ( [ gim ] ) is given by @xmath172    where the matrix @xmath173 is defined as @xmath174    establishing the equations for @xmath173 can be performed using the method of equations of motion for the @xmath175 operators . performing their decoupling due to hubbard @xcite , carrying out the fourier transformation and analytical continuation to the real frequency axis , and summing over @xmath25 and @xmath176 after ( [ hub ] ) we arrive to the net result 
   @xmath177the @xmath178 can be viewed in the matrix form ( [ hub ] ) with the following definition of a diagonal atomic green function@xmath179with @xmath111 being the total energies of the atom with @xmath25 electrons in @xmath169 the coefficients @xmath24 are the probabilities to find atom with @xmath25 electrons and were already discussed in connection to the formula ( [ leh ] ) for the atomic green function . 
   they are similar to the coefficients @xmath102 introduced within the sbmf method but now found from different set of equations . 
   these numbers are normalized to unity , @xmath180 and are expressed via diagonal elements of @xmath181 as follows : @xmath182their determination in principle assumes solving a non  linear set of equations while determining @xmath183 the mean number of electrons can be measured as follows : @xmath184 or as follows @xmath185 the numbers @xmath24 can be also used to find the averages @xmath186 @xmath187 in the way similar to what has been done in the sbmf approach .    , 
   ( b ) dependence of the spectral weight @xmath118 on concentration , and ( c ) density  density correlation function , @xmath188 versus filling , @xmath5 , in the two  band hubbard model in @xmath120 and for @xmath189 . 
   , title="fig:",scaledwidth=50.0% ] +    if we neglect by the hybridization @xmath35 in eq . 
   ( [ h1e ] ) , the probabilities @xmath24 become simply statistical weights : @xmath190we thus see that in principle there are several different ways to determine the coefficients @xmath24 , either via self  consistent determination ( scw ) , or using statistical formula ( [ stw ] ) , or taking them from sbmf equation ( [ guz ] ) , i.e. setting @xmath191 while still utilizing the functional dependence provided by the hubbard i method . to determine the best procedure 
   let us first consider limits of large and small @xmath7 s . 
   when @xmath192 , @xmath123 is reduced to @xmath193 i.e. the hubbard  i method reproduces the atomic limit . 
   setting @xmath194 gives @xmath195^{-1}$ ] , which is the correct band limit . 
   unfortunately , at half  filling this limit has a pathology connected to the instability towards mott transition at any interaction strength @xmath7 . to see this 
   , we consider a dynamical mean  field theory for the hubbard model 
   . using semicircular density of states , we obtain @xmath196^{-1}g_{at}(\omega ) $ ] and conclude that for any small @xmath7 the system opens a pathological gap in the spectrum . 
   clearly , using hubbard i only , the behavior of the green function at @xmath197 can not be reproduced and the quality of the numbers @xmath24 is at question . 
   this already emphasizes the importance of using the slave  boson treatment at small frequencies . 
   ultimately , making the comparisons with the qmc calculations is the best option in picking the most accurate procedure to compute the probabilities @xmath24 . to check the accuracy against the qmc we again consider the two  band hubbard model in @xmath120 symmetry as above . 
   the chemical potential@xmath144 mass renormalization and double occupancy are plotted versus filling in fig . 
   [ fighub1 ] . 
   all quantities here were computed with statistical weights after eq . 
   ( [ stw ] ) but we found a similar agreement while utilizing the self  consistent determination of @xmath24 after eq . 
   ( [ scw ] ) . 
   we first see that the hubbard i approximation does not give satisfactory agreement with the qmc data for @xmath198 because it misses the correct behavior at low frequencies . 
   the comparisons for @xmath199 plotted in fig . 
   [ fighub1](b ) surprisingly show a relatively good behavior . 
   however , the pathology of this approximation at half  filling would predict @xmath200 for any @xmath7 , which is a serious warning not to use it for extracting the quasiparticle weight . 
   [ fighub1](c ) shows @xmath119 as a function of average occupation @xmath5 . as this quantity 
   is directly related to the high  frequency expansion one may expect a better accuracy here . however , comparing fig . 
   [ fighub1](c ) and fig . 
   [ figsbmf](c ) , it is clear that the slave boson method gives more accurate double occupancy . 
   this is due to the fact that the density matrix obtained by the slave boson method is of higher quality than the one obtained from the hubbard i approximation . 
   the results of these comparisons suggest that the probabilities @xmath201 provided by the slave  boson method is a better way in determining the coefficients @xmath24 in the metallic region of parameters . 
   therefore it is preferable to use these numbers while establishing the equations for the unknown coefficients in the interpolational form ( [ rat ] ) . 
   however , the functional form ( [ atg ] ) of the hubbard i approximation with @xmath202 can still be used as it provides the positions of the poles @xmath28 and zeroes @xmath30 of the atomic green function necessary for the condition _ 
   d ) _ in the previous section . 
   this also ensures accurate high  frequency behavior of the interpolated self 
   energy since its moments expressed via multiple occupancies are directly related to @xmath203    interestingly , while more sophisticated qmc approach captures both the quasiparticle peak and the hubbard bands this is not the case for the slave  boson mean field method . to obtain the hubbard bands in this method fluctuations need to be computed , which would be very tedious in a general multiorbital situation . 
   however the slave  boson method delivers many parameters in a good agreement with the qmc results , and , hence , it can be used to give a functional dependence of the coefficients of the rational approximation . 
   by now the procedure to determine the coefficients is well established . 
   we use the sbmf method to determine @xmath204 @xmath205 as well as poles and zeroes of the atomic green function provided by the sbmf probabilities @xmath102 and by the bare atomic energy levels @xmath206 ( we omit the lifetime effects for simplicity)@xmath127 this generates a set of linear equations for coefficients @xmath207 , @xmath208 in the rational interpolation formula ( [ rat ] ) . in the present section 
   we show the trends our interpolative algorithm gives for the spectral functions in various regions of parameters as well as provide detailed comparisons for some values of doping for both imaginary and real axis spectral functions . the two  band hubbard model with semicircular density of states and 
   dmftself  consistency condition after ( [ scf ] ) is utilized in @xmath154 symmetry in all cases using the bandwidth @xmath209 and temperature @xmath210 . 
   [ figu2trend ] shows the behavior of the density of states @xmath211 for @xmath212 as a function of the chemical potential @xmath159 computed with respect to the particle  hole symmetry point @xmath213 and as a function of frequency @xmath214 the semicircular quasiparticle band is seen at the central part of the figure . 
   its bandwidth is only weakly renormalized by the interactions in this regime . 
   it is half 
   filled for @xmath215 ( i.e. when @xmath216 ) and gets fully emptied when chemical potential is shifted to negative values . 
   several weak satellites can be also seen on this figure which are due to atomic poles . 
   their spectral weight is extremely small in this case and any sizable lifetime effect ( which is not included while plotting this figure ) will smear these satellites out almost completely . while approaching fully emptied ( or fully filled situation ) the spectral weight of the hubbard bands disappears completely and only unrenormalized quasiparticle band remains . 
   it is clear that even without shifting the atomic poles to the complex axis , the numerical procedure of generating the self  energy is absolutely stable .     and frequency for the two  band hubbard model in @xmath120 and at @xmath217.,title="fig:",scaledwidth=50.0% ] +     and frequency for the two  band hubbard model in @xmath120  and at @xmath217 . 
   , title="fig:",scaledwidth=50.0% ] +     and frequency for the two  band hubbard model in @xmath120 and at @xmath189.,title="fig:",scaledwidth=50.0% ] +     and frequency for the two  band hubbard model in @xmath120  and at @xmath189 . 
   , title="fig:",scaledwidth=50.0% ] +    this trend can be directly compared with the simulations using a more accurate qmc  impurity solver . 
   we present this in fig . 
   [ figu2trendqmc ] for @xmath212 , which shows calculated density of states in the same region of parameters . 
   remarkably that again we can distinguish the renormalized quasiparticle band and very weak hubbard satellites . 
   the hubbard bands appear to be much more diffuse in this figure mainly due to the lifetime effects and partially due to maximum entropy method using for analytical continuation from imaginary to real axis . 
   otherwise the entire picture looks very much like the one on fig . 
   [ figu2trend ] , generated with much less computational effort . 
   [ figu4trend ] gives the same behavior of the density of states for the strongly correlated regime @xmath218 in this case the situation at integer filling is totally different as the system undergoes metal  insulator transition . 
   this is seen around the dopings levels with @xmath159 between 0 and -1 and between -3 and -5 where the wight of the quasiparticle band collapses while lower and upper hubbard bands acquire all spectral weight . in the remaining region of parameters 
   both strongly renormalized quasiparticle band and hubbard satellites remain . again , once full filling or full emptying is approached the quasiparticle bands restores its original bandwidth while the hubbard bands disappear . the qmc  result for the same region of parameters 
   is given in fig . 
   [ figu4trendqmc ] . 
   again we can distinguish the renormalized quasiparticle band and hubbard satellites as well as the areas of mott insulator and of strongly correlated metal . 
   the hubbard bands appear to be more sharp in this figure which signals on approaching the atomic limit . 
   we now turn to the comparison of the green functions and the self  energies obtained using the formulae ( [ imp ] ) , and ( [ rat ] ) respectively against the predictions of the quantum monte carlo method . 
   we will report our comparisons for the two  band hubbard model and sets of dopings corresponding to @xmath219 using the value of @xmath220 other tests for different degeneracies , doping levels and the interactions have been performed which display similar accuracy .     and 
   @xmath220,scaledwidth=65.0% ]    fig . 
   [ figgrn ] shows the comparison between the real and imaginary parts of the green function obtained by the interpolative method with the results of the qmc calculations . 
   as one can see almost complete agreement has been obtained for a wide regime of dopings . 
   the agreement gets less accurate once the half  filling is approached , but still very good giving an extraordinary computational speed of the given method compared to qmc 
   .     and @xmath220,title="fig:",scaledwidth=65.0% ] +    fig . 
   [ figsig ] shows similar comparison between the real and imaginary parts of the self 
   energies obtained by the interpolative and the qmc method . 
   we can see that the self 
   energies exhibit some noise which is intrinsic to stochastic qmc procedure . 
   the values of the self - energies near @xmath221 and @xmath222 are correctly captured with some residual discrepancies are attributed to slightly different chemical potentials used to reproduce given filling within every method . 
   the results at the imaginary axis show slightly underestimated slopes of the self 
   energies within the interpolative algorithm which is attributed to the underestimated values of @xmath4 obtained from the sbmf calculation . ultimately improving these numbers by inclusions of fluctuations beyond mean field 
   will further improve the comparisons . 
   however , even at the present stage of the accuracy all functional dependence brought by the sbmf method quantitatively captures the behavior of the self  energy seen from time consuming qmc simulation .      . and 
   @xmath220,title="fig:",scaledwidth=65.0% ] +    we also made detailed comparisons between calculated densities of states obtained at the real axis using the interpolative method and the qmcalgorithm . 
   the qmc  densities of states require an analytical continuation from the imaginary to real axis and were generated using the maximum entropy method . by itself 
   this procedure introduces some errors within the qmc especially at higher frequencies . in fig . 
   [ figdos ] , we show our results for the fillings corresponding to @xmath219 using the value of @xmath189 . 
   one can see the appearance of the quasiparticle band and two hubbard bands distanced by the value of @xmath7 . 
   it can be seen that the interpolative method remarkably reproduces the trend in shifting the hubbard bands upon changing the doping . 
   it automatically holds the distance between them to the value of @xmath7 while this is not always true in the quantum monte carlo method . despite this result , 
   the overall agreement between both methods is very satisfactory . 
   here we would like to discuss possible ways to further improve the accuracy of the method . 
   the inaccuracies are mainly seen in three different quantities : i ) the width of the hubbard bands , ii ) the mass renormalization @xmath223 which is borrowed from the sbmf method , and iii ) the number of electrons @xmath224 extracted from the interpolated impurity green function ( [ imp ] ) . the inaccuracy in the width of the hubbard band is mainly connected to neglecting the lifetime effect . 
   provided it is computed , this will shift the positions of atomic poles onto the complex plane which is in principle trivial to account for within our interpolative algorithm . to improve the accuracy of @xmath223 one 
   can , for example , work out a modified slave boson scheme which will account for the fluctuations around mean field solution . 
   the inaccuracy in @xmath224 is small in many regions of parameters and typically amounts to 510 per cent . 
   we can try to improve this agreement by requirement that @xmath224 obtained via interpolation matches with @xmath225 obtained by the sbmf method . 
   the latter agrees very well with the qmc for a wide region of parameters as it is evident from fig . [ figsbmf](a ) . in reality 
   , our analysis shows that in many cases the discrepancy in @xmath224 is connected with the overestimation of @xmath226 . 
   therefore , points ii ) and iii ) mentioned above are interrelated . 
   the requirement that @xmath227 can be enforced by adjusting the width of the quasiparticle band , and in many regions of parameters this is controlled by @xmath228 however , there are situations when the hubbard band appears in the vicinity of @xmath66 , and changing @xmath223 does not affect the bandwidth@xmath127 to gain a control in those cases it is better to replace the constraints @xmath229 by constraints of fixing the self  energies at two frequencies@xmath144 say @xmath230 and @xmath231 where @xmath232 is the frequency of the order of renormalized bandwidth@xmath127 we have found that this scheme brings mass renormalizations which are about 30% smaller than the sbmf ones , and the agreement with the qmc is significantly improved . 
   thus , inaccuracies ii ) and iii ) can be avoided with this very cheap trick . 
   however , we also would like to point out that the condition @xmath227 is essentially non  linear as the solution may not exist for all regions of parameters . 
   it is , for example , evident that in such points where @xmath224 is given by a symmetry ( as , e.g. , particle  hole symmetry point @xmath233 in the case considered above ) the mass renormalization does not affect the number of electrons . 
   as the philosophy of our approach is to get the best possible fit we are also open to implementing any kinds of _ ad hoc _ renormalizations constants . 
   one of such possibility could be the use of a quasiparticle residue 30 per cent smaller than @xmath234 as @xmath223 should go to unity when @xmath235 , the correction can , for example , be encoded into the formula @xmath236.$ ]    we finally would like to remark that the scheme defined by a set of linear equations for the coefficients ( [ co1])([co6 ] ) is absolutely robust as solutions exist for all regimes of parameters such as strength of the interaction , doping and degeneracy . in general , bringing any information on the self energy @xmath237 at some frequency point @xmath238 or its derivative @xmath239 would generate a linear relationship between the interpolation coefficients , thus keeping robustness of the method . on the other hand , fixing such relationships as numbers of electrons brings non  linearity to the problem which could lead to multiplicity or non  existence of the solutions . 
   it is also clear that by narrowing the regime of parameters , the accuracy of the interpolative algorithm can be systematically increased . 
   to summarize , this paper shows the possibility to interpolate the self 
   energies for a whole range of dopings , degeneracies and the interactions using a computationally efficient algorithm . 
   the parameters of the interpolation are obtained from a set of constraints in the slave boson mean field method combined with the functional form of the atomic green function . 
   the interpolative method reproduces all trends in remarkable agreement with such sophisticated and numerically accurate impurity solver as the qmc method . 
   we also obtain a very good quantitative agreement in a whole range of parameters for such quantities as mean level occupancies , spectral functions and self  energies . 
   some residual discrepancies remain which can be corrected provided better algorithms delivering the constraints will be utilized . 
   nevertheless , given the superior speed of the present approach , we have obtained a truly exceptional accuracy times efficiency of the proposed procedure . 
   the work was supported by nsf  dmr grants 0096462 , 02382188 , 0312478 , 0342290 and us doe grant no de 
   the authors also acknowledge the financial support from the computational material science network operated by us doe and from the ministry of education , science and sport of slovenia . 
   in the crystal field case we assume that @xmath29fold degenerate impurity level @xmath10 is split by a crystal field onto @xmath240 sublevels @xmath241 . 
   we assume that for each sublevel there is still some partial degeneracy @xmath242 so that @xmath243 in limiting case of @xmath21 degeneracy , @xmath244 , and in non  degenerate case , @xmath245 . 
   we need to discuss how a number of electrons @xmath25 can be accommodated over different sublevels @xmath246 . introducing numbers of electrons on each sublevel , @xmath247 
   , we obtain @xmath248 note the restrictions : @xmath249 @xmath250 and @xmath251 . in @xmath21 case , @xmath252 , and in non  degenerate case , @xmath253 , @xmath247 is either 0 or 1 . 
   total energy for the shell with @xmath25 electrons depends on particular configuration @xmath254 @xmath255 . 
   \label{slacrflev}\ ] ] many  body wave function is also characterized by a set of numbers @xmath254 , i.e. @xmath256 energy @xmath257 remains degenerate , which can be calculated as product of how many combinations exists to accommodate electrons in each sublevel , i.e. @xmath258 let us further introduce probabilities @xmath259 to find a shell in a given state with energy @xmath260 sum of all probabilities should be equal to @xmath261 , i.e.@xmath262    there are two green functions in gutzwiller method : impurity green function @xmath263 and quasiparticle green function @xmath264 where matrix coefficients @xmath265 represent generalized mass renormalizations parameters . all matrices are assumed to be diagonal and have diagonal elements enumerated as follows : @xmath266 each element in the green function is represented as follows@xmath267@xmath268and determines a mean number of electrons in each sublevel@xmath269the total mean number of electrons is thus : @xmath270 hybridization function @xmath271 is the matrix which is assumed to be diagonal , and it has diagonal elements enumerated as follows : @xmath272 mass renormalizations @xmath273 are determined in each sublevel . 
   the generalization of the non  linear equations ( [ guz ] ) has the form@xmath276 \psi _ { n_{1} ... n_{g}}+   \notag \\ & & \sum_{\alpha = 1}^{g}n_{\alpha } [ t\sigma _ { i\omega } \delta _ { \alpha } ( i\omega ) g_{g\alpha } ( i\omega ) ] b_{\alpha } \left [ r_{\alpha } l_{\alpha } \psi _ { n_{1} ... n_{\alpha } -1 ... n_{g}}+b_{\alpha } l_{\alpha } ^{2}\psi _ { n_{1} ... n_{\alpha } ... n_{g}}\right ] +   \notag \\ & & \sum_{\alpha = 1}^{g}(d_{\alpha } -n_{\alpha } ) [ t\sigma _ { i\omega } \delta _ { \alpha } ( i\omega ) g_{g\alpha } ( i\omega ) ] b_{\alpha } \left [ r_{\alpha } l_{\alpha } \psi _ { n_{1} ... n_{\alpha } + 1 ... n_{g}}+b_{\alpha } r_{\alpha } ^{2}\psi _ { n_{1} ... n_{\alpha } ... n_{g}}\right ] . 
   \label{slacrfscf}\end{aligned}\ ] ]            g. kotliar and s. y. savrasov , in _ new theoretical approaches to strongly correlated systems _ , edited by a. m. tsvelik , ( kluwer academic publishers , the netherlands , 2001 ) , p. 259 
   , ( available in cond 
   mat/020824 ) ; s. biermann , f. aryasetiawan , and a. georges , phys . 
   lett . * 90 * , 086402 ( 2003 ) . 
   k. held , i. a. nekrasov , g. keller , v. eyert , n. bluemer , a. k. mcmahan , r. t. scalettar , th . 
   pruschke , v. i. anisimov , and d. vollhardt , psi  k newsletter # * 56 * ( april 2003 ) , p. 65 ; a. i. lichtenstein , m. i. katsnelson , and g. kotliar , in _ electron correlations and materials properties _ 
   , ed . by a. gonis , n. kioussis and m. ciftan ( kluwer academic , plenum publishers , 2002 ) p. 428
summary: a rational representation for the self  energy is explored to interpolate the solution of the anderson impurity model in general orbitally degenerate case . 
   several constrains such as the friedel s sum rule , positions of the hubbard bands as well as the value of quasiparticle residue are used to establish the equations for the coefficients of the interpolation . 
   we employ two fast techniques , the slave  boson mean  field and the hubbard i approximations to determine the functional dependence of the coefficients on doping , degeneracy and the strength of the interaction . 
   the obtained spectral functions and self  energies are in good agreement with the results of numerically exact quantum monte carlo method .