Analysis of tensor approximation schemes for continuous functions

TL;DR

This paper investigates tensor approximation methods for continuous functions within Sobolev spaces, demonstrating that certain tensor formats can achieve dimension-robust approximation costs under specific conditions.

Contribution

It provides a detailed analysis of tensor approximation schemes for continuous functions, highlighting conditions for dimension-robustness in Tucker and tensor train formats.

Findings

Cost of tensor approximations is dimension-robust with proper weights.

Analysis applies to functions in isotropic Sobolev spaces.

Both Tucker and tensor train formats are effective under these conditions.

Abstract

In this article, we analyze tensor approximation schemes for continuous functions. We assume that the function to be approximated lies in an isotropic Sobolev space and discuss the cost when approximating this function in the continuous analogue of the Tucker tensor format or of the tensor train format. We especially show that the cost of both approximations are dimension-robust when the Sobolev space under consideration provides appropriate weights.

Equations

$$(A\otimes \cdots \otimes A){\mathcal{M}}_u={\mathcal{M}}_f\text{in}\Omega \times \cdots \times \Omega ,$$

$$f(\mathbf{x},\mathbf{y})\approx f_r(\mathbf{x},\mathbf{y}):=\sum _{\alpha =1}^r\sqrt{\lambda (\alpha )}\phi (\mathbf{x},\alpha )\psi (\alpha ,\mathbf{y}).$$

$${\mathcal{S}}_f:L^2({\Omega }_1)\to L^2({\Omega }_2),u\mapsto ({\mathcal{S}}_{f}u)(\mathbf{y}):={\int }_{{\Omega }_1}f(\mathbf{x},\mathbf{y})u(\mathbf{x})d\mathbf{x}.$$

$$\sqrt{\lambda(\alpha)}\psi(\alpha,\boldsymbol{y})=\big{(}\mathcal{S}_{f}\varphi(\alpha)\big{)}(\boldsymbol{y})\quad\text{and}\quad\sqrt{\lambda(\alpha)}\varphi(\boldsymbol{x},\alpha)=\big{(}\mathcal{S}_{f}^{\star}\psi(\alpha)\big{)}(\boldsymbol{x}),$$

$${\mathcal{S}}_f^{\star }:L^2({\Omega }_2)\to L^2({\Omega }_1),v\mapsto ({\mathcal{S}}_f^{\star }v)(\mathbf{x}):={\int }_{{\Omega }_2}f(\mathbf{x},\mathbf{y})v(\mathbf{y})d\mathbf{y}.$$

$${\mathcal{K}}_f\phi (\alpha )=\lambda (\alpha )\phi (\alpha )$$

$${\mathcal{K}}_f={\mathcal{S}}_f^{\star }{\mathcal{S}}_f:L^2({\Omega }_1)\to L^2({\Omega }_1),u\mapsto ({\mathcal{K}}_{f}u)(\mathbf{x}):={\int }_{{\Omega }_1}{k}_f(\mathbf{x},{\mathbf{x}}^{\prime })u({\mathbf{x}}^{\prime })d{\mathbf{x}}^{\prime }.$$

$$k_f(\mathbf{x},{\mathbf{x}}^{\prime })={\int }_{{\Omega }_2}f(\mathbf{x},\mathbf{y})f({\mathbf{x}}^{\prime },\mathbf{y})d\mathbf{y}\in L^2({\Omega }_1\times {\Omega }_1).$$

$$\lambda (1)\ge \lambda (2)\ge \cdots \ge \lambda (m)\stackrel{m\to \infty }{⟶}0$$

$${\tilde{\mathcal{K}}}_f\psi (\alpha )=\tilde{\lambda }(\alpha )\psi (\alpha )$$

$${\tilde{\mathcal{K}}}_f={\mathcal{S}}_f{\mathcal{S}}_f^{\star }:L^2({\Omega }_2)\to L^2({\Omega }_2),u\mapsto ({\tilde{\mathcal{K}}}_{f}u)(\mathbf{y}):={\int }_{{\Omega }_2}{\tilde{k}}_f(\mathbf{y},{\mathbf{y}}^{\prime })u({\mathbf{y}}^{\prime })d{\mathbf{y}}^{\prime }$$

$${\tilde{k}}_f(\mathbf{y},{\mathbf{y}}^{\prime })={\int }_{{\Omega }_1}f(\mathbf{x},\mathbf{y})f(\mathbf{x},{\mathbf{y}}^{\prime })d\mathbf{x}\in L^2({\Omega }_2\times {\Omega }_2).$$

$$H_{\mathrm{mix}}^{k,\ell }({\Omega }_1\times {\Omega }_2):=H^k({\Omega }_1)\otimes H^{\ell }({\Omega }_2),$$

$$\|f\|_{H_{\operatorname{mix}}^{k,\ell}(\Omega_{1}\times\Omega_{2})}:=\sqrt{\sum_{|\boldsymbol{\alpha}|\leq k}\sum_{|\boldsymbol{\beta}|\leq\ell}\bigg{\|}\frac{\partial^{|\boldsymbol{\alpha}|}\partial^{|\boldsymbol{\beta}|}}{\partial{\bf x}^{\boldsymbol{\alpha}}\partial{\bf y}^{\boldsymbol{\beta}}}f\bigg{\|}_{L^{2}(\Omega_{1}\times\Omega_{2})}^{2}}.$$

$$H^k({\Omega }_1\times {\Omega }_2)\subset H_{\mathrm{mix}}^{k,0}({\Omega }_1\times {\Omega }_2),H^k({\Omega }_1\times {\Omega }_2)\subset H_{\mathrm{mix}}^{0,k}({\Omega }_1\times {\Omega }_2).$$

$${\mathcal{S}}_f:L^2({\Omega }_1)\to H^k({\Omega }_2),{\mathcal{S}}_f^{\star }:L^2({\Omega }_2)\to H^k({\Omega }_1)$$

$$\big{\|}\mathcal{S}_{f}\big{\|}_{L^{2}(\Omega_{1})\to H^{k}(\Omega_{2})}\leq\|f\|_{H_{\operatorname{mix}}^{0,k}(\Omega_{1}\times\Omega_{2})},\quad\big{\|}\mathcal{S}_{f}^{\star}\big{\|}_{L^{2}(\Omega_{2})\to H^{k}(\Omega_{1})}\leq\|f\|_{H_{\operatorname{mix}}^{k,0}(\Omega_{1}\times\Omega_{2})}.$$

$$\displaystyle\big{\|}\mathcal{S}_{f}u\big{\|}_{H^{k}(\Omega_{2})}$$

$$=\underset{\Vert v{\Vert }_{H^{-k}({\Omega }_2)}=1}{\sup }(f,u\otimes v{)}_{L^2({\Omega }_1\times {\Omega }_2)}$$

$$\le \underset{\Vert v{\Vert }_{H^{-k}({\Omega }_2)}=1}{\sup }\Vert f{\Vert }_{H_{\mathrm{mix}}^{0,k}({\Omega }_1\times {\Omega }_2)}\Vert u\otimes v{\Vert }_{H_{\mathrm{mix}}^{0,-k}({\Omega }_1\times {\Omega }_2)}$$

$$\le \Vert f{\Vert }_{H_{\mathrm{mix}}^{0,k}({\Omega }_1\times {\Omega }_2)}\Vert u{\Vert }_{L^2({\Omega }_1)}.$$

$$\Vert \phi (\alpha ){\Vert }_{H^k({\Omega }_1)}\le \frac{1}{\sqrt{\lambda (\alpha )}}\Vert f{\Vert }_{H_{\mathrm{mix}}^{k,0}({\Omega }_1\times {\Omega }_2)},\Vert \psi (\alpha ){\Vert }_{H^k({\Omega }_2)}\le \frac{1}{\sqrt{\lambda (\alpha )}}\Vert f{\Vert }_{H_{\mathrm{mix}}^{0,k}({\Omega }_1\times {\Omega }_2)}.$$

$$\|\varphi(\alpha)\|_{H^{k}(\Omega_{1})}=\frac{1}{\sqrt{\lambda(\alpha)}}\big{\|}\mathcal{S}_{f}^{\star}\psi(\alpha)\big{\|}_{H^{k}(\Omega_{1})}\leq\frac{1}{\sqrt{\lambda(\alpha)}}\|f\|_{H_{\operatorname{mix}}^{k,0}(\Omega_{1}\times\Omega_{2})}\|\psi(\alpha)\|_{L^{2}(\Omega_{2})}.$$

$$\sum _{\alpha =1}^r\lambda (\alpha )\Vert \phi (\alpha ){\Vert }_{H^k({\Omega }_1)}^2\le r\Vert f{\Vert }_{H_{\mathrm{mix}}^{k,0}({\Omega }_1\times {\Omega }_2)}^2$$

$$\sum _{\alpha =1}^r\lambda (\alpha )\Vert \psi (\alpha ){\Vert }_{H^k({\Omega }_2)}^2\le r\Vert f{\Vert }_{H_{\mathrm{mix}}^{0,k}({\Omega }_1\times {\Omega }_2)}^2.$$

$$\Vert (I-P_r)w{\Vert }_{L^2({\Omega }_1)}\le c_k{r}^{-k/n_1}|w{|}_{H^k({\Omega }_1)}$$

$$\|f-(P_{r}\otimes I)f\|_{L^{2}(\Omega_{1}\times\Omega_{2})}^{2}=\operatorname{trace}\mathcal{K}_{f}-\operatorname{trace}\mathcal{K}_{f}^{\,r}=\sum_{\alpha=1}^{r}\big{(}\lambda(\alpha)-\lambda_{r}(\alpha)\big{)}+\sum_{\alpha=r+1}^{\infty}\lambda(\alpha).$$

$$\displaystyle\Bigg{\|}f-\sum_{\alpha=1}^{r}\sqrt{\lambda(\alpha)}\big{(}\varphi(\alpha)\otimes\psi(\alpha)\big{)}\Bigg{\|}_{L^{2}(\Omega_{1}\times\Omega_{2})}\!\!\!$$

$$f_{r}^{\operatorname{SVD}}=\sum_{\alpha=1}^{r}\sqrt{\lambda(\alpha)}\big{(}\varphi(\alpha)\otimes\psi(\alpha)\big{)}.$$

$$\Vert f-f_r^{\mathrm{SVD}}{\Vert }_{L^2({\Omega }_1\times {\Omega }_2)}=\sqrt{\sum _{\alpha =r+1}^{\infty }\lambda (\alpha )}\le c_k{r}^{-\frac{k}{\min \{n_1,n_2\}}}|f{|}_{H^k({\Omega }_1\times {\Omega }_2)}.$$

$$\lambda (\alpha )\lesssim {\alpha }^{-\frac{2k}{\min \{n_1,n_2\}}-1}\text{as α→∞}.$$

$$\sum _{\alpha =1}^{\infty }\lambda (\alpha )\Vert \phi (\alpha ){\Vert }_{H^k({\Omega }_1)}^2=\Vert f{\Vert }_{H_{\mathrm{mix}}^{k,0}({\Omega }_1\times {\Omega }_2)}^2$$

Full text

1. Introduction

The efficient approximate representation of multivariate functions is an important task in numerical analysis and scientific computing. In this article, we hence consider the approximation of functions which live on the product of $m$ bounded domains $\Omega_{1}\times\dots\times\Omega_{m}$, each of which satisfies $\Omega_{j}\subset\mathbb{R}^{n_{j}}$. Besides a sparse grid approximation of the function under consideration, being discussed in, e.g., [8, 17, 18, 50], one can also apply a low-rank approximation by means of a tensor approximation scheme, see, e.g., [15, 21, 22, 34, 35] and the references therein.

The low-rank approximation in the situation of the product of $m=2$ domains is well understood. It is related to the singular value decomposition and has been studied for arbitrary product domains in, e.g., [19, 20], see also [46, 47, 48] for the periodic case. However, the situation is not that clear for the product of $m>2$ domains, where one ends up with tensor decompositions. Such tensor decompositions are generalizations of the well known singular value decomposition and the corresponding low-rank matrix approximation methods of two dimensions to the higher-dimensional setting. There, besides the curse of dimension, we encounter – due to the non-existence of an Eckart-Young-Mirsky theorem – that the concepts of singular value decomposition and low-rank approximation can be generalized to higher dimensions in more than one way. Consequently, there exist many generalizations of the singular value decomposition of a function and of low-rank approximations to tensors. To this end, various schemes have been developed over the years in different areas of the sciences and have successfully been applied to various high-dimensional problems ranging from quantum mechanics and physics via biology and econometrics, computer graphics and signal processing to numerical analysis. Recently, tensor methods have even been recognized as special deep neural networks in machine learning and big data analysis [11, 28]. As tensor approximation schemes, we have, for example, matrix product states, DMRG, MERA, PEPS, CP, CANDECOMP, PARAFAC, Tucker, tensor train, tree tensor networks and hierarchical Tucker, to name a few. A mathematical introduction into tensor methods is given in the seminal book [21], while a survey on existing methods and their literature can be found in [16]. Also various software packages have been developed for an algebra of operators dealing with tensors.

Tensor methods are usually analyzed as low-rank approximations to a full discrete tensor of data with respect to the $\ell_{2}$-norm or Frobenius-norm. In this respective, they can be seen as compression methods which may avoid the curse of dimensionality, which is inherent in the full tensor representation. The main tool for studying tensor compression schemes is the so-called tensor-rank, compare [12, 13, 21]. Thus, instead of $\mathcal{O}(N^{n})$ storage, as less as $\mathcal{O}(nNr^{3})$ or even only $\mathcal{O}(nNr^{2})$ storage is needed, where $N$ denotes the number of data points in one coordinate direction, $n$ denotes the dimension of the tensor under consideration and $r$ denotes the respective tensor rank of the data. The cost complexity of the various algorithms working with sparse tensor representations is correspondingly reduced and working in a sparse tensor format allows to alleviate or to completely break the curse of dimension for suitable tensor data classes, i.e., for sufficiently small $r$.

However, the question where the tensor data stem from and the issue of the accuracy of the full tensor approximation, i.e., the discretization error of the full tensor itself and its relation to the error of a subsequent low-rank tensor approximation, is usually not adequately addressed.111We are only aware of [2, 3, 5, 39], where this question has been considered so far. Instead, only the approximation property of a low-rank tensor scheme with respect to the full tensor data is considered. But the former question is important since it clearly makes no sense to derive a tensor approximation with an error that is substantially smaller than the error which is already inherent in the full tensor data due to some discretization process for a continuous high-dimensional function which stems from some certain function class.

The approximation rates to continuous functions can be determined by a recursive use of the singular value decomposition, which is successively applied to convert the function into a specific continuous tensor format. We studied the singular value decomposition for arbitrary domains in [19, 20] and we now can apply these results to discuss approximation rates of continuous tensor formats. In the present article, given a function $f\in H^{k}(\Omega_{1}\times\dots\times\Omega_{m})$, we study the continuous analogues of the Tucker tensor decomposition and of the tensor train decomposition. We give bounds on the ranks required to ensure that the tensor decomposition admits a prescribed target accuracy. Especially, our analysis takes into account the influence of errors induced by truncating infinite expansions to finite ones. We therefore study an algorithm that computes the desired tensor expansion which is in contrast to the question of the smallest tensor-rank. We finally show that (isotropic) Sobolev spaces with dimension weights help to beat the curse of dimension when the number $m$ of product domains tends to infinity.

Besides the simple situation of $\Omega_{1}=\dots=\Omega_{m}=[0,1]$, which is usually considered in case of tensor decompositions, there are many more applications of our general setting. For example, non-Newtonian flow can be modeled by a coupled system which consists of the Navier Stokes equation for the flow in a three-dimensional geometry described by $\Omega_{1}$ and of the Fokker-Planck equation in a $3(\ell-1)$-dimensional configuration space $\Omega_{2}\times\dots\times\Omega_{\ell}$, consisting of $\ell-1$ spheres. Here, $\ell$ denotes the number of atoms in a chain-like molecule which constitutes the non-Newtonian behavior of the flow, for details see [4, 29, 31, 37]. Another example is homogenization. After unfolding [10], a two-scale homogenization problem gives raise to the product of the macroscopic physical domain and the periodic microscopic domain of the cell problem, see [32]. For multiple scales, several periodic microscopic domains appear which reflect the different separable scales, see e.g. [27]. Also the $m$-th moment of linear elliptic boundary value problems with random source terms, i.e. $Au(\omega)=f(\omega)$ in $\Omega$, are known to satisfy a deterministic partial differential equation on the $m$-fold product domain $\Omega\times\dots\times\Omega$. There, the solution’s $m$-th moment $\mathcal{M}_{u}$ is given by the equation

[TABLE]

see [40, 41]. This approach extends to boundary value problems with random diffusion and to random domains as well [9, 25]. Moreover, we find the product of several domains in quantum mechanics for e.g. the Schrödinger equation or the Langevin equation, where each domain is three-dimensional and corresponds to a single particle. Finally, we encounter it in uncertainty quantification, where one has the product of the physical domain $\Omega_{1}$ and of in general infinitely many intervals $\Omega_{2}=\Omega_{3}=\Omega_{4}=\dots$ for the random input parameter, which reflects its series expansion by the Karhunen-Lòeve decomposition or the Lévy-Ciesielski decomposition.

The remainder of this article is organized as follows: In Section 2, we give a short introduction to our results on the singular value decomposition, which are needed to derive the estimates for the continuous tensor decompositions. Then, in Section 3, we study the continuous Tucker tensor format, computed by means of the higher-oder singular value decomposition. Next, we study the continuous tensor train decomposition in Section 4, computed by means of a repeated use of a vector-valued singular value decomposition. Finally, Section 5 concludes with some final remarks.

Throughout this article, to avoid the repeated use of generic but unspecified constants, we denote by $C\lesssim D$ that $C$ is bounded by a multiple of $D$ independently of parameters which $C$ and $D$ may depend on. Obviously, $C\gtrsim D$ is defined as $D\lesssim C$, and $C\sim D$ as $C\lesssim D$ and $C\gtrsim D$. Moreover, given a Lipschitz-smooth domain $\Omega\subset\mathbb{R}^{n}$, $L^{2}(\Omega)$ means the space of square integrable functions on $\Omega$. For real numbers $k\geq 0$, the associated Sobolev space is denoted by $H^{k}(\Omega)$, where its norm $\|\cdot\|_{H^{k}(\Omega)}$ is defined in the standard way, compare [33, 45]. As usual, we have $H^{0}(\Omega)=L^{2}(\Omega)$. The seminorm in $H^{k}(\Omega)$ is denoted by $|\cdot|_{H^{k}(\Omega_{1})}$. Although not explicitly written, our subsequent analysis covers also the situation of $\Omega$ being not a domain but a (smooth) manifold.

2. Singular value decomposition

2.1. Definition and calculation

Let $\Omega_{1}\subset\mathbb{R}^{n_{1}}$ and $\Omega_{2}$ be Lipschitz-smooth domains. To represent functions $f\in L^{2}(\Omega_{1}\times\Omega_{2})$ on the tensor product domain $\Omega_{1}\times\Omega_{2}$ in an efficient way, we will consider low-rank approximations which separate the variables $\boldsymbol{x}\in\Omega_{1}$ and $\boldsymbol{y}\in\Omega_{2}$ in accordance with

[TABLE]

It is well known (see e.g. [30, 38, 42]) that the best possible representation (2.1) in the $L^{2}$-sense is given by the singular value decomposition, also called Karhunen-Lòeve expansion.222We refer the reader to [44] for a comprehensive historical overview on the singular value decomposition. Then, the coefficients $\sqrt{\lambda(\alpha)}\in\mathbb{R}$ are the singular values and the $\varphi(\alpha)\in L^{2}(\Omega_{1})$ and $\psi(\alpha)\in L^{2}(\Omega_{2})$ are the left and right ($L^{2}$-normalized) eigenfunctions of the integral operator

[TABLE]

This means that

[TABLE]

where

[TABLE]

is the adjoint of $\mathcal{S}_{f}$. Especially, the left and right eigenfunctions $\{\varphi(\alpha)\}_{\alpha=1}^{\infty}$ and $\{\psi(\alpha)\}_{\alpha=1}^{\infty}$ form orthonormal bases in $L^{2}(\Omega_{1})$ and $L^{2}(\Omega_{2})$, respectively.

In order to compute the singular value decomposition, we need to solve the eigenvalue problem

[TABLE]

for the integral operator

[TABLE]

Since $f\in L^{2}(\Omega_{1}\times\Omega_{2})$, the kernel

[TABLE]

is a symmetric Hilbert-Schmidt kernel. Hence, there exist countably many eigenvalues

[TABLE]

and the associated eigenfunctions $\{\varphi(\alpha)\}_{\alpha\in\mathbb{N}}$ constitute an orthonormal basis in $L^{2}(\Omega_{1})$.

Likewise, to obtain an orthonormal basis of $L^{2}(\Omega_{2})$, we can solve the eigenvalue problem

[TABLE]

for the integral operator

[TABLE]

with symmetric Hilbert-Schmidt kernel

[TABLE]

It holds $\lambda(\alpha)=\widetilde{\lambda}(\alpha)$ and the sequences $\{\varphi(\alpha)\}$ and $\{\psi(\alpha)\}$ are related by (2.2).

2.2. Regularity of the eigenfunctions

Now, we consider functions $f\in H^{k}(\Omega_{1}\times\Omega_{2})$. In the following, we collect results from [19, 20] concerning the singular value decomposition of such functions. We repeat the proof whenever needed for having explicit constants. To this end, we define the mixed Sobolev space $H_{\operatorname{mix}}^{k,\ell}(\Omega_{1}\times\Omega_{2})$ as a tensor product of Hilbert spaces

[TABLE]

which we equip with the usual cross norm

[TABLE]

Note that

[TABLE]

Lemma 2.1**.**

Assume that $f\in H^{k}(\Omega_{1}\times\Omega_{2})$ for some fixed $k\geq 0$. Then, the operators

[TABLE]

are continuous with

[TABLE]

Proof.

From $H^{k}(\Omega_{1}\times\Omega_{2})\subset H_{\operatorname{mix}}^{0,k}(\Omega_{1}\times\Omega_{2})$ it follows for $f\in H^{k}(\Omega_{1}\times\Omega_{2})$ that $f\in H_{\operatorname{mix}}^{0,k}(\Omega_{1}\times\Omega_{2})$. Therefore, the operator $\mathcal{S}_{f}:L^{2}(\Omega_{1})\to H^{k}(\Omega_{2})$ is continuous since

[TABLE]

Note that we have used here $H_{\operatorname{mix}}^{0,-k}(\Omega_{1}\times\Omega_{2})=L^{2}(\Omega_{1})\otimes H^{-k}(\Omega_{2})$. Proceeding likewise for $\mathcal{S}_{f}^{\star}:L^{2}(\Omega_{2})\to H^{k}(\Omega_{1})$ completes the proof. ∎

Lemma 2.2**.**

Assume that $f\in H^{k}(\Omega_{1}\times\Omega_{2})$ for some fixed $k\geq 0$. Then, it holds $\mathcal{S}_{f}\varphi(\alpha)\in H^{k}(\Omega_{2})$ and $\mathcal{S}_{f}^{\star}\psi(\alpha)\in H^{k}(\Omega_{1})$ for all $\alpha\in\mathbb{N}$ with

[TABLE]

Proof.

According to (2.2) and Lemma 2.1, we have

[TABLE]

This proves the first assertion. The second assertion follows by duality. ∎

As an immediate consequence of Lemma 2.2, we obtain

[TABLE]

and

[TABLE]

We will show later in Lemma 2.6 how to improve this estimate by sacrificing some regularity.

2.3. Truncation error

We next give estimates on the decay rate of the eigenvalues of the integral operator $\mathcal{K}_{f}=\mathcal{S}_{f}^{\star}\mathcal{S}_{f}$ with kernel (2.4). To this end, we exploit the smoothness in the function’s first variable and assume hence $f\in H_{\operatorname{mix}}^{k,0}(\Omega_{1}\times\Omega_{2})$. We introduce finite element spaces $U_{r}\subset L^{2}(\Omega_{1})$, which consist of $r$ discontinuous, piecewise polynomial functions of total degree $\lceil k\rceil$ on a quasi-uniform triangulation of $\Omega_{1}$ with mesh width $h_{r}\sim r^{-1/n_{1}}$. Then, given a function $w\in H^{k}(\Omega_{1})$, the $L^{2}$-orthogonal projection $P_{r}:L^{2}(\Omega_{1})\to U_{r}$ satisfies

[TABLE]

uniformly in $r$ due to the Bramble-Hilbert lemma, see e.g., [6, 7].

For the approximation of $f(\boldsymbol{x},\boldsymbol{y})$ in the first variable, i.e. \big{(}(P_{r}\otimes I)f\big{)}(\boldsymbol{x},\boldsymbol{y}), we obtain the following approximation result for the present choice of $U_{r}$, see [20] for the proof.

Lemma 2.3**.**

Assume that $f\in H^{k}(\Omega_{1}\times\Omega_{2})$ for some fixed $k\geq 0$. Let $\lambda(1)\geq\lambda(2)\geq\ldots\geq 0$ be the eigenvalues of the operator $\mathcal{K}_{f}=\mathcal{S}_{f}^{\star}\mathcal{S}_{f}$ and $\lambda_{r}(1)\geq\lambda_{r}(2)\geq\ldots\geq\lambda_{r}(r)\geq 0$ those of $\mathcal{K}_{f}^{\,r}:=P_{r}\mathcal{K}_{f}P_{r}$. Then, it holds

[TABLE]

By combining this lemma with the approximation estimate (2.6) and in view of $\lambda(\alpha)-\lambda_{r}(\alpha)\geq 0$ for all $\alpha\in\{1,2,\ldots,r\}$ according to the min-max theorem of Courant-Fischer, see [1] for example, we conclude that the truncation error of the singular value decomposition can be bounded by

[TABLE]

Since the eigenvalues of the integral operator $\mathcal{K}_{f}$ and its adjoint $\widetilde{\mathcal{K}}_{f}$ are the same, we can also exploit the smoothness of $f$ in the second coordinate by interchanging the roles of $\Omega_{1}$ and $\Omega_{2}$ in the above considerations. We thus obtain the following theorem:

Theorem 2.4**.**

Let $f\in H^{k}(\Omega_{1}\times\Omega_{2})$ for some fixed $k\geq 0$ and let

[TABLE]

Then, it holds

[TABLE]

Remark 2.5**.**

Theorem 2.4 implies that the eigenvalues $\{\lambda(\alpha)\}_{\alpha\in\mathbb{N}}$ in case of a function $f\in H^{k}(\Omega_{1}\times\Omega_{2})$ decay like

[TABLE]

Having the decay rate of the eigenvalues at hand, we are able to improve the result of Lemma 2.2 by sacrificing some regularity. Note that the proof of this result is based upon an argument from [43].

Lemma 2.6**.**

Assume that $f\in H^{k+\min\{n_{1},n_{2}\}}(\Omega_{1}\times\Omega_{2})$ for some fixed $k\geq 0$. Then, it holds

[TABLE]

and

[TABLE]

Proof.

Without loss of generality, we assume $n_{1}\leq n_{2}$. Then, since $f\in H^{k+n_{1}}(\Omega_{1}\times\Omega_{2})$, we conclude from (2.8) that

[TABLE]

where we used that $n_{1}=\min\{n_{1},n_{2}\}$. Moreover, by interpolating between $L^{2}(\Omega_{1})$ and $H^{k+n_{1}}(\Omega_{1})$, compare [33, 45] for example, we find

[TABLE]

that is

[TABLE]

As a consequence, we infer that

[TABLE]

with $\delta=\frac{n_{1}}{k+n_{1}}>0$. Therefore, it holds

[TABLE]

for any $\delta^{\prime}\in(0,\delta)$. Hence, the series

[TABLE]

converges for almost all $\boldsymbol{x}\in\Omega_{1}$, provided that $|\boldsymbol{\beta}|\leq k$. Likewise, the series

[TABLE]

converges for almost all $\boldsymbol{y}\in\Omega_{2}$. Thus, the series

[TABLE]

converges for almost all $\boldsymbol{x}\in\Omega_{1}$ and $\boldsymbol{y}\in\Omega_{2}$, provided that $|\boldsymbol{\beta}|\leq k$. Because of Egorov’s theorem, the pointwise absolute convergence almost everywhere implies uniform convergence. Hence, we can switch differentiation and summation to get

[TABLE]

Finally, we exploit the product structure of $L^{2}(\Omega_{1}\times\Omega_{2})$ and the orthonormality of $\{\psi(\alpha)\}_{\alpha\in\mathbb{N}}$ to derive the first assertion, i.e.,

[TABLE]

The second assertion follows in complete analogy. ∎

2.4. Vector-valued functions

In addition to the aforementioned results, we will also need the following result which concerns the approximation of vector-valued functions. Here and in the sequel, the vector-valued function $\boldsymbol{w}=[w(\alpha)]_{\alpha=1}^{m}$ is an element of $[L^{2}(\Omega)]^{m}$ and $[H^{k}(\Omega)]^{m}$ for some domain $\Omega\subset\mathbb{R}^{n}$, respectively, if the norms

[TABLE]

are finite. Likewise, the seminorm is defined in $[H^{k}(\Omega)]^{m}$.

Consider now a vector-valued function $\boldsymbol{w}\in[H^{k}(\Omega_{1})]^{m}$ of dimension $m$. Then, instead of (2.6), we find

[TABLE]

since $\boldsymbol{w}$ consists of $m$ components and we thus need $m$-times as many ansatz functions for our approximation argument. Hence, in case of a vector-valued function $\boldsymbol{f}\in[H_{\operatorname{mix}}^{k,0}(\Omega_{1}\times\Omega_{2})]^{m}\simeq[H^{k}(\Omega_{1})]^{m}\otimes L^{2}(\Omega_{2})$, we conclude by exploiting the smoothness in the first variable333Note that the kernel function of $\mathcal{S}_{f}^{\star}\mathcal{S}_{f}$ is matrix-valued while the kernel function of $\mathcal{S}_{f}\mathcal{S}_{f}^{\star}$ is scalar-valued. that the truncation error of the singular value decomposition can be estimated by

[TABLE]

Hence, the decay rate of the singular values is considerably reduced. Finally, we like to remark that Lemma 2.6 holds also in the vector case, i.e.,

[TABLE]

and

[TABLE]

provided that $\boldsymbol{f}$ has extra regularity in terms of $\boldsymbol{f}\in[H^{k+n_{1}}(\Omega_{1}\times\Omega_{2})]^{m}$. Here, analogously to above, $[H_{\operatorname{mix}}^{0,k}(\Omega_{1}\times\Omega_{2})]^{m}\simeq[L^{2}(\Omega_{1})]^{m}\otimes H^{k}(\Omega_{2})$.

After these preparations, we now introduce and analyze two types of continuous analogues of tensor formats, namely of the Tucker format [26, 49] and of the tensor train format [36, 34], and discuss their approximation properties for functions $f\in H^{k}(\Omega_{1}\times\dots\times\Omega_{m})$.

3. Tucker tensor format

3.1. Tucker decompostion

We shall consider from now on a product domain which consists of $m$ different domains $\Omega_{j}\subset\mathbb{R}^{n_{j}}$, $j=1,\ldots,m$. For given $f\in L^{2}(\Omega_{1}\times\dots\times\Omega_{m})$ and $j\in\{1,2,\ldots,m\}$, we apply the singular value decomposition to separate the variables $\boldsymbol{x}_{j}\in\Omega_{j}$ and $(\boldsymbol{x}_{1},\ldots,\boldsymbol{x}_{j-1},\boldsymbol{x}_{j+1},\ldots,\boldsymbol{x}_{m})\in\Omega_{1}\times\dots\times\Omega_{j-1}\times\Omega_{j+1}\times\dots\times\Omega_{m}$. We hence get

[TABLE]

where the left eigenfunctions $\{\varphi_{j}(\alpha_{j})\}_{\alpha_{j}\in\mathbb{N}}$ form an orthonormal basis in $L^{2}(\Omega_{j})$. Consequently, if we iterate over all $j\in\{1,2,\ldots,m\}$, this yields an orthonormal basis $\{\varphi_{1}(\alpha_{1})\otimes\cdots\otimes\varphi_{m}(\alpha_{m})\}_{\boldsymbol{\alpha}\in\mathbb{N}^{m}}$ of $L^{2}(\Omega_{1}\times\dots\times\Omega_{m})$, and we arrive at the representation

[TABLE]

Herein, the tensor \big{[}\omega(\boldsymbol{\alpha})\big{]}_{\boldsymbol{\alpha}\in\mathbb{N}^{m}} is the core tensor, where a single coefficient is given by

[TABLE]

3.2. Truncation error

If we intend to truncate the singular value decomposition (3.11) after $r_{j}$ terms such that the truncation error is bounded by $\varepsilon$, we have to choose

[TABLE]

according to Theorem 2.4. Doing so for all $j\in\{1,2,\ldots,m\}$, we obtain the approximation

[TABLE]

We have the following result on the Tucker decomposition:

Theorem 3.1**.**

Let $f\in H^{k}(\Omega_{1}\times\dots\times\Omega_{m})$ for some fixed $k>0$ and $0<\varepsilon<1$. If the ranks are chosen according to $r_{j}=\varepsilon^{-n_{j}/k}$ for all $j=1,\dots,m$. Then, the truncation error of the truncated Tucker decomposition is

[TABLE]

while the storage cost for the core tensor of $f_{r_{1},\ldots,r_{m}}^{\operatorname{TF}}$ are $\prod_{j=1}^{m}r_{j}=\varepsilon^{-(n_{1}+\dots+n_{m})/k}$.

Proof.

For the approximation of the core tensor, the sets of the univariate eigenfunctions $\{\varphi_{j}(\alpha_{j})\}_{\alpha_{j}=1}^{r_{j}}$ are used for all $j=1,\ldots,m$, cf. (3.12). Due to orthonormality, we find

[TABLE]

where we obtain $f_{\infty,\ldots,\infty}^{\operatorname{TF}}=f$ in case of $j=1$. Since

[TABLE]

for all $j\in\{1,2,\ldots,m\}$, we arrive with (3.13) and the summation over $j=1,\ldots,m$ at the desired error estimate. This completes the proof, since the estimate on the number of coefficients in the core tensor is obvious. ∎

3.3. Sobolev spaces with dimension weights

The cost of the core tensor of the Tucker decomposition exhibit the curse of dimension as the number $m$ of subdomains increases. This can be seen most simply for the example $n_{j}=n$. Then, the cost are $\varepsilon^{-nm/k}$, which expresses the curse of dimension as long as $k$ is not proportional to $m$. Nonetheless, in case of Sobolev spaces with dimension weights, the curse of dimension can be beaten.

For $f\in H^{k+n}(\Omega^{m})$, we shall discuss the situation $m\to\infty$ in more detail. To this end, we assume that all subdomains are identical to a single domain $\Omega\subset\mathbb{R}^{n}$ of dimension $n$ and note that the limit $m\to\infty$ only makes sense when weights are included in the underlying Sobolev spaces which ensure that higher dimensions become less important. For our proofs we choose as usual $m$ arbitrary but fixed and show the existence of $m$-independent constants in the convergence and cost estimates.

The Sobolev spaces $H_{\boldsymbol{\gamma}}^{k}(\Omega^{m})$ with dimension weights $\boldsymbol{\gamma}\in\mathbb{R}^{m}$ we consider are given by all functions $f\in H^{k}(\Omega^{m})$ such that

[TABLE]

The definition in (3.14) means that, given a function $f$ with norm $\|f\|_{H^{k}(\Omega^{m})}<\infty$, the partial derivatives with respect to $\boldsymbol{x}_{j}$ become less important as the dimension $j$ increases. Such functions appear for example in uncertainty quantification. Let be given a Karhunen-Loève expansion

[TABLE]

and insert it into a function $b:\mathbb{R}\to\mathbb{R}$ of finite smoothness $W^{k,\infty}(\mathbb{R})$. Then, the function b\big{(}u({\bf x},{\bf y})\big{)} satisfies (3.14) with respect to the ${\bf y}$-variable, where $\gamma_{j}=\sigma_{j}$. Hence, the solution of a given partial differential equation would satisfy a decay estimate similar to (3.14) whenever the stochastic field enters the partial differential equation through a non-smooth coefficient function $b$, compare [14, 23, 24] for example.

It turns out that algebraically decaying weights (3.15) are sufficient to beat the curse of dimension in case of the Tucker tensor decomposition.444In Theorem 3.2, no truncation of the dimension is applied, as it would be required in practice if the number $m$ of domains tends to infinity. Note that the dimension truncation is indeed here the same as for the tensor train decomposition later on, see also Theorem 4.3.

Theorem 3.2**.**

Given $\delta>0$, let $f\in H_{\boldsymbol{\gamma}}^{k}(\Omega^{m})$ for some fixed $k>0$ with weights (3.14) that decay like

[TABLE]

Then, for all $0<\varepsilon<1$, the error of the continuous Tucker decomposition with ranks

[TABLE]

is of order $\varepsilon$ while the storage cost for the core tensor of $f_{r_{1},\ldots,r_{m}}^{\operatorname{TF}}$ are bounded by $\varepsilon^{-n/k}$ independent of the dimension $m$.

Proof.

In view of Theorem 2.4 and (3.14), we deduce by choosing the ranks as in (3.13) that

[TABLE]

Therefore, we reach the desired over-all truncation error

[TABLE]

When the weights $\gamma_{j}$ decay as in (3.16), then the cost of the core tensor are

[TABLE]

with $\theta=(\delta^{\prime}-\delta)n/k>1$. Hence, the cost of the core tensor stay bounded independently of $m$ since

[TABLE]

∎

4. Tensor train format

4.1. Tensor train decomposition

For the discussion of the continuous tensor train decomposition, we should assume that the domains $\Omega_{j}\subset\mathbb{R}^{n_{j}}$, $j=1,\ldots,m$, are arranged in such a way that it holds $n_{1}\leq\dots\leq n_{m}$.555The considerations in this section are based upon [5]. Nonetheless, the results derived there are not correct. The authors did not consider the impact of the vector-valued singular value decomposition in a proper way, which indeed does result in the curse of dimension.

Now, consider $f\in H^{k}(\Omega_{1}\times\dots\times\Omega_{m})$ and separate the variables $\boldsymbol{x}_{1}\in\Omega_{1}$ and $(\boldsymbol{x}_{2},\ldots,\boldsymbol{x}_{m})\in\Omega_{2}\times\dots\times\Omega_{m}$ by the singular value decomposition

[TABLE]

Since

[TABLE]

we can separate $(\alpha_{1},\boldsymbol{x}_{2})\in\mathbb{N}\times\Omega_{2}$ from $(\boldsymbol{x}_{3},\ldots,\boldsymbol{x}_{m})\in\Omega_{3}\times\dots\times\Omega_{m}$ by means of a second singular value decomposition and arrive at

[TABLE]

By repeating the last step and successively separating $(\alpha_{j-1},\boldsymbol{x}_{j})\in\mathbb{N}\times\Omega_{j}$ from $(\boldsymbol{x}_{j+1},\ldots,\boldsymbol{x}_{m})\in\Omega_{j+1}\times\dots\times\Omega_{m}$ for $j=3,\ldots,m-1$ we finally arrive at the representation

[TABLE]

where

[TABLE]

In contrast to the Tucker format, we do not obtain a huge core tensor since each of the $m-1$ singular value decompositions of the tensor train decomposition removes the actual first spatial domain from the approximant. We just obtain a product of matrix-valued functions (except for the first and last factor which are vector-valued functions), each of which is related with a specific domain $\Omega_{j}$. This especially results in only $m-1$ sums in contrast to the $m$ sums for the Tucker format.

4.2. Truncation error

In practice, we truncate the singular value decomposition in step $j$ after $r_{j}$ terms, thus arriving at the representation

[TABLE]

One readily infers by using again Pythogoras’ theorem that the truncation error is bounded by

[TABLE]

see also [35]. Note that, for $j\geq 2$, the singular values $\{\lambda_{j}(\alpha)\}_{\alpha\in\mathbb{N}}$ in this estimate do not coincide with the singular values from the original continuous tensor train decomposition due to the truncation.

We next shall give bounds on the truncation error. In the $j$-th step of the algorithm, $j=2,3,\ldots,m-1$, one needs to approximate the vector-valued function

[TABLE]

by a vector-valued singular value decomposition. This means that we consider the singular value decomposition (2.9) for a vector-valued function in case of the domains $\Omega_{j}$ and $\Omega_{j+1}\times\dots\times\Omega_{m}$.

For $f\in H^{k+n_{m-1}}(\Omega_{1}\times\dots\times\Omega_{m})$, it holds $\boldsymbol{g}_{2}\in[H^{k+n_{m-1}}(\Omega_{2}\times\dots\times\Omega_{m})]^{r_{1}}$ and

[TABLE]

according to Lemma 2.6, precisely in its vectorized version (2.10). It follows $\boldsymbol{g}_{3}\in[H^{k+n_{m-1}}(\Omega_{3}\times\dots\times\Omega_{m})]^{r_{2}}$ and, again by (2.10),

[TABLE]

We hence conclude recursively $\boldsymbol{g}_{j}\in[H^{k+n_{m-1}}(\Omega_{j}\times\dots\times\Omega_{m})]^{r_{j-1}}$ and

[TABLE]

Estimate (4.19) shows that the $H^{k}$-seminorm of the vector-valued functions $\boldsymbol{g}_{j}$ stays bounded by $|f|_{H^{k}(\Omega_{1}\times\dots\times\Omega_{m})}$. But according to (2.9), we have in the $j$-th step only the truncation error estimate

[TABLE]

Hence, in view of (4.19), to achieve the target accuracy $\varepsilon$ per truncation, the truncation ranks need to be increased in accordance with

[TABLE]

We summarize our findings in the following theorem, which holds in this form also if the subdomains $\Omega_{j}\subset\mathbb{R}_{n_{j}}$ are not ordered in such a way that $n_{1}\leq\dots\leq n_{m}$.

Theorem 4.1**.**

Let $f\in H^{k+\max\{n_{1},\ldots,n_{m}\}}(\Omega_{1}\times\dots\times\Omega_{m})$ for some fixed $k>0$ and $0<\varepsilon<1$. Then, the over-all truncation error of the tensor train decomposition with truncation ranks (4.20) is

[TABLE]

The storage cost for $f_{r_{1},\ldots,r_{m-1}}^{\operatorname{TT}}$ are given by

[TABLE]

and hence are bounded by $\mathcal{O}(\varepsilon^{-(2m-1)\max\{n_{1},\ldots,n_{m-1}\}/k})$.

Remark 4.2**.**

If $n:=n_{1}=\dots=n_{m}$, then the cost of the tensor train decomposition are $\mathcal{O}(\varepsilon^{-(2m-1)n/k})$. Thus, the cost are quadratic compared to the cost of the Tucker decomposition. However, in practice, one performs $m/2$ forward steps and $m/2$ backward steps. This means one computes $m/2$ steps as described above to successively separate $\boldsymbol{x}_{1},\boldsymbol{x}_{2},\ldots,\boldsymbol{x}_{m/2}$ from the other variables. Then, one performs the algorithm in the opposite direction, i.e., one successively separates $\boldsymbol{x}_{m},\boldsymbol{x}_{m-1},\ldots,\boldsymbol{x}_{m/2+1}$ from the other variables. This way, the over-all cost are reduced to the order $\mathcal{O}(\varepsilon^{-mn/k})$.666If the spatial dimensions $n_{j}$, $j=1,\ldots,m$, of the subdomains are different, one can balance the number of forward and backward steps in a better way to reduce the cost further.

4.3. Sobolev spaces with dimension weights

Like for the Tucker decomposition, the cost of the tensor train decomposition suffer from the curse of dimension as the number $m$ of subdomains increases. We therefore discuss again appropriately Sobelev spaces with dimension weights, where we assume for reasons of simplicity that all subdomains are identical to a single domain $\Omega\subset\mathbb{R}^{n}$ of dimension $n$.

Theorem 4.3**.**

Given $\delta>0$, let $f\in H_{\boldsymbol{\gamma}}^{k+n}(\Omega^{m})$ for some fixed $k>0$ with weights (3.14) that decay like (3.15). For $0<\varepsilon<1$, choose the ranks successively in accordance with

[TABLE]

if $j\leq M$ and $r_{j}=0$ if $j>M$. Here, $M$ is given by

[TABLE]

Then, the error of the continuous tensor train decomposition is of order $\varepsilon$ while the storage cost of $f_{r_{1},\ldots,r_{m}}^{\operatorname{TT}}$ stay bounded by $M\exp(\varepsilon^{-n/k})^{2}$ independent of the dimension $m$.

Proof.

The combination of Theorem 2.4, (3.14) and (4.22) implies

[TABLE]

and

[TABLE]

Hence, as in the proof of Theorem 3.2, the approximation error of the continuous tensor train decomposition is bounded by a multiple of $\varepsilon$ independent of $m$.

Next, we observe for all $j\leq M$ that

[TABLE]

This recursively yields

[TABLE]

Hence, by using that $\theta=(\delta^{\prime}-\delta)n/k>1$, we obtain

[TABLE]

Therefore, the cost (4.21) are

[TABLE]

and, hence, are bounded independently of $m$ in view of (4.23). ∎

5. Discussion and conclusion

In the present article, we considered the continuous versions of the Tucker tensor format and of the tensor train format for the approximation of functions which live on an $m$-fold product of arbitrary subdomains. By considering (isotropic) Sobolev smoothness, we derived estimates on the ranks to be chosen in order to realize a prescribed target accuracy. These estimates exhibit the curse of dimension.

Both tensor formats have in common that always only the variable with respect to a single domain is separated from the other variables by means of the singular value decomposition. This enables cheaper storage schemes, while the influence of the over-all dimension of the product domain is reduced to a minimum.

We also examined the situation of Sobolev spaces with dimension weights. Having sufficiently fast decaying weights helps to beat the curse of dimension as the number of subdomains tends to infinity. It turned out that algebraically decaying weights are appropriate for both, the Tucker tensor format and the tensor train format.

We finally remark that we considered here only the ranks of the tensor decomposition in the continuous case, i.e., for functions and not for tensors of discrete data. Of course, an additional projection step onto suitable finite dimensional trial spaces on the individual domains would be necessary to arrive at a fully discrete approximation scheme that can really be used in computer simulations. This would impose a further error of discretization type which needs to be balanced with the truncation error of the particular continuous tensor format.

Acknowledgement

Michael Griebel was partially supported by the Sonderforschungsbereich 1060 The Mathematics of Emergent Effects funded by the Deutsche Forschungsgemeinschaft. Both authors like to thank Reinhold Schneider (Technische Universität Berlin) very much for fruitful discussions about tensor approximation.