The relation between integrable systems and orthogonal polynomials goes back to the seminar paper of Moser [1]}. There, the author made the remarkable observation that, given a family of orthogonal polynomials associated to a measure \(w\)  on the real axis, the related tridiagonal Jacobi matrix can be interpreted as a Lax matrix for the Toda equations, where the independent variable of equations plays the role of deformation parameters for the measure \(w\) . Starting from the nineties, this result and various generalizations (to orthogonal polynomials on the unit circle, to multiple orthogonal polynomials and many other types of special polynomials) played a central role in the field of random matrices, and the study of the equations satisfied by orthogonal polynomials found applications to 2D quantum gravity [2]} and to the computations of gap probabilities and partition functions associated to random models, both in the discrete and continuous setting (see, for instance, [3]} and [4]}).
