Search

Polynomial reconstruction on Cartesian grids is fundamental in many scientific and engineering applications, yet it is still an open problem how to construct for a finite subset K of Z D a lattice T ⊂ K so that multivariate polynomial interpolation on this lattice is unisolvent. In this work, we solve this open problem of poised lattice generation (PLG) via an interdisciplinary research of approximation theory, abstract algebra, and artificial intelligence. Specifically, we focus on the triangular lattices in approximation theory, study group actions of permutations upon triangular lattices, prove an isomorphism between the group of permutations and that of triangular lattices, and dynamically organize the state space of permutations so that a depth-first search of poised lattices has optimal efficiency. Based on this algorithm, we further develop the PLG finite difference method that retains the simplicity of Cartesian grids yet overcomes the disadvantage of legacy finite difference methods in handling irregular geometries. Results of various numerical tests demonstrate the effectiveness of our algorithm and the simplicity, flexibility, efficiency, and fourth-order accuracy of the PLG finite difference method. [ABSTRACT FROM AUTHOR]

In this study, we investigate an anisotropic weakly over-penalised symmetric interior penalty method for the Stokes equation on convex domains. Our approach is a simple discontinuous Galerkin method similar to the Crouzeix–Raviart finite element method. As our primary contribution, we show a new proof for the consistency term, which allows us to obtain an estimate of the anisotropic consistency error. The key idea of the proof is to apply the relation between the Raviart–Thomas finite element space and a discontinuous space. While inf-sup stable schemes of the discontinuous Galerkin method on shape-regular mesh partitions have been widely discussed, our results show that the Stokes element satisfies the inf-sup condition on anisotropic meshes. Furthermore, we provide an error estimate in an energy norm on anisotropic meshes. In numerical experiments, we compare calculation results for standard and anisotropic mesh partitions. [ABSTRACT FROM AUTHOR]

We consider the axisymmetric formulation of the equilibrium problem for a hot plasma in a tokamak. We adopt a non-overlapping mortar element approach, that couples C 0 piece-wise linear Lagrange finite elements in a region that does not contain the plasma and C 1 piece-wise cubic reduced Hsieh–Clough–Tocher finite elements elsewhere, to approximate the magnetic flux field on a triangular mesh of the poloidal tokamak section. The inclusion of ferromagnetic parts is simplified by assuming that they fit within the axisymmetric modeling and a new formulation of the Newton algorithm for the problem solution is stated, both in the static and quasi-static evolution cases. [ABSTRACT FROM AUTHOR]

Hermite interpolation is classically used to reconstruct smooth data when the function and its first order derivatives are available at certain nodes. If first order derivatives are not available, it is easy to set a system of equations imposing some regularity conditions at the data nodes in order to obtain them. This process leads to the construction of a Hermite spline. The problem of the described Hermite splines is that the accuracy is lost if the data contains singularities. The consequence is the appearance of oscillations, if there is a jump discontinuity in the function, that globally affects the accuracy of the spline, or the smearing of singularities, if the discontinuities are in the derivatives of the function. This paper is devoted to the construction and analysis of a new technique that allows for the computation of accurate first order derivatives of a function close to singularities using a Hermite spline. The idea is to correct the system of equations of the spline in order to attain the desired accuracy even close to the singularities. Once we have computed the first order derivatives with enough accuracy, a correction term is added to the Hermite spline in the intervals that contain a singularity. The aim is to reconstruct piecewise smooth functions with O (h 4) accuracy even close to the singularities. The process of adaption will require some knowledge about the position of the singularity and the jumps of the function and some of its derivatives at the singularity. The whole process can be used as a post-processing, where a correction term is added to the classical cubic Hermite spline. Proofs for the accuracy and regularity of the corrected spline and its derivatives are given. We also analyse the mechanism that eliminates the Gibbs phenomenon close to jump discontinuities in the function. The numerical experiments presented confirm the theoretical results obtained. [ABSTRACT FROM AUTHOR]

Finite elements with polynomial basis functions on the simplex with a symmetric distribution of nodes should have a unique polynomial representation. Unisolvence not only requires that the number of nodes equals the number of independent polynomials spanning a polynomial space of a given degree, but also that the Vandermonde matrix controlling their mapping to the Lagrange interpolating polynomials can be inverted. Here, a necessary condition for unisolvence is presented for polynomial spaces that have non-decreasing degrees when going from the edges and the various faces to the interior of the simplex. It leads to a proof of a conjecture on a necessary condition for unisolvence, requiring the node pattern to be the same as that of the regular simplex. [ABSTRACT FROM AUTHOR]

We propose an efficient mesh adaptive method for the numerical solution of time-dependent partial differential equations considered in the fixed space-time cylinder Ω × (0 , T) . We employ the space-time discontinuous Galerkin method which enables us to use different meshes at different time levels in a natural way. The mesh adaptive algorithm is based on control of the interpolation error in the L ∞ (0 , T ; L q (Ω)) -norm. The goal is to construct a sequence of conforming triangular meshes in such a way that the interpolation error bound is under a given tolerance and the number of degrees of freedom is minimal. The resulting grids consist of anisotropic mesh elements with varying polynomial approximation degrees with respect to space. We present a theoretical framework of this approach as well as several numerical examples demonstrating the accuracy, efficiency, and applicability of the method. [ABSTRACT FROM AUTHOR]

In one dimension, nodal locations that are distinct are necessary and sufficient to ensure that a unique polynomial interpolant exists for data provided at a set of nodes, i.e. that the set of nodes is unisolvent. In multiple dimensions however, unisolvency for a polynomial interpolant of degree p is not ensured even with nodal locations that are distinct and a set of n nodes, with n equal to the cardinality of a set of polynomial basis functions of at most degree p. In this paper a set of equations is derived for simplices of one to three dimensions with symmetrical nodal distributions to identify a combination of symmetry orbits that can provide a unisolvent set of nodes. The results suggest that there is a unique combination of symmetry orbits that can provide a unisolvent set of nodes for each degree of polynomial interpolant. Consequently, all other combinations of symmetry orbits cannot provide a unisolvent set of nodes for a degree p polynomial interpolant. This is verified numerically up to degree 10 for triangles and degree 7 for tetrahedra. The results suggest that the same is also true for higher-order polynomial interpolants. This significantly reduces the number of combination of symmetry orbits that needs to be considered. For example, for a tetrahedron with a degree seven interpolant, only one combination of symmetry orbits needs to be considered instead of the 161 different combinations of symmetry orbits that provide a set of nodes with n equal to the cardinality of the set of basis functions of at most degree seven. For a symmetrical nodal distribution in a simplex, the conditions presented are necessary but not sufficient to have a unisolvent set of nodes for polynomial interpolation. [ABSTRACT FROM AUTHOR]

The main purpose of this work is to provide an efficient scheme for constructing kernel-based reduced interpolation models. In the existing literature such problems are mainly addressed via the well-established knot insertion or knot removal schemes. Such iterative strategies are usually quite demanding from a computational point of view and our goal is to study an efficient implementation for data removal approaches, namely efficient reduced basis algorithm (ERBA). Focusing on kernel-based interpolation, the algorithm makes use of two iterative rules for removing data. The former, called ERBA-r, is based on classical residual evaluations. The latter, namely ERBA-p, is independent of the function values and relies on error bounds defined by the power function. In both cases, inspired by the so-called extended Rippa’s algorithm, our ERBA takes advantage of a fast implementation. [ABSTRACT FROM AUTHOR]

We propose an efficient, deterministic algorithm for constructing exponentially convergent deep neural network (DNN) approximations of multivariate, analytic maps f : [ - 1 , 1 ] K → R . We address in particular networks with the rectified linear unit (ReLU) activation function. Similar results and proofs apply for many other popular activation functions. The algorithm is based on collocating f in deterministic families of grid points with small Lebesgue constants, and by a-priori (i.e., “offline”) emulation of a spectral basis with DNNs to prescribed fidelity. Assuming availability of N function values of a possibly corrupted, numerical approximation f ˘ of f in [ - 1 , 1 ] K and a bound on ‖ f - f ˘ ‖ L ∞ ([ - 1 , 1 ] K) , we provide an explicit, computational construction of a ReLU DNN which attains accuracy ε (depending on N and ‖ f - f ˘ ‖ L ∞ ([ - 1 , 1 ] K) ) uniformly, with respect to the inputs. For analytic maps f : [ - 1 , 1 ] K → R , we prove exponential convergence of expression and generalization errors of the constructed ReLU DNNs. Specifically, for every target accuracy ε ∈ (0 , 1) , there exists N depending also on f such that the error of the construction algorithm with N evaluations of f ˘ as input in the norm L ∞ ([ - 1 , 1 ] K ; R) is smaller than ε up to an additive data-corruption bound ‖ f - f ˘ ‖ L ∞ ([ - 1 , 1 ] K) multiplied with a factor growing slowly with 1 / ε and the number of non-zero DNN weights grows polylogarithmically with respect to 1 / ε . The algorithmic construction of the ReLU DNNs which will realize the approximations, is explicit and deterministic in terms of the function values of f ˘ in tensorized Clenshaw–Curtis grids in [ - 1 , 1 ] K . We illustrate the proposed methodology by a constructive algorithm for (offline) computations of posterior expectations in Bayesian PDE inversion. [ABSTRACT FROM AUTHOR]

The data driven extrapolation requires the definition of a functional model depending on the available data and has the application scope of providing reliable predictions on the unknown dynamics. Since data might be scattered, we drive our attention towards kernel models that have the advantage of being meshfree. Precisely, the proposed numerical method makes use of the so-called Variably Scaled Kernels (VSKs), which are introduced to implement a feature augmentation-like strategy based on discrete data. Due to the possible uncertainty on the data and since we are interested in modelling the behaviour of the target functions, we seek for a regularized solution by ridge regression. Focusing on polyharmonic splines, we investigate their implementation in the VSK setting and we provide error bounds in Beppo–Levi spaces. The performances of the method are then tested on functions showing exponential or rational decay. Comparisons with Support Vector Regression (SVR) are also carried out and highlight that the proposed approach is effective, particularly since it does not require to train complex architecture constructions. [ABSTRACT FROM AUTHOR]

It is well understood that boundary conditions (BCs) may cause global radial basis function (RBF) methods to become unstable for hyperbolic conservation laws (CLs). Here we investigate this phenomenon and identify the strong enforcement of BCs as the mechanism triggering such stability issues. Based on this observation we propose a technique to weakly enforce BCs in RBF methods. In the case of hyperbolic CLs, this is achieved by carefully building RBF methods from the weak form of the CL, rather than the typically enforced strong form. Furthermore, we demonstrate that global RBF methods may violate conservation, yielding physically unreasonable solutions when the approximation does not take into account these considerations. Numerical experiments validate our theoretical results. [ABSTRACT FROM AUTHOR]

In this article we present a new adaptive algorithm for solving 2D interpolation problems of large scattered data sets through the radial basis function partition of unity method. Unlike other time-consuming schemes this adaptive method is able to efficiently deal with scattered data points with highly varying density in the domain. This target is obtained by decomposing the underlying domain in subdomains of variable size so as to guarantee a suitable number of points within each of them. The localization of such points is done by means of an efficient search procedure that depends on a partition of the domain in square cells. For each subdomain the adaptive process identifies a predefined neighborhood consisting of one or more levels of neighboring cells, which allows us to quickly find all the subdomain points. The algorithm is further devised for an optimal selection of the local shape parameters associated with radial basis function interpolants via leave-one-out cross validation and maximum likelihood estimation techniques. Numerical experiments show good performance of this adaptive algorithm on some test examples with different data distributions. The efficacy of our interpolation scheme is also pointed out by solving real world applications. [ABSTRACT FROM AUTHOR]

We propose a fast algorithm for the variable-order (VO) space-fractional advection-diffusion equations with nonlinear source terms on a finite domain. Due to the impact of the space-dependent the VO, the resulting coefficient matrices arising from the finite difference discretization of the fractional advection-diffusion equation are dense without Toeplitz-like structure. By the properties of the elements of coefficient matrices, we show that the off-diagonal blocks can be approximated by low-rank matrices. Then we present a fast algorithm based on the polynomial interpolation to approximate the coefficient matrices. The approximation can be constructed in O (k N) operations and requires O (k N) storage with N and k being the number of unknowns and the approximants, respectively. Moreover, the matrix-vector multiplication can be implemented in O (k N log N) complexity, which leads to a fast iterative solver for the resulting linear systems. The stability and convergence of the new scheme are also studied. Numerical tests are carried out to exemplify the accuracy and efficiency of the proposed method. [ABSTRACT FROM AUTHOR]

We introduce a new local meshfree method for the approximation of the Laplace–Beltrami operator on a smooth surface in R 3 . It is a direct method that uses radial basis functions augmented with multivariate polynomials. A key element of this method is that it does not need an explicit expression of the surface, which can be simply defined by a set of scattered nodes. Likewise, it does not require expressions for the surface normal vectors or for the curvature of the surface, which are approximated using explicit formulas derived in the paper. An additional advantage is that it is a local method and, hence, the matrix that approximates the Laplace–Beltrami operator is sparse, which translates into good scalability properties. The convergence, accuracy and other computational characteristics of the proposed method are studied numerically. Its performance is shown by solving two reaction–diffusion partial differential equations on surfaces; the Turing model for pattern formation, and the Schaeffer's model for electrical cardiac tissue behavior. [ABSTRACT FROM AUTHOR]

Including polynomials with small degree and stencil when designing very high order reconstructions is surely beneficial for their non oscillatory properties, but may bring loss of accuracy on smooth data unless special care is exerted. In this paper we address this issue with a new Central WENOZ (CWENOZ ) approach, in which the reconstruction polynomial is computed from a single set of non linear weights, but the linear weights of the polynomials with very low degree (compared to the final desired accuracy) are infinitesimal with respect to the grid size. After proving general results that guide the choice of the CWENOZ parameters, we study a concrete example of a reconstruction that blends polynomials of degree six, four and two, mimicking already published Adaptive Order WENO reconstructions (Arbogast et al. in SIAM J Numer Anal 56(3):1818-1947, 2018),(Balsara et al. in J Comput Phys 326:780-804, 2016). The novel reconstruction yields similar accuracy and oscillations with respect to the previous ones, but saves up to 20% computational time since it does not rely on a hierarchic approach and thus does not compute multiple sets of nonlinear weights in each cell. [ABSTRACT FROM AUTHOR]

Essentially non-oscillatory (ENO) and weighted ENO (WENO) methods on equidistant Cartesian grids are widely used to solve partial differential equations with discontinuous solutions. However, stable ENO/WENO methods on unstructured grids are less well studied. We propose a high-order ENO method based on radial basis function (RBF) to solve hyperbolic conservation laws on general two-dimensional grids. The radial basis function reconstruction offers a flexible way to deal with ill-conditioned cell constellations. We introduce a smoothness indicator based on RBFs and a stencil selection algorithm suitable for general meshes. Furthermore, we develop a stable method to evaluate the RBF reconstruction in the finite volume setting which circumvents the stagnation of the error and keeps the condition number of the reconstruction bounded. We conclude with several challenging numerical examples in two dimensions to show the robustness of the method. [ABSTRACT FROM AUTHOR]

In this paper we present a trivariate algorithm for fast computation of tetrahedral Shepard interpolants. Though the tetrahedral Shepard method achieves an approximation order better than classical Shepard formulas, it requires to detect suitable configurations of tetrahedra whose vertices are given by the set of data points. In doing that, we propose the use of a fast searching procedure based on the partitioning of domain and nodes in cubic blocks. This allows us to find the nearest neighbor points associated with each ball that need to be used in the 3D interpolation scheme. Numerical experiments show good performance of our interpolation algorithm. [ABSTRACT FROM AUTHOR]

We investigate adaptivity issues for the approximation of Poisson equations via radial basis function-based partition of unity collocation. The adaptive residual subsampling approach is performed with quasi-uniform node sequences leading to a flexible tool which however might suffer from numerical instability due to ill-conditioning of the collocation matrices. We thus develop a hybrid method which makes use of the so-called variably scaled kernels. The proposed algorithm numerically ensures the convergence of the adaptive procedure. [ABSTRACT FROM AUTHOR]

In this paper we present a method that uses radial basis functions to approximate the Laplace-Beltrami operator that allows to solve numerically diffusion (and reaction-diffusion) equations on smooth, closed surfaces embedded in R3. The novelty of the method is in a closed-form formula for the Laplace-Beltrami operator derived in the paper, which involve the normal vector and the curvature at a set of points on the surface of interest. An advantage of the proposed method is that it does not rely on the explicit knowledge of the surface, which can be simply defined by a set of scattered nodes. In that case, the surface is represented by a level set function from which we can compute the needed normal vectors and the curvature. The formula for the Laplace-Beltrami operator is exact for radial basis functions and it also depends on the first and second derivatives of these functions at the scattered nodes that define the surface. We analyze the converge of the method and we present numerical simulations that show its performance. We include an application that arises in cardiology. [ABSTRACT FROM AUTHOR]

This paper is devoted to the construction and analysis of an adapted and nonlinear multiresolution algorithm designed for interpolation or approximation of discontinuous univariate functions. The adaption attained allows to avoid numerical artifacts that appear when using linear algorithms and, at the same time, to obtain a high order of accuracy close to the singularities. It is known that linear algorithms are stable and convergent for smooth functions, but diffusion and Gibbs effect appear if the functions are piecewise continuous. Our aim is to develop an algorithm for function approximation with full accuracy that is capable to adapt to corners (kinks) and jump discontinuities, that uses a centered stencil and that does not use extrapolation. In order to reach this goal, we will need some information about the jumps in the function that we want to approximate and its derivatives. If this information is available, the algorithm is the most compact possible in the sense that the stencil is fixed and we do not need a stencil selection procedure as other algorithms do, such as ENO subcell resolution (ENO-SR). If the information about the jumps is not available, we will show a technique to approximate it. The algorithm is based on linear interpolation plus correction terms that provide the desired accuracy close to corners or jump discontinuities. [ABSTRACT FROM AUTHOR]

In this paper we use the idea of the weighted nonlocal Laplacian (Shi et al. in J Sci Comput, 2017) to deal with the constraints in the low dimensional manifold model (Osher et al. in SIAM J Imaging Sci, 2017). In the original LDMM, the constraints are enforced by the point integral method. The point integral method provides a correct way to deal with the constraints, however it is not very efficient due to the fact that the symmetry of the original Laplace-Beltrami operator is destroyed. WNLL provides another way to enforce the constraints in LDMM. In WNLL, the discretized system is symmetric and sparse and hence it can be solved very fast. Our experimental results show that the computational cost is reduced significantly with the help of WNLL. Moreover, the results in image inpainting and denoising are also better than the original LDMM and competitive with state-of-the-art methods. [ABSTRACT FROM AUTHOR]

In this paper, we are concerned with Jacobi polynomials Pn(α,β)(x) on the Bernstein ellipse with motivation mainly coming from recent studies of convergence rate of spectral interpolation. An explicit representation of Pn(α,β)(x) is derived in the variable of parametrization. This formula further allows us to show that the maximum value of Pn(α,β)(z) over the Bernstein ellipse is attained at one of the endpoints of the major axis if α+β≥-1. For the minimum value, we are able to show that for a large class of Gegenbauer polynomials (i.e., α=β), it is attained at two endpoints of the minor axis. These results particularly extend those previously known only for some special cases. Moreover, we obtain a more refined asymptotic estimate for Jacobi polynomials on the Bernstein ellipse. [ABSTRACT FROM AUTHOR]

We present a parallel algorithm for multivariate Radial Basis Function Partition of Unity Method (RBF-PUM) interpolation. The concurrent nature of the RBF-PUM enables designing parallel algorithms for dealing with a large number of scattered data-points in high space dimensions. To efficiently exploit this concurrency, our algorithm makes use of shared-memory parallel processors through the OpenCL standard. This efficiency is achieved by a parallel space partitioning strategy with linear computational time complexity with respect to the input and evaluation points. The speed of our approach allows for computationally more intensive construction of the interpolant. In fact, the RBF-PUM can be coupled with a cross-validation technique that searches for optimal values of the shape parameters associated with each local RBF interpolant, thus reducing the global interpolation error. The numerical experiments support our claims by illustrating the interpolation errors and the running times of our algorithm. [ABSTRACT FROM AUTHOR]

The partition of unity (PU) method, performed with local radial basis function (RBF) approximants, has been proved to be an effective tool for solving large scattered data interpolation problems. However, in order to achieve a good accuracy, the question about how many points we have to consider on each local subdomain, i.e. how large can be the local data sets, needs to be answered. Moreover, it is well-known that also the shape parameter affects the accuracy of the local RBF approximants and, as a consequence, of the PU interpolant. Thus here, both the shape parameter used to fit the local problems and the size of the associated linear systems are supposed to vary among the subdomains. They are selected by minimizing an a priori error estimate. As evident from extensive numerical experiments and applications provided in the paper, the proposed method turns out to be extremely accurate also when data with non-homogeneous density are considered. [ABSTRACT FROM AUTHOR]

Inspired by the nonlocal methods in image processing and the point integral method, we introduce a novel weighted nonlocal Laplacian method to compute a continuous interpolation function on a point cloud in high dimensional space. The numerical results in semi-supervised learning and image inpainting show that the weighted nonlocal Laplacian is a reliable and efficient interpolation method. In addition, it is fast and easy to implement. [ABSTRACT FROM AUTHOR]

We explore the use of radial basis functions (RBF) in the weighted essentially non-oscillatory (WENO) reconstruction process used to solve hyperbolic conservation laws, resulting in a numerical method of arbitrarily high order to solve problems with discontinuous solutions. Thanks to the mesh-less property of the RBFs, the method is suitable for non-uniform grids and mesh adaptation. We focus on multiquadric radial basis functions and propose a simple strategy to choose the shape parameter to control the balance between achievable accuracy and the numerical stability. We also develop an original smoothness indicator which is independent of the RBF for the WENO reconstruction step. Moreover, we introduce type I and type II RBF-WENO methods by computing specific linear weights. The RBF-WENO method is used to solve linear and nonlinear problems for both scalar and systems of conservation laws, including Burgers equation, the Buckley-Leverett equation, and the Euler equations. Numerical results confirm the performance of the proposed method. We finally consider an effective conservative adaptive algorithm that captures moving shocks and rapidly varying solutions well. Numerical results on moving grids are presented for both Burgers equation and the more complex Euler equations. [ABSTRACT FROM AUTHOR]

This paper is devoted to the construction and analysis of a new smoothness index for WENO interpolation capable of dealing with corner discontinuities. The new smoothness index presented is initially developed for the point-value framework of Harten's multiresolution. Even so, the ideas about how to extend the results to the cell-average framework are presented. The new smoothness index is inspired by the one proposed in Jiang and Shu (J Comput Phys 126(1):202-228, 1996). This index works very well for jump discontinuities as it was originally designed for the context of conservation laws in order to deal with problems that contain shocks and complicated fluid-structure interactions. Even so, it is easy to check that the mentioned index does not provide an appropriate performance for corner discontinuities. Our aim is to rise the order of accuracy of WENO interpolation near corner discontinuities. In order to do so, we will modify the original smoothness index proposed by Jiang and Shu such that the discontinuities in the first derivative of the function contribute effectively to the index. The modification proposed will produce a variation in the weights of WENO when dealing with a corner, that do not appear when using the smoothness indexes proposed by Jiang and Shu. The variation in the weights induced by the modification of the smoothness index will allow adaption to corner discontinuities, maintaining the adaption to jumps provided by the original smoothness index proposed by Jiang and Shu. The strategy proposed in Aràndiga et al. (SIAM J Numer Anal 49(2):893-915, 2011) can be adapted such that the accuracy is maintained near critical points at smooth zones. [ABSTRACT FROM AUTHOR]

A particular solution of a dynamical system is completely determined by its initial condition. When the omega limit set reduces to a point, the solution settles at steady state. The possible steady states of the system are completely determined by its parameters. However, with the same parameter set, it is possible that several steady states can originate from different initial conditions (multi-stability). In that case the outcome depends on the chosen initial condition. Therefore, it is important to assess the domain of attraction for each possible attractor. The algorithms presented here are general and robust enough so as to solve the problem of reconstructing the basin of attraction of each stable equilibrium point. In order to have a graphical representation of the separatrix manifolds, we focus on systems of two and three ordinary differential equations exhibiting bi- or tri-stability. For this purpose we have implemented several Matlab functions for the approximation of the points lying on the curves or on the surfaces determining the basins of attraction and for the reconstruction of such curves and surfaces. We approximate the latter with the implicit partition of unity method using radial basis functions as local approximants. Numerical results, obtained with a Matlab package made available to the scientific community, support our findings. [ABSTRACT FROM AUTHOR]