Smoothing out Green's function and correlation function integrals

Question

Often in quantum physics (for Green's functions, and response or correlation functions) we have integrals with terms such as $$ \frac{1}{i\eta+\varepsilon_{k+q}-\varepsilon_k} $$ which we must integrate over $k$ such that one state is above the Fermi surface and one is below it. For $q<2k_F$, this includes excitations arbitrarily close to $\varepsilon_{k+q}$ and $\varepsilon_k$ can be arbitrarily close.

I need to integrate $\frac{1}{i\eta+\varepsilon_{k+q}-\varepsilon_k}$ numerically over a grid of $k$-points for very small $\eta$. If I am not careful, the numerical error will be substantial. I cannot simply make $\eta$ large enough to broaden the Lorentzian peaks of the integrand because in my model $\eta$ is not a numerical parameter but rather has physical meaning.

Increasing the number of $k$-points in my integration grid helps, but is very costly and its benefit is still limited. What are some other useful ways to handle badly-behaved integrands like this?

Answer 1

Since the number of $q$-points with $q<2k_F$ is relatively small, I was able to evaluate most points ($q>2k_F$) with a rough $k$-point grid, and then the problematic few ($q<2k_F$) with a higher-resolution $k$-point grid. Although still a somewhat costly brute force approach, it works and is much more computationally feasible than evaluating all $q$-points on the high-resolution $k$-point grid.