{"id":2178,"date":"2025-08-13T22:07:54","date_gmt":"2025-08-13T22:07:54","guid":{"rendered":"https:\/\/www.ethanepperly.com\/?p=2178"},"modified":"2025-08-13T22:29:56","modified_gmt":"2025-08-13T22:29:56","slug":"vandermonde-matrices-are-merely-exponentially-ill-conditioned","status":"publish","type":"post","link":"https:\/\/www.ethanepperly.com\/index.php\/2025\/08\/13\/vandermonde-matrices-are-merely-exponentially-ill-conditioned\/","title":{"rendered":"Vandermonde Matrices are Merely Exponentially Ill-Conditioned"},"content":{"rendered":"\n

I am excited to share that my paper Does block size matter in randomized block Krylov low-rank approximation?<\/a><\/em> has recently been released on arXiv. In that paper, we study the randomized block<\/a> Krylov iteration<\/a> (RBKI) algorithm for low-rank approximation<\/a>. Existing<\/a> results<\/a> show<\/a> that RBKI is efficient at producing rank-\"k\" approximations with a large block size of \"k\" or a small block size of \"1\", but these results give poor results for intermediate block sizes \"1\ll b \ll k\". But often<\/a> these intermediate block sizes are the most efficient in practice. In our paper, we close this theoretical gap, showing RBKI is efficient for any block size \"1\le b \le k\". Check out the paper<\/a> for details!<\/p>\n\n\n\n

In our paper, the core technical challenge is understanding the condition number of a random block Krylov matrix of the form

  <\/span>   <\/span>\"\[K = \begin{bmatrix} G & AG & \cdots & A^{t-1}G \end{bmatrix},\]\"<\/p>where \"A\" is an \"n\times n\" real symmetric<\/a> positive semidefinite<\/a> matrix and \"G\" is an \"n\times (n/b)\" random Gaussian<\/a> matrix. Proving this result was not easy, and our proof required several ingredients. In this post, I want to talk about just one of them: Gautschi’s bound on the conditioning of Vandermonde matrices<\/a>. (Check out Gautschi’s original paper here<\/a><\/strong>.)<\/p>\n\n\n\n

Evaluating a Polynomial and Vandermonde Matrices<\/h2>\n\n\n\n

Let us begin with one the humblest but most important characters in mathematics, the univariate polynomial

  <\/span>   <\/span>\"\[p_a(u) = a_1 + a_2 u+ \cdots + a_t u^{t-1}.\]\"<\/p>Here, we have set the degree to be \"t-1\" and have written the polynomial \"p_a(\cdot)\" parametrically in terms of its vector of coefficients \"a = (a_1,\ldots,a_t)\". The polynomials of degree at most \"t-1\" form a linear space<\/a> of dimension \"t\", and the monomials \"1,u,\ldots,u^{t-1}\" provide a basis<\/a> for this space. In this post, we will permit the coefficients \"a \in \complex^t\" to be complex numbers.<\/p>\n\n\n\n

Given a polynomial, we often wish to evaluate it at a set of inputs. Specifically, let \"\lambda_1,\ldots,\lambda_s \in \complex\" be \"s\" (distinct<\/mark><\/em>) input locations<\/em>. If we evaluate \"p_a\" at each number, we obtain a list of (output) values, which we denote by \"f = (p_a(\lambda_1),\ldots,p_a(\lambda_s))\" of \"s\" (distinct) values, each of which given by the formula

  <\/span>   <\/span>\"\[p_a(\lambda_i) = \sum_{j=1}^t \lambda_i^{j-1} a_j.\]\"<\/p>Observe that the outputs \"f\" are a nonlinear<\/em> function of the input values \"\lambda\" but a linear<\/em> function of the coefficients \"a\". We may call the mapping \"a \mapsto f\" the coefficients-to-values map.<\/em><\/p>\n\n\n\n

Every linear transformation between vectors can<\/a> be realized as a matrix\u2013vector product, and the matrix for the coefficients-to-values map is called a Vandermonde matrix<\/a><\/em> \"V\". It is given by the formula

  <\/span>   <\/span>\"\[V= \begin{bmatrix} 1 & \lambda_1 & \lambda_1^2 & \cdots & \lambda_1^{t-1} \\ 1 & \lambda_2 & \lambda_2^2 & \cdots & \lambda_2^{t-1} \\ 1 & \lambda_3 & \lambda_3^2 & \cdots & \lambda_3^{t-1} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 1 & \lambda_s & \lambda_s^2 & \cdots & \lambda_s^{t-1} \end{bmatrix} \in \complex^{s\times t}.\]\"<\/p>The Vandermonde matrix defines the coefficients-to-values map in the sense that \"f = Va\".<\/p>\n\n\n\n

Interpolating by a Polynomial and Inverting a Vandermonde Matrix<\/h2>\n\n\n\n

Going forward, let us set \"s = t\" so the number of locations \"\lambda = (\lambda_1,\ldots,\lambda_t)\" equals the number of coefficients \"a = (a_1,\ldots,a_t)\". The Vandermonde matrix maps the vector of coefficients \"a\" to the vector of values \"f = (f_1,\ldots,f_t) = (p_a(\lambda_1),\ldots,p_a(\lambda_t))\". Its inverse \"V^{-1}\" maps a set of values \"f\" to a set of coefficients \"a\" defining a polynomial \"p_a\" which interpolates the values \"f\":

  <\/span>   <\/span>\"\[p_a(\lambda_i) = f_i \quad \text{for } i =1,\ldots,t .\]\"<\/p>More concisely, multiplying the inverse of a Vandermonde matrix solves the problem of polynomial interpolation<\/a><\/em>.<\/p>\n\n\n\n

To solve the polynomial interpolation problem with Vandermonde matrices, we can do the following. Given values \"f\", we first solve the linear system of equations \"Va = f\", obtaining a vector of coefficients \"a = V^{-1}f\". Then, define the interpolating polynomial

(1) <\/span>   <\/span>\"\[q(u) = a_1 + a_2 u + \cdots + a_t u^{t-1} \quad \text{with } a = V^{-1} f. \]\"<\/p>The polynomial \"q\" interpolates the values \"f\" at the locations \"\lambda_i\", \"q(\lambda_i) = p_a(\lambda_i) = f_i\".<\/p>\n\n\n\n

Lagrange Interpolation<\/h2>\n\n\n\n

Inverting a the Vandermonde matrix is one way to solve the polynomial interpolation problem, but the polynomial interpolation can also be solved directly. To do so, first notice that we can construct a special polynomial \"(u - \lambda_2)(u - \lambda_3)\cdots(u-\lambda_t)\" that is zero at the locations \"\lambda_2,\ldots,\lambda_t\" but nonzero at the first location \"\lambda_1\". (Remember that we have assumed that \"\lambda_1,\ldots,\lambda_t\" are distinct<\/mark><\/em>.) Further, by rescaling this polynomial to

  <\/span>   <\/span>\"\[\ell_1(u) = \frac{(u - \lambda_2)(u - \lambda_3)\cdots(u-\lambda_t)}{(\lambda_1 - \lambda_2)(\lambda_1 - \lambda_3)\cdots(\lambda_1-\lambda_t)},\]\"<\/p>we obtain a polynomial whose value at \"\lambda_1\" can be set to \"1\". Likewise, for each \"i\", we can define a similar polynomial

(2) <\/span>   <\/span>\"\[\ell_i(u) = \frac{(u - \lambda_1)\cdots(u-\lambda_{i-1}) (u-\lambda_{i+1})\cdots(u-\lambda_t)}{(\lambda_i - \lambda_1)\cdots(\lambda_i-\lambda_{i-1}) (\lambda_i-\lambda_{i+1})\cdots(\lambda_i-\lambda_t)}, \]\"<\/p>which is \"1\" at \"\lambda_i\" and zero at \"\lambda_j\" for \"j\ne i\". Using the Dirac delta symbol, we may write

  <\/span>   <\/span>\"\[\ell_i(\lambda_j) = \delta_{ij} = \begin{cases} 1, & i = j, \\ 0, & i\ne j. \end{cases}\]\"<\/p>The polynomials \"\ell_i\" are called the Lagrange polynomials<\/a><\/em> of the locations \"\lambda_1,\ldots,\lambda_q\". Below is an interactive illustration of the second Lagrange polynomial \"\ell_2\" associated with the points \"\lambda_i = i\" (with \"t = 5\").<\/p>\n\n\n\n