  This paper develops a method to carry out the large-$N$ asymptotic analysis
of a class of $N$-dimensional integrals arising in the context of the so-called
quantum separation of variables method. We push further ideas developed in the
context of random matrices of size $N$, but in the present problem, two scales
$1/N^{\alpha}$ and $1/N$ naturally occur. In our case, the equilibrium measure
is $N^{\alpha}$-dependent and characterised by means of the solution to a
$2\times 2$ Riemann--Hilbert problem, whose large-$N$ behavior is analysed in
detail. Combining these results with techniques of concentration of measures
and an asymptotic analysis of the Schwinger-Dyson equations at the
distributional level, we obtain the large-$N$ behavior of the free energy
explicitly up to $o(1)$. The use of distributional Schwinger-Dyson is a novelty
that allows us treating sufficiently differentiable interactions and the mixing
of scales $1/N^{\alpha}$ and $1/N$, thus waiving the analyticity assumptions
often used in random matrix theory.
