1. Introduction
The 1D Hubbard model [1] is one of the most important solvable models in non-perturbative quantum field theory [2]. It exhibits on-site Coulomb interaction and correlated hopping, which helps us to understand the mystery of high-Tc superconductivity. It is a paradigm of integrability in the strongly correlated systems.
In the past several decades, numerous approaches have been proposed to study the integrability and the physical properties of the 1D Hubbard model [3–12]. The Hubbard model with a periodic boundary condition was first exactly solved via the coordinate Bethe ansatz method [13, 14]. Shastry then constructed the corresponding R-matrix and the Lax matrix, and demonstrated the integrability of the 1D Hubbard model [15, 16]. The Hamiltonian of the conventional Hubbard model can be constructed by taking the derivation of the logarithm of the quantum transfer matrix at u = 0, {θm = 0}. Martins and his co-workers subsequently gave the solution of the conventional Hubbard model via the nested algebraic Bethe ansatz approach [17].
Our starting point is the construction of an extended 1D Hubbard model. We let all the inhomogeneous {θm} in the transfer matrix take the same nonzero value θ, i.e. u = θ, {θm = θ}. Then, the derivative of the logarithm of the quantum transfer matrix t(u) at u = θ gives another integrable Hamiltonian. This model depends on more free parameters. Compared to the conventional Hubbard model, the new model contains more possible nearest-neighbor interactions. Following the nested algebraic Bethe ansatz method, we solve the extended Hubbard model exactly. The T − Q relation and a set of Bethe ansatz equations (BAEs) are proposed.
This paper is organized as follows. In section 2 , we construct an integrable 1D extended Hubbard model. In section 3 , we formulate the nested algebraic Bethe ansatz for the extended Hubbard model and present our main results. Section 4 is devoted to the conclusion.
2. 1D extended Hubbard model
Let us recall the formulation of the integrability of the 1D Hubbard model [16]. The quantum R-matrix is given by [15],1 ) indeed satisfies the Yang–Baxter equation [18]:5 ) satisfies the RTT relation:
$\begin{eqnarray}\begin{array}{l}R(u,v)=\left(\begin{array}{cccccccccccccccc}{a}^{+} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & {b}^{+} & 0 & 0 & e & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & {b}^{+} & 0 & 0 & 0 & 0 & 0 & e & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & {c}^{+} & 0 & 0 & f & 0 & 0 & f & 0 & 0 & {d}^{+} & 0 & 0 & 0\\ 0 & e & 0 & 0 & {b}^{-} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & {a}^{-} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & f & 0 & 0 & {c}^{-} & 0 & 0 & {d}^{-} & 0 & 0 & f & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & {b}^{-} & 0 & 0 & 0 & 0 & 0 & e & 0 & 0\\ 0 & 0 & e & 0 & 0 & 0 & 0 & 0 & {b}^{-} & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & f & 0 & 0 & {d}^{-} & 0 & 0 & {c}^{-} & 0 & 0 & f & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & {a}^{-} & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & {b}^{-} & 0 & 0 & e & 0\\ 0 & 0 & 0 & {d}^{+} & 0 & 0 & f & 0 & 0 & f & 0 & 0 & {c}^{+} & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & e & 0 & 0 & 0 & 0 & 0 & {b}^{+} & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & e & 0 & 0 & {b}^{+} & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & {a}^{+}\end{array}\right),\end{array}\end{eqnarray}$
where, $\begin{eqnarray}\begin{array}{l}{a}^{\pm }(u,v)=\cos (u+v){\cos }^{2}(u-v)\cosh ({h}_{1}-{h}_{2})\pm {\cos }^{2}(u+v)\cos (u-v)\sinh ({h}_{1}-{h}_{2}),\\ {b}^{\pm }(u,v)=\sin (u-v)\cos (u-v)\cos (u+v)\cosh ({h}_{1}-{h}_{2})\pm \sin (u+v)\cos (u-v)\cos (u+v)\sinh ({h}_{1}-{h}_{2}),\\ {c}^{\pm }(u,v)={\sin }^{2}(u-v)\cos (u+v)\cosh ({h}_{1}-{h}_{2})\pm {\sin }^{2}(u+v)\cos (u-v)\sinh ({h}_{1}-{h}_{2}),\\ {d}^{\pm }(u,v)=\cos (u+v)\cosh ({h}_{1}-{h}_{2})\pm \cos (u-v)\sinh ({h}_{1}-{h}_{2}),\\ e(u,v)=\cos (u-v)\cos (u+v),\\ f(u,v)=\displaystyle \frac{\sin (u-v)}{\cosh ({h}_{1}+{h}_{2})}\cosh ({h}_{1}-{h}_{2})\cos (u+v),\end{array}\end{eqnarray}$
and functions h1 ≡ h(u), h2 ≡ h(v) are assumed to satisfy the constraint: $\begin{eqnarray}\begin{array}{l}\displaystyle \frac{U}{4}=\displaystyle \frac{\sinh 2h(u)}{\sin 2u}=\displaystyle \frac{\sinh 2h(v)}{\sin 2v}.\end{array}\end{eqnarray}$
Wadati proved that the R-matrix in ( $\begin{eqnarray}\begin{array}{l}{R}_{12}(u,v){R}_{13}(u,w){R}_{23}(v,w)={R}_{23}(v,w){R}_{13}(u,w){R}_{12}(u,v).\end{array}\end{eqnarray}$
We construct the monodromy matrix: $\begin{eqnarray}\begin{array}{l}{T}_{0}(u)={R}_{0N}(u,{\theta }_{N})\cdots {R}_{01}(u,{\theta }_{1}),\end{array}\end{eqnarray}$
where {θ1,…,θN} are inhomogeneous parameters. T0(u) in equation ( $\begin{eqnarray}\begin{array}{l}{R}_{12}(u,v){T}_{1}(u){T}_{2}(v)={T}_{2}(v){T}_{1}(u){R}_{12}(u,v).\end{array}\end{eqnarray}$
The transfer matrix is thus: $\begin{eqnarray}t(u)={{tr}}_{0}\{{T}_{0}(u)\},\end{eqnarray}$
which has the commutative property: $\begin{eqnarray}\begin{array}{l}[t(u),t(v)]=0.\end{array}\end{eqnarray}$
In the homogeneous limit {θm = θ}, the derivation of the logarithm of the transfer matrix at u = θ gives the following Hamiltonian:
$\begin{eqnarray}\begin{array}{rcl}H & = & {\left.\displaystyle \frac{\partial \mathrm{ln}t(u)}{\partial u}\right|}_{u=\theta ,\{{\theta }_{m}=\theta \}}+N\tan 2\theta \\ & = & \displaystyle \sum _{j=1}^{N}\left\{\displaystyle \frac{1}{2}(1+{\rm{sech}} \,2h(\theta ))\right.\,\left.({\sigma }_{j}^{+}{\sigma }_{j+1}^{-}+{\sigma }_{j}^{-}{\sigma }_{j+1}^{+}+{\tau }_{j}^{+}{\tau }_{j+1}^{-}+{\tau }_{j}^{-}{\tau }_{j+1}^{+})\right.\\ & & +\displaystyle \frac{1}{2}(1-{\rm{sech}} \,2h(\theta ))[{\sigma }_{j}^{z}{\sigma }_{j+1}^{z}({\tau }_{j}^{+}{\tau }_{j+1}^{-}+{\tau }_{j}^{-}{\tau }_{j+1}^{+})+({\sigma }_{j}^{+}{\sigma }_{j+1}^{-}+{\sigma }_{j}^{-}{\sigma }_{j+1}^{+}){\tau }_{j}^{z}{\tau }_{j+1}^{z}]\\ & & +\displaystyle \frac{U}{4}{\rm{sech}} \,2h(\theta )\left\{\displaystyle \frac{1}{2}({\cos }^{2}2\theta +1){\sigma }_{j}^{z}{\tau }_{j}^{z}+\displaystyle \frac{1}{4}({\cos }^{2}2\theta -1)({\sigma }_{j}^{z}{\tau }_{j+1}^{z}+{\sigma }_{j+1}^{z}{\tau }_{j}^{z})\right.\\ & & +\displaystyle \frac{1}{2}\sin 2\theta \cos 2\theta [({\sigma }_{j}^{z}+{\sigma }_{j+1}^{z})({\tau }_{j}^{-}{\tau }_{j+1}^{+}-{\tau }_{j}^{+}{\tau }_{j+1}^{-})+({\tau }_{j}^{z}+{\tau }_{j+1}^{z})({\sigma }_{j}^{-}{\sigma }_{j+1}^{+}-{\sigma }_{j}^{+}{\sigma }_{j+1}^{-})]\\ & & \left.\left.+{\sin }^{2}2\theta ({\sigma }_{j}^{+}{\sigma }_{j+1}^{-}-{\sigma }_{j}^{-}{\sigma }_{j+1}^{+})({\tau }_{j}^{+}{\tau }_{j+1}^{-}-{\tau }_{j}^{-}{\tau }_{j+1}^{+})\right\}\right\},\end{array}\end{eqnarray}$
where $\{{\sigma }_{j}^{\pm },{\sigma }_{j}^{z}\}$ and $\{{\tau }_{j}^{\pm },{\tau }_{j}^{z}\}$ are two commuting sets of Pauli matrices acting on site j. The periodic boundary condition implies ${\sigma }_{N+1}^{\pm }={\sigma }_{1}^{\pm },{\tau }_{N+1}^{\pm }={\tau }_{1}^{\pm }$.From the constraint in (3 ), one can obviously see that the function h(θ) is determined by θ and U. Therefore, the Hamiltonian depends on two independent parameters θ and U. The Hermitian condition of the Hamiltonian reads as follows:
$\begin{eqnarray}\begin{array}{l}\mathrm{Re}[\theta ]=0,\quad \mathrm{Im}[U]=0,\quad 16-{U}^{2}{\sinh }^{2}2| \theta | \gt 0.\end{array}\end{eqnarray}$
Moreover, in order to relate the coupled spin model in (9 ) to the Hubbard model, we have to perform the following inverse Jordan–Wigner transformation:9 ):
$\begin{eqnarray}\begin{array}{l}{\sigma }_{m}^{-}=\exp \left(-{\rm{i}}\pi \displaystyle \sum _{l=1}^{m-1}({n}_{l\uparrow }-1)\right){c}_{m\uparrow },\\ {\sigma }_{m}^{z}=2{n}_{m\uparrow }-1,\\ {\tau }_{m}^{-}=\exp \left(-{\rm{i}}\pi \displaystyle \sum _{l=1}^{m-1}({n}_{l\downarrow }-1)\right)\exp \left(-{\rm{i}}\pi \displaystyle \sum _{l=1}^{N}({n}_{l\uparrow }-1)\right){c}_{m\downarrow },\\ {\tau }_{m}^{z}=2{n}_{m\downarrow }-1,\end{array}\end{eqnarray}$
where cjσ and ${c}_{j\sigma }^{\dagger }$ are creation and annihilation fermion operators with spins (σ = ↑ , ↓ ) on site j, which satisfy anti-commutation relations$\{{c}_{i\sigma },{c}_{j\sigma ^{\prime} }\}=\{{c}_{i\sigma }^{\dagger },{c}_{j\sigma ^{\prime} }^{\dagger }\}=0$, $\{{c}_{i\sigma }^{\dagger },{c}_{j\sigma ^{\prime} }\}={\delta }_{i,j}{\delta }_{\sigma ,\sigma ^{\prime} }$ and ${n}_{j\sigma }={c}_{j\sigma }^{\dagger }{c}_{j\sigma }$ is the density operator. Using the inverse Jordan–Wigner transformation we can rewrite our Hamiltonian ( $\begin{eqnarray}\begin{array}{rcl}H & = & \displaystyle \sum _{j=1}^{N}{H}_{j,j+1},\\ {H}_{j,j+1} & = & \displaystyle \sum _{\sigma }{c}_{j,\sigma }^{\dagger }{c}_{j+1,\sigma }({\alpha }_{1}+{\alpha }_{3}({n}_{j,-\sigma }+{n}_{j+1,-\sigma })+{\alpha }_{4}{n}_{j,-\sigma }{n}_{j+1,-\sigma })\\ & & \quad +\displaystyle \sum _{\sigma }{c}_{j+1,\sigma }^{\dagger }{c}_{j,\sigma }({\alpha }_{2}+{\alpha }_{5}({n}_{j,-\sigma }+{n}_{j+1,-\sigma })+{\alpha }_{4}{n}_{j,-\sigma }{n}_{j+1,-\sigma })\\ & & \quad +{\alpha }_{6}\displaystyle \sum _{\sigma }({n}_{j,\sigma }{n}_{j+1,-\sigma }+{c}_{j,\sigma }^{\dagger }{c}_{j+1,-\sigma }^{\dagger }{c}_{j,-\sigma }{c}_{j+1,\sigma })+{\alpha }_{7}{n}_{j\uparrow }{n}_{j\downarrow }\\ & & \quad -{\alpha }_{6}({c}_{j\uparrow }^{\dagger }{c}_{j\downarrow }^{\dagger }{c}_{j+1\downarrow }{c}_{j+1\uparrow }+{c}_{j\uparrow }{c}_{j\downarrow }{c}_{j+1\downarrow }^{\dagger }{c}_{j+1\uparrow }^{\dagger })+{\alpha }_{8}\Space{0ex}{0.25ex}{0ex}({n}_{j}-\displaystyle \frac{1}{2}\Space{0ex}{0.25ex}{0ex}),\end{array}\end{eqnarray}$
where the parameters {α1,…,α8} are given by, $\begin{eqnarray}\begin{array}{l}{\alpha }_{1}=-1-\displaystyle \frac{U}{8}\sin 4\theta \,{\rm{sech}} \,2h(\theta ),\\ {\alpha }_{2}=-1+\displaystyle \frac{U}{8}\sin 4\theta \,{\rm{sech}} \,2h(\theta ),\\ {\alpha }_{3}=(1-{\rm{sech}} \,2h(\theta ))+\displaystyle \frac{U}{8}\sin 4\theta \,{\rm{sech}} \,2h(\theta ),\\ {\alpha }_{4}=-2(1-{\rm{sech}} \,2h(\theta )),\\ {\alpha }_{5}=(1-{\rm{sech}} \,2h(\theta ))-\displaystyle \frac{U}{8}\sin 4\theta \,{\rm{sech}} \,2h(\theta ),\\ {\alpha }_{6}=-\displaystyle \frac{U}{4}{\sin }^{2}2\theta \,{\rm{sech}} \,2h(\theta ),\\ {\alpha }_{7}=\displaystyle \frac{U}{2}\,{\rm{sech}} \,2h(\theta )({\cos }^{2}2\theta +1),\\ {\alpha }_{8}=-\displaystyle \frac{U}{2}{\cos }^{2}2\theta \,{\rm{sech}} \,2h(\theta ).\end{array}\end{eqnarray}$
The Hamiltonian (12 ) contains most of the possible nearest-neighbor interactions appearing in strongly correlated systems, e.g. the kinetic energy possessed by particles, the hopping terms that are also included in the conventional Hubbard model, the spin-spin interaction that is the familiar spin-exchange term of the Heisenberg XXX spin chain, and the pair hopping term that relates to the simultaneous hopping of two electrons from one site to a neighboring site.
The interaction intensities {α1,…,α8} all depend on θ and U. For finite θ and U, they are all of the same order of strength, which is clearly illustrated in figure 1.
Figure 1. Left: interaction intensity ∣αk∣ versus θ/i with U = 2.5. Right: interaction intensity ∣αk∣ versus U with θ = 0.5i. |
Compared to Alcaraz’s model [11], whose integrability has not been proved, the model we construct is integrable and Hermitian. Shiroishi presented two integrable Hamiltonians [19] that only depend on one free parameter. While, in this paper, we use a different R-matrix and construct a more general integrable Hamiltonian related to two free parameters θ and U.
The new Hamiltonian in (12 ) reduces to the conventional Hubbard model at θ = 0, namely:
$\begin{eqnarray}\begin{array}{l}{H}_{t}=-\displaystyle \sum _{j=1}^{N}({c}_{j,\sigma }^{\dagger }{c}_{j+1,\sigma }+{c}_{j+1,\sigma }^{\dagger }{c}_{j,\sigma })+U\displaystyle \sum _{j=1}^{N}\Space{0ex}{0.25ex}{0ex}({n}_{j\uparrow }-\displaystyle \frac{1}{2}\Space{0ex}{0.25ex}{0ex})\Space{0ex}{0.25ex}{0ex}({n}_{j\downarrow }-\displaystyle \frac{1}{2}\Space{0ex}{0.25ex}{0ex}).\end{array}\end{eqnarray}$
In conclusion, we construct a more general integrable Hamiltonian via the quantum inverse scattering method(QISM).
3. Exact diagonalization of the transfer matrix
In this section, we expect to diagonalize the transfer matrix and obtain the corresponding Bethe ansatz equations by following the procedure of the nested algebraic Bethe ansatz method [17, 20, 22]. We first represent the monodromy matrix (5 ) in the matrix form:
$\begin{eqnarray}\begin{array}{l}{T}_{0}(u)=\left(\begin{array}{cccc}B(u) & {B}_{1}(u) & {B}_{2}(u) & F(u)\\ {C}_{1}(u) & {A}_{11}(u) & {A}_{12}(u) & {E}_{1}(u)\\ {C}_{2}(u) & {A}_{21}(u) & {A}_{22}(u) & {E}_{2}(u)\\ {C}_{3}(u) & {C}_{4}(u) & {C}_{5}(u) & D(u)\end{array}\right).\end{array}\end{eqnarray}$
The transfer matrix can be expressed by, $\begin{eqnarray}\begin{array}{l}t(u)=B(u)+D(u)+{A}_{11}(u)+{A}_{22}(u).\end{array}\end{eqnarray}$
We introduce the local vacuum state at site j: $\begin{eqnarray}\begin{array}{l}| 0{\rangle }_{j}={\left(\displaystyle \genfrac{}{}{0em}{}{1}{0}\right)}_{\sigma }\otimes {\left(\displaystyle \genfrac{}{}{0em}{}{1}{0}\right)}_{\tau }.\end{array}\end{eqnarray}$
Then, the global vacuum is constructed as, $\begin{eqnarray}\begin{array}{l}| 0\rangle ={\otimes }_{j=1}^{N}| 0{\rangle }_{j}.\end{array}\end{eqnarray}$
The elements of the monodromy matrix T0(u) have the following effect on the reference state ∣0⟩: $\begin{eqnarray}\begin{array}{l}B(u)| 0\rangle =\displaystyle \prod _{m=1}^{N}{a}^{+}(u,{\theta }_{m})| 0\rangle ,\\ {A}_{k,k}(u)| 0\rangle =\displaystyle \prod _{m=1}^{N}{b}^{-}(u,{\theta }_{m})| 0\rangle ,\,\,k=1,2,\\ D(u)| 0\rangle =\displaystyle \prod _{m=1}^{N}{c}^{+}(u,{\theta }_{m})| 0\rangle ,\\ {B}_{k}(u)| 0\rangle \ne 0,\,\,{E}_{k}(u)| 0\rangle \ne 0,\,\,F(u)| 0\rangle \ne 0,\,\,k=1,2,\\ {A}_{12}(u)| 0\rangle ={A}_{21}(u)| 0\rangle =0,\\ {C}_{k}(u)| 0\rangle =0,\,\,k=1,\ldots ,5.\end{array}\end{eqnarray}$
One can see that the total number of particles is conserved and Bk(θ) is a creation operator. The eigenstate of t(u) can thus take the form:
$\begin{eqnarray}\begin{array}{l}| {\lambda }_{1},\ldots ,{\lambda }_{M}\rangle ={B}_{{a}_{1}}({\lambda }_{1}){B}_{{a}_{2}}({\lambda }_{2})\cdots {B}_{{a}_{M}}({\lambda }_{M}){{ \mathcal F }}^{{a}_{M}{a}_{M-1}\ldots {a}_{1}}| 0\rangle ,\end{array}\end{eqnarray}$
where {λ1,…,λM} is a set of Bethe roots and the repeated indices indicates the sum over the values 1 and 2, and $\{{{ \mathcal F }}^{{a}_{M}{a}_{M-1}\ldots {a}_{1}}\}$ are certain functions of {λj}.Before we go any further, let us introduce the following useful commutation relations:6 ). Here, the superscript represents the row and the subscript represents the column. The matrix r(u, v) in (23) is defined as,
$\begin{eqnarray}\begin{array}{rcl}B(u){B}_{k}(v) & = & -\displaystyle \frac{{a}^{-}(u,v)}{{b}^{-}(u,v)}{B}_{k}(v)B(u)+\displaystyle \frac{e(u,v)}{{b}^{-}(u,v)}{B}_{k}(u)B(v),\,\,k=1,2,\end{array}\end{eqnarray}$
$\begin{eqnarray}D(u){B}_{k}(v)=\displaystyle \frac{{b}^{+}(u,v)}{{c}^{+}(u,v)}{B}_{k}(v)D(u)+{\rm{other}}\,{\rm{terms}},\end{eqnarray}$
$\begin{eqnarray}{A}_{{ab}}(u){B}_{c}(v)=r(u,v{)}_{{bc}}^{{ed}}{B}_{e}(v){A}_{{ad}}(u)+\mathrm{other}\,\mathrm{terms},\end{eqnarray}$
which can be derived from the RTT relation ( $\begin{eqnarray}\begin{array}{l}r(u,v)=\displaystyle \frac{{a}^{-}(u,v)}{{b}^{-}(u,v)}\left(\begin{array}{cccc}1 & 0 & 0 & 0\\ 0 & \displaystyle \frac{U}{\tilde{u}-\tilde{v}+U} & \displaystyle \frac{-(\tilde{u}-\tilde{v})}{\tilde{u}-\tilde{v}+U} & 0\\ 0 & \displaystyle \frac{-(\tilde{u}-\tilde{v})}{\tilde{u}-\tilde{v}+U} & \displaystyle \frac{U}{\tilde{u}-\tilde{v}+U} & 0\\ 0 & 0 & 0 & 1\end{array}\right),\end{array}\end{eqnarray}$
with, $\begin{eqnarray}\begin{array}{l}\tilde{x}=\displaystyle \frac{\sin x}{\cos x}{{\rm{e}}}^{-2h(x)}-\displaystyle \frac{\cos x}{\sin x}{{\rm{e}}}^{2h(x)}.\end{array}\end{eqnarray}$
Using the commutation relations (23), we have: $\begin{eqnarray}\begin{array}{l}\displaystyle \sum _{a}{A}_{{aa}}(u){B}_{{a}_{1}}({\lambda }_{1}){B}_{{a}_{2}}({\lambda }_{2})\cdots {B}_{{a}_{M}}({\lambda }_{M})| 0\rangle \\ =\displaystyle \sum _{a}r(u,{\lambda }_{1}{)}_{{{aa}}_{1}}^{{e}_{1}{d}_{1}}{B}_{{e}_{1}}({\lambda }_{1}){A}_{{{ad}}_{1}}(u){B}_{{a}_{2}}({\lambda }_{2})\cdots {B}_{{a}_{M}}({\lambda }_{M})| 0\rangle +{\rm{u.t.}}\\ =\displaystyle \sum _{a}r(u,{\lambda }_{1}{)}_{{{aa}}_{1}}^{{e}_{1}{d}_{1}}r(u,{\lambda }_{2}{)}_{{d}_{1}{a}_{2}}^{{e}_{2}{d}_{2}}{B}_{{e}_{1}}({\lambda }_{1}){B}_{{e}_{2}}({\lambda }_{2}){A}_{{{ad}}_{2}}(u)\cdots {B}_{{a}_{M}}({\lambda }_{M})| 0\rangle +{\rm{u.t.}}\\ =\displaystyle \sum _{a}r(u,{\lambda }_{1}{)}_{{{aa}}_{1}}^{{e}_{1}{d}_{1}}r(u,{\lambda }_{2}{)}_{{d}_{1}{a}_{2}}^{{e}_{2}{d}_{2}}\cdots r(u,{\lambda }_{M}{)}_{{d}_{M-1}{a}_{M}}^{{e}_{M}a}{B}_{{e}_{1}}({\lambda }_{1}){B}_{{e}_{2}}({\lambda }_{2})\cdots {B}_{{e}_{M}}({\lambda }_{M}){A}_{{aa}}(u)| 0\rangle +{\rm{u.t.}}\\ =\displaystyle \sum _{a}[{\Pr }(u,{\lambda }_{1}){]}_{{{aa}}_{1}}^{{d}_{1}{e}_{1}}[{\Pr }(u,{\lambda }_{2}){]}_{{d}_{1}{a}_{2}}^{{d}_{2}{e}_{2}}\cdots [{\Pr }(u,{\lambda }_{M}){]}_{{d}_{M-1}{a}_{M}}^{{{ae}}_{M}}{B}_{{e}_{1}}({\lambda }_{1}){B}_{{e}_{2}}({\lambda }_{2})\cdots {B}_{{e}_{M}}({\lambda }_{M}){A}_{{aa}}(u)| 0\rangle +{\rm{u.t.}}\\ =\displaystyle \prod _{k\,=\,1}^{N}{a}^{+}(u,{\theta }_{k}){t}^{(1)}(u,\{{\lambda }_{j}\}{)}_{{a}_{1}{a}_{2}\ldots {a}_{M}}^{{e}_{1}{e}_{2}\ldots {e}_{M}}{B}_{{e}_{1}}({\lambda }_{1}){B}_{{e}_{2}}({\lambda }_{2})\cdots {B}_{{e}_{M}}({\lambda }_{M})| 0\rangle +{\rm{u.t.}},\end{array}\end{eqnarray}$
where u.t. denotes the unwanted terms and P is the permutation operator. Here, t(1)(u, {λj}) is the nested transfer matrix: $\begin{eqnarray}\begin{array}{rcl}{t}^{(1)}(u,\{{\lambda }_{j}\}) & = & {{tr}}_{0}\{{P}_{0,M}{r}_{M,0}(u,{\lambda }_{M})\cdots {P}_{\mathrm{0,1}}{r}_{\mathrm{1,0}}(u,{\lambda }_{1})\}\\ & = & \displaystyle \prod _{j=1}^{M}\displaystyle \frac{{a}^{-}(u,{\lambda }_{j})}{{b}^{-}(u,{\lambda }_{j})}{{tr}}_{0}\{{\tilde{R}}_{M0}(\tilde{u}-{\tilde{\lambda }}_{M}){\tilde{R}}_{M-10}(\tilde{u}-{\tilde{\lambda }}_{M-1})\cdots {\tilde{R}}_{10}(\tilde{u}-{\tilde{\lambda }}_{1})\},\end{array}\end{eqnarray}$
where, $\begin{eqnarray}\begin{array}{l}{\tilde{R}}_{n0}(u)=\left(\begin{array}{cccc}1 & 0 & 0 & 0\\ 0 & -\displaystyle \frac{u}{u+U} & \displaystyle \frac{U}{u+U} & 0\\ 0 & \displaystyle \frac{U}{u+U} & -\displaystyle \frac{u}{u+U} & 0\\ 0 & 0 & 0 & 1\end{array}\right).\end{array}\end{eqnarray}$
Applying the transfer matrix t(u) to the state ∣λ1,…,λM⟩ and using the commutation relations (21 )-(23 ) repeatedly, we obtain:27 ).
$\begin{eqnarray}\begin{array}{rcl}t(u)| {\lambda }_{1},\ldots ,{\lambda }_{M}\rangle & = & \left\{\displaystyle \prod _{n=1}^{N}{a}^{+}(u,{\theta }_{n})\displaystyle \prod _{j=1}^{M}\displaystyle \frac{-{a}^{-}(u,{\lambda }_{j})}{{b}^{-}(u,{\lambda }_{j})}+\displaystyle \prod _{n=1}^{N}{c}^{+}(u,{\theta }_{n})\displaystyle \prod _{j=1}^{M}\displaystyle \frac{{b}^{+}(u,{\lambda }_{j})}{{c}^{+}(u,{\lambda }_{j})}\right.\\ & & \left.+\displaystyle \prod _{n=1}^{N}{b}^{-}(u,{\theta }_{n}){{\rm{\Lambda }}}^{(1)}(u,\{{\lambda }_{j}\})\right\}| {\lambda }_{1},\ldots ,{\lambda }_{M}\rangle +{\rm{u}}.{\rm{t}}.,\end{array}\end{eqnarray}$
where Λ(1)(u, {λj}) is the eigenvalue of t(1)(u, {λj}) in (The function Λ(1)(u, {λj}) can be given by the algebraic Bethe ansatz method [22]:
$\begin{eqnarray}\begin{array}{l}{{\rm{\Lambda }}}^{(1)}(u,\{{\lambda }_{j}\})=\displaystyle \prod _{j=1}^{M}\displaystyle \frac{{a}^{-}(u,{\lambda }_{j})}{{b}^{-}(u,{\lambda }_{j})}\left[{\left(-1\right)}^{m}\displaystyle \prod _{k=1}^{m}\displaystyle \frac{\tilde{u}-{\mu }_{k}-U}{\tilde{u}-{\mu }_{k}}+{\left(-1\right)}^{M+m}\displaystyle \prod _{j=1}^{M}\displaystyle \frac{\tilde{u}-{\tilde{\lambda }}_{j}}{\tilde{u}-{\tilde{\lambda }}_{j}+U}\displaystyle \prod _{k=1}^{m}\displaystyle \frac{\tilde{u}-{\mu }_{k}+U}{\tilde{u}-{\mu }_{k}}\right],\end{array}\end{eqnarray}$
where m = 0,…,M and {μ1,…,μm} are the second set of Bethe roots.Define the following functions:30 ) into equation (29 ), the eigenvalue Λ(u) of the transfer matrix t(u) (7 ) can be parameterized as,
$\begin{eqnarray}\begin{array}{l}{z}_{+}(x)=\displaystyle \frac{\sin x}{\cos x}{{\rm{e}}}^{2h(x)},\,\,{z}_{-}(x)=\displaystyle \frac{\cos x}{\sin x}{{\rm{e}}}^{2h(x)}.\end{array}\end{eqnarray}$
Then, we can easily check the following useful relations: $\begin{eqnarray}\begin{array}{l}{z}_{\pm }(x+\pi )={z}_{\pm }(x),\,\,\,{z}_{\pm }(x+\tfrac{\pi }{2})=-\displaystyle \frac{1}{{z}_{\pm }(x)},\\ {z}_{+}(x)-\displaystyle \frac{1}{{z}_{+}(x)}+{z}_{-}(x)-\displaystyle \frac{1}{{z}_{-}(x)}=U.\end{array}\end{eqnarray}$
Substituting equation ( $\begin{eqnarray}\begin{array}{rcl}{\rm{\Lambda }}(u) & = & \displaystyle \prod _{l=1}^{N}{a}^{+}(u,{\theta }_{l})\displaystyle \prod _{j=1}^{M}\displaystyle \frac{-\cos u}{\sin u}\displaystyle \frac{{z}_{-}({\lambda }_{j})+{z}_{+}(u)}{{z}_{-}({\lambda }_{j})-{z}_{-}(u)}+\displaystyle \prod _{l=1}^{N}{c}^{+}(u,{\theta }_{l})\displaystyle \prod _{j=1}^{M}\displaystyle \frac{\cos u}{\sin u}\displaystyle \frac{{z}_{-}({\lambda }_{j})+1/{z}_{-}(u)}{{z}_{-}({\lambda }_{j})-1/{z}_{+}(u)}\\ & & +\displaystyle \prod _{l=1}^{N}{b}^{-}(u,{\theta }_{l})\left\{{\left(-1\right)}^{m}\displaystyle \prod _{j=1}^{M}\displaystyle \frac{\cos u}{\sin u}\displaystyle \frac{{z}_{-}({\lambda }_{j})+{z}_{+}(u)}{{z}_{-}({\lambda }_{j})-{z}_{-}(u)}\displaystyle \prod _{k=1}^{m}\displaystyle \frac{\tilde{u}-{\mu }_{k}-U}{\tilde{u}-{\mu }_{k}}\right.\\ & & \left.+{\left(-1\right)}^{M+m}\displaystyle \prod _{j=1}^{M}\displaystyle \frac{\cos u}{\sin u}\displaystyle \frac{{z}_{-}({\lambda }_{j})+1/{z}_{-}(u)}{{z}_{-}({\lambda }_{j})-1/{z}_{+}(u)}\displaystyle \prod _{k=1}^{m}\displaystyle \frac{\tilde{u}-{\mu }_{k}+U}{\tilde{u}-{\mu }_{k}}\right\},\end{array}\end{eqnarray}$
where $M,m\in {\mathbb{N}}$ and 0 ≤ m ≤ M ≤ 2N.We introduce the following short-hand notations:33 ) can be rewritten in a simpler form:29 ), the Bethe roots {λ1,…,λM} and {μ1,…,μm} should satisfy two sets of BAEs:9 ) in terms of the Bethe roots is:37 ) and (38 ) for the N = 2 case are shown in table 1. The energy spectrum given by Bethe roots is consistent with the ones from the exact diagonalization of the Hamiltonian.
$\begin{eqnarray}\begin{array}{l}{a}_{1}(u)=\displaystyle \prod _{l=1}^{N}{a}^{+}(u,{\theta }_{l}),\\ {a}_{2}(u)=\displaystyle \prod _{l=1}^{N}{b}^{-}(u,{\theta }_{l}),\\ {a}_{3}(u)=\displaystyle \prod _{l=1}^{N}{c}^{+}(u,{\theta }_{l}),\end{array}\end{eqnarray}$
and $\begin{eqnarray}\begin{array}{l}{Q}_{1}(u)=\displaystyle \prod _{j=1}^{M}\cos u\,\left({z}_{-}({\lambda }_{j})+{z}_{+}(u)\right),\\ {Q}_{2}(u)=\displaystyle \prod _{j=1}^{M}\sin u\,\left({z}_{-}({\lambda }_{j})-{z}_{-}(u)\right),\\ \tilde{Q}(\tilde{u})=\displaystyle \prod _{j=1}^{m}(\tilde{u}-{\mu }_{j}).\end{array}\end{eqnarray}$
Thus, the eigenvalue Λ(u) in ( $\begin{eqnarray}\begin{array}{rcl}{\rm{\Lambda }}(u) & = & {\left(-1\right)}^{M}{a}_{1}(u)\displaystyle \frac{{Q}_{1}(u)}{{Q}_{2}(u)}+{\left(-1\right)}^{m}{a}_{2}(u)\displaystyle \frac{{Q}_{1}(u)}{{Q}_{2}(u)}\displaystyle \frac{\tilde{Q}(\tilde{u}-U)}{\tilde{Q}(\tilde{u})}\\ & & +{\left(-1\right)}^{m}{a}_{2}(u)\displaystyle \frac{{Q}_{2}(u+\tfrac{\pi }{2})}{{Q}_{1}(u+\tfrac{\pi }{2})}\displaystyle \frac{\tilde{Q}(\tilde{u}+U)}{\tilde{Q}(\tilde{u})}\\ & & +{\left(-1\right)}^{M}{a}_{3}(u)\displaystyle \frac{{Q}_{2}(u+\tfrac{\pi }{2})}{{Q}_{1}(u+\tfrac{\pi }{2})}.\end{array}\end{eqnarray}$
To eliminate the unwanted terms in equation ( $\begin{eqnarray}\begin{array}{l}\displaystyle \prod _{k=1}^{N}\left[\displaystyle \frac{\sin {\theta }_{k}}{\cos {\theta }_{k}}\,\displaystyle \frac{1+{\nu }_{j}/{z}_{+}({\theta }_{k})}{1-{\nu }_{j}/{z}_{-}({\theta }_{k})}\right]\\ \,=\,{\left(-1\right)}^{m+1-M}\displaystyle \prod _{k=1}^{m}\displaystyle \frac{{\nu }_{j}^{-1}-{\nu }_{j}-{\bar{\mu }}_{k}-\tfrac{U}{2}}{{\nu }_{j}^{-1}-{\nu }_{j}-{\bar{\mu }}_{k}+\tfrac{U}{2}},\\ j=1,\ldots ,M,\end{array}\end{eqnarray}$
$\begin{eqnarray}\begin{array}{l}{\left(-1\right)}^{M+1}\displaystyle \prod _{k=1}^{M}\displaystyle \frac{{\bar{\mu }}_{j}-{\nu }_{k}^{-1}+{\nu }_{k}-\tfrac{U}{2}}{{\bar{\mu }}_{j}-{\nu }_{k}^{-1}+{\nu }_{k}+\tfrac{U}{2}}\\ \,=\,\displaystyle \displaystyle \prod _{l=1}^{m}\displaystyle \frac{{\bar{\mu }}_{j}-{\bar{\mu }}_{l}-U}{{\bar{\mu }}_{j}-{\bar{\mu }}_{l}+U},\,\,\,\,j=1,\ldots ,m,\end{array}\end{eqnarray}$
where, $\begin{eqnarray}\begin{array}{l}{\nu }_{j}={z}_{-}({\lambda }_{j}),\\ {\bar{\mu }}_{j}={\mu }_{j}+\displaystyle \frac{U}{2}.\end{array}\end{eqnarray}$
The eigenvalue of the Hamiltonian ( $\begin{eqnarray}\begin{array}{rcl}E & = & {\left.\displaystyle \frac{\partial \mathrm{ln}{\rm{\Lambda }}(u)}{\partial u}\right|}_{u=\theta ,\{{\theta }_{l}=\theta \}}+N\tan 2\theta \\ & = & \displaystyle \sum _{j=1}^{M}\displaystyle \frac{U\cos 2\theta \,{\rm{sech}} 2h(\theta )-2\left[{{\rm{e}}}^{2h(\theta )}{\nu }_{j}^{-1}+{{\rm{e}}}^{-2h(\theta )}{\nu }_{j}\right]}{-2\cos 2\theta +\sin 2\theta \left[-{{\rm{e}}}^{2h(\theta )}{\nu }_{j}^{-1}+{{\rm{e}}}^{-2h(\theta )}{\nu }_{j}\right]}\\ & & +\displaystyle \frac{U}{4}N{\cos }^{2}2\theta \,{\rm{sech}} 2h(\theta ).\end{array}\end{eqnarray}$
The numerical solutions of the BAEs (Table 1. The numerical solutions of the BAEs ( |
| ν1 | ν2 | ν3 | ν4 | ${\bar{\mu }}_{1}$ | ${\bar{\mu }}_{2}$ | ${\bar{\mu }}_{3}$ | ${\bar{\mu }}_{4}$ | E |
|---|---|---|---|---|---|---|---|---|
| −0.5786–0.8156i | −0.9863 + 0.1648i | − | − | 0.6508i | − | − | − | −4.0666 |
| −0.9020–0.4318i | − | − | − | − | − | − | − | −2.0000 |
| −0.9020–0.4318i | − | − | − | 0.8636i | − | − | − | −2.0000 |
| 0.2581 + 0.9661i | 0.4155-0.9096i | −0.9613–0.2755i | − | 0.4379i | − | − | − | −2.0000 |
| 0.2581 + 0.9661i | 0.4155-0.9096i | −0.9613–0.2755i | − | 2.0000i | 2.0000i | − | − | −2.0000 |
| 0.1126–0.9936i | −0.1126 + 0.9936i | − | − | − | − | − | − | -0.7327 |
| − | − | − | − | − | − | − | − | 0.7327 |
| 0.9757–0.2190i | − | − | − | − | − | − | − | 2.0000 |
| 0.9757–0.2190i | − | − | − | 0.4379i | − | − | − | 2.0000 |
| 0.5906–0.8070i | 0.2024 + 0.9793i | 0.7969–0.6041i | − | 0.8636i | − | − | − | 2.0000 |
| 0.5906–0.8070i | 0.7969–0.6041i | 0.2024 + 0.9793i | − | 2.0000i | −2.0000i | − | − | 2.0000 |
| 1.7146-0.4953i | 0.5383–0.1555i | − | − | 0.6508i | − | − | − | 4.0666 |
When θ = 0, our extended Hubbard model degenerates into the conventional one. As a consequence, the corresponding BAEs and the eigenvalue of the Hamiltonian reduce to,
$\begin{eqnarray}{\nu }_{j}^{N}={\left(-1\right)}^{m+1-M}\displaystyle \prod _{k=1}^{m}\displaystyle \frac{{\nu }_{j}^{-1}-{\nu }_{j}-{\bar{\mu }}_{k}-\tfrac{U}{2}}{{\nu }_{j}^{-1}-{\nu }_{j}-{\bar{\mu }}_{k}+\tfrac{U}{2}},\,\,\,\,j=1,\ldots ,M,\end{eqnarray}$
$\begin{eqnarray}\begin{array}{rcl} & & {\left(-1\right)}^{M+1}\displaystyle \prod _{k=1}^{M}\displaystyle \frac{{\bar{\mu }}_{j}-{\nu }_{k}^{-1}+{\nu }_{k}-\tfrac{U}{2}}{{\bar{\mu }}_{j}-{\nu }_{k}^{-1}+{\nu }_{k}+\tfrac{U}{2}}\\ & & \,=\displaystyle \prod _{l=1}^{m}\displaystyle \frac{{\bar{\mu }}_{j}-{\bar{\mu }}_{l}-U}{{\bar{\mu }}_{j}-{\bar{\mu }}_{l}+U},\,\,\,\,j=1,\ldots ,m,\end{array}\end{eqnarray}$
and $\begin{eqnarray}\begin{array}{rcl}E & = & \displaystyle \frac{U(N-2M)}{4}+\displaystyle \sum _{j=1}^{M}\left({\nu }_{j}+{\nu }_{j}^{-1}\right)\\ & = & \displaystyle \frac{U(N-2M)}{4}+2\displaystyle \sum _{j=1}^{M}\cos {k}_{j},\quad {\nu }_{j}={{\rm{e}}}^{{{ik}}_{j}}.\end{array}\end{eqnarray}$
4. Conclusion
In this paper, we study a 1D extended Hubbard model with a periodic boundary condition. We construct an integrable Hamiltonian (12 ) within the framework of the QISM. Compared with the conventional Hubbard model, the extended one contains more interaction terms. Using the nested algebraic Bethe ansatz method, the eigenvalue problem of the extended Hubbard model is solved by the homogeneous T − Q relation (36 ) and the associated BAEs (37 ) and (38 ). The numerical simulations imply that the solutions of the BAEs (37 ) and (38 ) indeed give the correct spectrum of the Hamiltonian. It should be noted that the T − Q relation (36 ) and BAEs (37 ) and (38 ) are constructed by selecting an all spin-up state as the vacuum state and they may not give the complete solutions. There also exists another T − Q relation with an all spin-down state being the vacuum. These two Bethe ansatz should give the complete set of eigenvalues and eigenstates of the Hamiltonian.
Furthermore, one can study the explicit form of the eigenstate in equation (20 ). In addition, based on our homogeneous BAEs, the thermodynamic properties of the model can also be studied via the well-known thermodynamic Bethe ansatz method [21].
Another interesting objective is to construct integrable extended Hubbard models with open boundary conditions. These models can be exactly solved via the off-diagonal Bethe ansatz method [22]. For open systems, we can study the thermodynamic limit of the model through the novel t − W scheme [23, 24].