Inner products in integrable Richardson-Gaudin models

TL;DR

This paper derives determinant formulas for inner products of eigenstates in Richardson-Gaudin models, linking Bethe rapidities and conserved charge eigenvalues, and connects different integrable model approaches.

Contribution

It introduces two classes of Gaudin-like determinant expressions for inner products, bridging Bethe rapidities and conserved charge eigenvalues in integrable models.

Findings

Derived determinant formulas for inner products.

Connected Bethe rapidities with conserved charge eigenvalues.

Expressed inner products as domain wall boundary partition functions.

Abstract

We present the inner products of eigenstates in integrable Richardson-Gaudin models from two different perspectives and derive two classes of Gaudin-like determinant expressions for such inner products. The requirement that one of the states is on-shell arises naturally by demanding that a state has a dual representation. By implicitly combining these different representations, inner products can be recast as domain wall boundary partition functions. The structure of all involved matrices in terms of Cauchy matrices is made explicit and used to show how one of the classes returns the Slavnov determinant formula. This framework provides a further connection between two different approaches for integrable models, one in which everything is expressed in terms of rapidities satisfying Bethe equations, and one in which everything is expressed in terms of the eigenvalues of conserved charges, satisfying quadratic equations.

Equations

$$[S^x(u),S^y(v)]$$

$$[S^y(u),S^z(v)]$$

$$[S^z(u),S^x(v)]$$

$$[S^{\kappa }(u),S^{\kappa }(v)]$$

$$S^{+}(u)=S^x(u)+i{S}^y(u),S^{-}(u)=S^x(u)-i{S}^y(u),$$

$$\left[S^{+}(u),S^{-}(v)\right]$$

$$\left[S^z(u),S^{\pm }(v)\right]$$

$${\mathbb{S}}^2(u)$$

$$=\frac{1}{2}\left(S^{+}(u)S^{-}(u)+S^{-}(u)S^{+}(u)\right)+S^z(u{)}^2.$$

$$[{\mathbb{S}}^2(u),{\mathbb{S}}^2(v)]=0,\forall u,v\in \mathbb{C},$$

$$|v_1\dots v_N⟩=S^{+}(v_1)\dots S^{+}(v_N)|0⟩=\prod _{a=1}^N{S}^{+}(v_a)|0⟩,$$

$${\mathbb{S}}^2(u)|v_1\dots v_N⟩=$$

$$+\sum _{a=1}^N\frac{2}{u-v_a}\left[F(v_a)-\sum _{b\ne a}^N\frac{1}{v_a-v_b}\right]|v_1\dots v_a\to u\dots v_N⟩,$$

$$F(v_a)-\sum _{b\ne a}^N\frac{1}{v_a-v_b}=0,\text{with}{S}^z(v_a)|0⟩=F(v_a)|0⟩,a=1\dots N,$$

$$[S_i^{+},S_j^{-}]=2{\delta }_{ij}{S}_i^z,[S_i^z,S_j^{\pm }]=\pm {\delta }_{ij}{S}_i^{\pm },$$

$$S^z(u)=-\frac{1}{g}-\sum _{i=1}^L\frac{S_i^z}{u-{ϵ}_i},S^{+}(u)=\sum _{i=1}^L\frac{S_i^{+}}{u-{ϵ}_i},S^{-}(u)=\sum _{i=1}^L\frac{S_i^{-}}{u-{ϵ}_i},$$

$${\mathbb{S}}^2(u)=\frac{2}{g}\sum _{i=1}^L\frac{R_i}{u-{ϵ}_i}+⟨0|{\mathbb{S}}^2(u)|0⟩,$$

$$R_i=\left(S_i^z+\frac{1}{2}\right)+g\sum _{j\ne i}^L\frac{1}{{ϵ}_i-{ϵ}_j}\left[\frac{1}{2}\left(S_i^{+}{S}_j^{-}+S_i^{-}{S}_j^{+}\right)+\left(S_i^z{S}_j^z-\frac{1}{4}\right)\right],i=1\dots L.$$

$$|v_1\dots v_N⟩=\prod _{a=1}^N\left(\sum _{i=1}^L\frac{S_i^{+}}{{ϵ}_i-v_a}\right)|\downarrow \cdots \downarrow ⟩,$$

$$\frac{1}{g}+\frac{1}{2}\sum _{i=1}^L\frac{1}{{ϵ}_i-v_a}-\sum _{b\ne a}^N\frac{1}{v_b-v_a}=0,a=1\dots N.$$

$$R_i|v_1\dots v_N⟩$$

$$=-\frac{g}{2}{\Lambda }_i(\{v_a\})|v_1\dots v_N⟩.$$

$$R_i^2=R_i-\frac{g}{2}\sum _{j\ne i}^L\frac{R_i-R_j}{{ϵ}_i-{ϵ}_j},i=1\dots L.$$

$${\Lambda }_i^2+\frac{2}{g}{\Lambda }_i-\sum _{j\ne i}^L\frac{{\Lambda }_i-{\Lambda }_j}{{ϵ}_i-{ϵ}_j}=0,i=1\dots L,$$

$${\Lambda }_i\equiv {\Lambda }_i(\{v_a\})=\sum _{a=1}^N\frac{1}{{ϵ}_i-v_a},i=1\dots L.$$

$$P(z)=\prod _{a=1}^N(z-v_a),$$

$$P^{\prime \prime }(z)+\left[\frac{2}{g}+\sum _{i=1}^L\frac{1}{{ϵ}_i-z}\right]P^{\prime }(z)-\left[\sum _{i=1}^L\frac{{\Lambda }_i}{{ϵ}_i-z}\right]P(z)=0.$$

$$⟨v_1\dots v_N⟩w_1\dots w_N=\frac{{\prod }_b{\prod }_{a\ne b}(v_a-w_b)}{{\prod }_{a<b}(v_b-v_a){\prod }_{b<a}(w_b-w_a)}\det S_N(\{v_a\},\{w_b\}),$$

$$S_N(\{v_a\},\{w_b\}{)}_{ab}=\frac{v_b-w_b}{v_a-w_b}\left(\sum _{i=1}^L\frac{1}{(v_a-{ϵ}_i)(w_b-{ϵ}_i)}-2\sum _{c\ne a}^N\frac{1}{(v_a-v_c)(w_b-v_c)}\right).$$

$$⟨v_1\dots v_N⟩w_1\dots w_N=(-1{)}^N{\left(\frac{g}{2}\right)}^{L-2N}\det J_L(\{v_a\},\{w_b\}),$$

$$J_L(\{v_a\},\{w_b\}{)}_{ij}=\{\begin{array}{ll}\frac{2}{g}+{\Lambda }_i(\{v_a\})+{\Lambda }_i(\{w_b\})-{\sum }_{k\ne i}^L\frac{1}{{ϵ}_i-{ϵ}_k}& \text{if}i=j\\ -\frac{1}{{ϵ}_i-{ϵ}_j}& \text{if}i\ne j\end{array},$$

$$J_L(\{v_a\},\{w_b\}{)}_{ij}=\{\begin{array}{ll}\frac{2}{g}+{\sum }_{a=1}^N\frac{1}{{ϵ}_i-v_a}+{\sum }_{b=1}^N\frac{1}{{ϵ}_i-w_b}-{\sum }_{k\ne i}^L\frac{1}{{ϵ}_i-{ϵ}_k}& \text{if}i=j\\ -\frac{1}{{ϵ}_i-{ϵ}_j}& \text{if}i\ne j\end{array},$$

$$⟨v_1\dots v_N⟩w_1\dots w_N=(-1{)}^N\det K_{2N}(\{v_a\},\{w_b\}),$$

$$K_{2N}(\{v_a\},\{w_b\}{)}_{ab}=\{\begin{array}{ll}\frac{2}{g}-{\sum }_{i=1}^L\frac{1}{x_a-{ϵ}_i}+{\sum }_{c\ne a}^{2N}\frac{1}{x_a-x_c}& \text{if}a=b\\ -\frac{1}{x_a-x_b}& \text{if}a\ne b\end{array}.$$

Full text

Inner products in integrable Richardson-Gaudin models

P.W. Claeys 1,2,3 *, D. Van Neck 2,3, S. De Baerdemacker 2,3,4

1 Institute for Theoretical Physics, University of Amsterdam, Science Park 904, 1098 XH Amsterdam, The Netherlands

2 Department of Physics and Astronomy, Ghent University, Proeftuinstraat 86, 9000 Ghent, Belgium

3 Center for Molecular Modeling, Ghent University, Technologiepark 903, 9052 Ghent, Belgium

4 Department of Inorganic and Physical Chemistry, Ghent University, Krijgslaan 281 (S3), 9000 Ghent, Belgium

• PieterW.Claeys@UGent.be

Abstract

We present the inner products of eigenstates in integrable Richardson-Gaudin models from two different perspectives and derive two classes of Gaudin-like determinant expressions for such inner products. The requirement that one of the states is on-shell arises naturally by demanding that a state has a dual representation. By implicitly combining these different representations, inner products can be recast as domain wall boundary partition functions. The structure of all involved matrices in terms of Cauchy matrices is made explicit and used to show how one of the classes returns the Slavnov determinant formula.

**Furthermore, this framework provides a further connection between two different approaches for integrable models, one in which everything is expressed in terms of rapidities satisfying Bethe equations, and one in which everything is expressed in terms of the eigenvalues of conserved charges, satisfying quadratic equations. **

Contents

1. 1 Introduction

2. 2 Richardson-Gaudin models

3. 2.1 Generalized Gaudin algebra

4. 2.2 Eigenvalue-based framework

5. 2.3 Moving between both frameworks

6. 3 Overview of the main results

7. 3.1 Inner products

8. 3.2 Gaudin determinant

9. 3.3 Izergin-Borchardt determinant

10. 4 Properties of Cauchy matrices

11. 5 From eigenvalue-based to Slavnov through dual states

12. 5.1 Dual states

13. 5.2 Inner products as DWPFs

14. 5.3 Reduction to Slavnov’s determinant

15. 6 Extension to hyperbolic models

16. 7 Conclusion

17. A Orthogonality of Bethe states

18. B Properties of Cauchy matrices

19. B.1 Inverse of Cauchy matrices

20. B.2 First order

21. B.3 Second order

22. B.4 Sylvester’s determinant identity

23. C Normalizations of the normal and the dual state

24. D Results for the hyperbolic model

1 Introduction

Integrable models have proven to be a powerful tool in the study of many-body physics. These models are often characterized by two key properties - they support a large amount of conserved quantities, and they can be solved exactly using Bethe ansatz methods [1, 2]. The eigenstates of integrable models are given by Bethe states, which are characterized by a set of variables, so-called rapidities or Bethe roots, which have to obey a set of coupled non-linear equations, the Bethe equations. The underlying framework of integrability then often allows for the exact calculation of physical observables from these eigenstates. For example correlation coefficients, one of the main means of extracting the underlying physics from a given wave function, can be efficiently calculated in these models [2]. More specifically, these correlation coefficients can often be reduced to (polynomial sums of) determinants. Since determinants can be efficiently evaluated numerically, such expressions have subsequently allowed for the investigation of various large-scale integrable systems in many different contexts, where e.g. their appearance in quantum quenches has attracted a lot of attention [3, 4]. Here, one of the key expressions is the well-known Slavnov formula for the inner product between two Bethe states, one with rapidities satisfying the Bethe equations (an on-shell state), and one with arbitrary rapidities (an off-shell state) [5]. Such determinant expressions provide a basic building block for the calculation of correlation coefficients from the Bethe states, which has allowed for massive simplifications in the calculations of correlation coefficients in these models [6, 7, 8, 9, 10]. Such inner products can also be interpreted as domain wall boundary partition functions (DWPF), which can similarly be expressed as determinants [11, 12, 13, 14].

In the following, we investigate the structure of these inner products in Richardson-Gaudin models [15, 16, 17, 5, 18, 7]. The existence of two distinct representations for each on-shell Bethe state allows such inner products to be recast as DWPFs, and we show here how the Slavnov determinant follows as a corollary. The demand that one of the two states is on-shell, a prerequisite for the existence of such determinant expressions, then follows naturally by requiring that a so-called dual representation of one of the two states exist, which is a known characteristic of on-shell states. In the (numerical) treatment of Richardson-Gaudin models, two separate approaches exist. In the first, on-shell Bethe states are obtained by solving the nonlinear coupled Richardson-Gaudin (or Bethe) equations for the rapidities, after which correlation coefficients are calculated from these rapidities using the Slavnov determinant. In the second approach, the so-called eigenvalue-based approach, eigenstates are characterized by solving for a set of variables determining the eigenvalues of the conserved charges for each eigenstate. The correlation coefficients can then also be calculated as DWPFs depending only on these new variables, keeping the rapidities implicit [19, 20, 21, 22]. In this work, the connection between these approaches is made explicit, highlighting the Cauchy structure of all involved matrices, which follows from the rational functions defining these models.

Recently, Richardson-Gaudin integrable models where no clear reference state exists have been under much attention [23, 24, 25, 26]. While complicating the structure of the Bethe equations, this only has minor influence on the eigenvalue-based equations. This also allows for a straightforward derivation of eigenvalue-based determinant expressions for inner products in these models. While it would be worthwhile to obtain similar connections as presented in this work for these models, we restrict ourselves to the better-known class of models containing a clear reference state.

The paper is organized as follows. In section 2, we give an overview of the theory of Richardson-Gaudin models, after which we give an overview of the main results of this work in section 3. The connection with previously known expressions is then made explicit in sections 4 and 5, where we first derive several crucial properties of Cauchy matrices, which allow for the use of dual representations and DWPFs. For the clarity of presentation, we restrict ourselves to the so-called rational model throughout the paper, and show in section 6 how this can be straightforwardly extended to hyperbolic models. Section 7 is then reserved for conclusions.

2 Richardson-Gaudin models

2.1 Generalized Gaudin algebra

In this section, we will give a brief overview of the theory of Richardson-Gaudin systems from the viewpoint of generalized Gaudin algebras, following the presentation from [27]. It is worth mentioning that the integrability of these models can also be derived within the framework of the Algebraic Bethe Ansatz through taking the so-called quasi-classical limit (see [28, 29, 18, 30, 8, 31, 32]).

As mentioned in the introduction, integrable models are characterized by two features going hand in hand - the existence of a large set of conserved quantities, and their exact solvability by Bethe ansatz. These both follow naturally from the concept of a Generalized Gaudin algebra (GGA). Such an algebra is defined by operators $S^{x}(u),S^{y}(u),S^{z}(u)$ satisfying the commutation relations

[TABLE]

with $u,v\in\mathbb{C}$. This is an infinite dimensional Lie algebra, highly reminiscent of the $su(2)$ algebra. Extending the similarity to $su(2)$, raising and lowering operators can be defined as

[TABLE]

which satisfy the commmutation relations

[TABLE]

Given such an algebra, a continuous family of mutually commuting operators111These provide the GGA equivalent of the transfer matrix in the Algebraic Bethe Ansatz [2, 33, 7]. can be defined as

[TABLE]

Note that although these resemble the Casimir operator of $su(2)$, they do not act as Casimir operators for the GGA since they do not commute with its generators. It follows from the given commutation relations of the GGA that

[TABLE]

so these operators generate a continuous set of commuting operators, leading to a continuous set of conserved charges.

Their exact solvability also follows from the definition of the GGA. Bethe ansatz states can be constructed as

[TABLE]

containing a product of raising operators acting on a vacuum reference state $\ket{0}$, which is defined as being both annihilated by the lowering operator $S^{-}(u)\ket{0}=0,\forall\ u\in\mathbb{C}$, and being an eigenstate of $S^{z}(u),\forall\ u\in\mathbb{C}$. This wave function then depends on a set of free parameters $\{v_{a}\}=\{v_{1}\dots v_{N}\}$, so-called rapidities or Bethe roots.

Using the commutation relations, the action of $\mathbb{S}^{2}(u)$ on such a Bethe state can be shown to be

[TABLE]

where the state $\ket{v_{1}\dots v_{a}\rightarrow u\dots v_{N}}$ is built as in Eq. (10), but with a single rapidity replaced by $u$, and $F(u)$ is defined as $S^{z}(u)\ket{0}=F(u)\ket{0}$. The action of $\mathbb{S}^{2}(u)$ gives rise to two contributions: one diagonal part and a set of $N$ off-diagonal (unwanted) states. So far, the Bethe state has been defined in terms of $N$ free rapidities. If these are chosen in such a way that the unwanted parts vanish, the resulting Bethe state is an eigenstate by construction. Indeed, if the rapidities satisfy the Richardson-Gaudin (or Bethe) equations

[TABLE]

the resulting Bethe state is an eigenstate.

A standard realization of this algebra exists within the context of spin-$1/2$ models, leading to both the Richardson [15, 16] and the Gaudin models [17]. For a set of $L$ $su(2)$-algebras $\{S_{i}^{z},S_{i}^{+},S_{i}^{-}\}$, labelled $i=1\dots L$ and satisfying

[TABLE]

a GGA can be realized by the following operators

[TABLE]

and $\ket{0}=\ket{\downarrow\dots\downarrow}$. Both $g$ and $\{\epsilon_{1}\dots\epsilon_{L}\}$ are free variables, with the latter playing the role of the inhomogeneities in the usual transfer matrix construction.

The resulting $\mathbb{S}^{2}(u)$ has single poles in $u=\epsilon_{i},\ i=1\dots L$ and can be expanded as

[TABLE]

where the residues lead to a set of commuting operators

[TABLE]

The common eigenstates of these operators now follow from the GGA as

[TABLE]

described by rapidities $\{v_{1}\dots v_{N}\}$ satisfying the equations

[TABLE]

The action of the conserved charges on these states can then be obtained from the eigenvalues of $\mathbb{S}^{2}(u)$ as

[TABLE]

Here, $\Lambda_{i}(\{v_{a}\})=\sum_{a=1}^{N}\frac{1}{\epsilon_{i}-v_{a}}$ encodes the dependence of the eigenvalues of the conserved charges on the rapidities, and is therefore also known as an eigenvalue-based variable [34].

2.2 Eigenvalue-based framework

Conventionally, solving integrable systems consists of first solving the Bethe equations (18) for the rapidities, which parametrize the eigenstates, and then extracting correlation coefficients from these rapidities using determinant expressions following from Slavnov’s determinant. An alternative approach avoids the explicit use of rapidities. Instead, it can be shown that the conserved charges (16) satisfy a set of quadratic equations [24, 35]

[TABLE]

Because these operators commute by construction, they can be diagonalized simultaneously, and this equation can be rewritten at the level of the eigenvalues. This leads to a set of quadratic equations for the eigenvalue-based variables

[TABLE]

with

[TABLE]

These equations are also known as the substituted or quadratic Bethe equations. Originally obtained by Babelon and Talalaev [19], these were later extended in a series of articles towards all known Richardson-Gaudin models [34, 36, 21, 22, 26]. They arose as a way of circumventing the singular behaviour of the regular Richardson-Gaudin equations, since they do not exhibit the singularities associated with the poles in the regular Bethe equations [37, 38, 39]. In many ways, this framework is similar to the use of Baxter’s $T-Q$-relation in the Algebraic Bethe Ansatz in immediately solving for the eigenvalues of the conserved charges [40].

Remarkably, it is possible to obtain expressions for inner products between two Bethe states which only depend explicitly on these new variables, further cementing their usefulness, so it is not always necessary to extract the rapidities from these variables [12, 13, 21, 22, 14].

2.3 Moving between both frameworks

Nevertheless, not all possible correlation coefficients or matrix elements can be expressed in terms of the eigenvalue-based variables, so it is necessary to obtain an efficient way of moving between representations (rapidity-based and eigenvalue-based). From a numerical point of view these different representations presents a trade-off. On the one hand, compared to the Richardson-Gaudin equations the eigenvalue-based equations can be solved in a much more straightforward way, but it is then not trivial to obtain the rapidities from the eigenvalue-based variables. On the other hand, while directly solving the Bethe equations (18) for the rapidities is a more involved task, the eigenvalue-based variables can afterwards be immediately obtained by simply plugging the rapidities into Eq. (22).

One way of (numerically) extracting the rapidities from the eigenvalue-based variables and inverting Eq. (22) is by defining a polynomial with the rapidities as roots

[TABLE]

which satisfies the ordinary differential equation [41, 36, 35]

[TABLE]

Once the variables $\Lambda_{i}$ are known, this differential equation is well-defined, and efficient algorithms exist for finding the roots of a polynomial satisfying such a differential equation [36]. Remarkably, the differential equation (24) can be derived starting from either the Bethe equations (18) or the eigenvalue-based equations (21), and can be considered the connection between both frameworks [35].

This has profound consequences, as recently pointed out in [35]. Since the eigenvalue-based equations follow from operator identities, they are necessarily complete, and the equivalence between the Bethe equations and the eigenvalue-based equations then implies the completeness of the Bethe equations and the associated Bethe ansatz states.

Given two different ways of solving the Bethe equations, the differential equation (24) then poses another way of solving for the eigenstates, combining both approaches. Instead of first solving for the eigenvalue-based variables and afterwards extracting the rapidities, the expansion of the eigenvalue-based variables in terms of the rapidities can be introduced, and the full equation can be solved for the rapidities in a self-consistent way, as done in [42]. The correspondence between Bethe equations and differential equations has similarly led to a plethora of techniques for solving the Bethe equations [42, 43, 44, 45, 46, 47, 48, 49].

3 Overview of the main results

3.1 Inner products

The main object of interest in this work is the inner product of Bethe states. In this section, a series of determinant expressions for such inner products will be presented without proof, after which the mathematical identities and reasoning underlying these results will be presented in later sections.

Suppose we have two sets of rapidities, one ($\{v_{a}\}=\{v_{1}\dots v_{N}\}$) obtained by solving the Richardson-Gaudin equations, and another ($\{w_{b}\}=\{w_{1}\dots w_{N}\}$) given by arbitrary complex variables. If the rapidities satisfy the Richardson-Gaudin equations, the resulting Bethe state is also known as an on-shell state, whereas a Bethe state defined by arbitrary variables is known as an off-shell state. The inner product between an on-shell and an off-shell state is then famously given by Slavnov’s determinant formula [5, 2, 50, 7]

[TABLE]

with $S_{N}(\{v_{a}\},\{w_{b}\})$ an $N\times N$ matrix defined as

[TABLE]

The first main results of this paper are the alternative expressions

[TABLE]

with $J_{L}(\{v_{a}\},\{w_{b}\})$ an $L\times L$ matrix defined as

[TABLE]

where the diagonal elements depend on the eigenvalue-based variables as defined in Eq. (22), and the dependence on the rapidities can be made explicit as

[TABLE]

and the related and equivalent determinant expression

[TABLE]

with $K_{2N}(\{v_{a}\},\{w_{b}\})$ an $2N\times 2N$ matrix defined in terms of $\{x_{1}\dots x_{2N}\}=\{v_{1}\dots v_{N}\}\cup\{w_{1}\dots w_{N}\}$ as

[TABLE]

These two expressions are clearly related by exchanging the role of the rapidities and the inhomogeneities, a feature which will be further detailed in following sections. Furthermore, the matrix (28) is well-suited for the eigenvalue-based framework, since it only depends on the rapidities through the eigenvalue-based variables in the diagonal elements. Because of this specific structure, it is not necessary to explicitly know the rapidities in order to evaluate the matrix elements. Note that the diagonal elements are reminiscent of the eigenvalue-based equations (21), a feature which will be further commented on in later sections. Similarly, the matrix (30) only depends on the inhomogeneities through the diagonal elements, in a form similar to the Bethe equations (18). The orthogonality of two different eigenstates can then easily be shown by exploiting the similarity of the diagonal elements to the eigenvalue-based/Bethe equations, as shown in Appendix A.

These matrices also exhibit the same structure as the Gaudin matrix for the normalization of an on-shell state [17], hence the denomination of Gaudin-like matrices. For this kind of matrices, the off-diagonal elements only depend on the difference of two rapidities (inhomogeneities). The diagonal elements then contain two summations, one over all but one rapidities (inhomogeneities), and one over all inhomogeneities (rapidities). Remarkably, these Gaudin-like structures are not limited to Richardson-Gaudin models and have been observed in general integrable models such as the XXZ spin chain [51, 52, 53, 54] and the Lieb-Liniger model [55, 56].

The second main result of this paper is the identification of the Cauchy matrix as the fundamental building block underneath all matrix expressions, being directly responsible for the Gaudin-like structure. In the following sections, it will also be shown that the Slavnov determinant (25) can be derived starting from Eq. (30), similarly exposing its structure in terms of Cauchy matrices. We conclude this section with two well-known specific cases of inner products which have special use in calculations, leading to the Gaudin determinant and the Izergin-Borchardt determinant, and show how they fit within this scheme.

3.2 Gaudin determinant

If the two sets of rapidities coincide, Slavnov’s determinant expression (25) reduces to the Gaudin determinant for the normalization of Bethe states [17], given by

[TABLE]

with $G_{N}(v_{1}\dots v_{N})$ an $N\times N$ matrix defined as

[TABLE]

The Gaudin matrix is also obtained by taking the limit of coinciding rapidities for Eq. (30). Using elementary row and column operations, the determinant of the $2N\times 2N$ matrix from Eq. (30) can in this limit be reduced to that of the $N\times N$ Gaudin matrix.

Interestingly, the eigenvalue-based expression for the inner product (29) leads to an alternative expression for the normalization as

[TABLE]

with $J_{L}(\{v_{a}\})$ an $L\times L$ matrix defined as

[TABLE]

It has already been mentioned how the structure of this matrix resembles that of the presented determinant expressions (29) and (30). This similarity and the related use of both determinant expressions can be even further established, since the Gaudin matrix is identical to the Jacobian of the Bethe equations (18), which seems to be a common property of integrable systems [57, 58, 59, 60, 61, 62], while the alternative matrix (35) is exactly the Jacobian of the eigenvalue-based equations (21).

3.3 Izergin-Borchardt determinant

Another frequently encountered case is where the off-shell state simplifies to a simple product state $\ket{i_{1}\dots i_{N}}=\prod_{b=1}^{N}S^{\dagger}_{i_{b}}\ket{\downarrow\dots\downarrow}$, where Slavnov’s determinant expression results in an Izergin-Borchardt determinant [63, 64, 65]

[TABLE]

which contains the determinant of an $N\times N$ matrix with matrix elements $\frac{1}{(\epsilon_{i_{b}}-v_{a})^{2}}$. The alternative matrix following from Eq. (29) corresponding to the Izergin-Borchardt determinant then reads

[TABLE]

with $J_{N}(\{v_{a}\},\{i_{b}\})$ an $N\times N$ matrix defined as

[TABLE]

while Eq. (30) gives rise to

[TABLE]

with $K_{N}(\{v_{a}\},\{i_{b}\})$ an $N\times N$ matrix defined as

[TABLE]

In fact, these results were originally obtained in [12], and served as the building blocks for this work. We should stress that the matrices in Eqs. (29) and (35) also appeared in [12] where it was shown that inner products and normalizations are proportional to the determinants of these respective matrices. We show here that the proportionality factor can be immediately evaluated as a single constant, independent of both the rapidities and the inhomogeneities of the model. This will be proven in section 5.

An overview of all relevant matrices and their counterparts in both frameworks is given in Figure 1.

4 Properties of Cauchy matrices

In order to derive the previously-presented determinant expressions, the Cauchy structure of all matrices needs to be made clear. Several properties of Cauchy matrices will then be used in order to navigate between different determinant expressions, which will be presented in this section. Necessary proofs have been moved to Appendix B for clarity of presentation.

To set the stage, assume we have two sets of variables $\{\epsilon_{1}\dots\epsilon_{N}\}=\{\epsilon_{i}\}$ and $\{x_{1}\dots x_{N}\}=\{x_{\alpha}\}$. Given two such sets, an $N\times N$ Cauchy matrix $C$ is defined by the matrix elements

[TABLE]

The inverse of this matrix is related to the transpose of this matrix through two diagonal matrices [66], which are defined in terms of two polynomials $p(x)=\prod_{i}(x-\epsilon_{i})$ and $q(x)=\prod_{\alpha}(x-x_{\alpha})$ as

[TABLE]

such that

[TABLE]

Permanents of Cauchy matrices are ubiquitous in the theory of Richardson-Gaudin integrability, and an important result by Borchardt [63] showed how the permanent of a Cauchy matrix can be evaluated as a ratio of determinants given by

[TABLE]

with $*$ the Hadamard product, defined as $(A*B)_{ij}=A_{ij}B_{ij}$. The denominator $\det\left[C\right]$ can be explicitly evaluated as

[TABLE]

However, instead of directly evaluating these determinants, it is also possible to rewrite this as

[TABLE]

as noted in [14]. Because of the known structure of the inverse of a Cauchy matrix, these matrices can be explicitly calculated as

[TABLE]

with the $J$-matrices given by (see Appendix B)

[TABLE]

and the $D$-matrices are those arising in the equation for the inverse (44). Now multiple determinant expressions for the permanent can be used, since

[TABLE]

Because of the product structure of the $J$-matrices, the following identity holds

[TABLE]

Note that so far all involved matrices were $N\times N$ matrices defined in terms of two sets of variables $\{\epsilon_{1}\dots\epsilon_{N}\}$ and $\{x_{1}\dots x_{N}\}$. Remarkably, this final expression (52) can immediately be generalized towards matrices of different dimensions, defined in terms of sets with a different number of variables, playing an important role in connecting the inner products to DWPFs. Given two such sets $\{\epsilon_{1}\dots\epsilon_{L}\}$ and $\{x_{1}\dots x_{N}\}$, with $L\neq N$, the following also holds

[TABLE]

with $J_{x}$ and $N\times N$ matrix and $J_{\epsilon}$ an $L\times L$ matrix defined as in Eq. (49), but now in terms of sets with a different number of variables. This result can be obtained by generalizing Eq. (52) towards non-square matrices and applying Sylvester’s determinant identity

[TABLE]

with $A$ an arbitrary $N\times L$ matrix and $B$ an arbitrary $L\times N$ matrix, as commented on in Appendix B.4.

Returning to $N\times N$ matrices, one final property of Cauchy matrices which will be key in relating the eigenvalue-based determinant to Slavnov’s determinant, is that

[TABLE]

as also proven in Appendix B. This can be seen as a higher-order extension of the previous equality

[TABLE]

and leads to the determinant identity

[TABLE]

5 From eigenvalue-based to Slavnov through dual states

5.1 Dual states

The crucial element in all eigenvalue-based expressions is the (implicit) existence of dual states for on-shell Bethe states. Any eigenstate of the integrable models under study can be constructed in two different ways: either by creating excitations on top of a vacuum, or by annihilating excitations from a dual vacuum state [12, 21]. The former approach was used in Section 2, where eigenstates were constructed as

[TABLE]

with rapidities $\{v_{1}\dots v_{N}\}$ satisfying the equations

[TABLE]

However, the latter approach states that eigenstates also have a dual representation given by

[TABLE]

with the rapidities of this dual state satisfying

[TABLE]

Note that these equations are related to the previous ones by changing the sign of $g$, which can be seen as a consequence of the spin-flip symmetry of the conserved charges (16). Because total spin-projection $S^{z}=\sum_{i=1}^{L}S_{i}^{z}$ is a symmetry of the system, this also implies that the same state will in general be described by a different number of rapidities in both representations.

Two such states are eigenstates of a given integrable Hamiltonian by construction, and these can be made to represent the same eigenstate by demanding that their eigenvalues coincide, leading to

[TABLE]

While the correspondence between the rapidities of both representations is not intuitive, the correspondence between the eigenvalue-based variables is simple and given by

[TABLE]

A special case is the $g\to\infty$ limit, where both sets of variables are equal apart from a number of diverging rapidities going to infinity. In this limit, the conserved charges have an additional $su(2)$ total spin symmetry, which is reflected in $\lim_{v\to\infty}S^{+}(v)\sim\sum_{i}S_{i}^{+}$, as was also observed in [11].

These two different representations of a single state obviously have different normalizations, but it is shown in Appendix C that the ratio of these normalizations is simply given by

[TABLE]

This ratio is directly responsible for the appearance of the prefactor in Eq. (29), as will now be shown.

5.2 Inner products as DWPFs

Since the existence of dual states is guaranteed for on-shell Bethe states, we can always express the inner product in terms of the dual state. The overlap between a (possibly off-shell) state $\ket{w_{1}\dots w_{N}}$ and an on-shell state $\ket{v_{1}\dots v_{N}}$ can now be written as

[TABLE]

The overlap between a dual state and a normal state is simply given by

[TABLE]

which can be interpreted as the overlap of a Bethe state defined by $L$ rapidities $\{x_{\alpha}\}=\{x_{1}\dots x_{L}\}=\{v^{\prime}\}\cup\{w\}$ with the dual vacuum $\ket{\uparrow\dots\uparrow}$,

[TABLE]

which is also known as a domain wall boundary partition function (DWPF), as introduced by Korepin [67]. This has a graphical interpretation and is illustrated in Figure 2.

For the models at hand, the expression (68) is exactly the definition of the permanent of a matrix

[TABLE]

with $C$ the Cauchy matrix defined as

[TABLE]

Using the results on Cauchy matrices from the previous section (51), this permanent can be rewritten as a determinant of

[TABLE]

with $J_{\epsilon}$ an $L\times L$ matrix defined as

[TABLE]

We can now simply reintroduce the rapidities $\{x_{1}\dots x_{L}\}=\{v^{\prime}\}\cup\{w\}$ in order to write

[TABLE]

where, because of the correspondence between the eigenvalues of the original state and the dual state (62), we can write this entirely in terms of the original rapidities to obtain

[TABLE]

resulting in the proposed determinant expression (29) for the inner product if we take into account the ratio of normalizations for the original state and the dual state (64). Note that the existence of the dual state was necessary to derive these determinant expressions, but is only implicit in the final results due to the demand that one of the two states in the inner product must be on-shell.

At this point, the first proposed expression (29) has been derived, which is expressed in the eigenvalue-based variables. In order to derive the equivalent expression (30) for the rapidities, the crucial observation is that the matrix (74) has the structure $(\mathbb{1}_{L}+\cdots)$ after multiplying all matrix elements with $g/2$. This leads to the exact same structure obtained when applying Sylvester’s determinant identity in Eq. (53). The involved matrices are now defined in terms of two sets of variables of unequal size, being the $L$ inhomogeneities $\{\epsilon_{1}\dots\epsilon_{L}\}$ and the $2N$ combined rapidities $\{v_{1}\dots v_{N}\}\cup\{w_{1}\dots w_{N}\}$. Applying the relation (53) then connects the determinant of the $L\times L$ matrix $(\mathbb{1}_{L}+\cdots)$ to that of an $2N\times 2N$ matrix $(\mathbb{1}_{2N}+\cdots)$, equalling the matrix proposed in (30) after correcting for the missing prefactors of $(g/2)$ in the matrix elements by changing the prefactor of the determinant.

Afterwards, the related Gaudin and Izergin-Borchardt determinants follow immediately by taking the appropriate limits.

5.3 Reduction to Slavnov’s determinant

Starting from the eigenvalue-based inner product, it is possible to rederive Slavnov’s determinant expression, obtaining other determinant expressions in the process and shedding some light on the structure of all involved identities.

Firstly, by making using of the Bethe equations (18) the well-known Slavnov’s determinant expression (25) can be straightforwardly rewritten as

[TABLE]

[TABLE]

Note that, by rewriting it in this way, the only explicit dependence on the inhomogeneities is through $\sum_{i}\frac{1}{w_{b}-\epsilon_{i}}$, bringing to mind the eigenvalue-based determinants and Eq. (30). The prefactor is the inverse of the determinant of the Cauchy matrix $U$ defined by

[TABLE]

Defining a diagonal matrix $\Lambda$ in order to absorb the dependency of this matrix on the inhomogeneities as

[TABLE]

the matrix in Slavnov’s determinant can be decomposed as

[TABLE]

This can be summarized in

[TABLE]

We will now work towards this expression starting from the eigenvalue-based expressions using the special properties of Cauchy matrices. The inner product as given in Eq. (30) depends on a $2N\times 2N$ matrix, which can be interpreted as a $2\times 2$ block matrix of $N\times N$ matrices as

[TABLE]

with

[TABLE]

and the diagonal matrices are given by

[TABLE]

The dependence of $J_{w}$ on the inhomogeneities can be absorbed in the same diagonal matrix $\Lambda$ as defined for Slavnov’s determinant and the off-diagonal matrices are then determined by the same Cauchy matrix defined previously. Using the Richardson-Gaudin equations (18) for $\{v_{a}\}$ in the diagonal elements of $J_{v}$ then results in

[TABLE]

By rewriting the matrix in this way, all submatrices exhibit the structures previously introduced since $J_{v}\sim U^{-1}(U*U)$ and $J_{w}\sim(U*U)U^{-1}$. One final step now consists of relating the determinant of a $2N\times 2N$ matrix to that of an $N\times N$ one. For this, it is possible to evaluate the determinant of a $2\times 2$ block matrix as

[TABLE]

Applying this to the equation for the inner product yields

[TABLE]

where we have used $\det(J_{v})=\det(U*U)/\det(U)$ and Eq. (55) to evaluate the inverse of $J_{v}$ as

[TABLE]

This corresponds exactly to Slavnov’s determinant expression. Alternative ways of calculating the determinant of the block matrix would result in alternative $N\times N$ matrices, but their structure is more involved than that of the Slavnov determinant. Again, the on-shell requirement for one of the two states is crucial, arising here due to the implied existence of a dual state.

In conclusion, there are three ways of evaluating the inner product - by calculating the determinants of $L\times L$, $2N\times 2N$, or $N\times N$ matrices. The structure of all these matrices is intricately related to Cauchy matrices, where the $2N\times 2N$ matrix can be reduced to the $N\times N$ Slavnov determinant, while the $L\times L$ matrix follows from the eigenvalue-based framework. Such results were also obtained in [11] and [68] for the rational six-vertex model, connecting DWPFs with Slavnov determinants, however without invoking the dual representation of Bethe states. Alternatively, similar results were obtained by Kitanine et al. for the integrable XXX Heisenberg chain through both a separation of variables and an Algebraic Bethe Ansatz approach [69].

6 Extension to hyperbolic models

The presentation thus far was focused on the rational (XXX) Richardson-Gaudin models, because of their clear-cut connection with Cauchy matrices. However, it is possible to extend all results to the hyperbolic (XXZ) Richardson-Gaudin case, for which we will only provide an overview and refer to Appendix D for a detailed derivation. Such models are a generalization of the rational (or XXX) model treated so far, often arising in the context of $p+ip$ superconductivity, where they are known to exhibit a richer phase diagram [70, 71, 72, 32, 73].

Conserved charges and Bethe ansatz

The class of XXZ Richardson-Gaudin models are defined by a similar set of commuting operators222We follow the presentation from [73, 23], which differs from the usual one in the antisymmetry in the last term. However, this allows for a clearer presentation and can be connected to the fact that these models can also be constructed from a non-skew-symmetric $r$-matrix [74, 75, 76]. given by

[TABLE]

The common eigenstates of these operators are given by

[TABLE]

with eigenvalues

[TABLE]

provided the rapidities $\{v_{1}\dots v_{N}\}$ satisfy the equations

[TABLE]

Eigenvalue-based equations

The conserved charges of these models again satisfy a set of similar quadratic equations [24]

[TABLE]

Remarkably, the eigenvalues of the conserved charges are determined by eigenvalue-based variables [21] defined in the same way as

[TABLE]

where the eigenvalue-based variables now follow from the quadratic equations

[TABLE]

The connection between this quadratic equation and the original Bethe equations can then be made through the differential equation

[TABLE]

Inner products

The normalizations of the original and the dual state are now related through (for $L-2N>0$)

[TABLE]

This is a more involved expression compared to the rational model (64), which can be seen as a consequence of the various symmetries and dualities present in this model - for specific values $g^{-1}\in\mathbb{Z}$ the dual wave function vanishes because some of the dual rapidities diverge [77]. Similar terms also arise in the investigation of the phase diagram of these models [70, 72, 78]. The origin of this prefactor is shown in Appendix D and will be further discussed after the presentation of the inner products.

From this, it is possible to obtain a similar set of determinant expressions for inner products in the hyperbolic model, and we find (for on-shell $\{v_{a}\}$)

[TABLE]

with $J_{L}(\{v_{a}\},\{w_{b}\})$ an $L\times L$ matrix given by

[TABLE]

The diagonal elements of these matrices again reflect the structure of the eigenvalue-based equations, which can also be used to show the orthogonality of different eigenstates. When calculating the normalization of an on-shell state, this reduces to the Jacobian of the eigenvalue-based equations (94), where all derivatives are now taken w.r.t. $\epsilon_{i}\Lambda_{i}$ instead of $\Lambda_{i}$. Alternatively, there is again a direct equivalence with the determinant of a $2N\times 2N$ matrix defined in terms of $\{x_{1}\dots x_{2N}\}=\{v_{1}\dots v_{N}\}\cup\{w_{1}\dots w_{N}\}$ as

[TABLE]

with

[TABLE]

where the diagonal elements now clearly reflect the Bethe equations (91).

The Slavnov determinant for the inner product of an on-shell and an off-shell state in the XXZ model [32] then follows from evaluating this matrix as a block matrix, where all steps in the derivation are fundamentally the same as for the rational XXX model. The resulting expression is given by

[TABLE]

with

[TABLE]

Discussion

As mentioned before, the prefactor in Eq. (96) reflects the various symmetries and dualities present in the hyperbolic model. Since it arises from the ratio of the normalizations of the Bethe state and its dual representation, it can be connected to the behaviour of the dual rapidities. In fact, it is well-established that (dual) rapidities in the XXZ model are allowed to diverge or collapse to zero at specific values of $g^{-1}\in\mathbb{Z}$, and how this reflects additional symmetries in the model at these points [70, 72, 78, 77]. For diverging dual rapidities, the normalization of the dual state will vanish, resulting in a diverging prefactor in Eq. (97). However, this is not a fundamental problem, since this is controllable via a proper normalization of the dual state. It can also easily be checked that this leads to a vanishing determinant, resulting in a well-conditioned finite value of the inner product.

This behaviour can be better understood by making the connection between the dual rapidities and the original rapidities explicit. For the hyperbolic model, the correspondence between the two sets of rapidities is given by [21] as

[TABLE]

There are now two cases where the dual state vanishes. First, the case where some of the $v_{a}=0$ is known to occur only at the so-called Read-Green points [70, 72, 78]. At these points $g^{-1}=-p$, with $p=1\dots N$ the number of vanishing rapidities. The dual rapidities are then given by the $N-p$ non-zero rapidities of the original state supplemented with $L-2N+p$ diverging rapidities. This can easily be checked as

[TABLE]

where the diverging rapidities do not contribute to the summation and the zero rapidities result in a term $p/\epsilon_{i}$. This hampers a straightforward evaluation of all involved determinants, not because of the prefactor, but because it necessitates taking the appropriate limit where multiple rapidities coincide. This is a common difficulty in all proposed determinant expressions, requiring a more in-depth analysis falling outside the scope of the current work.

Secondly, the other case where the prefactor might pose trouble is when $g^{-1}=+p$, with $p=0\dots L-2N-1$. Here, all rapidities of the original state are known to be finite and non-zero, and the dual rapidities are given by the $N$ rapidities of the original state, supplemented with $L-2N-p$ diverging dual rapidities and $p$ dual rapidities collapsing to zero. It can again be easily checked that this leads to the correct number of dual rapidities $N+(L-2N-p)+p=L-N$ and

[TABLE]

In this case, Eqs. (100) and (101) can be evaluated without difficulties since no rapidities coincide.

The diverging prefactor hence reflects the divergence of a subset of dual rapidities at specific and predictable values of $g^{-1}$ where the XXZ model exhibits an additional symmetry [77], allowing us to relate the dual rapidities to the original rapidities.

7 Conclusion

Integrable Richardson-Gaudin models can be investigated in two distinct ways - either by solving the Bethe equations for the rapidities, or by solving a set of quadratic Bethe equations for the eigenvalues of the conserved charges. Similarly, inner products of eigenstates of these integrable models can be expressed in two distinct ways, as determinants of matrices depending on either the rapidities or the eigenvalues. As a direct consequence, the normalization of such eigenstates can be calculated as the determinant of either the Jacobian of the Bethe equations or the Jacobian of the quadratic Bethe equations.

In this work, we investigated the connection between both approaches and presented two classes of determinant expressions for the inner product between on-shell and off-shell states in these integrable models. Slavnov’s determinant is then obtained as a special case of one of the two classes of determinants. The crucial element in this connection is the existence of two distinct representations for each eigenstate, allowing inner products to be recast in domain wall boundary partition functions, and the structure of all involved matrices in terms of Cauchy matrices, which has been made explicit here.

Acknowledgements

P.W.C. acknowledges support from a Ph.D. fellowship and a travel grant for a long stay abroad at the University of Amsterdam from the Research Foundation Flanders (FWO Vlaanderen). P.W.C. thanks J.-S. Caux for valuable discussions and for his hospitality at the University of Amsterdam.

Appendix A Orthogonality of Bethe states

Starting from Eq. (28) or Eq. (30), the orthogonality of different Bethe states can easily be shown. Suppose we have two non-equal on-shell states $\{v_{a}\}$ and $\{w_{b}\}$. Then it is straightforward to show that the rows of (29) are linearly dependent.

[TABLE]

which vanish since these are exactly the eigenvalue-based equations satisfied by on-shell states (21). The linear dependence of the rows implies that the determinant of this matrix also vanishes, proving the orthogonality of different eigenstates

The orthogonality also follows from the similarity of the diagonal elements of (30) to the Bethe equations (18). Taking $\{x_{1}\dots x_{N}\}=\{v_{1}\dots v_{N}\}$ and $\{x_{N+1}\dots x_{2N}\}=\{w_{1}\dots w_{N}\}$, this results in

[TABLE]

where we have identified $w_{b}$ and $w_{N+b}$ in the second line. Since these are exactly the Bethe equations (18), the rows are again linearly dependent, leading to a vanishing determinant and orthogonal eigenstates.

Appendix B Properties of Cauchy matrices

B.1 Inverse of Cauchy matrices

In this Appendix, proofs will be provided for the properties of Cauchy matrices used throughout the main text. For this, the crucial element is the explicit expression for the inverse of a Cauchy matrix. Starting from a Cauchy matrix defined by two sets of variables $\{\epsilon_{1}\dots\epsilon_{N}\}$ and $\{x_{1}\dots x_{N}\}$ as

[TABLE]

the inverse is given by (see e.g. [66])

[TABLE]

with $p(x)=\prod_{i}(x-\epsilon_{i})$ and $q(x)=\prod_{\alpha}(x-x_{\alpha})$. Explicitly writing out the action of the inverse on the Cauchy matrix results in the identity

[TABLE]

which is a necessary result for the following. Note that this equality also holds if the sets $\{\epsilon_{i}\}$ contains less variables than the set $\{x_{\alpha}\}$, since this can be obtained by simply sending an appropriate amount of $\epsilon_{i}$ to infinity.

B.2 First order

In the notation of Section 4, the matrix elements of $(C*C)C^{-1}$ can be explicitly calculated. The off-diagonal elements ($i\neq j$) result in

[TABLE]

Here, the properties of the inverse as written out in Eq. (110) have been used twice in order to evaluate the summations, first for a vanishing summation in Eq. (112) to Eq. (113), making use of $i\neq j$, and later in order to evaluate the non-zero summation when moving from Eq. (113) to Eq. (114). The diagonal elements can similarly be calculated as

[TABLE]

This can be simplified by performing a partial fraction decomposition in $\epsilon_{i}$, which has single poles in $x_{\alpha}$ and $\epsilon_{j}$ $(j\neq i)$. This results in

[TABLE]

where again only the properties of the inverse of the Cauchy matrix as in Eq. (110) have been used to evaluate the single summation. This then leads to

[TABLE]

with

[TABLE]

and $D_{\epsilon}$ a diagonal matrix with $(D_{\epsilon})_{ii}=q(\epsilon_{i})/p^{\prime}(\epsilon_{i})$. Equivalently, changing the order of the product leads to

[TABLE]

with

[TABLE]

and $D_{x}$ a diagonal matrix with $(D_{x})_{\alpha\alpha}={p(x_{\alpha})}/{q^{\prime}(x_{\alpha})}$.

B.3 Second order

The first-order relations (51) can be extended to higher-order Hadamard products, where the relevant expression for our purpose is

[TABLE]

By using the definitions of $J_{x}$ and $J_{\epsilon}$, relating the inverse of $C$ to its transposed (44), and multiplying with $C*C$, this can again be expressed purely in terms of Cauchy matrices and diagonal matrices as

[TABLE]

So, written out explicitly, the equality to be proven is

[TABLE]

where the previously obtained expression for $(C*C)C^{-1}$ has been inserted. The second term on the right-hand side can be expanded as

[TABLE]

simplifying this summation to

[TABLE]

where again only the properties of the inverse of the Cauchy matrix have been used. Plugging this into the original equation and multiplying with $(\epsilon_{i}-x_{\alpha})^{2}$, the equality to be proven is

[TABLE]

This can be done by doing a partial fraction expansion of the left-hand side in $\epsilon_{i}$, which has single poles in $\epsilon_{i}=x_{\beta},\ \forall\beta$ and $\epsilon_{i}=\epsilon_{j},\ j\neq i$. This calculation is highly reminiscent of the one used in the previous identities, resulting in

[TABLE]

thus proving the identity.

B.4 Sylvester’s determinant identity

In Section 4, it was already mentioned how the factorization of the $J$-matrices leads to the determinant identity

[TABLE]

with $J_{x}$ and $J_{\epsilon}$ two $N\times N$ matrices defined according to Eqs. (119) and (121) in terms of two sets of variables $\{\epsilon_{1}\dots\epsilon_{N}\}$ and $\{x_{1}\dots x_{N}\}$. This can now be generalized towards an $L\times L$ matrix $J_{\epsilon}$ and an $N\times N$ matrix $J_{x}$ defined in terms of differently-sized sets $\{\epsilon_{1}\dots\epsilon_{L}\}$ and $\{x_{1}\dots x_{N}\}$, leading to

[TABLE]

This can be done in two ways. The first is by generalizing the factorization properties (118) and (120) towards non-square Cauchy matrices, where $C^{-1}$ is no longer the inverse, but can still be explicitly defined through Eq. (109). The calculation here is the exact same one as for square matrices, since all involved identities hold for differently-sized sets of variables.

Another way is by starting from two sets of equal size, e.g. $L$ if $L>N$, and taking the limit where the variables $\{x_{N+1}\dots x_{L}\}$ go to infinity as $x_{N+1}\ll x_{N+2}\ll\dots\ll x_{L}$. These then drop out of the summation in the diagonal elements, their off-diagonal elements vanish, and the matrix identity for equal sizes reduces to

[TABLE]

Appendix C Normalizations of the normal and the dual state

From the identities presented in Appendix B, it is possible to derive the ratio of the normalizations of a dual Bethe state and a regular Bethe state. Since these are both eigenstates, we already know that they describe the same state, but with different normalizations. We now calculate the ratio of the normalization of these states by taking the overlap of both states with an (arbitrary) reference state, and show how it leads to our proposed expression (64).

To reiterate, a Bethe ansatz eigenstate for the Richardson-Gaudin models can be written in two representations as

[TABLE]

with the (dual) rapidities satisfying the Richardson-Gaudin equations

[TABLE]

describing the same eigenstate if these variables are coupled through

[TABLE]

The crucial element in this proof is the existence of eigenvalue-based expressions for the inner product of a Bethe state and a reference state [12]. Defining a reference state from a set of occupied levels can again be done in two ways

[TABLE]

where both states are already normalized. Then the overlap of the original Bethe state with this reference state is given by

[TABLE]

with $J_{N}(\{v_{a}\},\{i_{occ}\})$ an $N\times N$ matrix defined as

[TABLE]

Alternatively, the overlap of the dual state with the same reference state can be written as

[TABLE]

This final expression is the complex conjugate of the inner product of a regular Bethe state defined by the dual rapidities and a different reference state. Because all eigenvalue-based variables are real, these inner products are always real, and we can use the same determinant expressions to write

[TABLE]

with $J_{L-N}(\{v_{a}^{\prime}\},\{i|i\notin\{i_{occ}\}\})$ an $(L-N)\times(L-N)$ matrix defined as

[TABLE]

Because this only depends on the dual rapidities through the eigenvalue-based variables in the diagonal elements, the correspondence (135) can be used to express this determinant in the rapidities of the original state as

[TABLE]

Now Sylvester’s determinant identity can be used, as explained in Section 4, to exchange the role of rapidities and inhomogeneities and obtain

[TABLE]

with $K_{N}(\{v_{a}\},\{i|i\notin\{i_{occ}\}\})$ an $N\times N$ matrix defined as

[TABLE]

Here, it is important to note the similarity between the diagonal elements of this matrix and the Richardson-Gaudin equations. By slightly rewriting the Richardson-Gaudin equations, these diagonal elements can be rewritten as

[TABLE]

resulting in

[TABLE]

Both $g/2$ and the minus sign in the diagonal elements have been absorbed in a prefactor, since multiplying all rows with $-1$ and taking the transpose of the matrix returns the original matrix but with the sign of the diagonal elements exchanged. The roles of both variables can again be exchanged in this final expression to obtain the original equality

[TABLE]

Since the set of reference states forms a complete basis and this equality holds for all such states, this yields the proposed expression (64) as

[TABLE]

Appendix D Results for the hyperbolic model

Applying the matrix determinant lemma

All derivations for the hyperbolic model are essentially the same as those for the rational model, except for one additional result. From the matrix determinant lemma, it can be shown using the same arguments as in Appendix B that

[TABLE]

for two sets of equal size $\{\epsilon_{1}\dots\epsilon_{L}\}$ and $\{x_{1}\dots x_{L}\}$ and arbitrary $G\in\mathbb{C}$, where we again have $J_{\epsilon}$ an $L\times L$ as an matrix defined as

[TABLE]

and $J_{x}$ also an $L\times L$ matrix defined as

[TABLE]

In this expression $1/\epsilon$ and $1/x$ are shorthand for $L\times L$ diagonal matrices with diagonal elements $1/\epsilon_{i}$ and $1/x_{\alpha}$ respectively. This follows from the version of the matrix determinant lemma relating

[TABLE]

with $A$ an invertible $n\times n$ matrix and $U,V$ are $n\times m$ matrices. This reduces to Sylvester’s determinant identity in the special case $A=\mathbb{1}$, and here it is applied by setting $A=1/\epsilon$, $U=C*C$ and $V^{T}=C^{-1}$. Since $A$ is a diagonal matrix, both its inverse and its determinant can be explicitly calculated, and all resulting expressions are highly similar to those presented in Appendix B, and will hence not be repeated here.

The limit where the two sets do not contain an equal amount of variables, leading to two matrices of unequal dimensions, can again be obtained by sending the appropriate amount of variables to infinity, similar to the techniques used in [22]. As will be shown explicitly, starting from two sets containing $L$ variables and taking the limit where $L-N$ variables $x_{\alpha}$ go subsequently to infinity, each diverging variable adds a prefactor due to its appearance in the diagonal elements, giving rise to

[TABLE]

The origin of the prefactors follows from the diagonal elements of the right-hand side and will be derived in full detail, since these are not as trivial as in the rational case. Firstly, starting from Eq. (149), the limit where a single variable is sent to infinity can be obtained. Selecting $x_{L}\to\infty$ while all other variables remain finite, both sides of Eq. (149) can be evaluated. In the left-hand side, the only dependence on $x_{L}$ is in the diagonal elements, which reduce to

[TABLE]

where the limit $x_{L}\to\infty$ simply removes a single variable from the summation. In order to evaluate the right-hand side, the term $x_{L}$ from the prefactor can be absorbed in the determinant by multiplying row and column $L$ of the matrix with $x_{L}^{1/2}$. Taking the limit to infinity, the off-diagonal elements in this row and column will vanish, since these are given by $\pm x_{L}^{1/2}/(x_{L}-x_{\alpha})\sim x_{L}^{-1/2}\to 0$. The diagonal element in row/column $L$ is then given by

[TABLE]

while the other diagonal elements will be given by

[TABLE]

where the term containing $x_{L}$ again drops out of the summation. Since the off-diagonal elements vanish, the diagonal element will simply contribute a prefactor in the determinant, and the evaluation of the determinant of the $L\times L$ matrix reduces to that of a similar $(L-1)\times(L-1)$ matrix with an additional prefactor $(G+1)/2-1$. The effect of multiply diverging variables can be obtained by subsequently sending $x_{L-1}$ to infinity and again absorbing the prefactor $x_{L-1}$ in the matrix. The off-diagonal elements in the row/column $(L-1)$ will again vanish, and the related diagonal element contributes a prefactor

[TABLE]

This can be repeated until $N$ variables $\{x_{1}\dots x_{N}\}$ remain, resulting in the proposed expression

[TABLE]

with $J_{\epsilon}$ an $L\times L$ matrix again defined as

[TABLE]

and $J_{x}$ an $N\times N$ matrix as

[TABLE]

Normalizations of the normal and the dual state

The derivation from Appendix C can now be repeated for the hyperbolic model, where we will pay special attention to the points where it deviates from the rational model. In the hyperbolic case the Bethe ansatz eigenstate and its dual representation are given by

[TABLE]

with the rapidities satisfying the Richardson-Gaudin equations

[TABLE]

where the correspondence between both is given by[21]

[TABLE]

A reference state can again be introduced as

[TABLE]

and the overlap of the original Bethe state with this reference state is given by[21]

[TABLE]

with $J_{N}(\{v_{a}\},\{i_{occ}\})$ an $N\times N$ matrix defined as

[TABLE]

This is constructed in the exact same way as for the rational model, with only an additional prefactor because of the terms $\sqrt{\epsilon_{i}}$ in the Bethe state. Similarly, the overlap of the dual state with the same reference state is given by

[TABLE]

with $J_{L-N}(\{v_{a}^{\prime}\},\{i|i\notin\{i_{occ}\}\})$ an $(L-N)\times(L-N)$ matrix defined as

[TABLE]

The eigenvalue-based variables can again be identified in the diagonal elements, and the correspondence (165) then leads to

[TABLE]

where $1/\epsilon$ here stands for the diagonal matrix with matrix elements $\epsilon_{i},i\in\{i_{occ}\}$. This is where the matrix determinant lemma (149) comes into play, this time exchanging the role of rapidities and inhomogeneities as

[TABLE]

with $K_{N}(\{v_{a}\},\{i|i\notin\{i_{occ}\}\})$ an $N\times N$ matrix defined as

[TABLE]

and $1/v$ shorthand for a diagonal matrix with matrix elements $1/v_{a}$. Remarkably, the diagonal elements again resemble the Richardson-Gaudin equations (163), although now those for the hyperbolic model, and can be rewritten as

[TABLE]

leading to

[TABLE]

As for the rational model, the minus sign in the diagonal elements has been absorbed in the prefactor. As a final step, the roles of both variables can be exchanged by using Sylvester’s determinant identity, in order to obtain

[TABLE]

or, equivalently,

[TABLE]