Mott physics beyond Brinkman-Rice scenario

TL;DR

This paper enhances the Brinkman-Rice approach to the Mott transition by incorporating high-energy virtual processes via a variational Schrieffer-Wolff transformation, achieving results close to dynamical mean-field theory in a Bethe lattice.

Contribution

It introduces a novel method combining Brinkman-Rice and Schrieffer-Wolff transformation to better capture high-energy processes in Mott transitions.

Findings

The method accurately reproduces the Mott transition in the Bethe lattice.

Results closely match those from dynamical mean-field theory.

The approach offers potential for studying more complex models.

Abstract

The main flaw of the well-known Brinkman-Rice description, obtained through the Gutzwiller approximation, of the paramagnetic Mott transition in the Hubbard model is in neglecting high-energy virtual processes that generate for instance the antiferromagnetic exchange $J\sim t^2/U$. Here we propose a way to capture those processes by combining the Brinkman-Rice approach with a variational Schrieffer-Wolff transformation, and apply this method to study the single-band metal-to-insulator transition in a Bethe lattice with infinite coordination number, where the Gutzwiller approximation becomes exact. We indeed find for the Mott transition a description very close to the real one provided by dynamical mean-field theory; an encouraging result in view of possible applications to more involved models.

Equations

$$H=-\frac{t}{\sqrt{z}}\sum_{\langle{\bf i}{\bf j}\rangle}T_{{\bf i}{\bf j}}+\frac{\displaystyle U}{\displaystyle 2}\sum_{{\bf i}}\big{(}n_{{\bf i}\uparrow}+n_{{\bf i}\downarrow}-1\big{)}^{2},$$

$$\begin{array}{r}{H}_{\mathrm{int}}=\frac{U}{2}\sum _{\mathbf{i}}\sum _n(n-1{)}^2{P}_{\mathbf{i}}(n),\end{array}$$

$$\begin{split}\tilde{T}_{{\bf i}{\bf j}}&\equiv\Big{(}P_{{\bf i}}(2)P_{{\bf j}}(0)+P_{{\bf i}}(0)P_{{\bf j}}(2)\Big{)}T_{{\bf i}{\bf j}}\Big{(}P_{{\bf i}}(1)P_{{\bf j}}(1)\Big{)},\\ \tilde{T}_{{\bf i}{\bf j}}^{\dagger}&\equiv\Big{(}P_{{\bf i}}(1)P_{{\bf j}}(1)\Big{)}T_{{\bf i}{\bf j}}\Big{(}P_{{\bf i}}(2)P_{{\bf j}}(0)+P_{{\bf i}}(0)P_{{\bf j}}(2)\Big{)}.\end{split}$$

$$\mathcal{U}(ϵ)=\exp \left(\frac{ϵ}{\sqrt{z}}S\right),$$

$$S=\sum_{\langle{\bf i}{\bf j}\rangle}S_{{\bf i}{\bf j}}=\sum_{\langle{\bf i}{\bf j}\rangle}\big{(}\tilde{T}_{{\bf i}{\bf j}}-\tilde{T}_{{\bf i}{\bf j}}^{\dagger}\big{)}.$$

$$\begin{split}\mathcal{H}=&\mathcal{U}(\epsilon)^{\dagger}H\mathcal{U}(\epsilon)\simeq H-\frac{\epsilon}{\sqrt{z}}[S,H]\\ &+\frac{\epsilon^{2}}{2z}\big{[}S,[S,H]\big{]}-\frac{\epsilon^{3}}{6z^{3/2}}\Big{[}S,\big{[}S,[S,H]\big{]}\!\Big{]}\,,\end{split}$$

$$\begin{split}\mathcal{H}&\simeq H+\epsilon\bigg{(}\frac{t}{z}\sum_{{\bf i}{\bf j}{\bf k}{\bf l}}[S_{{\bf i}{\bf j}},T_{{\bf k}{\bf l}}]+\frac{U}{\sqrt{z}}\sum_{{\bf i}{\bf j}}\mathbb{T}_{{\bf i}{\bf j}}\bigg{)}\\ &-\frac{\epsilon^{2}}{2}\bigg{(}\!\frac{t}{z^{3/2}}\!\sum_{\begin{subarray}{c}{\bf i}{\bf j}{\bf k}\\ {\bf l}{\bf m}{\bf n}\end{subarray}}\!\big{[}S_{{\bf i}{\bf j}},[S_{{\bf k}{\bf l}},T_{{\bf m}{\bf n}}]\big{]}+\frac{U}{z}\sum_{{\bf i}{\bf j}{\bf k}{\bf l}}[S_{{\bf i}{\bf j}},\mathbb{T}_{{\bf k}{\bf l}}]\!\bigg{)}\\ &+\frac{\epsilon^{3}}{6}\bigg{(}\frac{t}{z^{2}}\!\!\sum_{\begin{subarray}{c}{\bf i}{\bf j}{\bf k}{\bf l}\\ {\bf m}{\bf n}{\bf p}{\bf q}\end{subarray}}\!\Big{[}S_{{\bf i}{\bf j}},\big{[}S_{{\bf k}{\bf l}},[S_{{\bf m}{\bf n}},T_{{\bf p}{\bf q}}]\big{]}\!\Big{]}\\ &\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ +\frac{U}{z^{3/2}}\sum_{\begin{subarray}{c}{\bf i}{\bf j}{\bf k}\\ {\bf l}{\bf m}{\bf n}\end{subarray}}\big{[}S_{{\bf i}{\bf j}},[S_{{\bf k}{\bf l}},\mathbb{T}_{{\bf m}{\bf n}}]\big{]}\bigg{)}\,,\end{split}$$

$$\Big{[}\sum_{{\bf i}{\bf j}}S_{{\bf i}{\bf j}},H_{\rm int}\Big{]}=-U\sum_{{\bf i}{\bf j}}\mathbb{T}_{{\bf i}{\bf j}}.$$

$$\mid {\psi }_G⟩\equiv \prod _{\mathbf{i}}{\mathcal{P}}_{\mathbf{i}}\mid {\psi }_0⟩.$$

$$\mathcal{P}_{{\bf i}}=\sqrt{2}\Big{(}\sin\theta\ \!P_{{\bf i}}(0)+\cos\theta\ \!P_{{\bf i}}(1)+\sin\theta\ \!P_{{\bf i}}(2)\Big{)}\,,$$

$$\mid \Psi ⟩=\mathcal{U}(ϵ)\mid {\psi }_G⟩,$$

$$8T_{0}N\mathcal{F}(ϵ,\theta )=⟨\Psi \mid H\mid \Psi ⟩=⟨{\psi }_G\mid \mathcal{H}\mid {\psi }_G⟩,$$

$$\begin{split}d\equiv&\frac{1}{N}\,\langle\psi_{G}\mid\mathcal{U}(\epsilon)^{\dagger}\,\Big{(}\sum_{{\bf i}}\,P_{{\bf i}}(2)\Big{)}\,\mathcal{U}(\epsilon)\mid\psi_{G}\rangle\,,\end{split}$$

$$\begin{split}\mathcal{F}(\epsilon,\theta)&\simeq\frac{1}{8}\Bigg{[}-\Big{(}1-\frac{\epsilon u}{2\tau}\Big{)}\sin^{2}2\theta+2u\Big{(}1-\cos 2\theta\Big{)}-8\epsilon\,\tau\,\bigg{(}1-\frac{\epsilon u}{2\tau}\bigg{)}\,\cos 2\theta\\ &-\frac{3\,\epsilon}{8\,\tau}\Bigg{(}1-\frac{\epsilon u}{4\tau}\Bigg{)}\sin^{2}2\theta\,\cos 2\theta+\epsilon^{2}\Bigg{(}\frac{9}{4}-\frac{\epsilon u}{2\tau}\Bigg{)}\sin^{2}2\theta+\frac{32\,\epsilon^{3}\tau}{3}\cos 2\theta\\ &-\frac{\epsilon^{2}}{128\,\tau^{2}}\bigg{(}1-\frac{\epsilon u}{6\tau}\bigg{)}\bigg{(}-\frac{7}{2}\,\sin^{4}2\theta+9\,\sin^{2}2\theta\,\cos^{2}2\theta\bigg{)}+\frac{155\,\epsilon^{3}}{192\,\tau}\sin^{2}2\theta\,\cos 2\theta\Bigg{]}\,,\end{split}$$

$$u=\frac{U}{8T_0},\tau =\frac{t}{8T_0}=\frac{3\pi }{64}.$$

$$\begin{split}d(\epsilon,\theta&)=\frac{1}{8}\Bigg{[}2-\bigg{(}2-4\epsilon^{2}\bigg{)}\cos 2\theta\\ &+\frac{\epsilon}{2\tau}\big{(}1-\epsilon^{2}\big{)}\sin^{2}2\theta+\frac{3\,\epsilon^{2}}{32\,\tau^{2}}\sin^{2}2\theta\,\cos 2\theta\\ &+\frac{\epsilon^{3}}{768\tau^{3}}\bigg{(}-\frac{7}{2}\sin^{4}2\theta+9\sin^{2}2\theta\,\cos^{2}2\theta\bigg{)}\Bigg{]}\,.\end{split}$$

$$\partial \mathcal{F}/\partial ϵ=\partial \mathcal{F}/\partial \theta =0,$$

$$\epsilon_{\rm ins}=\frac{\tau}{u}-4\Big{(}\frac{\tau}{u}\Big{)}^{3}+\mathcal{O}\left(\frac{\tau^{5}}{u^{5}}\right),$$

$$E_{\rm ins}=-8T_{0}\Big{(}\frac{\tau^{2}}{2u}-\frac{4\tau^{4}}{3u^{3}}\Big{)}+\mathcal{O}\Big{(}\frac{\tau^{6}}{u^{5}}\Big{)}\simeq-\frac{t^{2}}{2U}+\frac{4t^{4}}{3U^{3}}.$$

$$d_{\rm ins}=\frac{1}{2}\Big{(}\frac{\tau}{u}\Big{)}^{2}-4\Big{(}\frac{\tau}{u}\Big{)}^{4}+\mathcal{O}\Big{(}\frac{\tau^{6}}{u^{6}}\Big{)},$$

Full text

A metal to insulator transition driven by the electron-electron repulsion was envisioned by Mott more than fifty years ago Mott (1961). Since then, the underlying physics of this phenomenon has been studied by large variety of quantum many-body tools in models for strongly correlated systems Brinkman and Rice (1970); Georges et al. (1996); Imada et al. (1998); Capello et al. (2005).

One of the earliest microscopic descriptions of the Mott localisation is owned to Brinkman and Rice Brinkman and Rice (1970), and obtained through the Gutzwiller approximation applied to the half-filled Hubbard model. In their scenario the transition to the insulating state occurs when the hopping is fully hampered by repulsion, i.e. its expectation value in the variational wavefunction strictly vanishes. This result is elegant in many ways. It is fully analytical and provides a very intuitive and physically transparent, almost classical, interpretation of the Mott phenomenon.

Nonetheless, this description, frequently called Brinkman–Rice transition, has a severe drawback: the expectation value of the hopping cannot be zero, and it is so in the Gutzwiller approximation only because there is a complete, static and dynamic, locking of charge degrees of freedom. In reality, dynamical charge fluctuations do play a role even deep in the Mott phase, and in particular they mediate the antiferromagnetic spin-exchange, as clear by the large $U$ mapping onto the Heisenberg model that can be formally derived through the Schrieffer-Wolff transformation Chao et al. (1977).

Since the result of Brinkman and Rice, a variety of quantum many-body tools have been constructed that are generically able to sensibly capture those dynamical processes Georges et al. (1996); Capello et al. (2005); Schollwöck (2005); Toschi et al. (2007). In particular, when applied to the Hubbard model, they provide satisfying descriptions of the Mott transition, though relying on heavy numerical computation. Nowadays, the scenario provided by dynamical mean-field theory (DMFT) Georges et al. (1996), which becomes exact in infinite dimensions Metzner and Vollhardt (1989), has become an invaluable benchmark to compare with.

In the present work we revisit the problem of the Brinkman–Rice transition, and complement it with the inclusion of the dynamical processes in a semi-analytic manner. In order to achieve this goal, we construct a method that combines the Gutzwiller’s variational approach with a variational Schrieffer-Wolff transformation. As a case study, we apply our technique to the half-filled Hubbard model in the paramagnetic phase on the infinitely coordinated Bethe lattice. The energy functional to be minimised can be obtained fully analytically. Its minimisation leads to a significantly improved description of the Mott transition as compared to the standard Brinkman–Rice scenario, and much closer to the exact DMFT one Georges et al. (1996). The improvement is in particular highlighted in: (i) a sizeable lowering of the critical interaction strength for a transition; (ii) a lower value of the insulator energy that includes a non-zero expectation value of the hopping $\sim-t^{2}/U$; and (iii) a proper balance of kinetic and potential energies at the transition.

The starting point of our analysis is the half-filled single-band Hubbard model on the infinitely coordinated Bethe lattice,

[TABLE]

where ${z\to\infty}$ is the coordination number of the Bethe lattice. ${T_{{\bf i}{\bf j}}=T_{{\bf j}{\bf i}}\equiv\sum_{\sigma}(c_{{\bf i}\sigma}^{\dagger}c_{{\bf j}\sigma}^{\phantom{\dagger}}+c_{{\bf j}\sigma}^{\dagger}c_{{\bf i}\sigma}^{\phantom{\dagger}})}$ is the hermitian hopping operator between neighbouring sites $i$ and $j$, and $n_{{\bf i}\sigma}=c_{{\bf i}\sigma}^{\dagger}c_{{\bf i}\sigma}^{\phantom{\dagger}}$ the local density of spin ${\sigma=\uparrow,\downarrow}$ electrons. We rewrite the interaction, last term in Eq. (1), as

[TABLE]

where $P_{{\bf i}}(n)$ is the projector at site ${\bf i}$ onto the subspace with $n$ electrons.

In order to construct a partial Schrieffer-Wolff transformation Chao et al. (1977) that accounts for not complete projection of double occupancies, we separately define components of the hopping operator, $T_{{\bf i}{\bf j}}$ projected on the right or on the left onto the configurations where both sites ${\bf i}$ and ${\bf j}$ are singly occupied,

[TABLE]

Their sum is gathered under the form of the new operator $\mathbb{T}_{{\bf i}{\bf j}}\equiv\tilde{T}_{{\bf i}{\bf j}}+\tilde{T}_{{\bf i}{\bf j}}^{\dagger}$, while the remaining part of the hopping operator under $\mathcal{T}_{{\bf i}{\bf j}}\equiv T_{{\bf i}{\bf j}}-\mathbb{T}_{{\bf i}{\bf j}}$. We construct the partial Schrieffer-Wolff transformation

[TABLE]

through the anti-hermitian operator

[TABLE]

The transformed Hamiltonian reads,

[TABLE]

where $\epsilon$ is variationally determined so as to minimise the energy. In the following, we shall assume that for any value of $U$ the optimal $\epsilon$ is small enough to safely neglect higher order terms, $\mathcal{O}(\epsilon^{4})$ in Eq. (6). A posteriori, we shall check the validity of such assumption. We rewrite the transformed Hamiltonian in a more useful form,

[TABLE]

where we made use of the following equality

[TABLE]

The transformed low energy Hamiltonian (7) is then analysed by a variational approach. Specifically, the ground state of $\mathcal{H}$ is approximated by a variational Gutzwiller wave function $|\psi_{G}\rangle$ constructed from the uncorrelated Fermi sea $|\psi_{0}\rangle$ through

[TABLE]

$\mathcal{P}_{{\bf i}}$ is a linear operator that, in the presence of particle-hole symmetry, can be parametrised as

[TABLE]

where $\theta$ is a variational parameter bounded by ${\theta\in\{0,\pi/4\}}$, where $\theta=\pi/4$ corresponds to the uncorrelated (metallic) state, whereas $\theta=0$ projects out of $|\psi_{0}\rangle$ all configurations with doubly occupied and empty sites. In other words, the actual variational wavefunction for the ground state of original Hamiltonian $H$ is

[TABLE]

and depends both on $\theta$ and $\epsilon$. The variational energy functional per lattice site (where $N$ is the total number of sites), $\mathcal{F}(\epsilon,\theta)$ can be now obtained as the expectation value in the Gutzwiller wave function (9) of the Hamiltonian (7), which can be analytically computed in the infinitely coordinated Bethe lattice,

[TABLE]

in reduced units of $8T_{0}$. Here $-T_{0}$ is the hopping energy per site of $|\psi_{0}\rangle$, which in a Bethe lattice reads ${T_{0}=8t/3\pi}$.

Already at this point, qualitative differences with respect to the standard Brinkman–Rice transition emerge clearly. In our approach $\theta=0$ providing vanishing double occupancies in $|\psi_{G}\rangle$ does not yield the same for the actual wavefunction $|\Psi\rangle$. Explicitly, the density of doubly occupied sites, $d$, can be calculated as

[TABLE]

which generically does not provide $d\!=\!0$ even if ${\theta\!=\!0}$.

In order to evaluate (12) as well as (13) we apply Wick’s theorem. Each expectation value resulting from this procedure can be conveniently visualised by a diagram with nodes denoting sites and edges being averages of the inter-site single particle density matrix. Sum of diagrams with the same number $x$ of nodes we shall shortly denote as an $x$-vertex. We checked that a satisfying accuracy is obtained by keeping all $x$-vertices up to $x=4$.

The resulting energy functional reads

[TABLE]

where interaction strength, $U$ and hopping amplitude, $t$ are rescaled as

[TABLE]

The first line in Eq. (14) includes the 2-vertex contribution, the second line is the 3-vertex one, and finally the third is the 4-vertex correction. Additionally, the expectation value of the double occupancy $d$ reads

[TABLE]

From (14) and (16) we can easily recover the results of the standard Gutzwiller approximation applied to the Hubbard model by setting $\epsilon=0$. In this case the Brinkman–Rice transition takes place for $u_{\rm BR}=1$ and the insulating state is characterised by $d=0$.

We search for minima of the functional $\mathcal{F}$ with increasing $u$ by standard methods. Namely, for each $u$ we look for the pair of variables $\{\epsilon,\theta\}$ satisfying

[TABLE]

under the condition that the Hessian is positive definite. We start observing that for $u=0$, the minimum of the functional (14) is correctly determined by $\epsilon=0$ and ${\theta=\pi/4}$ that correspond to fully uncorrelated metal. For interaction strength roughly up to $u\simeq 0.4$, the optimised energy is almost coincident with that obtained either by Gutzwiller approximation or by DMFT.

For stronger correlations, $u\gtrsim 0.4$, the Gutzwiller approximation starts to deviate appreciably with respect to DMFT, while our variational energy remains quite close. In Fig. 1(a) we plot the total energy, as well as separately kinetic and potential energies, of the minimum of functional Eq. (14), as compared with DMFT Weber et al. (2012), and with the sole Gutzwiller approximation, for which we just show the total energy.

Following Brinkman and Rice Brinkman and Rice (1970), we associate the Mott insulating state with $\theta=0$, which is always a saddle point of the functional (14). However, this saddle point becomes minimum only when metal becomes unstable; the metal to insulator transition is thus continuous and occurs at a critical interaction, $u_{c}\simeq 0.822$, which is sizeably lower than the Brinkman-Rice value, $u_{\rm BR}=1$, and quite close to DMFT, $u_{\rm c2;DMFT}\simeq 0.854$. In Fig. 1(b) we show the values of the variational parameters $\epsilon$ and $\theta$ on the both sides of the transition. In the same figure, we also plot the average double occupancy $d$ (from Eq. (16)), which is non-zero in the insulating phase and decreases almost linearly in the metallic state.

In the insulating phase, $u\gtrsim u_{c}$ when the optimal $\theta=0$, we can analytically calculate several quantities. For instance, the saddle point value of $\epsilon$ can be obtained in power series of $\tau/u$:

[TABLE]

whereas the energy per site is

[TABLE]

Additionally, the average double occupancy in powers of $\tau/u$ reads

[TABLE]

which is indeed finite.

Let us now compare more in detail the above results with the exact DMFT ones  Georges et al. (1996). Alike DMFT, we find continuous metal-insulator transition for quite similar critical interaction $U$. However, in our case there is no coexistence region of the insulating and metallic solutions, which in DMFT spreads over significant region ($u_{\rm c1}$ and $u_{\rm c2}$ obtained by DMFT are marked in Fig. 1(b)). In spite of this deficiency, we do find an energy balance across the transition close to DMFT, and quite different from the Gutzwiller approximation. Indeed, hopping and potential energies have kinks at the transition, though the total energy is smooth, and the insulator energy at the leading order scales as $\sim-t^{2}/2U$ (cf. Fig.1(a)).

In summary, we have analysed a very simple variational wavefunction for a correlated system that consists of a Gutzwiller wavefunction combined with a variational Schrieffer-Wolff transformation. We have benchmarked this wavefunction against the exact DMFT results  Georges et al. (1996) for the paramagnetic Mott transition in the half-filled single-band Hubbard model on a Bethe lattice with infinite coordination number. Although there are obviously differences with exact results, nevertheless our variational wavefunction provides a description of the Mott transition much closer to reality than the Brinkman-Rice scenario. More importantly, our wavefunction is able to portray a Mott insulator where charge fluctuations are not completely suppressed as in the Brinkman-Rice scenario, and which therefore has a non-zero expectation value of the hopping. This variational technique might open new possibilities to access Mott physics or related phenomena in more realistic models with minimal computational effort.

Acknowledgements. We are grateful to Adriano Amaricci for providing us data by DMFT. MMW acknowledges support from Polish Ministry of Science and Higher Education under the “Mobilność Plus” program, Agreement No. 1265/MOB/IV/2015/0. MF acknowledges support from European Union under the H2020 Framework Programme, ERC Advanced Grant No. 692670 “FIRSTORM”.