Bivariate Hermite interpolation by a limiting case of the cross approximation algorithm

Abstract

Computing low-rank approximations of a given function is a key step for implementing efficiently numerous algorithms in various fields, including the discretisation of non-local integral operators and Isogeometric Analysis. The adaptive cross approximation (ACA) algorithm is an efficient method requiring few computational resources introduced by Bebendorf. We introduce in the present paper the new paradigm of approximating the given function by a piecewise low-rank function with $C^1$-regularity. The proposed approximation is based on the ACA algorithm and our main contribution is the extension of the interpolation property characterising this algorithm to Hermite interpolation. Therefore, we introduce a new method for low-rank Hermite interpolation using a limiting case of the ACA algorithm. The proposed method has full approximation order. We then propose a piecewise low-rank approximation with adaptive refinement using either the ACA algorithm or our new method to compute each piece. We finally compare the results obtained for the two methods.

Introduction

In the present paper we focus on the approximation of a given bivariate function with a sum of products of univariate functions. This so-called low-rank approximation has received a lot of attention, see [1] and references therein, and [2]. A first approach to perform this approximation uses the singular value decomposition [3], [4]. This method has clear theoretical foundations, but is computationally expensive. A method reducing drastically the resources needed, the adaptive cross approximation (ACA) algorithm, has been introduced by Bebendorf in [2].

To estimate the error made by this approximation, the author focused on the case of integral operators with kernels satisfying a reasonable asymptotic decay. Schneider then showed in [5] that, under a particular choice of the points, the error can be controlled by the infimum of the error in the $L^{\infty }$-norm over all the functions having the expected rank. One can find in [6] a review of low-rank approximations of matrices, including the ACA algorithm for the case where the matrix arises from the sampling of a function on a grid.

Such approximations are notably used to implement efficient discretisations of integral operators. Indeed, a specificity of these non-local operators is that the resulting matrices are dense. Nevertheless, under some reasonable asymptotic assumptions on the kernel, a data sparse approximation of these matrices follows from their decomposition into blocks with low-rank [7], [8], [9]. Moreover, the low-rank approximation of the blocks is a direct consequence of the low-rank approximation of the considered kernel.

More recently, in the field of Isogeometric Analysis [10], the authors of [11] have used low-rank approximations to improve significantly the implementation of the discretisation of partial differential equations. Indeed, to operate isogeometric discretisations of PDEs, one mainly uses tensor spline bases mapped to the physical domain by a global parameterisation. The authors of [11] then used the tensor structure to decompose the Galerkin matrix into Kronecker products of matrix factors with small dimensions. Similarly to the case of integral equations, the key step is a low-rank approximation of the kernel (density of the integral).

An algorithm to approximate a bivariate function by low-rank splines using the ACA algorithm has then been proposed in [12]. The decomposition of the matrix obtained in [11] is even crucial when considering higher dimensions to avoid prohibitive computational costs [13], [14]. The rank of the aforementioned kernel depending directly on the rank of the parameterisation, another way to obtain low-rank kernels is to directly construct low-rank parameterisations. To this end, an optimisation approach has successfully been applied in [15] to construct low-rank parameterisations of physical domains delimited by four given curves.

On the other hand, interpolation methods, such as Coons patches [16], can be applied to construct parameterisations. An interpolation algorithm generating low-rank parameterisations has then been proposed in [17]. In a second paper, the authors have generalised this method in order to interpolate multiple curves [18]. They also showed that this method is equivalent to the ACA algorithm, this latter having the advantage of choosing adaptively the points. Indeed, the ACA algorithm performs a low-rank approximation of a given function using the restrictions of this function to some $x$-lines and $y$-lines, permitting it to be equivalently used as an interpolation method.

Contrary to the error estimates of the ACA algorithm presented in [2], [5], the authors of [18] focus on local estimations and show that this method has full approximation order. This motivates the approach taken in the present paper: we approximate hierarchically a given function by a piecewise low-rank function. The main bottleneck is to obtain $C^1$-continuity of the result, leading us to low-rank Hermite interpolation. We note that we could obtain Hermite interpolation using blending functions [19], but the result may have high rank. Therefore, we focus in the first sections of this paper on the development of a new method for low-rank Hermite interpolation.

The proposed method results from a limiting case of the ACA algorithm, renamed as the cross algorithm (CA) in what follows to underline that the points are fixed in our approach. A closed-form solution is first proposed. We then present an efficient implementation of this closed-form solution using an extension of the CA algorithm. We also prove that this new method has full approximation order. We finally show the efficiency of the piecewise low-rank approximation by numerical examples in the last section.