Non-existence of the Luttinger-Ward functional and misleading convergence of skeleton diagrammatic series for Hubbard-like models

The Luttinger-Ward functional $\Phi[\mathbf{G}]$, which expresses the thermodynamic grand potential in terms of the interacting single-particle Green's function $\mathbf{G}$, is found to be ill-defined for fermionic models with the Hubbard on-site interaction. In particular, we show that the self-energy $\mathbf{\Sigma}[\mathbf{G}] \propto \delta\Phi[\mathbf{G}]/\delta \mathbf{G}$ is not a single-valued functional of $\mathbf{G}$: in addition to the physical solution for $\mathbf{\Sigma}[\mathbf{G}]$, there exists at least one qualitatively distinct unphysical branch. This result is demonstrated for several models: the Hubbard atom, the Anderson impurity model, and the full two-dimensional Hubbard model. Despite this pathology, the skeleton Feynman diagrammatic series for $\mathbf{\Sigma}$ in terms of $\mathbf{G}$ is found to converge at least for moderately low temperatures. However, at strong interactions, its convergence is to the unphysical branch. This reveals a new scenario of breaking down of diagrammatic expansions. In contrast, the bare series in terms of the non-interacting Green's function $\mathbf{G}_0$ converges to the correct physical branch of $\mathbf{\Sigma}$ in all cases currently accessible by diagrammatic Monte Carlo. Besides their conceptual importance, these observations have important implications for techniques based on the explicit summation of diagrammatic series.