Vanishing off-diagonal elements of Green's functions

[szankerhqs]

Calculating Green's function's for interacting clusters with off-diagonal hopping we run into strange behavior of CTHYB. For a hybridization function without off-diagonal elements the off-diagonal elements of Green's functions after solving are zero although the local Hamiltonian supports off-diagonal hopping.

A simple example is the 'complex_Gtau_ED.py' test. Setting the off-diagonal elements in the hybridization to zero this behavior occurs.
Changing the original Hamiltonian in the test (first two orbitals are the impurity)
H_mat = np.array([[-0.2 , -1.4 , 0.5 , 0.4 ], [-1.4 , -0.3 , 0.4 , 0.5 ], [ 0.5 , 0.4 , 0.1 , 0.0 ], [ 0.4 , 0.5 , 0.0 , 0.0 ]])
to
H_mat = np.array([[-0.2 , -1.4 , 0.5 , 0.4 ], [-1.4 , -0.3 , 0.4 , 0.5 ], [ 0.5 , 0.4 , 0.1 , 0.0 ], [ 0.4 , 0.5 , 0.0 , 0.0 ]])
and solving for the impurity Green's function, the off diagonal elements S.G_iw['ud'][0,1] vanish (see code below for working example).

import numpy as np
from pytriqs.operators import *
from pytriqs.applications.impurity_solvers.cthyb import *
from pytriqs.gf.local import *
from pytriqs.archive import HDFArchive
from pytriqs.utility.h5diff import h5diff

# set solver parameters
p = {}
p["random_seed"] = 123 
p["length_cycle"] = 100
p["n_warmup_cycles"] = 10000
p["n_cycles"] = 50000

# uncomment to include orbital flips in global moves
gm = {}
gm["swap_orbs"]  = {("ud",0) : ("ud",1), ("ud",1) : ("ud",0)}
p["move_global"] = gm
p["move_global_prob"] = 0.6

U = 0.5
H_int = U*n("ud",0)*n("ud",1)
#
# hybridization with off-diagonal elements ~ 0.4
#
# the Hamiltonian of the system
# the first two orbitals are the impurity, the other two are the bath
H_mat = np.array([[-0.2  , -1.4 , 0.5  , 0.04 ],
                  [-1.4  , -0.3 , 0.04 , 0.5 ],
                  [ 0.5  , 0.04 , 0.1  , 0.0 ],
                  [ 0.04 , 0.5  , 0.0  , 0.0 ]])
corr_dim = 2
beta = 10.0

G0_iw = GfImFreq(beta=beta,indices=range(len(H_mat)),n_points=101)
G0_iw << inverse(iOmega_n - H_mat)


S = Solver(beta=beta,gf_struct={"ud":range(corr_dim)},n_tau=203,n_iw=101)
S.G0_iw << G0_iw[:corr_dim,:corr_dim]
S.solve(h_int=H_int,**p)


#
# hybridization without off-diagonal elements
#

# the Hamiltonian of the system
# the first two orbitals are the impurity, the other two are the bath
H_mat = np.array([[-0.2  ,-1.4  , 0.5 ,  0.0 ],
                  [-1.4  ,-0.3  , 0.0 ,  0.5 ],
                  [ 0.5  , 0.0  , 0.1 ,  0.0 ],
                  [ 0.0  , 0.5  , 0.0 ,  0.0 ]])

G0_iw = GfImFreq(beta=beta,indices=range(len(H_mat)),n_points=101)
G0_iw << inverse(iOmega_n - H_mat)


S2 = Solver(beta=beta,gf_struct={"ud":range(corr_dim)},n_tau=203,n_iw=101)
S2.G0_iw << G0_iw[:corr_dim,:corr_dim]
S2.solve(h_int=H_int,**p)

# plot the off-diagonal elements of G_iw.
from pytriqs.plot.mpl_interface import oplot, plt
fig = plt.figure
oplot(S.G_iw['ud'][0,1], '-o', label=r'$\Delta$ with off-diag', mode='R')
oplot(S2.G_iw['ud'][0,1], '-o', label=r'$\Delta$ without off-diag', mode='R')

[krivenko]

What you describe is a fundamental limitation of the hybridization expansion QMC algorithm in general. Or, at least, of the sampling/measurement method most implementations use.

Please, have a look at the equation (9) of our paper: https://arxiv.org/pdf/1507.00175.pdf
This is the accumulation formula for G(\tau).

Non-vanishing contributions to G_{\alpha\beta}(\tau), where \alpha\neq\beta, involve a matrix element M_{ij} such that \alpha=\alpha_i and \beta=\alpha'_j.
If the hybridization function is diagonal, no elements with said property can exists in matrix M. Therefore, there is no way for an off-diagonal element of G_{\alpha\beta}(\tau) to obtain a finite contribution.

@HugoStrand Should we add this question to the FAQ?