Efficient sum-of-exponentials approximations for the heat kernel and their applications

In this paper, we show that efficient separated sum-of-exponentials approximations can be constructed for the heat kernel in any dimension. In one space dimension, the heat kernel admits an approximation involving a number of terms that is of the order $O(\log(\frac{T}{\delta}) (\log(\frac{1}{\epsilon})+\log\log(\frac{T}{\delta})))$ for any x\in\bbR and $\delta \leq t \leq T$, where $\epsilon$ is the desired precision. In all higher dimensions, the corresponding heat kernel admits an approximation involving only $O(\log^2(\frac{T}{\delta}))$ terms for fixed accuracy $\epsilon$. These approximations can be used to accelerate integral equation-based methods for boundary value problems governed by the heat equation in complex geometry. The resulting algorithms are nearly optimal. For $N_S$ points in the spatial discretization and $N_T$ time steps, the cost is $O(N_S N_T \log^2 \frac{T}{\delta})$ in terms of both memory and CPU time for fixed accuracy $\epsilon$. The algorithms can be parallelized in a straightforward manner. Several numerical examples are presented to illustrate the accuracy and stability of these approximations.Comment: 23 pages, 5 figures, 3 table

Efficient high-order integral equation methods for the heat equation

Efficient high-order integral equation methods have been developed for solving the boundary value problems of the heat equation with complex geometries in two and three dimensions. First of all, the classical heat potential theory is applied to convert such problems to Volterra integral equations of the second kind via the heat layer potentials. Some advantages of the integral formulation as compared with standard finite difference and finite element methods include reduction of the dimension of the problem by one, high order accuracy, unconditional stability, insensitivity to different geometries, and elimination of truncating the computational domain and the need of artificial boundary conditions for exterior problems. However, the heat layer potentials contains convolution integrals in both space and time whose direct evaluation requires O(NS2NT2) work and O(NSNT) storage, where NS is the total number of discretization points in the spatial boundary and NT is the total number of time steps. This is excessively expensive even for problems of modest size, especially for three-dimensional problems. In order to evaluate the heat layer potentials accurately and efficiently, they are split into two parts - the local part which contains the temporal integration from t - δ to t and the history part which contains the temporal integration from 0 to t - δ. For the local part, Product integration is applied on the temporal integral to convert it to a sum of several spatial convolution integrals where the so-called local kernels have logarithmic singularity in two dimensions and 1r singularity in three dimensions. These weakly singular integrals are discretized via high-order quadratures and the resulting discrete summations can then be evaluated via fast algorithms such as the fast multipole method and its descendants. For the history part, efficient separated sum-of-exponentials approximations can be constructed for the heat kernel in any dimension. Specifically, in one space dimension, the heat kernel admits an approximation involving a number of terms that is of the order O(log(T δ )(log(1 ) + log log(T δ ))) for any x Î R and δ ≤ t ≤ T, where E is the desired precision. In all higher dimensions, the corresponding heat kernel admits an approximation involving only O (log2(Tδ )) terms for fixed accuracy E. These approximations can be used to accelerate the evaluation of the history part of the heat layer potentials for stationary geometries. For two-dimensional problems with complex stationary geometries, the sum-of-exponentials approximation is used for the heat kernel and all local and history kernels are compressed only once. The resulting algorithm is very efficient with quasilinear complexity in both space and time for both interior and exterior problems. For two-dimensional problems with complex moving geometries, the spectral Fourier approximation is applied for the heat kernel and NUFFT is used to speed up the evaluation of the history part of the heat potentials. The complexity of the algorithm is again quasilinear in both space and time, albeit only for the interior problem. For three-dimensional problems, the sum-of-exponentials approximation is applied to speed up the evaluation of the history part. The singular surface integrals in the local kernels are treated with a spectrally accurate integrator. The algorithm is applicable for both interior and exterior problems and has quasilinear complexity with respect to the temporal variable. All these algorithms can be parallelized in a straightforward manner and their performance is demonstrated with extensive numerical experiments