Matrices over a Hilbert space and their low-rank cross approximation2


Abstract

Motivated by applications in reduced-order modeling (ROM) of parametric partial differential equations, we investigate the algebraic properties of Bochner matrices—matrices with entries in an abstract Hilbert space. Low-rank cross approximation is extended to Bochner matrices and its approximation guarantees are derived. The high-dimensional nature of the entries is shown to manifest itself in maximum-volume bounds, making them smaller than in the classical setting. An analogue of adaptive cross approximation is proposed and validated as a non-intrusive ROM method in numerical experiments with parametric nonlinear Stokes equations.

vector-valued matrices, low-rank approximation, cross approximation, maximum volume, reduced-order modeling, parametric partial differential equations

15A69, 65D40, 65F55, 65N99

1 Introduction↩︎

The kind of entries a matrix has is determined by an application from which it arises. The most common are real and complex matrices, appearing universally across scientific computing, signal processing, and data science. Matrices over finite fields [1] are used in coding theory [2] and cryptography [3].

Matrices over commutative rings [4] have also been extensively studied. Matrices of polynomials [5], Laurent polynomials [6], and analytic functions [7], [8] appear in multichannel broadband signal processing [9]. Matrices over the convolution ring of vectors [10][12] are well-suited to represent visual data [13].

Matrices over the non-commutative ring of quaternions [14] are used in signal and image processing [15][17], while matrices of linear operators appear in the context of partial differential equations (PDEs) [18], [19].

The motivating application behind the present work, however, produces matrices whose entries form neither a field nor a ring, but rather reside in a Hilbert space.

1.1 Reduced-order modeling of parametric PDEs↩︎

Let \(\Omega \subset \mathbb{R}^{N}\) be open and bounded, and consider a Dirichlet boundary-value problem for the Poisson equation in \(\Omega\) with right-hand side \(f \in \mathrm{L}^2(\Omega)\) and varying coefficient \(\kappa \in \mathrm{L}^\infty(\Omega)\) such that \(\kappa \geq 1\) almost everywhere: \[-\mathrm{div}(\kappa \nabla u) = f~\text{ in }~\Omega, \qquad u = 0~\text{ on }~\partial\Omega.\] This problem admits a unique weak solution in the Sobolev space \(\mathrm{H}^{1}_0(\Omega)\) [20]. Suppose now that \(f = f_{\alpha,\beta}\) and \(\kappa = \kappa_{\alpha,\beta}\) depend on parameters \(\alpha,\beta \in [0,1]\); then there is a well-defined solution map \((\alpha, \beta) \mapsto u_{\alpha,\beta}\). Upon introducing grids \(\{ \alpha_i \}_{i = 1}^{m}\) and \(\{ \beta_j \}_{j = 1}^{n}\), we obtain an \(m \times n\) matrix with entries in \(\mathrm{H}^{1}_0(\Omega)\).

In reduced-order modeling (ROM), one often seeks to construct a data-driven surrogate model to approximate a quantity of interest (QoI) associated with the solution of a parametric PDE, allowing for rapid evaluation without querying the PDE solver. This essentially reduces to approximating a function \([0,1]^2 \to \mathbb{R}\) from samples, a problem with numerous established approaches. When the QoI is discretized into a matrix, cross approximation [21], [22] provides a surrogate model by inspecting a small number of rows and columns and solving the PDE only a small number of times [23].

Our central thesis is that cross approximation can be extended to matrices over a Hilbert space and used as a non-intrusive ROM method to approximate the solution map itself, rather than merely an associated scalar QoI.

1.2 Contributions and outline↩︎

We formally introduce our objects of study, Bochner matrices, in 2 and then analyze them through four distinct algebraic lenses: as linear operators (3), as quasimatrices (4), as matrices (5), and as third-order tensors (6). Following this progression, we identify how the algebraic properties of Bochner matrices emerge from their structure, providing the necessary foundation for cross approximation.

We investigate the theoretical properties of Bochner cross approximation in 7, deriving approximation guarantees based on (i) the operator norm of the pseudoinverse of submatrices of singular factors and (ii) the maximum-volume property. Notably, we prove that the maximum-volume error bounds are smaller for Bochner matrices than for classical matrices of numbers (7.3).

8 presents an extension of adaptive cross approximation [22], [24] to Bochner matrices as a greedy algorithm for index selection. We validate our algorithm, and the framework of Bochner cross approximation as a whole, as a non-intrusive ROM method in numerical experiments with parametric nonlinear Stokes equations.

Some concluding remarks are collected in 9.

1.3 Notation↩︎

We write vectors as \(\mathsf{a,b}\) and matrices as \(\mathsf{A,B}\), denoting the \(n \times n\) identity matrix by \(\mathsf{I}_n\). The transpose is written as \(\mathsf{A}^\intercal\), the complex conjugate as \(\overline{\mathsf{A}}\), and the conjugate transpose as \(\mathsf{A}^\ast\). We let \([n] = \{ 1, \ldots, n \}\). By convention, the inner product \(\langle \cdot, \cdot \rangle_{\mathrm{\mathrm{H}}}\) of a Hilbert space \(\mathrm{H}\) is linear in the second argument.

2 Bochner matrices↩︎

Consider a Hilbert space \(\mathrm{H}\) over a field \(\mathbb{F}\in \{ \mathbb{R}, \mathbb{C}\}\). Let \(d \in \mathbb{N}\) and \(n_1, \ldots, n_d \in \mathbb{N}\). We say that an \(n_1 \times \cdots \times n_d\) array with entries in \(\mathrm{H}\) is a Bochner tensor of order \(d\), which we refer to as a Bochner matrix when \(d = 2\) and a Bochner vector when \(d = 1\). Equivalently, a Bochner tensor is a map \[\boldsymbol{{\mathscr{A}}} : [n_1] \times \cdots \times [n_d] \to \mathrm{H}, \quad \mathfrak{i} = (i_1, \ldots, i_d) \mapsto \boldsymbol{{\mathscr{A}}}(i_1, \ldots, i_d) = {\boldsymbol{a}}_{\mathfrak{i}}.\] Tensors over the field \(\mathbb{F}\) will be referred to simply as “tensors,” or as “classical tensors” for additional contrast.

Under entrywise addition and scalar multiplication, the set \(\mathrm{H}^{n_1 \times \cdots \times n_d}\) of Bochner tensors is a vector space. We endow it with the family of \(\ell_p(\mathrm{H})\) norms: \[\| \boldsymbol{{\mathscr{A}}} \|_{\mathrm{\ell_p(\mathrm{H})}} = \begin{cases} \left( \sum_{\mathfrak{i} \in [n_1] \times \cdots \times [n_d]} \| {\boldsymbol{a}}_{\mathfrak{i}} \|_{\mathrm{\mathrm{H}}}^p \right)^{1/p}, & 1 \leq p < \infty, \\ \max_{\mathfrak{i} \in [n_1] \times \cdots \times [n_d]} \| {\boldsymbol{a}}_{\mathfrak{i}} \|_{\mathrm{\mathrm{H}}}, & p = \infty. \end{cases}\] These norms are pairwise equivalent, satisfying the bounds \[\| \boldsymbol{{\mathscr{A}}} \|_{\mathrm{\ell_q(\mathrm{H})}} \leq \| \boldsymbol{{\mathscr{A}}} \|_{\mathrm{\ell_{p}(\mathrm{H})}} \leq (n_1 \cdots n_d)^{\frac{1}{p} - \frac{1}{q}} \| \boldsymbol{{\mathscr{A}}} \|_{\mathrm{\ell_q(\mathrm{H})}}, \quad p \leq q,\] and satisfy Hölder’s inequality: \[|\langle \boldsymbol{{\mathscr{A}}}, \boldsymbol{{\mathscr{B}}} \rangle_{\mathrm{\ell_2(\mathrm{H})}}| \leq \| \boldsymbol{{\mathscr{A}}} \|_{\mathrm{\ell_p(\mathrm{H})}} \| \boldsymbol{{\mathscr{B}}} \|_{\mathrm{\ell_{p_\ast}(\mathrm{H})}}, \quad p_{\ast} = \tfrac{p}{p-1},\] where \(\langle \boldsymbol{{\mathscr{A}}}, \boldsymbol{{\mathscr{B}}} \rangle_{\mathrm{\ell_2(\mathrm{H})}} = \sum_{\mathfrak{i}} \langle {\boldsymbol{a}}_{\mathfrak{i}}, {\boldsymbol{b}}_{\mathfrak{i}} \rangle_{\mathrm{\mathrm{H}}}\) is the inner product inducing the \(\ell_2(\mathrm{H})\) norm.

A basic Cauchy-sequence argument shows that the vector space of Bochner tensors equipped with the \(\ell_p(\mathrm{H})\) norm is complete, and so is itself a Hilbert space when \(p = 2\). Let us remark in passing that Bochner tensors can be identified with functions from Bochner–Lebesgue spaces with respect to the discrete counting measure.

3 Bochner matrices as linear operators↩︎

An \(m \times n\) Bochner matrix \(\boldsymbol{\mathsf{A}}\) can be multiplied from the right by classical vectors in \(\mathbb{F}^n\). This induces a linear operator \(\mathsf{x} \mapsto \boldsymbol{\mathsf{A}} \mathsf{x}\) between Hilbert spaces \(\mathbb{F}^{n}\) and \(\mathrm{H}^{m}\). The operator norm of \(\boldsymbol{\mathsf{A}}\) is equivalent to the \(\ell_p(\mathrm{H})\) norms, satisfying the bounds \(\| \boldsymbol{\mathsf{A}} \|_{\mathrm{\ell_{\infty}(\mathrm{H})}} \leq \| \boldsymbol{\mathsf{A}} \|_{\mathrm{\ell_2 \to \ell_2(\mathrm{H})}} \leq \| \boldsymbol{\mathsf{A}} \|_{\mathrm{\ell_1(\mathrm{H})}}\), and thus the operator is bounded. Furthermore, it has a finite rank not exceeding \(n\), which we denote by \(\mathrm{rank} ( \boldsymbol{\mathsf{A}} )\) and that satisfies the following basic property.

Lemma 1. Let \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{m \times n}\) and \(\mathsf{B} \in \mathbb{F}^{n \times k}\). Then \[\mathrm{rank} ( \boldsymbol{\mathsf{A}} \mathsf{B} ) \leq \min\{ \mathrm{rank} ( \boldsymbol{\mathsf{A}} ), \mathrm{rank} ( \mathsf{B} ) \}, \quad \mathrm{rank} ( \boldsymbol{\mathsf{A}} \mathsf{B} ) = \begin{cases} \mathrm{rank} ( \boldsymbol{\mathsf{A}} ), & \mathrm{rank} ( \mathsf{B} ) = n, \\ \mathrm{rank} ( \mathsf{B} ), & \mathrm{rank} ( \boldsymbol{\mathsf{A}} ) = n. \end{cases}\]

Proof. The inequality is standard. The equalities follow directly from Sylvester’s rank inequality, which applies as the operator of \(\boldsymbol{\mathsf{A}}\) has a finite-dimensional domain. ◻

As a bounded operator, \(\boldsymbol{\mathsf{A}}\) admits a bounded adjoint \(\boldsymbol{\mathsf{A}}^{\ast} : \mathrm{H}^{m} \to \mathbb{F}^{n}\) given by \[\boldsymbol{\mathsf{A}}^\ast \boldsymbol{\mathsf{y}} = \begin{bmatrix} \langle \boldsymbol{\mathsf{A}}(:,1), \boldsymbol{\mathsf{y}} \rangle_{\mathrm{\ell_2(\mathrm{H})}} & \cdots & \langle \boldsymbol{\mathsf{A}}(:,n), \boldsymbol{\mathsf{y}} \rangle_{\mathrm{\ell_2(\mathrm{H})}} \end{bmatrix}^\intercal.\] We say that \(\boldsymbol{\mathsf{A}}\) has orthonormal columns when its induced operator is a partial isometry, i.e., \(\boldsymbol{\mathsf{A}}^{\ast} \boldsymbol{\mathsf{A}} = \mathsf{I}_n\), where we adopt the convention that the adjoint \(\boldsymbol{\mathsf{A}}^{\ast}\) can be applied to Bochner matrices columnwise. Since the range of the operator induced by \(\boldsymbol{\mathsf{A}}\) is finite-dimensional (and hence closed), there exists a unique bounded pseudoinverse \(\boldsymbol{\mathsf{A}}^{\dagger} : \mathrm{H}^{m} \to \mathbb{F}^{n}\) defined by the standard Moore–Penrose equations [25].

Since \(\boldsymbol{\mathsf{A}}\) corresponds to a finite-rank (and hence compact) operator, it admits a singular value decomposition (SVD) with \(\mathrm{rank} ( \boldsymbol{\mathsf{A}} )\) uniquely defined nonzero singular values [26]. We formulate this classical result in terms of Bochner matrices.

Theorem 1. Let \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{m \times n}\) be nonzero and \(r = \mathrm{rank} ( \boldsymbol{\mathsf{A}} )\). There exists a unique \(\mathsf{\Sigma} = \mathrm{diag}(\sigma_1, \ldots, \sigma_r) \in \mathbb{R}^{r \times r}\) with \(\sigma_1 \geq \cdots \geq \sigma_r > 0\), and there exist \(\boldsymbol{\mathsf{U}} \in \mathrm{H}^{m \times r}\) and \(\mathsf{V} \in \mathbb{F}^{n \times r}\) with orthonormal columns such that \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{U}} \mathsf{\Sigma} \mathsf{V}^\ast\). Furthermore, the pseudoinverse is given by \(\boldsymbol{\mathsf{A}}^{\dagger} = \mathsf{V} \mathsf{\Sigma}^{-1} \boldsymbol{\mathsf{U}}^{\ast}\).

The optimality of the truncated SVD for low-rank approximation of compact operators is well-established in the context of symmetrically normed operator ideals [26]. Finite-rank operators belong to all symmetrically normed operator ideals, and hence the Eckart–Young–Mirsky theorem holds for Bochner matrices. A more direct proof of this result can be obtained via a reduction to classical matrices. In the formulation of the theorem below, we set \(\sigma_j(\boldsymbol{\mathsf{A}}) = 0\) for all \(j > \mathrm{rank} ( \boldsymbol{\mathsf{A}} )\).

Theorem 2. Let \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{m \times n}\). Let \(\phi\) be a symmetric gauge function3 inducing a unitarily invariant norm \(\| \boldsymbol{\mathsf{A}} \|_{\mathrm{\phi}} = \phi(\sigma_1(\boldsymbol{\mathsf{A}}), \ldots, \sigma_n(\boldsymbol{\mathsf{A}}))\). For every \(1 \leq k < \mathrm{rank} ( \boldsymbol{\mathsf{A}} )\), the best approximation of \(\boldsymbol{\mathsf{A}}\) by a Bochner matrix of rank not exceeding \(k\) is achieved by the \(k\)-truncated SVD, denoted \(\lfloor \boldsymbol{\mathsf{A}} \rfloor_{k} = \boldsymbol{\mathsf{U}}_k \mathsf{\Sigma}_k \mathsf{V}_k^\ast\): \[\| \boldsymbol{\mathsf{A}} - \lfloor \boldsymbol{\mathsf{A}} \rfloor_{k} \|_{\mathrm{\phi}} = \min \{ \| \boldsymbol{\mathsf{A}} - \boldsymbol{\mathsf{B}} \|_{\mathrm{\phi}} : \mathrm{rank} ( \boldsymbol{\mathsf{B}} ) \leq k \}.\]

Clearly, the \(\ell_2(\mathrm{H})\) norm and the operator \(\ell_2 \to \ell_2(\mathrm{H})\) norms are unitarily invariant. Let us also note that every such norm \(\| \cdot \|_{\mathrm{\phi}}\) is equivalent to the \(\ell_p(\mathrm{H})\) norms, a property which we state without proof for brevity.

4 Bochner matrices as quasimatrices↩︎

A quasimatrix [27][29] is a generalized matrix with columns that are functions in \(\mathrm{L}^2(0,1)\). Adopting a broader view, the columns can lie in an abstract Hilbert space, which allows us to treat a Bochner matrix \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{m \times n}\) as a quasimatrix with columns in \(\mathrm{H}^{m}\). Consequently, \(\boldsymbol{\mathsf{A}}\) has exactly \(\mathrm{rank} ( \boldsymbol{\mathsf{A}} )\) linearly independent columns and admits an interpolative decomposition.

Let \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{m \times n}\) be nonzero and \(r = \mathrm{rank} ( \boldsymbol{\mathsf{A}} )\). There exists an index set \({J} \subseteq [n]\) of cardinality \(r\) together with \(\mathsf{B} \in \mathbb{F}^{r \times n}\) such that \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{A}}(:, {J}) \mathsf{B}\).

As the columns of \(\boldsymbol{\mathsf{A}}\) belong to a Hilbert space, they can be orthogonalized via the Gram–Schmidt procedure, yielding a thin QR decomposition. The QR decomposition also provides a practical way to compute the SVD of a Bochner matrix, since it reduces the original problem to computing the SVD of a classical triangular matrix; the \(n \times n\) Gram matrix of \(\boldsymbol{\mathsf{A}}\) could also be used to approach the computation of the SVD.

Theorem 3. Let \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{m \times n}\) satisfy \(\mathrm{rank} ( \boldsymbol{\mathsf{A}} ) = n\). Then there exists a unique \(\boldsymbol{\mathsf{Q}} \in \mathrm{H}^{m \times n}\) with orthonormal columns, and a unique upper triangular \(\mathsf{R} \in \mathbb{F}^{n \times n}\) with strictly positive diagonal entries, such that \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{Q}} \mathsf{R}\).

Corollary 1. Let \(\boldsymbol{\mathsf{A}} \neq 0\) and \(r = \mathrm{rank} ( \boldsymbol{\mathsf{A}} ) < n\). There exist \(\boldsymbol{\mathsf{Q}} \in \mathrm{H}^{m \times r}\) with orthonormal columns, an upper triangular \(\mathsf{R}_1 \in \mathbb{F}^{r \times r}\) with strictly positive diagonal entries, \(\mathsf{R}_2 \in \mathbb{F}^{r \times (n-r)}\), and an \(n \times n\) permutation matrix \(\mathsf{\Pi}\) such that \(\boldsymbol{\mathsf{A}}\mathsf{\Pi} = \boldsymbol{\mathsf{Q}} [ \mathsf{R}_1~\mathsf{R}_2 ]\).

Proof. Apply 3 to \(\boldsymbol{\mathsf{A}}(:, {J})\) from [proposition:1s95interpol]. ◻

The columnwise structure of Bochner matrices permits their low-rank approximation via approximate interpolative decompositions as in [proposition:1s95interpol]. Given selected column indices \({J}\), the approximation can be obtained as \(\boldsymbol{\mathsf{A}} \approx \boldsymbol{\mathsf{A}}(:, {J}) \boldsymbol{\mathsf{A}}(:, {J})^{\dagger} \boldsymbol{\mathsf{A}}\).

Let \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{m \times n}\) and \({J} \subseteq [n]\), and denote \(\boldsymbol{\mathsf{C}} = \boldsymbol{\mathsf{A}}(:, {J})\). Then \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{C}}^{\dagger} \boldsymbol{\mathsf{A}}\) if and only if \(\mathrm{rank} ( \boldsymbol{\mathsf{C}} ) = \mathrm{rank} ( \boldsymbol{\mathsf{A}} )\).

Proof. Follows from 1 and the fact that \(\boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{C}}^{\dagger}\) is an orthogonal projection operator onto the column span of \(\boldsymbol{\mathsf{C}}\). ◻

When the goal is to select a small number of columns, the column subset selection problem aims to identify the optimal ones and can be approached via rank-revealing QR decompositions or volume-based sampling [30]—relying only on the inner-product structure of the space of columns.

The rows of a quasimatrix remain undefined unless its columns possess additional structure, e.g., its column-functions having continuous representatives in \(\mathrm{L}^2(0,1)\) [28] or lying in a reproducing kernel Hilbert space [29]. Bochner matrices exhibit a different enabling structure.

5 Bochner matrices as matrices↩︎

Finally, we acknowledge that Bochner matrices are two-dimensional arrays, and hence can be transposed. However, the ranks of \(\boldsymbol{\mathsf{A}}\) and \(\boldsymbol{\mathsf{A}}^{\intercal}\) are, in general, distinct. An example demonstrating this divergence is an \(m \times n\) Bochner matrix whose entries constitute an orthonormal system in \(\mathrm{H}\): it has \(\mathrm{rank} ( \boldsymbol{\mathsf{A}} ) = n\) and \(\mathrm{rank} ( \boldsymbol{\mathsf{A}}^\intercal ) = m\). Therefore, we distinguish between the column rank and row rank of a Bochner matrix: \[\mathrm{rank_{c}} ( \boldsymbol{\mathsf{A}} ) = \mathrm{rank} ( \boldsymbol{\mathsf{A}} ), \quad \mathrm{rank_{r}} ( \boldsymbol{\mathsf{A}} ) = \mathrm{rank} ( \boldsymbol{\mathsf{A}}^\intercal ).\] All of the decompositions discussed so far apply equally well to transposed Bochner matrices, though it is important to keep in mind that the like decompositions of \(\boldsymbol{\mathsf{A}}\) and \(\boldsymbol{\mathsf{A}}^\intercal\) are not connected in any immediate way. Another property of classical matrices that fails for Bochner matrices is the commutativity of transposition with taking the adjoint or pseudoinverse, because the latter are defined strictly as linear operators.

Since we now allow ourselves to access the rows, we can multiply Bochner matrices on the left by classical matrices. While 1 describes the row rank of \(\mathsf{B} \boldsymbol{\mathsf{A}}\) via transposition, its column rank needs to be treated separately.

Lemma 2. Let \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{m \times n}\) and \(\mathsf{B} \in \mathbb{F}^{k \times m}\). If \(\mathrm{rank} ( \mathsf{B} ) = m\) then \(\mathrm{rank} ( \mathsf{B} \boldsymbol{\mathsf{A}} ) = \mathrm{rank} ( \boldsymbol{\mathsf{A}} )\), otherwise \(\mathrm{rank} ( \mathsf{B} \boldsymbol{\mathsf{A}} ) \leq \mathrm{rank} ( \boldsymbol{\mathsf{A}} )\).

Proof. Applying the rank-nullity theorem to the linear operator \(\mathsf{B} \boldsymbol{\mathsf{A}} : \mathbb{F}^n \to \mathrm{H}^k\) and representing its kernel as a direct sum \(\ker(\mathsf{B} \boldsymbol{\mathsf{A}}) = \ker(\boldsymbol{\mathsf{A}}) \oplus (\ker(\boldsymbol{\mathsf{A}})^\perp \cap \ker(\mathsf{B} \boldsymbol{\mathsf{A}}))\), we obtain the following identity: \[\mathrm{rank} ( \mathsf{B} \boldsymbol{\mathsf{A}} ) = \mathrm{rank} ( \boldsymbol{\mathsf{A}} ) - \dim(\ker(\boldsymbol{\mathsf{A}})^\perp \cap \ker(\mathsf{B} \boldsymbol{\mathsf{A}})).\] The rank inequality is its immediate consequence. To prove the equality, it is sufficient to show that \(\ker(\boldsymbol{\mathsf{A}}) = \ker(\mathsf{B} \boldsymbol{\mathsf{A}})\). Let \(\boldsymbol{\mathsf{z}} \in \mathrm{H}^m\) satisfy \(\mathsf{B} \boldsymbol{\mathsf{z}} = 0\). For any \({\boldsymbol{h}} \in \mathrm{H}\), we have \[\begin{bmatrix} 0 \\ \vdots \\ 0 \end{bmatrix} = \begin{bmatrix} \langle \sum_{j = 1}^{m} b_{1,j} {\boldsymbol{z}}_j, {\boldsymbol{h}} \rangle_{\mathrm{\mathrm{H}}} \\ \vdots \\ \langle \sum_{j = 1}^{m} b_{k,j} {\boldsymbol{z}}_j, {\boldsymbol{h}} \rangle_{\mathrm{\mathrm{H}}} \end{bmatrix} = \mathsf{B} \begin{bmatrix} \langle {\boldsymbol{z}}_1, {\boldsymbol{h}} \rangle_{\mathrm{\mathrm{H}}} \\ \vdots \\ \langle {\boldsymbol{z}}_m, {\boldsymbol{h}} \rangle_{\mathrm{\mathrm{H}}} \end{bmatrix}.\] Since \(\mathsf{B}\) has full column rank by assumption, \(\langle {\boldsymbol{z}}_j, {\boldsymbol{h}} \rangle_{\mathrm{\mathrm{H}}} = 0\) for each \(j\). But as \({\boldsymbol{h}}\) is arbitrary, we conclude that \(\boldsymbol{\mathsf{z}} = 0\) and the kernels coincide. ◻

Furthermore, the access to rows provides us with a different way to construct the approximate interpolative decomposition, generalizing [proposition:cca95exact].

Let \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{m \times n}\), \({I} \subseteq [m]\), and \({J} \subseteq [n]\). Denote \(\boldsymbol{\mathsf{C}} = \boldsymbol{\mathsf{A}}(:, {J})\), \(\boldsymbol{\mathsf{R}} = \boldsymbol{\mathsf{A}}({I}, :)\), and \(\boldsymbol{\mathsf{G}} = \boldsymbol{\mathsf{A}}({I}, {J})\). Then \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^{\dagger} \boldsymbol{\mathsf{R}}\) if and only if \(\mathrm{rank} ( \boldsymbol{\mathsf{G}} ) = \mathrm{rank} ( \boldsymbol{\mathsf{A}} )\).

Proof. If \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}}\) then \(\mathrm{rank} ( \boldsymbol{\mathsf{A}} ) \leq \mathrm{rank} ( \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} ) \leq \mathrm{rank} ( \boldsymbol{\mathsf{G}} ) \leq \mathrm{rank} ( \boldsymbol{\mathsf{A}} )\). Conversely, we have \(\mathrm{rank} ( \boldsymbol{\mathsf{A}} ) = \mathrm{rank} ( \boldsymbol{\mathsf{G}} ) \leq \mathrm{rank} ( \boldsymbol{\mathsf{C}} ) \leq \mathrm{rank} ( \boldsymbol{\mathsf{A}} )\). By [proposition:cca95exact], \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{C}}^\dagger\boldsymbol{\mathsf{A}}\) and thus \(\boldsymbol{\mathsf{R}} = \boldsymbol{\mathsf{G}} \boldsymbol{\mathsf{C}}^\dagger\boldsymbol{\mathsf{A}}\). On the other hand, \(\mathrm{rank} ( \boldsymbol{\mathsf{A}} ) = \mathrm{rank} ( \boldsymbol{\mathsf{G}} ) \leq \mathrm{rank} ( \boldsymbol{\mathsf{R}} ) \leq \mathrm{rank} ( \boldsymbol{\mathsf{A}} )\), so that \(\boldsymbol{\mathsf{R}} = \boldsymbol{\mathsf{G}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}}\) by [proposition:cca95exact]. Let \(\boldsymbol{\mathsf{C}} = \boldsymbol{\mathsf{C}}(:, {K}) \mathsf{B}\) be an interpolative decomposition from [proposition:1s95interpol], and note that \(\mathsf{B}\) has full row rank. We then have \(\boldsymbol{\mathsf{G}} = \boldsymbol{\mathsf{G}}(:, {K}) \mathsf{B}\), and \(\boldsymbol{\mathsf{G}}(:, {K})\) has full column rank by 1. To finish the proof, \[\boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} = \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{G}} \boldsymbol{\mathsf{C}}^\dagger\boldsymbol{\mathsf{A}} = \mathsf{B}^\dagger\boldsymbol{\mathsf{G}}(:, {K})^\dagger\boldsymbol{\mathsf{G}}(:, {K}) \mathsf{B} \mathsf{B}^\dagger\boldsymbol{\mathsf{C}}(:, {K})^\dagger\boldsymbol{\mathsf{A}} = \mathsf{B}^\dagger\boldsymbol{\mathsf{C}}(:, {K})^\dagger\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{C}}^\dagger\boldsymbol{\mathsf{A}}.\] ◻

The interpolative decomposition in [proposition:cgr95exact] is written in exactly the same form as the skeleton decomposition of classical matrices [21], i.e., \(\mathsf{A}(:, {J}) \mathsf{A}({I}, {J})^{\dagger} \mathsf{A}({I}, :)\). Despite this syntactic similarity, there are two crucial differences. First, the classical skeleton decomposition is a product of three factors, whereas the decomposition in [proposition:cgr95exact] fundamentally collapses into two factors because the pseudoinverse operator cannot stand alone as a purely algebraic middle factor. Second, the skeleton decomposition of classical matrices is a CUR decomposition [31]; that is, every column is expressed as a linear combination of selected columns and every row is expressed as a linear combination of selected rows. In contrast, the decomposition in [proposition:cgr95exact] is interpolative only with respect to the columns and fails to leverage the row rank.

In fact, CUR decompositions in the sense described above are structurally impossible for Bochner matrices: they would contain two Bochner-matrix factors, yet those cannot be multiplied within our framework. In this regard, consider a verbatim extension of the projection-based classical CUR decomposition \(\mathsf{C} \mathsf{C}^{\dagger} \mathsf{A} \mathsf{R}^{\dagger} \mathsf{R}\) [32], which we interpret not as a three-factor product with \(\mathsf{C}^{\dagger} \mathsf{A} \mathsf{R}^{\dagger}\) in the middle, but as the simultaneous application of column and row orthogonal projectors. For Bochner matrices, we can write this down as \((\boldsymbol{\mathsf{R}}^\intercal(\boldsymbol{\mathsf{R}}^\intercal)^\dagger(\boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{C}}^\dagger\boldsymbol{\mathsf{A}})^\intercal)^\intercal\) or \(\boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{C}}^\dagger(\boldsymbol{\mathsf{R}}^\intercal(\boldsymbol{\mathsf{R}}^\intercal)^\dagger\boldsymbol{\mathsf{A}}^\intercal)^\intercal\). Consider a simple \(2 \times 2\) example where the entries of \(\boldsymbol{\mathsf{A}}\) form an orthonormal system in \(\mathrm{H}\), and let \({I} = {J} = \{ 1 \}\). Then the two projection-based CUR extensions evaluate to \[(\boldsymbol{\mathsf{R}}^\intercal(\boldsymbol{\mathsf{R}}^\intercal)^\dagger(\boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{C}}^\dagger\boldsymbol{\mathsf{A}})^\intercal)^\intercal= \frac{1}{2} \begin{bmatrix} \boldsymbol{\mathsf{A}}(1,1) & \boldsymbol{\mathsf{A}}(1,2) \\ 0 & 0 \end{bmatrix}, \quad \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{C}}^\dagger(\boldsymbol{\mathsf{R}}^\intercal(\boldsymbol{\mathsf{R}}^\intercal)^\dagger\boldsymbol{\mathsf{A}}^\intercal)^\intercal= \frac{1}{2} \begin{bmatrix} \boldsymbol{\mathsf{A}}(1,1) & 0 \\ \boldsymbol{\mathsf{A}}(2,1) & 0 \end{bmatrix}.\] They do not agree, the column and row projections overriding each other, and neither is strictly a CUR decomposition. Note also that the subspaces we selected to project on are dominant left singular subspaces of \(\boldsymbol{\mathsf{A}}\) and \(\boldsymbol{\mathsf{A}}^\intercal\).

6 Bochner matrices as third-order tensors↩︎

As demonstrated, a CUR-style application of the one-sided interpolative decompositions from [proposition:cca95exact] [proposition:cgr95exact] does not yield a simultaneous reduction of both the column and row ranks. Instead, we will extend classical two-sided interpolative decompositions [33] into the Bochner setting.

Let \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{m \times n}\) be nonzero with \(\rho = \mathrm{rank_{r}} ( \boldsymbol{\mathsf{A}} )\), \(\kappa = \mathrm{rank_{c}} ( \boldsymbol{\mathsf{A}} )\). There exist index sets \({I} \subseteq [m]\) of cardinality \(\rho\) and \({J} \subseteq [n]\) of cardinality \(\kappa\), together with \(\mathsf{B} \in \mathbb{F}^{m \times \rho}\) and \(\mathsf{C} \in \mathbb{F}^{\kappa \times n}\), such that \(\boldsymbol{\mathsf{A}} = \mathsf{B} \boldsymbol{\mathsf{A}}({I}, {J}) \mathsf{C}\).

Proof. Applying [proposition:1s95interpol] to \(\boldsymbol{\mathsf{A}}\) and \(\boldsymbol{\mathsf{A}}^\intercal\), we get \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{A}}(:, {J}) \mathsf{C} = \mathsf{B} \boldsymbol{\mathsf{A}}({I}, :)\). ◻

The structure of the decomposition in [proposition:2s95interpol], having a Bochner matrix as its inner core and classical matrices as its outer factors, strictly mirrors the Tucker decomposition of classical tensors [34]. Indeed, in the finite-dimensional special case of \(\mathrm{H}= \mathbb{F}^{k}\) this reduces to the classical Tucker2 decomposition; and in the general case, the space of Bochner matrices \(\mathrm{H}^{m \times n}\) is isometrically isomorphic to the tensor-product Hilbert space \(\mathrm{H}\mathbin{\otimes}\mathbb{F}^{m} \mathbin{\otimes}\mathbb{F}^{n}\) [35]. Motivated by these correspondences, we call every decomposition of the form \(\boldsymbol{\mathsf{A}} = \mathsf{B} \boldsymbol{\mathsf{G}} \mathsf{D}\) a Tucker decomposition and refer to it as a cross decomposition when its core \(\boldsymbol{\mathsf{G}}\) is a submatrix of \(\boldsymbol{\mathsf{A}}\). Accordingly, we define the Tucker rank of a Bochner matrix as a tuple \(\mathrm{rank_{T}} ( \boldsymbol{\mathsf{A}} ) = (\mathrm{rank_{r}} ( \boldsymbol{\mathsf{A}} ), \mathrm{rank_{c}} ( \boldsymbol{\mathsf{A}} ))\).

The following theorem extends [proposition:cgr95exact] to cross decompositions and connects them to the CUR-style application of one-sided interpolative decompositions. This generalizes properties of classical matrices and tensors derived in [31], [36], [37] to Bochner matrices.

Theorem 4. Let \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{m \times n}\), \({I} \subseteq [m]\), and \({J} \subseteq [n]\). Denote \(\boldsymbol{\mathsf{C}} = \boldsymbol{\mathsf{A}}(:, {J})\), \(\boldsymbol{\mathsf{R}} = \boldsymbol{\mathsf{A}}({I},:)\), and \(\boldsymbol{\mathsf{G}} = \boldsymbol{\mathsf{A}}({I}, {J})\). Then the following are equivalent:

  1. \(\mathrm{rank_{c}} ( \boldsymbol{\mathsf{C}} ) = \mathrm{rank_{c}} ( \boldsymbol{\mathsf{A}} )\) and \(\mathrm{rank_{r}} ( \boldsymbol{\mathsf{R}} ) = \mathrm{rank_{r}} ( \boldsymbol{\mathsf{A}} )\);

  2. \(\mathrm{rank_{T}} ( \boldsymbol{\mathsf{G}} ) = \mathrm{rank_{T}} ( \boldsymbol{\mathsf{A}} )\);

  3. \(\boldsymbol{\mathsf{A}} = (\boldsymbol{\mathsf{R}}^\intercal(\boldsymbol{\mathsf{R}}^\intercal)^\dagger(\boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{C}}^\dagger\boldsymbol{\mathsf{A}})^\intercal)^\intercal= \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{C}}^\dagger(\boldsymbol{\mathsf{R}}^\intercal(\boldsymbol{\mathsf{R}}^\intercal)^\dagger\boldsymbol{\mathsf{A}}^\intercal)^\intercal\);

  4. \(\boldsymbol{\mathsf{A}} = ((\boldsymbol{\mathsf{G}}^\intercal)^\dagger\boldsymbol{\mathsf{C}}^\intercal)^\intercal\boldsymbol{\mathsf{G}} (\boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}})\).

Proof. \(\mathbf{(1 \Rightarrow 2).}\) [proposition:cca95exact] gives \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{C}}^\dagger\boldsymbol{\mathsf{A}}\), and thus \(\boldsymbol{\mathsf{R}} = \boldsymbol{\mathsf{G}} \boldsymbol{\mathsf{C}}^\dagger\boldsymbol{\mathsf{A}}\). It follows that \(\mathrm{rank_{c}} ( \boldsymbol{\mathsf{G}} ) = \mathrm{rank_{c}} ( \boldsymbol{\mathsf{R}} )\) since \(\boldsymbol{\mathsf{G}}\) is a submatrix of \(\boldsymbol{\mathsf{R}}\). Similarly, \(\mathrm{rank_{r}} ( \boldsymbol{\mathsf{G}} ) = \mathrm{rank_{r}} ( \boldsymbol{\mathsf{C}} )\), and there exists \({K} \subseteq [m]\) together with \(\mathsf{X}\) of full column rank such that \(\boldsymbol{\mathsf{C}} = \mathsf{X} \boldsymbol{\mathsf{G}}({K}, :)\). By 2, we get \(\mathrm{rank_{c}} ( \boldsymbol{\mathsf{C}} ) = \mathrm{rank_{c}} ( \boldsymbol{\mathsf{G}}({K}, :) ) \leq \mathrm{rank_{c}} ( \boldsymbol{\mathsf{G}} ) \leq \mathrm{rank_{c}} ( \boldsymbol{\mathsf{A}} ) = \mathrm{rank_{c}} ( \boldsymbol{\mathsf{C}} )\). By analogy, \(\mathrm{rank_{r}} ( \boldsymbol{\mathsf{G}} ) = \mathrm{rank_{r}} ( \boldsymbol{\mathsf{A}} )\).

\(\mathbf{(2 \Rightarrow 1).}\) We have \(\mathrm{rank_{c}} ( \boldsymbol{\mathsf{A}} ) = \mathrm{rank_{c}} ( \boldsymbol{\mathsf{G}} ) \leq \mathrm{rank_{c}} ( \boldsymbol{\mathsf{C}} ) \leq \mathrm{rank_{c}} ( \boldsymbol{\mathsf{A}} )\) as \(\boldsymbol{\mathsf{G}}\) is a submatrix of \(\boldsymbol{\mathsf{C}}\). Similarly, \(\mathrm{rank_{r}} ( \boldsymbol{\mathsf{R}} ) \leq \mathrm{rank_{r}} ( \boldsymbol{\mathsf{A}} )\).

\(\mathbf{(1 \Rightarrow 3).}\) By [proposition:cca95exact], \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{C}}^\dagger\boldsymbol{\mathsf{A}} = (\boldsymbol{\mathsf{R}}^\intercal(\boldsymbol{\mathsf{R}}^\intercal)^\dagger\boldsymbol{\mathsf{A}}^\intercal)^\intercal= \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{C}}^\dagger(\boldsymbol{\mathsf{R}}^\intercal(\boldsymbol{\mathsf{R}}^\intercal)^\dagger\boldsymbol{\mathsf{A}}^\intercal)^\intercal\). Repeat the argument for \(\boldsymbol{\mathsf{A}}^\intercal\) to get the second representation.

\(\mathbf{(3 \Rightarrow 1).}\) By 1, \(\mathrm{rank_{c}} ( \boldsymbol{\mathsf{A}} ) \leq \mathrm{rank_{c}} ( \boldsymbol{\mathsf{C}} ) \leq \mathrm{rank_{c}} ( \boldsymbol{\mathsf{A}} )\), and similarly for \(\boldsymbol{\mathsf{R}}^\intercal\).

\(\mathbf{(2 \Rightarrow 4).}\) [proposition:cgr95exact] yields \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} = (\boldsymbol{\mathsf{R}}^\intercal(\boldsymbol{\mathsf{G}}^\intercal)^\dagger\boldsymbol{\mathsf{C}}^\intercal)^\intercal\), whence we obtain \(\boldsymbol{\mathsf{C}}^\intercal= \boldsymbol{\mathsf{G}}^\intercal(\boldsymbol{\mathsf{G}}^\intercal)^\dagger\boldsymbol{\mathsf{C}}^\intercal\).

\(\mathbf{(4 \Rightarrow 2).}\) Note that \(\boldsymbol{\mathsf{C}} = ((\boldsymbol{\mathsf{G}}^\intercal)^\dagger\boldsymbol{\mathsf{C}}^\intercal)^\intercal\boldsymbol{\mathsf{G}} (\boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{G}}) = ((\boldsymbol{\mathsf{G}}^\intercal)^\dagger\boldsymbol{\mathsf{C}}^\intercal)^\intercal\boldsymbol{\mathsf{G}}\) by the definition of pseudoinverse. Then \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}}\), and [proposition:cgr95exact] gives \(\mathrm{rank_{c}} ( \boldsymbol{\mathsf{G}} ) = \mathrm{rank_{c}} ( \boldsymbol{\mathsf{A}} )\). The same argument gives \(\mathrm{rank_{r}} ( \boldsymbol{\mathsf{G}} ) = \mathrm{rank_{r}} ( \boldsymbol{\mathsf{A}} )\) when applied to \(\boldsymbol{\mathsf{A}}^\intercal\). ◻

The interpretation of Bochner matrices as third-order tensors paves the way for their low-rank approximation by means of higher-order SVD (HOSVD) [38]. However, as the example at the end of 5 shows, unlike the case of classical tensors, low-rank truncation for each unfolding (i.e., \(\boldsymbol{\mathsf{A}}\) and \(\boldsymbol{\mathsf{A}}^\intercal\)) cannot be performed via column-space projection: the two projections do not commute and interfere with one another. Instead, we must apply classical orthogonal projection matrices on the right, thereby projecting in the row space. Such construction is perfectly aligned with the definition of mode unfoldings in abstract tensor-product Hilbert spaces [35].

Theorem 5. Let \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{m \times n}\) be nonzero, and let \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{U}} \mathsf{\Sigma} \mathsf{V}^\ast\) and \(\boldsymbol{\mathsf{A}}^\intercal= \boldsymbol{\mathsf{P}} \mathsf{\Lambda} \mathsf{Q}^\ast\) be its SVDs. For every pair of \(1 \leq \rho \leq \mathrm{rank_{r}} ( \boldsymbol{\mathsf{A}} )\) and \(1 \leq \kappa \leq \mathrm{rank_{c}} ( \boldsymbol{\mathsf{A}} )\), denote by \(\lfloor \boldsymbol{\mathsf{A}} \rfloor_{\kappa} = \boldsymbol{\mathsf{U}}_{\kappa} \mathsf{\Sigma}_{\kappa} \mathsf{V}_{\kappa}^\ast\) and \(\lfloor \boldsymbol{\mathsf{A}}^\intercal \rfloor_{\rho} = \boldsymbol{\mathsf{P}}_{\rho} \mathsf{\Lambda}_{\rho} \mathsf{Q}_{\rho}^\ast\) the corresponding truncated SVDs. Then \[\begin{align} \| \boldsymbol{\mathsf{A}} - \overline{\mathsf{Q}}_\rho \mathsf{Q}_\rho^\intercal\boldsymbol{\mathsf{A}} \mathsf{V}_\kappa \mathsf{V}_\kappa^\ast \|_{\mathrm{\ell_2(\mathrm{H})}} &\leq \sqrt{\| \mathsf{\Sigma} - \mathsf{\Sigma}_\kappa \|_{\mathrm{\ell_2}}^2 + \| \mathsf{\Lambda} - \mathsf{\Lambda}_\rho \|_{\mathrm{\ell_2}}^2} \\ &\leq \sqrt{2} \min \Big\{ \| \boldsymbol{\mathsf{A}} - \boldsymbol{\mathsf{B}} \|_{\mathrm{\ell_2(\mathrm{H})}} : \mathrm{rank_{r}} ( \boldsymbol{\mathsf{B}} ) \leq \rho,~\mathrm{rank_{c}} ( \boldsymbol{\mathsf{B}} ) \leq \kappa \Big\}. \end{align}\]

Proof. The Pythagorean theorem gives \[\begin{align} \| \boldsymbol{\mathsf{A}} - \overline{\mathsf{Q}}_\rho \mathsf{Q}_\rho^\intercal\boldsymbol{\mathsf{A}} \mathsf{V}_\kappa \mathsf{V}_\kappa^\ast \|_{\mathrm{\ell_2(\mathrm{H})}}^2 &= \| \boldsymbol{\mathsf{A}} (\mathsf{I} - \mathsf{V}_\kappa \mathsf{V}_\kappa^\ast) \|_{\mathrm{\ell_2(\mathrm{H})}}^2 + \| (\mathsf{I} - \overline{\mathsf{Q}}_\rho \mathsf{Q}_\rho^\intercal) \boldsymbol{\mathsf{A}} \mathsf{V}_\kappa \mathsf{V}_\kappa^\ast \|_{\mathrm{\ell_2(\mathrm{H})}}^2 \\ &\leq \| \boldsymbol{\mathsf{A}} (\mathsf{I} - \mathsf{V}_\kappa \mathsf{V}_\kappa^\ast) \|_{\mathrm{\ell_2(\mathrm{H})}}^2 + \| \boldsymbol{\mathsf{A}}^\intercal(\mathsf{I} - \mathsf{Q}_\rho \mathsf{Q}_\rho^\ast) \|_{\mathrm{\ell_2(\mathrm{H})}}^2 \\ &= \| \mathsf{\Sigma} - \mathsf{\Sigma}_\kappa \|_{\mathrm{\ell_2}}^2 + \| \mathsf{\Lambda} - \mathsf{\Lambda}_\rho \|_{\mathrm{\ell_2}}^2. \end{align}\] A standard corollary of the Eckart–Young–Mirsky theorem (2) is that the set of Bochner matrices with bounded Tucker rank is closed, and therefore the minimal distance is achieved at some \(\boldsymbol{\mathsf{B}}_{\circ}\). Invoking 2 again, we obtain \[\| \mathsf{\Sigma} - \mathsf{\Sigma}_\kappa \|_{\mathrm{\ell_2}}^2 \leq \| \boldsymbol{\mathsf{A}} - \boldsymbol{\mathsf{B}}_{\circ} \|_{\mathrm{\ell_2(\mathrm{H})}}^2, \quad \| \mathsf{\Lambda} - \mathsf{\Lambda}_\rho \|_{\mathrm{\ell_2}}^2 \leq \| \boldsymbol{\mathsf{A}}^\intercal- \boldsymbol{\mathsf{B}}_{\circ}^\intercal \|_{\mathrm{\ell_2(\mathrm{H})}}^2.\] ◻

The proof is a direct adaptation of the HOSVD argument for classical tensors [38], and it extends to tensors in abstract tensor-product Hilbert spaces because their mode unfoldings are Hilbert–Schmidt operators (even though [35] is restricted to classical tensors). A similar quasioptimality bound can be shown for the sequentially truncated HOSVD of Bochner matrices (see [39] and [35]).

7 Cross approximation of Bochner matrices↩︎

In the preceding sections, we established that Bochner matrices can be approximated with low Tucker rank via the truncated HOSVD (5) and derived necessary and sufficient conditions for their cross decompositions to be exact (4). However, our motivating application introduced in 1.1 relies on approximate cross decompositions, and therefore this section investigates the theoretical properties of cross approximation: \[\mathrm{cross}(\boldsymbol{\mathsf{B}}, {I}, {J} ) = ((\boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\intercal)^\dagger\boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}}^\intercal)^\intercal\boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}} (\boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}}),\] where \(\boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} = \boldsymbol{\mathsf{B}}(:, {J})\), \(\boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}} = \boldsymbol{\mathsf{B}}({I}, {J})\), \(\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} = \boldsymbol{\mathsf{B}}({I}, :)\) are the selected submatrices. We will drop the subscripts for the submatrices of \(\boldsymbol{\mathsf{A}}\) and use regular font for the submatrices of classical matrices.

7.1 Interpolation property↩︎

When the conditions of 4 are not met and cross approximation fails to recover the entire Bochner matrix, it still interpolates entries residing in the selected rows and columns. The following theorem formally extends the interpolation property of classical cross approximation to Bochner matrices. For instance, it shows that every cross approximant is a cross decomposition of itself.

Theorem 6. Let \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{m \times n}\) and \(\boldsymbol{\mathsf{X}} = \mathrm{cross}(\boldsymbol{\mathsf{A}}, {I}, {J} )\). Then

  1. \(\boldsymbol{\mathsf{X}}({I}, {J}) = \boldsymbol{\mathsf{G}}\) and \(\mathrm{rank_{T}} ( \boldsymbol{\mathsf{X}} ) = \mathrm{rank_{T}} ( \boldsymbol{\mathsf{G}} )\);

  2. \(\boldsymbol{\mathsf{X}}({I},:) = \boldsymbol{\mathsf{R}}\) if and only if \(\mathrm{rank_{c}} ( \boldsymbol{\mathsf{G}} ) = \mathrm{rank_{c}} ( \boldsymbol{\mathsf{R}} )\);

  3. \(\boldsymbol{\mathsf{X}}(:,{J}) = \boldsymbol{\mathsf{C}}\) if and only if \(\mathrm{rank_{r}} ( \boldsymbol{\mathsf{G}} ) = \mathrm{rank_{r}} ( \boldsymbol{\mathsf{C}} )\);

  4. \(\boldsymbol{\mathsf{X}} = \boldsymbol{\mathsf{A}}\) if and only if \(\mathrm{rank_{T}} ( \boldsymbol{\mathsf{G}} ) = \mathrm{rank_{T}} ( \boldsymbol{\mathsf{A}} )\).

Proof. \(\mathbf{1.}\) By the defining Moore–Penrose equations of pseudoinverse: \[\boldsymbol{\mathsf{X}}({I}, {J}) = ((\boldsymbol{\mathsf{G}}^\intercal)^\dagger\boldsymbol{\mathsf{G}}^\intercal)^\intercal\boldsymbol{\mathsf{G}} (\boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{G}}) = ((\boldsymbol{\mathsf{G}}^\intercal)^\dagger\boldsymbol{\mathsf{G}}^\intercal)^\intercal\boldsymbol{\mathsf{G}} = (\boldsymbol{\mathsf{G}}^\intercal(\boldsymbol{\mathsf{G}}^\intercal)^\dagger\boldsymbol{\mathsf{G}}^\intercal)^\intercal= \boldsymbol{\mathsf{G}}.\] By 1 2, \(\mathrm{rank_{T}} ( \boldsymbol{\mathsf{X}} )\) is upper bounded by \(\mathrm{rank_{T}} ( \boldsymbol{\mathsf{G}} )\). As \(\boldsymbol{\mathsf{G}}\) is a submatrix of \(\boldsymbol{\mathsf{X}}\), the reverse inequality holds as well.

\(\mathbf{2.}\) By definition, \(\boldsymbol{\mathsf{X}}({I},:) = \boldsymbol{\mathsf{G}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}}\). [proposition:cca95exact] guarantees that \(\boldsymbol{\mathsf{X}}({I},:) = \boldsymbol{\mathsf{R}}\) if and only if \(\mathrm{rank_{c}} ( \boldsymbol{\mathsf{G}} ) = \mathrm{rank_{c}} ( \boldsymbol{\mathsf{R}} )\).

\(\mathbf{3.}\) Apply the same argument to \(\boldsymbol{\mathsf{X}}^\intercal\).

\(\mathbf{4.}\) Follows from 4. ◻

7.2 Error bound for one-sided interpolation↩︎

As an intermediate step towards the error bound for Bochner cross approximation, we derive an error bound for the approximate one-sided interpolation ([proposition:cgr95exact]). The proof mirrors the argument used for cross approximation of classical perturbed low-rank matrices [40]. We begin with a technical lemma.

Lemma 3. Let \(\boldsymbol{\mathsf{B}} \in \mathrm{H}^{m \times n}\) be nonzero and satisfy \(\mathrm{rank_{c}} ( \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}} ) = \mathrm{rank_{c}} ( \boldsymbol{\mathsf{B}} )\), and let \(\boldsymbol{\mathsf{B}} = \boldsymbol{\mathsf{U}} \mathsf{\Sigma} \mathsf{V}^\ast\) be its SVD. Then \[\| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2(\mathrm{H})}} = \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}}, \quad \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} \|_{\mathrm{\ell_2 \to \ell_2}} = \| \mathsf{C}_{\mathsf{V}^\ast}^\dagger \|_{\mathrm{\ell_2 \to \ell_2}}.\]

Proof. By definition, \(\boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}} = \boldsymbol{\mathsf{U}}({I}, :) \mathsf{\Sigma} \mathsf{V}({J}, :)^\ast = \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}} \mathsf{\Sigma} \mathsf{C}_{\mathsf{V}^\ast}\); likewise, \(\boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} = \boldsymbol{\mathsf{U}} \mathsf{\Sigma} \mathsf{C}_{\mathsf{V}^\ast}\) and \(\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} = \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}} \mathsf{\Sigma} \mathsf{V}^\ast\). The rank assumption guarantees, via 1, that \(\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}\) has full column rank and \(\mathsf{C}_{\mathsf{V}^\ast}\) has full row rank. Therefore, \(\boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger= \mathsf{C}_{\mathsf{V}^\ast}^\dagger\mathsf{\Sigma}^{-1} \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger\) and thus \[\boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger= \boldsymbol{\mathsf{U}} \mathsf{\Sigma} \mathsf{C}_{\mathsf{V}^\ast} \mathsf{C}_{\mathsf{V}^\ast}^\dagger\mathsf{\Sigma}^{-1} \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger= \boldsymbol{\mathsf{U}} \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger, \quad \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} = \mathsf{C}_{\mathsf{V}^\ast}^\dagger\mathsf{\Sigma}^{-1} \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}} \mathsf{\Sigma} \mathsf{V}^\ast = \mathsf{C}_{\mathsf{V}^\ast}^\dagger\mathsf{V}^\ast.\] It remains to note that \(\boldsymbol{\mathsf{U}}\) and \(\mathsf{V}^\ast\), being partial isometries, do not alter the respective operator norms. ◻

Theorem 7. Let \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{B}} + \boldsymbol{\mathsf{E}} \in \mathrm{H}^{m \times n}\) with \(\boldsymbol{\mathsf{B}} \neq 0\). Let \(\boldsymbol{\mathsf{B}} = \boldsymbol{\mathsf{U}} \mathsf{\Sigma} \mathsf{V}^\ast\) be an SVD, and let \(\mathrm{rank_{c}} ( \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}} ) = \mathrm{rank_{c}} ( \boldsymbol{\mathsf{B}} )\). Let \(\| \cdot \|_{\mathrm{\phi}}\) be a unitarily invariant norm induced by a symmetric gauge function \(\phi\). Denote \(\eta_r = \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}}\) and \(\eta_c = \| \mathsf{C}_{\mathsf{V}^\ast}^\dagger \|_{\mathrm{\ell_2 \to \ell_2}}\). Then \[\begin{align} \| \boldsymbol{\mathsf{B}} &- \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\phi}} \leq \eta_r \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \eta_c \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + 3 \eta_r \eta_c \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} \\ &+ \| \boldsymbol{\mathsf{G}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}} \Big( \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \big(\eta_r \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \eta_c \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \eta_r \eta_c \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}}\big) \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} \Big). \end{align}\]

Proof. The rank assumption together with [proposition:cgr95exact] yields \[\| \boldsymbol{\mathsf{B}} - \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\phi}} \leq \| \boldsymbol{\mathsf{B}} - \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} \|_{\mathrm{\phi}} + \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} - \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\phi}} = \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} - \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\phi}},\] and therefore \[\begin{align} \| \boldsymbol{\mathsf{B}} - \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\phi}} &\leq \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} - \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} \|_{\mathrm{\phi}} \\ &+ \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} - \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\phi}} + \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} - \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\phi}}. \end{align}\] The standard properties of unitarily invariant norms of operators [26] provide bounds for the last two terms: \[\begin{align} \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} - \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\phi}} &\leq \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2(\mathrm{H})}} \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}}, \\ \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} - \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\phi}} &\leq \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} \| \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_2 \to \ell_2}}. \end{align}\] To bound the first term, we refer to [proposition:cgr95exact] again and rely on the equalities \(\boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} = \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}\) and \(\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} = \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}}\). When combined with the defining equations of pseudoinverses, this leads to the equality \[\label{eq:cgr95error95part} \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} - \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} \|_{\mathrm{\phi}} = \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} - \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger(\boldsymbol{\mathsf{G}} - \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}}) \boldsymbol{\mathsf{G}}^\dagger(\boldsymbol{\mathsf{G}} - \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}}) \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} \|_{\mathrm{\phi}},\tag{1}\] which we bound in turn using the triangle inequality: \[\begin{align} \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} - \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} \|_{\mathrm{\phi}} &\leq \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2(\mathrm{H})}} \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}} \|_{\mathrm{\ell_2 \to \ell_2}} \\ &\cdot \big( \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \| \boldsymbol{\mathsf{G}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{G}} \|_{\mathrm{\phi}} + \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} \big). \end{align}\] Since \(\boldsymbol{\mathsf{G}} \boldsymbol{\mathsf{G}}^\dagger\) and \(\boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{G}}\) are orthogonal projection operators in \(\mathrm{H}^{|{I}|}\) and \(\mathbb{F}^{|{J}|}\), respectively, the expression in parentheses is bounded by4 \((3 + \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} \| \boldsymbol{\mathsf{G}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}}) \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}}\). To finish the proof, we need to bound the operator norms of \(\boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\), \(\boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{B}}}\), \(\boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}^\dagger\), \(\boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}}\). The first two bounds follow from 3. For the third bound, consider \[\begin{align} \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2(\mathrm{H})}} &\leq \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2(\mathrm{H})}} \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2(\mathrm{H})}} \\ &= \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}} \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}} \boldsymbol{\mathsf{G}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2(\mathrm{H})}} \\ &\leq \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}} \big(1 + \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \boldsymbol{\mathsf{G}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2(\mathrm{H})}} \big) \\ &\leq \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}} \big(1 + \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} \| \boldsymbol{\mathsf{G}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}} \big), \end{align}\] where we use \(\| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\ell_2 \to \ell_2(\mathrm{H})}} \leq \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}}\) to simplify the overall bound. Likewise, \[\| \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_2 \to \ell_2}} \leq \| \mathsf{C}_{\mathsf{V}^\ast}^\dagger \|_{\mathrm{\ell_2 \to \ell_2}} \big(1 + \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} \| \boldsymbol{\mathsf{G}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}} \big) + \| \boldsymbol{\mathsf{G}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}} \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}}.\] Combining all the intermediate bounds results in the final bound. ◻

7 quantifies how well the approximate one-sided interpolation recovers an underlying Bochner matrix \(\boldsymbol{\mathsf{B}}\) with low column rank; an error bound for \(\boldsymbol{\mathsf{A}}\) itself follows by triangle inequality. Notably, the low-rank assumption is not imposed explicitly; rather, we replace it with a rank-equality condition, essentially demanding that the selected index sets be sufficiently rich for the particular \(\boldsymbol{\mathsf{B}}\).

The decomposition produced by 7 possesses an inner dimension \(|{J}|\) and satisfies a rank bound \(\mathrm{rank_{c}} ( \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} ) \leq \mathrm{rank_{c}} ( \boldsymbol{\mathsf{G}} ) \leq |{J}|\). For a fixed selection of indices, we can further reduce the column rank by truncating \(\boldsymbol{\mathsf{G}}\) to lower column rank prior to taking its pseudoinverse. Recall that \(\lfloor \boldsymbol{\mathsf{G}} \rfloor_{k}\) denotes the SVD of \(\boldsymbol{\mathsf{G}}\) truncated to \(k\) nonzero singular values. Similarly, let us denote by \(\lfloor \boldsymbol{\mathsf{G}} \rfloor_{\tau}\) the result of discarding all singular values below \(\tau > 0\); we distinguish between the two based on whether the subscript is an integer or a real number. With \(\lfloor \boldsymbol{\mathsf{G}} \rfloor_{\tau}\) used in place of \(\boldsymbol{\mathsf{G}}\), the second-order error term can be controlled better, overcoming the stability issues due to the possible ill-conditioning of \(\boldsymbol{\mathsf{G}}\)—which is common when the indices are generously oversampled. However, the interpolation properties break down once the pseudoinverse is regularized.

Corollary 2. Let \(\tau > 0\). In the setting of 7, \[\begin{align} \| \boldsymbol{\mathsf{B}} - \boldsymbol{\mathsf{C}} \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\tau}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\phi}} &\leq \eta_r \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \eta_c \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \eta_r \eta_c (2 \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}} - \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\tau} \|_{\mathrm{\phi}}) \\ &+ \frac{1}{\tau} \Big( \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \big(\eta_r \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \eta_c \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \eta_r \eta_c \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}}\big) \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} \Big). \end{align}\]

Proof. The proof of 7 extends almost verbatim. The only step that requires special attention is 1 . Expanding its right-hand side, we note that \(\boldsymbol{\mathsf{G}} \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\tau}^\dagger\boldsymbol{\mathsf{G}} = \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\tau}^\dagger\) and, subsequently, rely on \(\boldsymbol{\mathsf{G}} \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\tau}^\dagger\) and \(\lfloor \boldsymbol{\mathsf{G}} \rfloor_{\tau}^\dagger\boldsymbol{\mathsf{G}}\) also being orthogonal projection operators due to the SVD truncation. ◻

A similar corollary holds if we truncate the rank rather than apply thresholding to the singular values. Below, we bound the norm of the pseudoinverse of the truncated \(\boldsymbol{\mathsf{G}}\) in terms of \(\boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}\) and \(\boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}}\), separating their contributions.

Corollary 3. In the setting of 7, let \(1 \leq k \leq \mathrm{rank_{c}} ( \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}} )\) be such that \(\gamma = \sigma_{k}(\boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}}) - \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\ell_2 \to \ell_2(\mathrm{H})}} > 0\). Then \[\begin{align} \| \boldsymbol{\mathsf{B}} - \boldsymbol{\mathsf{C}} \lfloor \boldsymbol{\mathsf{G}} \rfloor_{k}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\phi}} &\leq \eta_r \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \eta_c \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \eta_r \eta_c (2 \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}} - \lfloor \boldsymbol{\mathsf{G}} \rfloor_{k} \|_{\mathrm{\phi}}) \\ &+ \frac{1}{\gamma} \Big( \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \big(\eta_r \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \eta_c \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} + \eta_r \eta_c \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}}\big) \| \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\phi}} \Big). \end{align}\]

Proof. By Weyl’s inequalities, which apply to Bochner matrices, \(\sigma_k(\boldsymbol{\mathsf{G}}) \geq \gamma > 0\), and thus the rank truncation \(\lfloor \boldsymbol{\mathsf{G}} \rfloor_{k}\) is well-defined. The proof remains the same, except for the bound \(\| \lfloor \boldsymbol{\mathsf{G}} \rfloor_{k}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}} = \sigma_k(\boldsymbol{\mathsf{G}})^{-1} \leq \gamma^{-1}\). ◻

7.3 Well-conditioned submatrices↩︎

The approximation error bound in 7 is governed by the operator norms of the pseudoinverses of submatrices of singular factors. The existence of well-conditioned submatrices of classical \(n \times \kappa\) matrices with orthonormal columns is a long-standing question. It was conjectured in [21] that there always exists a \(\kappa \times \kappa\) submatrix for which the operator norm of its inverse is bounded by \(\sqrt{n}\); see [41], [42] for recent progress in this direction.

Maximum-volume submatrices do not resolve this conjecture, but are nonetheless well-conditioned [21], [43][45]. The volume is defined as \(\nu(\mathsf{A}) = \sqrt{\det(\mathsf{A}^\ast \mathsf{A})}\), and if a \(|{J}| \times \kappa\) (\(|{J}| \geq \kappa\)) submatrix has the largest volume among all \(|{J}| \times \kappa\) submatrices then the operator norm of its pseudoinverse is upper bounded by \[\label{eq:maxvol95classical} \sqrt{1 + \frac{(n - |{J}|)\kappa}{|{J}| - \kappa + 1}}.\tag{2}\] We generalize this to Bochner matrices and define their volume as \(\nu(\boldsymbol{\mathsf{A}}) = \sqrt{\det(\boldsymbol{\mathsf{A}}^\ast \boldsymbol{\mathsf{A}})}\), stressing that \(\nu(\boldsymbol{\mathsf{A}}^\intercal) \neq \nu(\boldsymbol{\mathsf{A}})\) in general. In the following statements, we also denote by \(\boldsymbol{\mathsf{A}}_{-i}\) the Bochner matrix with its \(i\)th row \(\boldsymbol{\mathsf{a}}_i^\intercal\) removed.

Lemma 4. Let \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{N \times \kappa}\) have \(\mathrm{rank_{c}} ( \boldsymbol{\mathsf{A}} ) = \kappa\). Let \(\mathsf{\Gamma} = \boldsymbol{\mathsf{A}}^\ast \boldsymbol{\mathsf{A}}\) and \(\mathsf{\Gamma}_i = (\boldsymbol{\mathsf{a}}_{i}^\intercal)^\ast \boldsymbol{\mathsf{a}}_{i}^\intercal\) for every \(1 \leq i \leq N\). Then \[\frac{\nu^2(\boldsymbol{\mathsf{A}}_{-i})}{\nu^2(\boldsymbol{\mathsf{A}})} = \det(\mathsf{I}_{\kappa} - \mathsf{S}_i), \quad \mathsf{S}_i = \mathsf{\Gamma}^{-\tfrac{1}{2}} \mathsf{\Gamma}_i \mathsf{\Gamma}^{-\tfrac{1}{2}}.\] Let \(\Delta_{i}(\boldsymbol{\mathsf{A}}) = \det(\mathsf{I}_{\kappa} - \mathsf{S}_i) - (1 - \mathrm{tr}(\mathsf{S}_i))\) and \(\Delta(\boldsymbol{\mathsf{A}}) = \sum_{i = 1}^{N} \Delta_{i}(\boldsymbol{\mathsf{A}})\). Then \[\sum_{i = 1}^{N} \frac{\nu^2(\boldsymbol{\mathsf{A}}_{-i})}{\nu^2(\boldsymbol{\mathsf{A}})} = N - \kappa + \Delta(\boldsymbol{\mathsf{A}}).\]

Proof. By assumption, \(\mathsf{\Gamma} \in \mathbb{F}^{\kappa \times \kappa}\) has full rank and \(\nu(\boldsymbol{\mathsf{A}}) \neq 0\); consequently, all \(\mathsf{S}_i\) are well-defined positive semidefinite matrices. Furthermore, \(\mathsf{\Gamma} = \sum_{i = 1}^{N} \mathsf{\Gamma}_i\) and \(\sum_{i = 1}^{N} \mathsf{S}_i = \mathsf{I}_{\kappa}\). We expand the definition of \(\nu^2(\boldsymbol{\mathsf{A}}_{-i})\) to prove the first equality: \[\nu^2(\boldsymbol{\mathsf{A}}_{-i}) = \det(\boldsymbol{\mathsf{A}}_{-i}^\ast \boldsymbol{\mathsf{A}}_{-i}) = \det(\mathsf{\Gamma} - \mathsf{\Gamma}_i) = \nu^2(\boldsymbol{\mathsf{A}}) \det(\mathsf{I}_{\kappa} - \mathsf{S}_i).\] To obtain the second equality, observe that \[\sum_{i = 1}^{N} (1 - \mathrm{tr}(\mathsf{S}_i)) = N - \mathrm{tr}\left(\sum_{i = 1}^N \mathsf{S}_i\right) = N - \mathrm{tr}(\mathsf{I}_{\kappa}) = N - \kappa.\] ◻

We call \(\Delta_i(\boldsymbol{\mathsf{A}})\) the intra-row rank excess of \(\boldsymbol{\mathsf{A}}\) and \(\Delta(\boldsymbol{\mathsf{A}})\) its inter-row rank excess. Since the matrices \(\mathsf{S}_i\) are positive semidefinite and satisfy \(\sum_{i = 1}^{N} \mathsf{S}_i = \mathsf{I}_{\kappa}\), their eigenvalues lie in \([0, 1]\). Thus, by the Weierstrass product inequality, we have \(\Delta_i(\boldsymbol{\mathsf{A}}) \geq 0\); moreover, \(\Delta_i(\boldsymbol{\mathsf{A}})\) vanishes if and only if \(\mathrm{rank} ( \mathsf{\Gamma}_i ) \leq 1\). Based on the SVD of \(\boldsymbol{\mathsf{a}}_i^\intercal\), we can establish a simple yet fundamental relationship: \[\mathrm{rank} ( \mathsf{\Gamma}_i ) = \mathrm{rank_{c}} ( \boldsymbol{\mathsf{a}}_i^\intercal ).\] Therefore, the inter-row rank excess is always zero for classical matrices.

Lemma 5. In the setting of 4, let \(i_\ast\) be such that \(\mathrm{rank_{c}} ( \boldsymbol{\mathsf{A}}_{-i_\ast} ) = \kappa\), and denote \(\mathsf{\Gamma}_{-i_\ast} = \boldsymbol{\mathsf{A}}_{-i_\ast}^\ast \boldsymbol{\mathsf{A}}_{-i_\ast}\). Then \[\frac{\nu^2(\boldsymbol{\mathsf{A}})}{\nu^2(\boldsymbol{\mathsf{A}}_{-i_\ast})} = \det\Big(\mathsf{I}_{\kappa} + \mathsf{\Gamma}_{-i_\ast}^{-\tfrac{1}{2}} \mathsf{\Gamma}_{i_\ast} \mathsf{\Gamma}_{-i_\ast}^{-\tfrac{1}{2}}\Big) \geq 1 + \| \boldsymbol{\mathsf{a}}_{i_\ast}^\intercal\boldsymbol{\mathsf{A}}_{-i_\ast}^{\dagger} \|_{\mathrm{\ell_2(\mathrm{H}) \to \mathrm{H}}}^2.\]

Proof. The equality follows from the definition of volume. To prove the inequality, consider vectors \(\mathsf{x} \in \mathbb{F}^{\kappa}\) of unit \(\ell_2\) norm and let \(\boldsymbol{\mathsf{A}}_{-i_\ast} = \boldsymbol{\mathsf{P}} \mathsf{\Lambda} \mathsf{Q}^\ast\) be an SVD. Then \[\begin{align} \Big\| \mathsf{\Gamma}_{-i_\ast}^{-\tfrac{1}{2}} \mathsf{\Gamma}_{i_\ast} \mathsf{\Gamma}_{-i_\ast}^{-\tfrac{1}{2}} \Big\|_{\mathrm{\ell_2 \to \ell_2}} &= \max_{\mathsf{x}} \Big| \mathsf{x}^\ast \mathsf{\Gamma}_{-i_\ast}^{-\tfrac{1}{2}} \mathsf{\Gamma}_{i_\ast} \mathsf{\Gamma}_{-i_\ast}^{-\tfrac{1}{2}} \mathsf{x} \Big| = \max_{\mathsf{x}} \Big\| \boldsymbol{\mathsf{a}}_{i_\ast}^\intercal\mathsf{\Gamma}_{-i_\ast}^{-\tfrac{1}{2}} \mathsf{x} \Big\|_{\mathrm{\mathrm{H}}}^2 = \Big\| \boldsymbol{\mathsf{a}}_{i_\ast}^\intercal\mathsf{\Gamma}_{-i_\ast}^{-\tfrac{1}{2}} \Big\|_{\mathrm{\ell_2 \to \mathrm{H}}}^2 \\ &= \Big\| \boldsymbol{\mathsf{a}}_{i_\ast}^\intercal\mathsf{Q} \mathsf{\Lambda}^{-1} \Big\|_{\mathrm{\ell_2 \to \mathrm{H}}}^2 = \Big\| \boldsymbol{\mathsf{a}}_{i_\ast}^\intercal\mathsf{Q} \mathsf{\Lambda}^{-1} \boldsymbol{\mathsf{P}}^\ast \Big\|_{\mathrm{\ell_2(\mathrm{H}) \to \mathrm{H}}}^2 = \| \boldsymbol{\mathsf{a}}_{i_\ast}^\intercal\boldsymbol{\mathsf{A}}_{-i_\ast}^{\dagger} \|_{\mathrm{\ell_2(\mathrm{H}) \to \mathrm{H}}}^2, \end{align}\] where we use the unitary invariance of operator norms. The inequality follows since the remaining eigenvalues are nonnegative. ◻

Combining 4 5, we can bound the volume of \(\boldsymbol{\mathsf{A}}\) as \[\Big( 1 + \| \boldsymbol{\mathsf{a}}_{i_\ast}^\intercal\boldsymbol{\mathsf{A}}_{-i_\ast}^{\dagger} \|_{\mathrm{\ell_2(\mathrm{H}) \to \mathrm{H}}}^2 \Big) \nu^2(\boldsymbol{\mathsf{A}}_{-i_\ast}) \leq \nu^2(\boldsymbol{\mathsf{A}}) = \frac{N}{N - \kappa + \Delta(\boldsymbol{\mathsf{A}})} \bigg(\frac{1}{N} \sum_{i = 1}^{N} \nu^2(\boldsymbol{\mathsf{A}}_{-i})\bigg).\] Suppose that the squared volume of the submatrix \(\boldsymbol{\mathsf{A}}_{-i_\ast}\) is at least the average squared volume across all \((N-1) \times \kappa\) submatrices of \(\boldsymbol{\mathsf{A}}\). Then this volume bound yields \[\label{eq:maxvol95lsq95bound} \| \boldsymbol{\mathsf{a}}_{i_\ast}^\intercal\boldsymbol{\mathsf{A}}_{-i_\ast}^{\dagger} \|_{\mathrm{\ell_2(\mathrm{H}) \to \mathrm{H}}}^2 \leq \frac{\kappa - \Delta(\boldsymbol{\mathsf{A}})}{N - \kappa + \Delta(\boldsymbol{\mathsf{A}})},\tag{3}\] which extends the classical maximum-volume least-squares bound to Bochner matrices and relies on a relaxed above-average volume condition. We are finally ready to derive an operator-norm bound for the pseudoinverse of a maximum-volume submatrix.

Theorem 8. Let \(\boldsymbol{\mathsf{U}} \in \mathrm{H}^{m \times \kappa}\) have orthonormal columns. Let indices \({I} \subseteq [m]\) be such that \(\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}\) has the largest volume among all \(|{I}| \times \kappa\) submatrices of \(\boldsymbol{\mathsf{U}}\). If \(\nu(\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}) > 0\), \[\| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}} \leq \sqrt{1 + \sum_{i \not\in {I}} \frac{\kappa - \Delta(\tilde{\boldsymbol{\mathsf{R}}}_{i})}{|{I}| + 1 - \kappa + \Delta(\tilde{\boldsymbol{\mathsf{R}}}_{i})}}, \quad \tilde{\boldsymbol{\mathsf{R}}}_{i} = \begin{bmatrix} \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}} \\ \boldsymbol{\mathsf{u}}_{i}^\intercal \end{bmatrix}.\]

Proof. Using the standard properties of the operator norm, we have \[\begin{align} \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}}^2 &= \| \boldsymbol{\mathsf{U}} \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2(\mathrm{H})}}^2 = \| (\boldsymbol{\mathsf{U}} \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger)^\ast (\boldsymbol{\mathsf{U}} \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger) \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2(\mathrm{H})}} \\ &= \bigg\| (\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}} \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger)^\ast (\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}} \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger) + \sum_{i \not\in {I}} (\boldsymbol{\mathsf{u}}_{i}^\intercal\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger)^\ast (\boldsymbol{\mathsf{u}}_{i}^\intercal\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger) \bigg\|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2(\mathrm{H})}} \\ &= \bigg\| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}} \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger+ \sum_{i \not\in {I}} (\boldsymbol{\mathsf{u}}_{i}^\intercal\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger)^\ast (\boldsymbol{\mathsf{u}}_{i}^\intercal\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger) \bigg\|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2(\mathrm{H})}}\\ &\leq 1 + \sum_{i \not\in {I}} \| (\boldsymbol{\mathsf{u}}_{i}^\intercal\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger)^\ast (\boldsymbol{\mathsf{u}}_{i}^\intercal\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger) \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2(\mathrm{H})}} \\ &= 1 + \sum_{i \not\in {I}} \| \boldsymbol{\mathsf{u}}_{i}^\intercal\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \mathrm{H}}}^2, \end{align}\] since \(\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}} \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger\) is an orthogonal projector. We then bound each term in the sum using 3 ; this is allowed because the positive volume of \(\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}\) implies \(\mathrm{rank_{c}} ( \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}} ) = \kappa\). ◻

For classical matrices, 8 recovers the standard well-conditioning bound for maximum-volume submatrices. Meanwhile, the presence of inter-row rank excesses in the bound shows that maximum-volume submatrices are even better conditioned in the Bochner setting; indeed, the right-hand side in 3 is a monotonically decreasing function of the inter-row rank excess.

To illustrate the conditioning gains, consider \(\boldsymbol{\mathsf{Q}} \in \mathrm{H}^{m \times \kappa}\) whose entries are mutually orthonormal in \(\mathrm{H}\) and set \(\boldsymbol{\mathsf{U}} = \boldsymbol{\mathsf{Q}} / \sqrt{m}\). Due to symmetry, all \(|{I}| \times \kappa\) submatrices have the same volume, and it holds for every row-augmentation in 8 that \[\Delta(\tilde{\boldsymbol{\mathsf{R}}}_{i}) = (|{I}| + 1) \left( 1 - \frac{1}{|{I}| + 1} \right)^{\kappa} - (|{I}| + 1 - \kappa),\] which we can substitute into the bound, simplifying to \[\| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}} \leq \sqrt{1 + (m - |{I}|) \Big[ \Big(1 + \frac{1}{|{I}|}\Big)^\kappa - 1 \Big]} < \sqrt{1 + (m - |{I}|) \big( e^{\kappa / |{I}|} - 1 \big)}.\] It is instructive to compare this bound with the bound for classical maximum-volume submatrices 2 ; see 1. When a square submatrix is selected, classical matrices enjoy the bound \(\sqrt{1 + m \kappa}\), whereas the bound for \(\boldsymbol{\mathsf{U}}\) is \(\sqrt{1 + (e-1)m}\) and does not grow with \(\kappa\). A standard way to regularize classical submatrices is to sample twice as many rows as there are columns, improving the bound to \(\sqrt{1 + m}\). For the specific \(\boldsymbol{\mathsf{U}}\), such oversampling leads to an even smaller bound \(\sqrt{1 + (e^{1/2}-1)m}\). Unlike classical matrices, the rows of Bochner matrices can be undersampled, while still delivering well-conditioned submatrices. For example, the selection of \(\kappa / \log(\kappa)\) rows results in a bound that is similar to the case of classical square submatrices.

Table 1: Comparison of operator-norm bounds for the pseudoinverse of maximum-volume submatrices of classical and Bochner matrices with orthonormal columns.
\(|\Index{I}|\) Classical matrices Bochner matrix \(\Batrix{U}\)
\(\kappa / \log(\kappa)\) \(\sqrt{1 + (m - \kappa / \log(\kappa)) (\kappa - 1)}\)
\(\kappa\) \(\sqrt{1 + (m - \kappa) \kappa}\) \(\sqrt{1 + (m - \kappa) (e - 1)}\)
\(c\kappa,~c \geq 1\) \(\sqrt{1 + (m - c\kappa) \frac{\kappa}{(c-1)\kappa + 1}}\) \(\sqrt{1 + (m - c\kappa) (e^{1/c} - 1)}\)

Let us note that 8 holds for any \(\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}\) that has the largest volume among all submatrices that disagree with it in exactly one row. Relaxed further, it is sufficient for \(\boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}\) to have above-average volume within each individual augmentation \(\tilde{\boldsymbol{\mathsf{R}}}_{i}\).

7.4 Entrywise error bound for one-sided interpolation↩︎

Another important class of approximation guarantees for classical cross approximation are entrywise error bounds for when the submatrix at the intersection of the selected columns and rows itself exhibits maximum-volume properties. Following [45], we define the projective volume of a Bochner matrix as \(\nu_{\kappa}(\boldsymbol{\mathsf{A}}) = \prod_{s = 1}^{\kappa} \sigma_s(\boldsymbol{\mathsf{A}})\) and truncate the rank of the pseudoinverse for regularization. We also introduce a projective inter-row rank excess \(\Delta_{\kappa}(\boldsymbol{\mathsf{A}}) = \sup_{\mathsf{Q}_{\kappa}} \Delta(\boldsymbol{\mathsf{A}} \mathsf{Q}_{\kappa})\), where \(\mathsf{Q}_\kappa\) are (non-unique) truncated right singular factors of \(\boldsymbol{\mathsf{A}}\). This quantity is well-defined when \(\mathrm{rank_{c}} ( \boldsymbol{\mathsf{A}} ) \geq \kappa\).

Theorem 9. Let \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{m \times n}\) and \(\kappa \in \mathbb{N}\). Let \(\boldsymbol{\mathsf{G}}\) have the maximum \(\kappa\)-projective volume among all \(|{I}| \times |{J}|\) submatrices of \(\boldsymbol{\mathsf{A}}\). If \(\nu_{\kappa}(\boldsymbol{\mathsf{G}}) > 0\), \[\| \boldsymbol{\mathsf{A}} - \boldsymbol{\mathsf{C}} \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_{\infty}(\mathrm{H})}} \leq \sqrt{1 + \frac{\kappa - \Delta_\kappa'}{|{I}| + 1 - \kappa + \Delta_\kappa'}} \sqrt{1 + \frac{\kappa}{|{J}| + 1 - \kappa}} \sigma_{\kappa + 1}(\boldsymbol{\mathsf{A}})\] with \(\Delta_\kappa' = \min_{i \not\in {I}} \min_{1 \leq j \leq n} \Delta_{\kappa}\big( \boldsymbol{\mathsf{A}}({I} \cup \{i\}, {J} \cup \{j\}) \big)\).

Proof. For fixed indices \(1 \leq i \leq m\) and \(1 \leq j \leq n\), denote \[{\boldsymbol{a}} = \boldsymbol{\mathsf{A}}(i, j), \quad \boldsymbol{\mathsf{c}}^\intercal= \boldsymbol{\mathsf{C}}(i, :), \quad \boldsymbol{\mathsf{r}} = \boldsymbol{\mathsf{R}}(:, j), \quad {\boldsymbol{e}} = \boldsymbol{\mathsf{a}} - \boldsymbol{\mathsf{c}}^\intercal\lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{r}}.\] Our goal is to bound \(\| {\boldsymbol{e}} \|_{\mathrm{\mathrm{H}}}\) depending on whether \(i\) and \(j\) fall into the index sets \({I}\) and \({J}\), respectively. We shall assume that \({\boldsymbol{e}} \neq 0\), for otherwise the bounds are trivial.

First, let \(i \in {I}\) and \(j \in {J}\). Then \[\begin{align} \| {\boldsymbol{e}} \|_{\mathrm{\mathrm{H}}} &\leq \| \boldsymbol{\mathsf{G}} - \boldsymbol{\mathsf{G}} \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{G}} \|_{\mathrm{\ell_{\infty}(\mathrm{H})}} = \| \boldsymbol{\mathsf{G}} - \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa} \|_{\mathrm{\ell_{\infty}(\mathrm{H})}} \\ &\leq \| \boldsymbol{\mathsf{G}} - \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa} \|_{\mathrm{\ell_2 \to \ell_2(\mathrm{H})}} = \sigma_{\kappa + 1}(\boldsymbol{\mathsf{G}}) \leq \sigma_{\kappa + 1}(\boldsymbol{\mathsf{A}}), \end{align}\] where we use standard properties of the pseudoinverse and singular values.

Second, let \(i \not\in {I}\) and \(j \not\in {J}\). Construct an augmented submatrix \[\boldsymbol{\mathsf{M}} = \begin{bmatrix} \boldsymbol{\mathsf{G}} & \boldsymbol{\mathsf{r}} \\ \boldsymbol{\mathsf{c}}^\intercal& \boldsymbol{\mathsf{a}} \end{bmatrix} \in \mathrm{H}^{(|{I}| + 1) \times (|{J}| + 1)}\] and let \(\lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa} = \boldsymbol{\mathsf{U}}_{\kappa} \mathsf{\Sigma}_{\kappa} \mathsf{V}_{\kappa}^\ast\) be a truncated SVD. Consider an auxiliary matrix \[\mathsf{T} = \begin{bmatrix} \mathsf{V}_{\kappa} & -\lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{r}} \\ 0 & 1 \end{bmatrix} \in \mathbb{F}^{(|{J}| + 1) \times (\kappa + 1)}\] and multiply \(\boldsymbol{\mathsf{M}}\) with it on the right: \[\boldsymbol{\mathsf{K}} = \boldsymbol{\mathsf{M}} \mathsf{T} = \begin{bmatrix} \boldsymbol{\mathsf{U}}_{\kappa} \mathsf{\Sigma}_{\kappa} & \boldsymbol{\mathsf{r}}^\perp \\ \boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa} & {\boldsymbol{e}} \end{bmatrix} \in \mathrm{H}^{(|{I}| + 1) \times (\kappa + 1)}, \quad \boldsymbol{\mathsf{r}}^\perp = \boldsymbol{\mathsf{r}} - \boldsymbol{\mathsf{U}}_{\kappa} \boldsymbol{\mathsf{U}}_{\kappa}^\ast \boldsymbol{\mathsf{r}} \in \mathrm{H}^{|{I}|}.\] By the definition of volume, we get \[\nu^2(\boldsymbol{\mathsf{K}}) = \det(\boldsymbol{\mathsf{K}}^\ast \boldsymbol{\mathsf{K}}) = \det \left( \begin{bmatrix} \mathsf{\Sigma}_{\kappa}^2 + (\boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa})^\ast \boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa} & (\boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa})^\ast {\boldsymbol{e}} \\ {\boldsymbol{e}}^\ast \boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa} & \| \boldsymbol{\mathsf{r}}^\perp \|_{\mathrm{\ell_2(\mathrm{H})}}^2 + \| {\boldsymbol{e}} \|_{\mathrm{\mathrm{H}}}^2 \end{bmatrix} \right).\] Since \({\boldsymbol{e}} \neq 0\) by our assumption, we can invoke Schur’s complement: \[\nu^2(\boldsymbol{\mathsf{K}}) = (\| \boldsymbol{\mathsf{r}}^\perp \|_{\mathrm{\ell_2(\mathrm{H})}}^2 + \| {\boldsymbol{e}} \|_{\mathrm{\mathrm{H}}}^2) \det\bigg( \mathsf{\Sigma}_{\kappa}^2 + (\boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa})^\ast \boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa} - \frac{(\boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa})^\ast {\boldsymbol{e}} {\boldsymbol{e}}^\ast \boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa}}{\| \boldsymbol{\mathsf{r}}^\perp \|_{\mathrm{\ell_2(\mathrm{H})}}^2 + \| {\boldsymbol{e}} \|_{\mathrm{\mathrm{H}}}^2} \bigg).\] The matrix \[(\boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa})^\ast \boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa} - \frac{(\boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa})^\ast {\boldsymbol{e}} {\boldsymbol{e}}^\ast \boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa}}{\| {\boldsymbol{e}} \|_{\mathrm{\mathrm{H}}}^2} = (\boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa} - {\boldsymbol{e}} {\boldsymbol{e}}^\dagger\boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa})^\ast (\boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa} - {\boldsymbol{e}} {\boldsymbol{e}}^\dagger\boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa})\] is a Gram matrix and hence PSD; therefore, \[(\boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa})^\ast \boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa} - \frac{(\boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa})^\ast {\boldsymbol{e}} {\boldsymbol{e}}^\ast \boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa}}{\| \boldsymbol{\mathsf{r}}^\perp \|_{\mathrm{\ell_2(\mathrm{H})}}^2 + \| {\boldsymbol{e}} \|_{\mathrm{\mathrm{H}}}^2} \succcurlyeq (\boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa})^\ast \boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa} - \frac{(\boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa})^\ast {\boldsymbol{e}} {\boldsymbol{e}}^\ast \boldsymbol{\mathsf{c}}^\intercal\mathsf{V}_{\kappa}}{\| {\boldsymbol{e}} \|_{\mathrm{\mathrm{H}}}^2} \succcurlyeq 0,\] and the monotonicity of determinant for PSD matrices implies \[\nu^2(\boldsymbol{\mathsf{K}}) \geq (\| \boldsymbol{\mathsf{r}}^\perp \|_{\mathrm{\ell_2(\mathrm{H})}}^2 + \| {\boldsymbol{e}} \|_{\mathrm{\mathrm{H}}}^2) \det(\mathsf{\Sigma}_{\kappa}^2) \geq \| {\boldsymbol{e}} \|_{\mathrm{\mathrm{H}}}^2 \det(\mathsf{\Sigma}_{\kappa}^2).\] Let us relate the volumes of \(\boldsymbol{\mathsf{K}}\) and \(\boldsymbol{\mathsf{M}}\). Using Schur’s complement, we can show that the auxiliary matrix \(\mathsf{T}\) satisfies \[\nu_{\kappa + 1}^2(\mathsf{T}) = \nu^2(\mathsf{T}) = \det(\mathsf{T}^\ast \mathsf{T}) = \det\left(\begin{bmatrix} \mathsf{I}_{\kappa} & -\mathsf{\Sigma}_\kappa^{-1} \boldsymbol{\mathsf{U}}_\kappa^\ast \boldsymbol{\mathsf{r}} \\ (-\mathsf{\Sigma}_\kappa^{-1} \boldsymbol{\mathsf{U}}_\kappa^\ast \boldsymbol{\mathsf{r}})^\ast & 1 + \| \mathsf{\Sigma}_\kappa^{-1} \boldsymbol{\mathsf{U}}_\kappa^\ast \boldsymbol{\mathsf{r}} \|_{\mathrm{\ell_2}}^2 \end{bmatrix}\right) = 1.\] Let \(\boldsymbol{\mathsf{M}} = \boldsymbol{\mathsf{P}} \mathsf{\Lambda} \mathsf{Q}^\ast\) be an SVD, then by the log-majorization of singular values of a matrix product [46] we obtain \[\nu(\boldsymbol{\mathsf{K}}) = \nu(\boldsymbol{\mathsf{P}} \boldsymbol{\mathsf{P}}^\ast \boldsymbol{\mathsf{M}} \mathsf{T}) = \nu(\boldsymbol{\mathsf{P}}^\ast \boldsymbol{\mathsf{M}} \mathsf{T}) = \nu_{\kappa + 1}(\boldsymbol{\mathsf{P}}^\ast \boldsymbol{\mathsf{M}} \mathsf{T}) \leq \nu_{\kappa + 1}(\boldsymbol{\mathsf{P}}^\ast \boldsymbol{\mathsf{M}}) \nu_{\kappa + 1}(\mathsf{T}) = \nu_{\kappa + 1}(\boldsymbol{\mathsf{M}}),\] and therefore \[\| {\boldsymbol{e}} \|_{\mathrm{\mathrm{H}}} \leq \frac{\nu(\boldsymbol{\mathsf{K}})}{\nu_{\kappa}(\boldsymbol{\mathsf{G}})} \leq \frac{\nu_{\kappa + 1}(\boldsymbol{\mathsf{M}})}{\nu_{\kappa}(\boldsymbol{\mathsf{G}})} = \frac{\nu_{\kappa}(\boldsymbol{\mathsf{M}})}{\nu_{\kappa}(\boldsymbol{\mathsf{G}})} \sigma_{\kappa + 1}(\boldsymbol{\mathsf{M}}) \leq \frac{\nu_{\kappa}(\boldsymbol{\mathsf{M}})}{\nu_{\kappa}(\boldsymbol{\mathsf{G}})} \sigma_{\kappa + 1}(\boldsymbol{\mathsf{A}}).\] To bound the volume ratio, we shall expand \(\nu_{\kappa}^2(\boldsymbol{\mathsf{M}})\) into a sum of squared volumes of submatrices with removed columns and rows based on the fundamental identities \[\nu_{\kappa}(\boldsymbol{\mathsf{M}}) = \nu(\boldsymbol{\mathsf{P}}_{\kappa}^\ast \boldsymbol{\mathsf{M}}) = \nu(\boldsymbol{\mathsf{M}} \mathsf{Q}_{\kappa}) = \max_{\boldsymbol{\mathsf{\tilde{P}}}^\ast \boldsymbol{\mathsf{\tilde{P}}} = \mathsf{I}_{\kappa}} \nu(\boldsymbol{\mathsf{\tilde{P}}}^\ast \boldsymbol{\mathsf{M}}) = \max_{\mathsf{\tilde{Q}}^\ast \mathsf{\tilde{Q}} = \mathsf{I}_{\kappa}} \nu(\boldsymbol{\mathsf{M}} \mathsf{\tilde{Q}}).\] Denote by \(\boldsymbol{\mathsf{M}}_{-s}\) the submatrix with \(s\)th row removed, by \(\boldsymbol{\mathsf{M}}_{:, -t}\) the submatrix with \(t\)th column removed, and by \(\boldsymbol{\mathsf{M}}_{-s,-t}\) the submatrix with both \(s\)th row and \(t\)th column removed. Applying 4 to \(\boldsymbol{\mathsf{M}} \mathsf{Q}_\kappa\), which is allowed since \(\nu_{\kappa}(\boldsymbol{\mathsf{G}}) > 0\), we obtain \[\sum_{s = 1}^{|{I}| + 1} \nu^2(\boldsymbol{\mathsf{M}}_{-s} \mathsf{Q}_\kappa) = \Big(|{I}| + 1 - \kappa + \Delta(\boldsymbol{\mathsf{M}} \mathsf{Q}_\kappa) \Big) \nu_{\kappa}^2(\boldsymbol{\mathsf{M}}).\] Next, we bound each term via \(\nu(\boldsymbol{\mathsf{M}}_{-s} \mathsf{Q}_\kappa) \leq \nu_{\kappa}(\boldsymbol{\mathsf{M}}_{-s})\). If \(\nu_{\kappa}(\boldsymbol{\mathsf{M}}_{-s}) = 0\) then we move on; otherwise, let \(\boldsymbol{\mathsf{P}}_{s,\kappa}\) be the truncated left singular factor of \(\boldsymbol{\mathsf{M}}_{-s}\). Then we use the equality \(\nu_{\kappa}(\boldsymbol{\mathsf{M}}_{-s}) = \nu(\boldsymbol{\mathsf{P}}_{s,\kappa}^\ast \boldsymbol{\mathsf{M}}_{-s})\) and apply 4 to \((\boldsymbol{\mathsf{P}}_{s,\kappa}^\ast \boldsymbol{\mathsf{M}}_{-s})^\intercal\). Because this is a classical matrix, its inter-row excess is zero, and hence \[\sum_{t = 1}^{|{J}| + 1} \nu^2(\boldsymbol{\mathsf{P}}_{s,\kappa}^\ast \boldsymbol{\mathsf{M}}_{-s,-t}) = \Big(|{J}| + 1 - \kappa \Big) \nu_{\kappa}^2(\boldsymbol{\mathsf{M}}_{-s}).\] Combining the two sums yields \[\begin{align} \sum_{s = 1}^{|{I}| + 1} \sum_{t = 1}^{|{J}| + 1} \nu_{\kappa}^2(\boldsymbol{\mathsf{M}}_{-s,-t}) &\geq \sum_{s = 1}^{|{I}| + 1} \sum_{t = 1}^{|{J}| + 1} \nu^2(\boldsymbol{\mathsf{P}}_{s,\kappa}^\ast \boldsymbol{\mathsf{M}}_{-s,-t}) = \sum_{s = 1}^{|{I}| + 1} \Big(|{J}| + 1 - \kappa \Big) \nu_{\kappa}^2(\boldsymbol{\mathsf{M}}_{-s}) \\ &\geq \Big(|{J}| + 1 - \kappa \Big) \sum_{s = 1}^{|{I}| + 1} \nu^2(\boldsymbol{\mathsf{M}}_{-s} \mathsf{Q}_\kappa) \\ &= \Big(|{J}| + 1 - \kappa \Big) \Big(|{I}| + 1 - \kappa + \Delta(\boldsymbol{\mathsf{M}} \mathsf{Q}_\kappa) \Big) \nu_{\kappa}^2(\boldsymbol{\mathsf{M}}), \end{align}\] and it remains to recall that \(\nu_{\kappa}(\boldsymbol{\mathsf{M}}_{-s,-t}) \leq \nu_{\kappa}(\boldsymbol{\mathsf{G}})\) and maximize over the SVDs of \(\boldsymbol{\mathsf{M}}\).

Third, let \(i \in {I}\) and \(j \not\in {J}\). Then we augment as \(\boldsymbol{\mathsf{M}} = [\boldsymbol{\mathsf{G}}~\boldsymbol{\mathsf{r}}]\) and use the same auxiliary matrix \(\mathsf{T}\) to form \(\boldsymbol{\mathsf{K}}\). We can then readily show that \[\nu^2(\boldsymbol{\mathsf{K}}) = \| \boldsymbol{\mathsf{r}}^\perp \|_{\mathrm{\ell_2(\mathrm{H})}}^2 \det(\mathsf{\Sigma}_{\kappa}^2) \geq \| {\boldsymbol{e}} \|_{\mathrm{\mathrm{H}}}^2 \det(\mathsf{\Sigma}_{\kappa}^2),\] since \({\boldsymbol{e}}\) is now a component of \(\boldsymbol{\mathsf{r}}^\perp\), proceed to the same bound \(\| {\boldsymbol{e}} \|_{\mathrm{\mathrm{H}}} \leq \tfrac{\nu_{\kappa}(\boldsymbol{\mathsf{M}})}{\nu_{\kappa}(\boldsymbol{\mathsf{G}})} \sigma_{\kappa + 1}(\boldsymbol{\mathsf{A}})\), and skip the row-removal step in volume expansions to get \[\| {\boldsymbol{e}} \|_{\mathrm{\mathrm{H}}} \leq \sqrt{1 + \frac{\kappa}{|{J}| + 1 - \kappa}} \sigma_{\kappa + 1}(\boldsymbol{\mathsf{A}}).\]

Fourth, let \(i \not\in {I}\) and \(j \in {J}\). Then we augment as \(\boldsymbol{\mathsf{M}} = [\boldsymbol{\mathsf{G}}^\intercal~\boldsymbol{\mathsf{c}}]^\intercal\in \mathrm{H}^{(|{I}| + 1) \times |{J}|}\) and introduce an auxiliary matrix \(\mathsf{T} = [\mathsf{V}_{\kappa}~\mathsf{e}_\omega - \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{r}}] \in \mathbb{F}^{|{J}| \times (\kappa + 1)}\), where \(\omega\) is the local position of \(j\) within \({J}\). This modification of \(\mathsf{T}\) ensures that \(\boldsymbol{\mathsf{K}}\) has exactly the same block structure as in the second considered case. An identical argument yields the same bound \(\| {\boldsymbol{e}} \|_{\mathrm{\mathrm{H}}} \leq \tfrac{\nu_{\kappa}(\boldsymbol{\mathsf{M}})}{\nu_{\kappa}(\boldsymbol{\mathsf{G}})} \sigma_{\kappa + 1}(\boldsymbol{\mathsf{A}})\), where now \[\det(\mathsf{T}^\ast \mathsf{T}) = \| \mathsf{e}_\omega - \mathsf{V}_{\kappa} \mathsf{V}_{\kappa}^\ast \mathsf{e}_\omega \|_{\mathrm{\ell_2}}^2 \leq 1.\] Applying 4 we obtain \[\sum_{s = 1}^{|{I}| + 1} \nu_{\kappa}^2(\boldsymbol{\mathsf{M}}_{-s}) \geq \sum_{s = 1}^{|{I}| + 1} \nu^2(\boldsymbol{\mathsf{M}}_{-s} \mathsf{Q}_\kappa) = \Big(|{I}| + 1 - \kappa + \Delta(\boldsymbol{\mathsf{M}} \mathsf{Q}_\kappa) \Big) \nu_{\kappa}^2(\boldsymbol{\mathsf{M}})\] and can directly bound \(\nu_{\kappa}(\boldsymbol{\mathsf{M}}_{-s}) \leq \nu_{\kappa}(\boldsymbol{\mathsf{G}})\) since \(\boldsymbol{\mathsf{M}}_{-s}\) is of size \(|{I}| \times |{J}|\). ◻

A standard remark is in order that 9 holds when the projective volume of \(\boldsymbol{\mathsf{G}}\) is maximum only among the submatrices that differ from it in at most one row and at most one column.

7.5 Error bounds for cross approximation↩︎

The preparatory work on the approximation guarantees for one-sided interpolation can now be streamlined to derive error bounds for cross approximation. We begin with a bound in the \(\ell_2(\mathrm{H})\) norm, the only unitarily invariant norm that is also invariant to transposition.

Lemma 6. Let \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{B}} + \boldsymbol{\mathsf{E}} \in \mathrm{H}^{m \times n}\) and suppose that \[\| \boldsymbol{\mathsf{B}} - \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_2(\mathrm{H})}} \leq \epsilon, \quad \| \boldsymbol{\mathsf{B}}^\intercal- \boldsymbol{\mathsf{R}}^\intercal(\boldsymbol{\mathsf{G}}^\intercal)^\dagger\boldsymbol{\mathsf{C}}^\intercal \|_{\mathrm{\ell_2(\mathrm{H})}} \leq \delta.\] Then \[\| \boldsymbol{\mathsf{B}} - \mathrm{cross}(\boldsymbol{\mathsf{A}}, {I}, {J} ) \|_{\mathrm{\ell_2(\mathrm{H})}} \leq \epsilon + \| \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_2 \to \ell_2}} (\delta + \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\ell_2(\mathrm{H})}})\]

Proof. We expand the error by triangle inequality to obtain \[\begin{align} \| \boldsymbol{\mathsf{B}} - \mathrm{cross}(\boldsymbol{\mathsf{A}}, {I}, {J} ) \|_{\mathrm{\ell_2(\mathrm{H})}} &\leq \| \boldsymbol{\mathsf{B}} - \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_2(\mathrm{H})}} + \| \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} - \mathrm{cross}(\boldsymbol{\mathsf{A}}, {I}, {J} ) \|_{\mathrm{\ell_2(\mathrm{H})}} \\ &\leq \epsilon + \| \boldsymbol{\mathsf{C}} - \big( (\boldsymbol{\mathsf{G}}^\intercal)^\dagger\boldsymbol{\mathsf{C}}^\intercal\big)^\intercal\boldsymbol{\mathsf{G}} \|_{\mathrm{\ell_2(\mathrm{H})}} \| \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_2 \to \ell_2}} \\ &= \epsilon + \| \boldsymbol{\mathsf{C}}^\intercal- \boldsymbol{\mathsf{G}}^\intercal(\boldsymbol{\mathsf{G}}^\intercal)^\dagger\boldsymbol{\mathsf{C}}^\intercal \|_{\mathrm{\ell_2(\mathrm{H})}} \| \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_2 \to \ell_2}} \\ &\leq \epsilon + \big( \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{B}}}^\intercal- \boldsymbol{\mathsf{G}}^\intercal(\boldsymbol{\mathsf{G}}^\intercal)^\dagger\boldsymbol{\mathsf{C}}^\intercal \|_{\mathrm{\ell_2(\mathrm{H})}} + \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\ell_2(\mathrm{H})}} \big) \| \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_2 \to \ell_2}} \\ &\leq \epsilon + \big( \| \boldsymbol{\mathsf{B}}^\intercal- \boldsymbol{\mathsf{R}}^\intercal(\boldsymbol{\mathsf{G}}^\intercal)^\dagger\boldsymbol{\mathsf{C}}^\intercal \|_{\mathrm{\ell_2(\mathrm{H})}} + \| \boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{E}}} \|_{\mathrm{\ell_2(\mathrm{H})}} \big) \| \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_2 \to \ell_2}}. \end{align}\] ◻

The following theorem builds upon 7, and we simplify its error bound for the sake of presentation. The proof relies on the sequential processing of a Bochner matrix and its transpose, so a symmetric bound holds too. The effects induced by the regularization of two pseudoinverses can be incorporated in the cross-approximation bound by means of 2 3, though we omit such extension for brevity.

Theorem 10. Let \(\boldsymbol{\mathsf{A}} = \boldsymbol{\mathsf{B}} + \boldsymbol{\mathsf{E}} \in \mathrm{H}^{m \times n}\) with \(\boldsymbol{\mathsf{B}} \neq 0\), and let \(\boldsymbol{\mathsf{B}} = \boldsymbol{\mathsf{U}} \mathsf{\Sigma} \mathsf{V}^\ast\) and \(\boldsymbol{\mathsf{B}}^\intercal= \boldsymbol{\mathsf{P}} \mathsf{\Lambda} \mathsf{Q}^\ast\) be SVDs. Denote \[\begin{gather} \eta_r = \| \boldsymbol{\mathsf{R}}_{\boldsymbol{\mathsf{U}}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}}, \quad \eta_c = \| \mathsf{C}_{\mathsf{V}^\ast}^\dagger \|_{\mathrm{\ell_2 \to \ell_2}}, \quad \xi_r = \| \mathsf{R}_{\mathsf{Q}}^\dagger \|_{\mathrm{\ell_2 \to \ell_2}}, \quad \xi_c = \| (\boldsymbol{\mathsf{C}}_{\boldsymbol{\mathsf{P}}^\intercal}^\intercal)^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}}, \\ \tilde{\eta} = \eta_r + \eta_c + 3 \eta_r \eta_c, \quad \tilde{\xi} = \xi_c + \xi_r + 3 \xi_c \xi_r. \end{gather}\] Let \(\mathrm{rank_{T}} ( \boldsymbol{\mathsf{G}}_{\boldsymbol{\mathsf{B}}} ) = \mathrm{rank_{T}} ( \boldsymbol{\mathsf{B}} )\), then \[\begin{align} \| \boldsymbol{\mathsf{B}} - \mathrm{cross}(\boldsymbol{\mathsf{A}}, {I}, {J} ) \|_{\mathrm{\ell_2(\mathrm{H})}} &\leq \tilde{\eta} (2 + \tilde{\xi}) \| \boldsymbol{\mathsf{E}} \|_{\mathrm{\ell_2(\mathrm{H})}} \\ &+ (1 + \tilde{\eta})(2 + \tilde{\xi}) \Big( \| \boldsymbol{\mathsf{G}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}} + \| (\boldsymbol{\mathsf{G}}^\intercal)^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}} \Big) \| \boldsymbol{\mathsf{E}} \|_{\mathrm{\ell_2(\mathrm{H})}}^2. \end{align}\]

Proof. Coarsening of the bound in 7 gives \[\| \boldsymbol{\mathsf{B}} - \boldsymbol{\mathsf{C}} \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_2(\mathrm{H})}} \leq \tilde{\eta} \| \boldsymbol{\mathsf{E}} \|_{\mathrm{\ell_2(\mathrm{H})}} + \| \boldsymbol{\mathsf{G}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}} (1 + \tilde{\eta}) \| \boldsymbol{\mathsf{E}} \|_{\mathrm{\ell_2(\mathrm{H})}}^2.\] By the same token, we bound \[\| \boldsymbol{\mathsf{B}}^\intercal- \boldsymbol{\mathsf{R}}^\intercal(\boldsymbol{\mathsf{G}}^\intercal)^\dagger\boldsymbol{\mathsf{C}}^\intercal \|_{\mathrm{\ell_2(\mathrm{H})}} \leq \tilde{\xi} \| \boldsymbol{\mathsf{E}} \|_{\mathrm{\ell_2(\mathrm{H})}} + \| (\boldsymbol{\mathsf{G}}^\intercal)^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}} (1 + \tilde{\xi}) \| \boldsymbol{\mathsf{E}} \|_{\mathrm{\ell_2(\mathrm{H})}}^2.\] To finalize the proof via 6, we need a bound on \(\| \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_2 \to \ell_2}}\), which we take from the end of the proof of 7: \[\| \boldsymbol{\mathsf{G}}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_2 \to \ell_2}} \leq \eta_c + \| \boldsymbol{\mathsf{G}}^\dagger \|_{\mathrm{\ell_2(\mathrm{H}) \to \ell_2}} (1 + \eta_c) \| \boldsymbol{\mathsf{E}} \|_{\mathrm{\ell_2(\mathrm{H})}}.\] We substitute the bounds in 6 and massage the expression using \(\eta_c \leq \tilde{\eta}\). ◻

To extend the maximum-volume bound from 9 to cross approximation, we need to assume that the projective volumes of both \(\boldsymbol{\mathsf{G}}\) and \(\boldsymbol{\mathsf{G}}^\intercal\) are simultaneously maximal across the submatrices of \(\boldsymbol{\mathsf{A}}\) and \(\boldsymbol{\mathsf{A}}^\intercal\), respectively. However, this is a stringent assumption, because submatrices maximizing both volumes at the same time do not exist for Bochner matrices in general. Therefore, we rely instead on a quasi-maximum volume assumption and require that \(\theta \nu_\kappa(\boldsymbol{\mathsf{G}})\) be larger than or equal to the maximum projective volume for \(\theta \geq 1\).

Theorem 11. Let \(\boldsymbol{\mathsf{A}} \in \mathrm{H}^{m \times n}\). Let \(\rho,\kappa \in \mathbb{N}\) and \(\theta_c, \theta_r \geq 1\). Suppose \(\boldsymbol{\mathsf{G}}\) has the \(\theta_c\)-quasi-maximum \(\kappa\)-projective volume among all \(|{I}| \times |{J}|\) submatrices of \(\boldsymbol{\mathsf{A}}\) and \(\boldsymbol{\mathsf{G}}^\intercal\) has the \(\theta_r\)-quasi-maximum \(\rho\)-projective volume among all \(|{J}| \times |{I}|\) submatrices of \(\boldsymbol{\mathsf{C}}^\intercal\). Denote \(\boldsymbol{\mathsf{X}} = (\lfloor \boldsymbol{\mathsf{G}}^\intercal \rfloor_{\rho}^\dagger\boldsymbol{\mathsf{C}}^\intercal)^\intercal\boldsymbol{\mathsf{G}} (\lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{R}})\). If \(\nu_{\kappa}(\boldsymbol{\mathsf{G}}) > 0\) and \(\nu_{\rho}(\boldsymbol{\mathsf{G}}^\intercal) > 0\), \[\begin{align} \| \boldsymbol{\mathsf{A}} - \boldsymbol{\mathsf{X}} \|_{\mathrm{\ell_{\infty}(\mathrm{H})}} &\leq \theta_c \sqrt{1 + \frac{\kappa - \Delta_\kappa'}{|{I}| + 1 - \kappa + \Delta_\kappa'}} \sqrt{1 + \frac{\kappa}{|{J}| + 1 - \kappa}} \sigma_{\kappa + 1}(\boldsymbol{\mathsf{A}}) \\ &+ \theta_r \theta_c \sqrt{1 + \frac{\rho}{|{I}| + 1 - \rho}} \sqrt{1 - \frac{1}{\theta_c^2} + \frac{\kappa}{|{J}| + 1 - \kappa}} \sigma_{\rho + 1}(\boldsymbol{\mathsf{C}}^\intercal). \end{align}\] with \(\Delta_\kappa' = \min_{i \not\in {I}} \min_{1 \leq j \leq n} \Delta_{\kappa}\big( \boldsymbol{\mathsf{A}}({I} \cup \{i\}, {J} \cup \{j\}) \big)\).

Proof. By triangle inequality, we have \[\| \boldsymbol{\mathsf{A}} - \boldsymbol{\mathsf{X}} \|_{\mathrm{\ell_{\infty}(\mathrm{H})}} \leq \| \boldsymbol{\mathsf{A}} - \boldsymbol{\mathsf{C}} \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_{\infty}(\mathrm{H})}} + \| \boldsymbol{\mathsf{C}} \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{R}} - \boldsymbol{\mathsf{X}} \|_{\mathrm{\ell_{\infty}(\mathrm{H})}}.\] We bound the first term via a quasi-maximum-volume variant of 9, the proof of which is a direct modification of the original proof: \[\| \boldsymbol{\mathsf{A}} - \boldsymbol{\mathsf{C}} \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_{\infty}(\mathrm{H})}} \leq \theta_c \sqrt{1 + \frac{\kappa - \Delta_\kappa'}{|{I}| + 1 - \kappa + \Delta_\kappa'}} \sqrt{1 + \frac{\kappa}{|{J}| + 1 - \kappa}} \sigma_{\kappa + 1}(\boldsymbol{\mathsf{A}}).\] Suppose that the maximum error in the second term is attained at the \((i,j)\)th entry and denote \(\boldsymbol{\mathsf{c}}^\intercal= \boldsymbol{\mathsf{C}}(i,:)\) and \(\boldsymbol{\mathsf{r}} = \boldsymbol{\mathsf{R}}(:, j)\) as previously: \[\begin{align} \| \boldsymbol{\mathsf{C}} \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{R}} - \boldsymbol{\mathsf{X}} \|_{\mathrm{\ell_{\infty}(\mathrm{H})}} &= \| (\boldsymbol{\mathsf{C}} - (\lfloor \boldsymbol{\mathsf{G}}^\intercal \rfloor_{\rho}^\dagger\boldsymbol{\mathsf{C}}^\intercal)^\intercal\boldsymbol{\mathsf{G}}) \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{R}} \|_{\mathrm{\ell_{\infty}(\mathrm{H})}} \\ &= \| (\boldsymbol{\mathsf{c}} - \boldsymbol{\mathsf{G}}^\intercal\lfloor \boldsymbol{\mathsf{G}}^\intercal \rfloor_{\rho}^\dagger\boldsymbol{\mathsf{c}})^\intercal\lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{r}} \|_{\mathrm{\mathrm{H}}} \\ &\leq \| \boldsymbol{\mathsf{c}} - \boldsymbol{\mathsf{G}}^\intercal\lfloor \boldsymbol{\mathsf{G}}^\intercal \rfloor_{\rho}^\dagger\boldsymbol{\mathsf{c}} \|_{\mathrm{\ell_2(\mathrm{H})}} \| \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{r}} \|_{\mathrm{\ell_2}}. \end{align}\] As in the proof of 9, we consider different cases based on the position of \((i,j)\) relative to \(({I}, {J})\). When \(i \in {I}\) and \(\boldsymbol{\mathsf{c}} = \boldsymbol{\mathsf{G}}^\intercal\mathsf{e}_{\omega}\), \[\begin{align} \| \boldsymbol{\mathsf{c}} - \boldsymbol{\mathsf{G}}^\intercal\lfloor \boldsymbol{\mathsf{G}}^\intercal \rfloor_{\rho}^\dagger\boldsymbol{\mathsf{c}} \|_{\mathrm{\ell_2(\mathrm{H})}} &= \| (\boldsymbol{\mathsf{G}}^\intercal- \boldsymbol{\mathsf{G}}^\intercal\lfloor \boldsymbol{\mathsf{G}}^\intercal \rfloor_{\rho}^\dagger\boldsymbol{\mathsf{G}}^\intercal) \mathsf{e}_{\omega} \|_{\mathrm{\ell_2(\mathrm{H})}} \\ &\leq \| \boldsymbol{\mathsf{G}}^\intercal- \lfloor \boldsymbol{\mathsf{G}}^\intercal \rfloor_{\rho} \|_{\mathrm{\ell_2 \to \ell_2(\mathrm{H})}} = \sigma_{\rho + 1}(\boldsymbol{\mathsf{G}}^\intercal) \leq \sigma_{\rho + 1}(\boldsymbol{\mathsf{C}}^\intercal). \end{align}\] When \(i \not\in {I}\), we repeat the argument from the third part of the proof of 9. Namely, we form an augmentation \(\boldsymbol{\mathsf{M}} = [\boldsymbol{\mathsf{G}}^\intercal~\boldsymbol{\mathsf{c}}]\), take an auxiliary matrix \(\mathsf{T}\), multiply \(\boldsymbol{\mathsf{K}} = \boldsymbol{\mathsf{M}} \mathsf{T}\), obtain an equality \[\nu(\boldsymbol{\mathsf{K}}) = \nu_{\rho}(\boldsymbol{\mathsf{G}}^{\intercal}) \| \boldsymbol{\mathsf{c}} - \boldsymbol{\mathsf{G}}^\intercal\lfloor \boldsymbol{\mathsf{G}}^\intercal \rfloor_{\rho}^\dagger\boldsymbol{\mathsf{c}} \|_{\mathrm{\ell_2(\mathrm{H})}},\] and derive the bound \[\| \boldsymbol{\mathsf{c}} - \boldsymbol{\mathsf{G}}^\intercal\lfloor \boldsymbol{\mathsf{G}}^\intercal \rfloor_{\rho}^\dagger\boldsymbol{\mathsf{c}} \|_{\mathrm{\ell_2(\mathrm{H})}} = \frac{\nu(\boldsymbol{\mathsf{K}})}{\nu_{\rho}(\boldsymbol{\mathsf{G}}^{\intercal})} \leq \frac{\nu_{\rho + 1}(\boldsymbol{\mathsf{M}})}{\nu_{\rho}(\boldsymbol{\mathsf{G}}^{\intercal})} \leq \theta_r \sqrt{1 + \frac{\rho}{|{I}| + 1 - \rho}} \sigma_{\rho + 1}(\boldsymbol{\mathsf{C}}^\intercal).\] When \(j \in {J}\) and \(\boldsymbol{\mathsf{r}} = \boldsymbol{\mathsf{G}} \mathsf{e}_{\omega}\), we have \[\| \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{r}} \|_{\mathrm{\ell_2}} = \| \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{G}} \mathsf{e}_{\omega} \|_{\mathrm{\ell_2}} \leq \| \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{G}} \|_{\mathrm{\ell_2 \to \ell_2}} \leq 1.\] When \(j \not\in {J}\), consider \(\boldsymbol{\mathsf{M}} = [\boldsymbol{\mathsf{G}}~\boldsymbol{\mathsf{r}}]\), and let \(\lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa} = \boldsymbol{\mathsf{U}}_\kappa \mathsf{\Sigma}_\kappa \mathsf{V}_\kappa^\ast\) be an SVD. Then, by the matrix determinant lemma, \[\begin{align} \nu_\kappa^2(\boldsymbol{\mathsf{M}}) \geq \nu^2(\boldsymbol{\mathsf{U}}_\kappa^\ast \boldsymbol{\mathsf{M}}) &= \det\Big( \mathsf{\Sigma}_\kappa^2 + (\boldsymbol{\mathsf{U}}_\kappa^\ast \boldsymbol{\mathsf{r}}) (\boldsymbol{\mathsf{U}}_\kappa^\ast \boldsymbol{\mathsf{r}})^\ast \Big) \\ &= \nu_\kappa^2(\boldsymbol{\mathsf{G}}) (1 + \| \mathsf{\Sigma}_\kappa^{-1} \boldsymbol{\mathsf{U}}_\kappa^\ast \boldsymbol{\mathsf{r}} \|_{\mathrm{\ell_2}}^2) = \nu_\kappa^2(\boldsymbol{\mathsf{G}}) (1 + \| \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{r}} \|_{\mathrm{\ell_2}}^2), \end{align}\] and consequently \[\| \lfloor \boldsymbol{\mathsf{G}} \rfloor_{\kappa}^\dagger\boldsymbol{\mathsf{r}} \|_{\mathrm{\ell_2}} \leq \sqrt{\frac{\nu_\kappa^2(\boldsymbol{\mathsf{M}})}{\nu_\kappa^2(\boldsymbol{\mathsf{G}})} - 1} \leq \sqrt{\theta_c^2 + \frac{\theta_c^2 \kappa}{|{J}| + 1 - \kappa} - 1} = \theta_c \sqrt{1 - \frac{1}{\theta_c^2} + \frac{\kappa}{|{J}| + 1 - \kappa}}.\] ◻

The bound in 11 is asymmetric, and a similar one holds with the role of rows and columns reversed. In addition, the statement still holds when \(\boldsymbol{\mathsf{G}}\) optimizes the volume among submatrices of \(\boldsymbol{\mathsf{A}}\) that differ from it in at most one row and at most one column, and \(\boldsymbol{\mathsf{G}}^\intercal\) optimizes the volume among submatrices of \(\boldsymbol{\mathsf{C}}^\intercal\) that differ from it in at most one column.

8 Numerical experiments↩︎

8.1 Adaptive Bochner cross approximation↩︎

For classical matrices, adaptive cross approximation (ACA, [22], [24]) is a widely adopted approach to index selection for cross approximation, whereby columns and rows are picked one by one. At its every iteration, the position of an entry of large absolute value is chosen in the current residual. In full pivoting, the largest entry in the whole residual is selected, but then every entry needs to be inspected. In practice, rook pivoting delivers a good balance between stability and efficiency and is de facto the standard way to select indices in ACA [47], [48]. [alg:rook] is a direct extension of rook pivoting to Bochner matrices.

\(i_* \gets \mathrm{argmax}_{i \in [m]} \| \boldsymbol{\mathsf{A}}(i, j_*) \|_{\mathrm{\mathrm{H}}}\)

The ACA algorithm maintains a decomposition \(\mathsf{A} = \mathsf{B} + \mathsf{E}\), starting with \(\mathsf{B} = 0\) and increasing its rank by one at every iteration by adding a rank-one cross approximation of \(\mathsf{E}\) to it. A fundamental property of ACA is additivity; that is, if \(\mathsf{B} = \mathrm{cross}(\mathsf{A}, {I}, {J} )\) with invertible submatrix \(\mathsf{A}({I}, {J})\) and \(\mathsf{E}(i,j) \neq 0\) then \[\mathsf{B} + \mathrm{cross}(\mathsf{E}, i,j ) = \mathrm{cross}(\mathsf{A}, {I} \cup \{ i \}, {J} \cup \{ j \} ).\] Therefore, ACA iteratively builds up cross approximation based on rank-one updates.

However, this property ceases to hold for Bochner matrices. Recall the example from 5 of a \(2 \times 2\) Bochner matrix whose entries are mutually orthonormal: \[\boldsymbol{\mathsf{B}} = \mathrm{cross}(\boldsymbol{\mathsf{A}}, 1,1 ) = \begin{bmatrix} \boldsymbol{\mathsf{A}}(1,1) & 0 \\ 0 & 0 \end{bmatrix}, \quad \boldsymbol{\mathsf{B}} + \mathrm{cross}(\boldsymbol{\mathsf{E}}, 2,2 ) = \begin{bmatrix} \boldsymbol{\mathsf{A}}(1,1) & 0 \\ 0 & \boldsymbol{\mathsf{A}}(2,2) \end{bmatrix},\] whereas \(\mathrm{cross}(\boldsymbol{\mathsf{A}}, \{1,2\}, \{1,2\} ) = \boldsymbol{\mathsf{A}}\). Rank-one updates do not add up to a higher-rank cross approximation, requiring a proper recomputation. We take this into account in our adaptive Bochner cross (ABC) approximation algorithm ([alg:abc]).

\({I} \gets \varnothing,~{J} \gets \varnothing,~\boldsymbol{\mathsf{B}} \gets 0\)

Note the presence of the inner loop, in which the residual undergoes fast approximate rank-one corrections. When \(n_{\mathrm{sub}} = 1\), the ABC algorithm extends an ACA-like algorithm developed in [36] from classical third-order tensors to Bochner matrices (neglecting one mode). With \(n_{\mathrm{sub}} > 1\), we attempt to garner more indices at reduced computational cost, making the algorithm more flexible.

Unlike classical ACA, the ABC algorithm can select the same indices multiple times. This follows from the weaker interpolation properties of Bochner cross approximation (6), rooted in the difference between column and row ranks.

8.2 Reduced-order modeling: Nonlinear Stokes equations↩︎

To validate the ABC algorithm as a non-intrusive ROM method for parametric PDEs, we consider the two-dimensional Stokes equations in \(\Omega \subset \mathbb{R}^2\), which describe the stationary flow of a fluid with dominating viscous forces. Let \(\mathsf{u} : \Omega \to \mathbb{R}^2\) be the velocity and \(p : \Omega \to \mathbb{R}\) be the pressure. When there are no external forces, the Stokes equations read \[\label{eq:stokes} -\mathrm{div}(2 \nu \epsilon(\mathsf{u}) - p \mathsf{I}_2) = 0,~\mathrm{div}(\mathsf{u}) = 0~\text{ in }~\Omega,\tag{4}\] where \(\epsilon(\mathsf{u}) = \tfrac{1}{2} \nabla \mathsf{u} + \tfrac{1}{2} (\nabla \mathsf{u})^\intercal\) and \(\nu\) is the viscosity, which we model as \[\nu_{\alpha,\beta}(\mathsf{u}) = \nu_0 + \alpha \| \nabla \mathsf{u} \|_{\mathrm{\ell_2}}^{\beta},\] making the equations nonlinear and parametric. We partition the boundary as \(\partial\Omega = \partial\Omega_D \cup \partial\Omega_N\) and impose the following boundary conditions: \[\label{eq:stokes95bcs} \mathsf{u} = 0~\text{ on }~\partial\Omega_D, \qquad -(2 \nu \epsilon(\mathsf{u}) - p \mathsf{I}_2) \mathsf{n} = p_0 \mathsf{n}~\text{ on }~\partial\Omega_N,\tag{5}\] where \(\mathsf{n} : \partial\Omega_N \to \mathbb{R}^2\) is the outward normal and \(p_0 : \partial\Omega_N \to \mathbb{R}\).

This boundary-value problem is a standard example used to illustrate the FEniCS5 library. We take the domain, mesh, and boundary conditions from [49] to make the problem 4 5 concrete and use FEniCS to solve it using Taylor–Hood finite elements with 6495 degrees of freedom. The mesh and the solution for \(\nu_0 = 0.01\) and \(\alpha = \beta = 0\) are shown in 1.

a

b

c

Figure 1: Finite-element solution to 4 5 obtained with FEniCS as in [49]..

a

b

c

d

Figure 2: Cross approximation of the solution map of 4 5 with ABC over 50 random seeds. Left: median and 90th percentile of the relative approximation error. Right: median of the row ranks (solid) and column ranks (dashed)..

We apply ABC to approximate the solution map of 4 5 on a \(100 \times 100\) uniform grid in the parameter space \((\alpha,\beta) \in [0,1] \times [0,1]\). To define the corresponding Bochner matrix, we equip the 6495-dimensional space of finite-element solutions with the \(\mathrm{H}^1(\Omega; \mathbb{R}^3)\) inner product. In the experiments, we fix the product \(n_{\mathrm{iter}} n_{\mathrm{sub}} = 16\) and vary \(n_{\mathrm{sub}} \in \{ 1, 2, 4, 8, 16 \}\), plotting the statistics over 50 random seeds in 2.

When \(n_{\mathrm{rook}} = 0\), the column selection is completely random and the approximation error stagnates, indicating that although the algorithm continues to blindly expand the index sets, the chosen subspaces fail to capture the structure of the solution map. Conversely, when \(n_{\mathrm{rook}} = 1\), both rows and columns are chosen adaptively, leading to decreasing approximation errors and fewer selected indices. As expected, the errors decay faster per iteration for smaller values of \(n_{\mathrm{sub}}\), because the algorithm relies less frequently on the approximate residual. Finally, we observe that fewer rows are selected than columns, which signifies that the dependence on \(\beta\) is stronger; indeed, \(\beta\) controls the degree of nonlinearity.

To assess the accuracy of ABC approximation, we compare it with the quasioptimal HOSVD approximation of the same Tucker rank. The results in 2 for ABC with \(n_{\mathrm{rook}} = 1\) and \(n_{\mathrm{sub}} = 1\) show that it achieves errors that are about twice the HOSVD errors, despite querying the PDE solver only a small number of times.

Table 2: Relative \(\ell_2(\Hilb^1)\) approximation errors of the solution map of [eq:stokes] [eq:stokes95bcs] obtained with ABC and HOSVD for the same Tucker rank.
\(n_{\mathrm{iter}}\) 4 8 12 16
\((|\Index{I}|, |\Index{J}|)\) \((3,3)\) \((6,7)\) \((7,11)\) \((9,15)\)
ABC \(2.2 \times 10^{-2}\) \(1.6 \times 10^{-4}\) \(1.8 \times 10^{-5}\) \(2.6 \times 10^{-6}\)
HOSVD \(1.0 \times 10^{-2}\) \(7.3 \times 10^{-5}\) \(8.3 \times 10^{-6}\) \(1.3 \times 10^{-6}\)

9 Remarks↩︎

Unrestricted Hilbert space. The developed framework of Bochner matrices and their cross approximation does not impose any constraints on the underlying Hilbert space, relying exclusively on its inner-product structure. In contrast, the framework of quasimatrices [28], [29] and quasitensors [50], [51] requires a reproducing kernel Hilbert space or assumes that the entries are continuous functions. Similarly, the infinite-dimensional extension of tubal tensors in [12] requires separability.

Structured sampling. The formalism of Bochner matrices and tensors provides an abstraction for data models where the naturally accessible unit is a fiber or slice, rather than a single entry of a classical tensor. Tensor completion with such structured sampling was studied in [52] and applied to tensorial parametric ROM of dynamical systems in [53][55].

Cross approximation of Bochner tensors. The present work establishes the theoretical foundation of cross approximation for Bochner matrices. In a follow-up paper [56], we extend this framework to the Tucker-cross approximation of Bochner tensors, focusing primarily on its applications to ROM of parametric PDEs.

Acknowledgements↩︎

I thank Vladimir Kazeev and Maxim Olshanskii for the discussions of tensorial parametric ROM, which have lead to this article.

References↩︎

[1]
D. Hachenberger and D. Jungnickel, Topics in Galois Fields, Springer, 2020, https://doi.org/10.1007/978-3-030-60806-4.
[2]
D. Silva, F. R. Kschischang, and R. Kotter, Communication over finite-field matrix channels, IEEE Trans Inf Theory, 56 (2010), pp. 1296–1305, https://doi.org/10.1109/TIT.2009.2039167.
[3]
T. Kleinjung, K. Aoki, J. Franke, A. K. Lenstra, E. Thomé, J. W. Bos, P. Gaudry, A. Kruppa, P. L. Montgomery, D. A. Osvik, et al., Factorization of a 768-bit RSA modulus, in CRYPTO 2010, 2010, pp. 333–350, https://doi.org/10.1007/978-3-642-14623-7_18.
[4]
W. C. Brown, Matrices over Commutative Rings, CRC Press, 1992. isbn: 9780824787554.
[5]
J. A. Foster, J. G. McWhirter, M. R. Davies, and J. A. Chambers, An algorithm for calculating the QR and singular value decompositions of polynomial matrices, IEEE Trans Signal Process, 58 (2009), pp. 1263–1274, https://doi.org/10.1109/TSP.2009.2034325.
[6]
D. Cescato and H. Bölcskei, QR decomposition of Laurent polynomial matrices sampled on the unit circle, IEEE Trans Inf Theory, 56 (2010), pp. 4754–4761, https://doi.org/10.1109/TIT.2010.2054454.
[7]
A. Bunse-Gerstner, R. Byers, V. Mehrmann, and N. K. Nichols, Numerical computation of an analytic singular value decomposition of a matrix valued function, Numer Math, 60 (1991), pp. 1–39, https://doi.org/10.1007/bf01385712.
[8]
S. Weiss, I. K. Proudler, G. Barbarino, J. Pestana, and J. G. McWhirter, Properties and structure of the analytic singular value decomposition, IEEE Trans Signal Process, 72 (2024), pp. 2260–2275, https://doi.org/10.1109/TSP.2024.3387726.
[9]
V. W. Neo, S. Redif, J. G. McWhirter, J. Pestana, I. K. Proudler, S. Weiss, and P. A. Naylor, Polynomial eigenvalue decomposition for multichannel broadband signal processing: a mathematical technique offering new insights and solutions, IEEE Signal Process Mag, 40 (2023), pp. 18–37, https://doi.org/10.1109/MSP.2023.3269200.
[10]
M. E. Kilmer and C. D. Martin, Factorization strategies for third-order tensors, Linear Algebra Appl, 435 (2011), pp. 641–658, https://doi.org/10.1016/j.laa.2010.09.020.
[11]
M. E. Kilmer, K. Braman, N. Hao, and R. C. Hoover, Third-order tensors as operators on matrices: A theoretical and computational framework with applications in imaging, SIAM J Matrix Anal Appl, 34 (2013), pp. 148–172, https://doi.org/10.1137/110837711.
[12]
U. Mor and H. Avron, Quasitubal tensor algebra over separable Hilbert spaces, 2025, https://doi.org/10.48550/arXiv.2504.16231.
[13]
R. Dian and S. Li, Hyperspectral image super-resolution via subspace-based low tensor multi-rank regularization, IEEE Trans Image Process, 28 (2019), pp. 5135–5146, https://doi.org/10.1109/TIP.2019.2916734.
[14]
F. Zhang, Quaternions and matrices of quaternions, Linear Algebra Appl, 251 (1997), pp. 21–57, https://doi.org/10.1016/0024-3795(95)00543-9.
[15]
S. Miron, N. Le Bihan, and J. I. Mars, Quaternion-MUSIC for vector-sensor array processing, IEEE Trans Signal Process, 54 (2006), pp. 1218–1229, https://doi.org/10.1109/TSP.2006.870630.
[16]
Z. Jia, M. K. Ng, and G.-J. Song, Robust quaternion matrix completion with applications to image inpainting, Numer Linear Algebra Appl, 26 (2019), p. e2245, https://doi.org/10.1002/nla.2245.
[17]
S. Miron, J. Flamant, N. Le Bihan, P. Chainais, and D. Brie, Quaternions in signal and image processing: A comprehensive and objective overview, IEEE Signal Process Mag, 40 (2023), pp. 26–40, https://doi.org/10.1109/MSP.2023.3278071.
[18]
C. Tretter, Spectral Theory of Block Operator Matrices and Applications, World Scientific, 2008, https://doi.org/10.1142/p493.
[19]
G. Barbarino, C. Garoni, and S. Serra-Capizzano, Block generalized locally Toeplitz sequences: theory and applications in the unidimensional case, Electron Trans Numer Anal, 53 (2020), pp. 28–112, https://doi.org/10.1553/etna_vol53s28.
[20]
H. Brezis, Functional Analysis, Sobolev Spaces and Partial Differential Equations, Springer, 2011, https://doi.org/10.1007/978-0-387-70914-7.
[21]
S. A. Goreinov, E. E. Tyrtyshnikov, and N. L. Zamarashkin, A theory of pseudoskeleton approximations, Linear Algebra Appl, 261 (1997), pp. 1–21, https://doi.org/10.1016/S0024-3795(96)00301-1.
[22]
M. Bebendorf, Approximation of boundary element matrices, Numer Math, 86 (2000), pp. 565–589, https://doi.org/10.1007/PL00005410.
[23]
J. Ballani and L. Grasedyck, Hierarchical tensor approximation of output quantities of parameter-dependent PDEs, SIAM/ASA J Uncertain Quantif, 3 (2015), pp. 852–872, https://doi.org/10.1137/140960980.
[24]
M. Bebendorf and S. Rjasanow, Adaptive low-rank approximation of collocation matrices, Computing, 70 (2003), pp. 1–24, https://doi.org/10.1007/s00607-002-1469-6.
[25]
A. Ben-Israel and T. N. Greville, Generalized Inverses: Theory and Applications, Springer, 2 ed., 2003, https://doi.org/10.1007/b97366.
[26]
I. Gohberg and M. Krein, Introduction to the Theory of Linear Nonselfadjoint Operators in Hilbert Space, AMS, 1969. isbn: 978-1-4704-4436-5.
[27]
L. N. Trefethen, Householder triangularization of a quasimatrix, IMA J Numer Anal, 30 (2010), pp. 887–897, https://doi.org/10.1093/imanum/drp018.
[28]
A. Townsend and L. N. Trefethen, Continuous analogues of matrix factorizations, Proc R Soc A: Math Phys Eng Sci, 471 (2015), 20140585, https://doi.org/10.1098/rspa.2014.0585.
[29]
P. F. Shustin and H. Avron, Semi-infinite linear regression and its applications, SIAM J Matrix Anal Appl, 43 (2022), pp. 479–511, https://doi.org/10.1137/21M1411950.
[30]
D. Kressner, T. Ni, and A. Uschmajew, On the approximation of vector-valued functions by volume sampling, J Complex, 86 (2025), p. 101887, https://doi.org/10.1016/j.jco.2024.101887.
[31]
K. Hamm and L. Huang, Perspectives on CUR decompositions, Appl Comput Harmon Anal, 48 (2020), pp. 1088–1099, https://doi.org/10.1016/j.acha.2019.08.006.
[32]
M. W. Mahoney and P. Drineas, CUR matrix decompositions for improved data analysis, Proc Natl Acad Sci, 106 (2009), pp. 697–702, https://doi.org/10.1073/pnas.0803205106.
[33]
S. Voronin and P.-G. Martinsson, Efficient algorithms for CUR and interpolative matrix decompositions, Adv Comput Math, 43 (2017), pp. 495–516, https://doi.org/10.1007/s10444-016-9494-8.
[34]
G. Ballard and T. G. Kolda, Tensor Decompositions for Data Science, CUP, 2025, https://doi.org/10.1017/9781009471664.
[35]
W. Hackbusch, Tensor Spaces and Numerical Tensor Calculus, Springer, 2 ed., 2010, https://doi.org/10.1007/978-3-030-35554-8.
[36]
C. F. Caiafa and A. Cichocki, Generalizing the column–row matrix decomposition to multi-way arrays, Linear Algebra Appl, 433 (2010), pp. 557–573, https://doi.org/10.1016/j.laa.2010.03.020.
[37]
H. Cai, K. Hamm, L. Huang, and D. Needell, Mode-wise tensor decompositions: Multi-dimensional generalizations of CUR decompositions, J Mach Learn Res, 22 (2021), pp. 1–36, https://doi.org/10.48550/arXiv.2103.11037.
[38]
L. De Lathauwer, B. De Moor, and J. Vandewalle, A multilinear singular value decomposition, SIAM J Matrix Anal Appl, 21 (2000), pp. 1253–1278, https://doi.org/10.1137/s0895479896305696.
[39]
N. Vannieuwenhoven, R. Vandebril, and K. Meerbergen, A new truncation strategy for the higher-order singular value decomposition, SIAM J Sci Comput, 34 (2012), pp. A1027–A1052, https://doi.org/10.1137/110836067.
[40]
K. Hamm and L. Huang, Perturbations of CUR decompositions, SIAM J Matrix Anal Appl, 42 (2021), pp. 351–375, https://doi.org/10.1137/19M128394X.
[41]
Y. Nesterenko, About subspaces the most deviating from the coordinate ones, 2025, https://doi.org/10.48550/arXiv.2511.02387.
[42]
R. Sengupta and M. Pautov, On the submatrices with the best-bounded inverses, 2026, https://doi.org/10.48550/arXiv.2604.05944.
[43]
F. De Hoog and R. Mattheij, Subset selection for matrices, Linear Algebra Appl, 422 (2007), pp. 349–359, https://doi.org/10.1016/j.laa.2006.08.034.
[44]
A. Mikhalev and I. V. Oseledets, Rectangular maximum-volume submatrices and their applications, Linear Algebra Appl, 538 (2018), pp. 187–211, https://doi.org/10.1016/j.laa.2017.10.014.
[45]
A. Osinsky and N. L. Zamarashkin, Pseudo-skeleton approximations with better accuracy estimates, Linear Algebra Appl, 537 (2018), pp. 221–249, https://doi.org/10.1016/j.laa.2017.09.032.
[46]
R. A. Horn and C. R. Johnson, Topics in Matrix Analysis, CUP, 1994. isbn: https://www.cambridge.org/us/universitypress/subjects/mathematics/algebra/topics-matrix-analysis.
[47]
S. Dolgov and D. Savostyanov, Parallel cross interpolation for high-precision calculation of high-dimensional integrals, Comput Phys Commun, 246 (2020), 106869, https://doi.org/10.1016/j.cpc.2019.106869.
[48]
T. Shi, D. Hayes, and J.-M. Qiu, Distributed memory parallel adaptive tensor-train cross approximation, SIAM J Sci Comput, (2026), pp. S486–S512, https://doi.org/10.1137/24M1677587.
[49]
M. E. Rognes and G. N. Wells, Fenics course. Lecture 3: Static nonlinear PDEs. https://pub.fenicsproject.org/course/lectures/2017-nordic-phdcourse, 2017. Accessed 16/11/2025.
[50]
B. W. Larsen, T. G. Kolda, A. R. Zhang, and A. H. Williams, Tensor decomposition meets RKHS: Efficient algorithms for smooth and misaligned data, 2024, https://doi.org/10.48550/arXiv.2408.05677.
[51]
R. Han, P. Shi, and A. R. Zhang, Guaranteed functional tensor singular value decomposition, J Am Stat Assoc, 119 (2024), pp. 995–1007, https://doi.org/10.1080/01621459.2022.2153689.
[52]
M. Sørensen, S. Hendrikx, and L. De Lathauwer, Multilinear singular value decomposition–based completion with fibers observed in a single mode, SIAM J Matrix Anal Appl, 46 (2025), pp. 1061–1090, https://doi.org/10.1137/23M1622830.
[53]
A. V. Mamonov and M. A. Olshanskii, Interpolatory tensorial reduced order models for parametric dynamical systems, Comput Methods Appl Mech Eng, 397 (2022), 115122, https://doi.org/10.1016/j.cma.2022.115122.
[54]
A. V. Mamonov and M. A. Olshanskii, Model order reduction of parametric dynamical systems by slice sampling tensor completion, 2024, https://doi.org/10.48550/arXiv.2411.07151.
[55]
A. V. Mamonov and M. A. Olshanskii, Tensorial parametric model order reduction of nonlinear dynamical systems, SIAM J Sci Comput, 46 (2024), pp. A1850–A1878, https://doi.org/10.1137/23M1553789.
[56]
S. Budzinskiy, V. Kazeev, and M. Olshanskii, Low-rank cross approximation of function-valued tensors for reduced-order modeling of parametric PDEs, 2025, https://doi.org/10.48550/arXiv.2511.22650.

  1. Faculty of Mathematics, University of Vienna, Kolingasse 14-16, 1090 Vienna, Austria ().↩︎

  2. This research was funded by the Disruptive Innovation – Early Career Seed Money funding program of the Austrian Academy of Sciences (ÖAW) and the Austrian Science Fund (FWF).↩︎

  3. To simplify submultiplicative norm bounds, we assume throughout the text our symmetric gauge functions to be normalized as \(\phi(1, 0, \ldots, 0) = 1\). For instance, such normalization is used in [26].↩︎

  4. Owing to the normalization of symmetric gauge functions, as adopted throughout the text.↩︎

  5. https://fenicsproject.org/↩︎