A Matlab package computing simultaneous Gaussian quadrature rules for multiple orthogonal polynomials

The aim of this paper is to describe a Matlab package for computing the simultaneous Gaussian quadrature rules associated with a variety of multiple orthogonal polynomials. Multiple orthogonal polynomials can be considered as a generalization of classical orthogonal polynomials, satisfying orthogonality constraints with respect to r different measures, with r >= 1 . Moreover, they satisfy ( r + 2)-term recurrence relations. In this manuscript, without loss of generality, r is considered equal to 2 . The so-called simultaneous Gaussian quadrature rules associated with multiple orthogonal polynomials can be computed by solving a banded lower Hessenberg eigenvalue problem. Unfortunately, computing the eigendecomposition of such a matrix turns out to be strongly ill-conditioned and the Matlab function balance.m does not improve the condition of the eigenvalue problem. Therefore, most procedures for computing simultaneous Gaussian quadrature rules are implemented with variable precision arithmetic. Here, we propose a Matlab package that allows to reliably compute the simultaneous Gaussian quadrature rules in floating point arithmetic. It makes use of a variant of a new balancing procedure, recently developed by the authors of the present manuscript, that drastically reduces the condition of the Hessenberg eigenvalue problem.