Solving a Single Impurity Anderson model

In this section we use the methods in EDIpack2.0 to solve a simple example of Anderson quantum impurity problem.

Looking forward for a DMFT application, here we consider the Bethe lattice DOS \(\rho(x)=\frac{1}{2D}\sqrt{D^2-x^2}\) and build the corresponding non-interacting Green’s function \(G_0(z) = \int_{\mathbb{R}}\rho(\epsilon)\left[ z -\epsilon \right]^{-1}\).

We construct a discretized bath by fitting such function on the Matsubara frequencies \(G_0(i\omega_n)\) using the methods in \chi^2 Fit. Finally we input this bath into the ed_solve() solver of EDIpack2.0 in presence of local interaction on the impurity.

The initialization of the code is:

program lancED
   USE EDIPACK2
   USE SCIFOR
   implicit none
   ! Bethe lattice half-bandwidth = energy unit
   real(8),parameter                       :: D=1d0
   ! Bath size and Nso=Nspin*Norb (here =1)
   integer                                 :: Nb,Nso
   ! User bath, allocatable see below
   real(8),allocatable                     :: Bath(:)
   ! Non-interacting Bethe lattice Green's function (naming
   ! convention will be clear in the following section)
   complex(8),allocatable,dimension(:,:,:) :: Weiss

   !> READ THE input using EDIpack procedure:
   call ed_read_input('inputED.conf')

where we load both the EDIpack2.0 and SciFortran libraries through their main module edipack2 and scifor. We also define some local variables and proceed with reading the input file "inputED.conf" using the EDIpack2.0 function ed_read_input().

Next we construct the non-interacting Green’s function, using procedures available in SciFortran:

!> Get the Bethe lattice non-interacting Matsubara GF as a guess for the bath
allocate(Weiss(Nso,Nso,Lmats))
call bethe_guess_g0(Weiss(1,1,:),D,beta,hloc=0d0)

Then we initialize the solver. This step requires the user to pass the user bath as a rank-1 double precision array of a given size. The correct size of the bath array is evaluated internally by the EDIpack2.0 code through the function ed_get_bath_dimension().

!> Init solver:
Nb=ed_get_bath_dimension()
allocate(bath(Nb))
call ed_init_solver(bath)

Upon initialization the bath is guessed from a flat distribution centered around zero and with half-width ed_hw_bath. Here we update the bath optimizing it against the non-interacting Bethe lattice Green’s function:

!> Fit the bath against G0 guess: the outcome is a bath discretizing the Bethe DOS.
call ed_chi2_fitgf(Weiss,bath,ispin=1,iorb=1)

We are now ready to solve the quantum impurity problem for a given set of parameters specified in the input file (see below)

!> Solve SIAM with this given bath
call ed_solve(bath)

Here is a snapshot of the results obtained for \(U=1.0, 10.0\).

../_images/01_anderson_fig.svg

In the top panel we report the impurity spectral functions \(-\Im G^{im}(\omega)/\pi\) compared to the Bethe density of states (filled curve). In the bottom panel we show the real part of the impurity real-axis self-energy functions \(\Sigma(\omega)\) near the Fermi level \(\omega=0\). The linear fit \(y = A\omega\) gives a direct estimate of the derivatives \(A\simeq \tfrac{\partial\Re\Sigma}{\partial\omega}_{|_{\omega\rightarrow 0}}\) and thus of the quasi-particle renormalization constants \(Z=\left( 1 - A \right)^{-1}\) as reported in the legend.

As a direct comparison we report also the values of \(Z\) estimated from the Matsubara axis using the relation \(\frac{\Im\Sigma(i\omega_n}{\omega_n}_{|_{\omega_n\rightarrow 0}}= \frac{1}{\pi}\int_{\mathbb R}d\epsilon \frac{\Re\Sigma(\epsilon)}{\epsilon^2}= \frac{\partial\Re\Sigma}{\partial\omega}_{|_{\omega\rightarrow 0}}\).

We obtained: \(Z=0.74\) and \(Z=0.002\), respectively, for \(U=1.0\) and \(U=10.0\).

Here is an example of input file used in actual calculations:

NORB=1                                        !Number of impurity orbitals (max 5).
NBATH=9                                       !Number of bath sites
ULOC=1d0
ED_MODE=normal                                !Flag to set ED type
BATH_TYPE=normal                              !flag to set bath type
BETA=1000.000000000                           !Inverse temperature, at T=0 is used as a IR cut-off.
ED_VERBOSE=3                                  !Verbosity level: 0=almost nothing --> 5:all. Really: all
LMATS=4096                                    !Number of Matsubara frequencies.
LFIT=2048
EPS=1.000000000E-02                           !Broadening on the real-axis.
CG_FTOL=1.000000000E-10                       !Conjugate-Gradient tolerance.
CG_NITER=2048                                 !Max. number of Conjugate-Gradient iterations.
ED_TWIN=T
LANC_NGFITER=500                              !Number of Lanczos iteration in GF determination. Number of momenta.
ED_HW_BATH=1.000000000                        !half-bandwidth for the bath initialization: flat in -ed_hw_bath:ed_hw_bath