Showing papers in "Ima Journal of Numerical Analysis in 2021"

Optimal-order error estimates of finite element approximations to variable-order time-fractional diffusion equations without regularity assumptions of the true solutions

Abstract: We study a fully discretized finite element approximation to variable-order Caputo and Riemann–Liouville time-fractional diffusion equations (tFDEs) in multiple space dimensions, which model solute transport in heterogeneous porous media and related applications. We prove error estimates for the proposed methods, which are discretized on an equidistant or graded temporal partition predetermined by the behavior of the variable order at the initial time, only under the regularity assumptions of the variable order, coefficients and the source term but without any regularity assumption of the true solutions. Roughly, we prove that the finite element approximations to variable-order Caputo tFDEs have optimal-order convergence rates on a uniform temporal partition. In contrast the finite element approximations to variable-order Riemann–Liouville tFDEs discretized on a uniform temporal partition achieve an optimal-order convergence rate if $\\alpha (0)=\\alpha ^{\\prime}(0) = 0$ but a suboptimal-order convergence rate if $\\alpha (0)>0$. In the latter case, optimal-order convergence rate can be proved by employing the graded temporal partition. We conduct numerical experiments to investigate the performance of the numerical methods and to verify the mathematical analysis.

Blow-up of error estimates in time-fractional initial-boundary value problems

Abstract: Time-fractional initial-boundary value problems of the form $D_t^\\alpha u-p \\varDelta u +cu=f$ are considered, where $D_t^\\alpha u$ is a Caputo fractional derivative of order $\\alpha \\in (0,1)$ and the spatial domain lies in $\\mathbb{R}^d$ for some $d\\in \\{1,2,3\\}$. As $\\alpha \\to 1^-$ we prove that the solution $u$ converges, uniformly on the space-time domain, to the solution of the classical parabolic initial-boundary value problem where $D_t^\\alpha u$ is replaced by $\\partial u/\\partial t$. Nevertheless, most of the rigorous analyses of numerical methods for this time-fractional problem have error bounds that blow up as $\\alpha \\to 1^-$, as we demonstrate. We show that in some cases these analyses can be modified to obtain robust error bounds that do not blow up as $\\alpha \\to 1^-$.

Accurately Computing the Log-Sum-Exp and Softmax Functions

Abstract: Evaluating the log-sum-exp function or the softmax function is a key step in many modern data science algorithms, notably in inference and classification. Because of the exponentials that these functions contain, the evaluation is prone to overflow and underflow, especially in low precision arithmetic. Software implementations commonly use alternative formulas that avoid overflow and reduce the chance of harmful underflow, employing a shift or another rewriting. Although mathematically equivalent, these variants behave differently in floating-point arithmetic ew{and shifting can introduce subtractive cancellation}. We give rounding error analyses of different evaluation algorithms and interpret the error bounds using condition numbers for the functions. We conclude, based on the analysis and numerical experiments, that the shifted formulas are of similar accuracy to the unshifted ones, so can safely be used, but that a division-free variant of softmax can suffer from loss of accuracy.

Strong error analysis for stochastic gradient descent optimization algorithms

Abstract: Stochastic gradient descent (SGD) type optimization schemes are fundamental ingredients in a large number of machine learning based algorithms In particular, SGD type optimization schemes are frequently employed in applications involving natural language processing, object and face recognition, fraud detection, computational advertisement, and numerical approximations of partial differential equations In mathematical convergence results for SGD type optimization schemes there are usually two types of error criteria studied in the scientific literature, that is, the error in the strong sense and the error with respect to the objective function In applications one is often not only interested in the size of the error with respect to the objective function but also in the size of the error with respect to a test function which is possibly different from the objective function The analysis of the size of this error is the subject of this article In particular, the main result of this article proves under suitable assumptions that the size of this error decays at the same speed as in the special case where the test function coincides with the objective function

The Euler-Maruyama Scheme for SDEs with Irregular Drift: Convergence Rates via Reduction to a Quadrature Problem

Abstract: We study the strong convergence order of the Euler-Maruyama scheme for scalar stochastic differential equations with additive noise and irregular drift. We provide a general framework for the error analysis by reducing it to a weighted quadrature problem for irregular functions of Brownian motion. Assuming Sobolev-Slobodeckij-type regularity of order $\kappa \in (0,1)$ for the non-smooth part of the drift, our analysis of the quadrature problem yields the convergence order $\min\{3/4,(1+\kappa)/2\}-\epsilon$ for the equidistant Euler-Maruyama scheme (for arbitrarily small $\epsilon>0$). The cut-off of the convergence order at $3/4$ can be overcome by using a suitable non-equidistant discretization, which yields the strong convergence order of $(1+\kappa)/2-\epsilon$ for the corresponding Euler-Maruyama scheme.

Anderson acceleration for contractive and noncontractive operators

Abstract: A one-step analysis of Anderson acceleration with general algorithmic depths is presented. The resulting residual bounds within both contractive and noncontractive settings reveal the balance between the contributions from the higher and lower order terms, which are both dependent on the success of the optimization problem solved at each step of the algorithm. The new residual bounds show the additional terms introduced by the extrapolation produce terms that are of a higher order than was previously understood. In the contractive setting, these bounds sharpen previous convergence and acceleration results. The bounds rely on sufficient linear independence of the differences between consecutive residuals, rather than assumptions on the boundedness of the optimization coefficients, allowing the introduction of a theoretically sound safeguarding strategy. Several numerical tests illustrate the analysis primarily in the noncontractive setting, and demonstrate the use of the method, the safeguarding strategy, and theory-based guidance on dynamic selection of the algorithmic depth, on a p-Laplace equation, a nonlinear Helmholtz equation, and the steady Navier-Stokes equations with high Reynolds number in three spatial dimensions.

Stability of Galerkin discretizations of a mixed space–time variational formulation of parabolic evolution equations

Abstract: We analyze Galerkin discretizations of a new well-posed mixed space–time variational formulation of parabolic partial differential equations. For suitable pairs of finite element trial spaces, the resulting Galerkin operators are shown to be uniformly stable. The method is compared to two related space–time discretization methods introduced by Andreev (2013, Stability of sparse space-time finite element discretizations of linear parabolic evolution equations. IMA J. Numer. Anal., 33, 242–260) and by Steinbach (2015, Space-time finite element methods for parabolic problems. Comput. Methods Appl. Math., 15, 551–566).

Adaptive cubic regularization methods with dynamic inexact Hessian information and applications to finite-sum minimization

Abstract: We consider the Adaptive Regularization with Cubics approach for solving nonconvex optimization problems and propose a new variant based on inexact Hessian information chosen dynamically. The theoretical analysis of the proposed procedure is given. The key property of ARC framework, constituted by optimal worst-case function/derivative evaluation bounds for first- and second-order critical point, is guaranteed. Application to large-scale finite-sum minimization based on subsampled Hessian is discussed and analyzed in both a deterministic and probabilistic manner and equipped with numerical experiments on synthetic and real datasets.

Simulation of McKean Vlasov SDEs with super linear growth

Abstract: We present two fully probabilistic Euler schemes, one explicit and one implicit, for the simulation of McKean-Vlasov Stochastic Differential Equations (MV-SDEs) with drifts of super-linear growth and random initial condition We provide a pathwise propagation of chaos result and show strong convergence for both schemes on the consequent particle system The explicit scheme attains the standard $1/2$ rate in stepsize From a technical point of view, we successfully use stopping times to prove the convergence of the implicit method although we avoid them altogether for the explicit one The combination of particle interactions and random initial condition makes the proofs technically more involved Numerical tests recover the theoretical convergence rates and illustrate a computational complexity advantage of the explicit over the implicit scheme A comparative analysis is carried out on a stylized non-Lipschitz MV-SDE and a mean-field model for FitzHugh-Nagumo neurons We provide numerical tests illustrating "particle corruption" effect where one single particle diverging can "corrupt" the whole particle system Moreover, the more particles in the system the more likely this divergence is to occur

Convergence in Hölder norms with applications to Monte Carlo methods in infinite dimensions

Abstract: We show that if a sequence of piecewise affine linear processes converges in the strong sense with a positive rate to a stochastic process that is strongly Holder continuous in time, then this sequence converges in the strong sense even with respect to much stronger Holder norms and the convergence rate is essentially reduced by the Holder exponent. Our first application hereof establishes pathwise convergence rates for spectral Galerkin approximations of stochastic partial differential equations. Our second application derives strong convergence rates of multilevel Monte Carlo approximations of expectations of Banach-space-valued stochastic processes.

Domain decomposition preconditioners for high-order discretisations of the heterogeneous Helmholtz equation

Abstract: We consider one-level additive Schwarz domain decomposition preconditioners for the Helmholtz equation with variable coefficients (modelling wave propagation in heterogeneous media), subject to boundary conditions that include wave scattering problems. Absorption is included as a parameter in the problem. This problem is discretised using $H^1$-conforming nodal finite elements of fixed local degree $p$ on meshes with diameter $h = h(k)$, chosen so that the error remains bounded with increasing $k$. The action of the one-level preconditioner consists of the parallel solution of problems on subdomains (which can be of general geometry), each equipped with an impedance boundary condition. We prove rigorous estimates on the norm and field of values of the left- or right-preconditioned matrix that show explicitly how the absorption, the heterogeneity in the coefficients and the dependence on the degree enter the estimates. These estimates prove rigorously that, with enough absorption and for $k$ large enough, GMRES is guaranteed to converge in a number of iterations that is independent of $k,p,$ and the coefficients. The theoretical threshold for $k$ to be large enough depends on $p$ and on the local variation of coefficients in subdomains (and not globally). Extensive numerical experiments are given for both the absorptive and the propagative cases; in the latter case we investigate examples both when the coefficients are nontrapping and when they are trapping. These experiments (i) support our theory in terms of dependence on polynomial degree and the coefficients; (ii) support the sharpness of our field of values estimates in terms of the level of absorption required.

A unified theory for continuous-in-time evolving finite element space approximations to partial differential equations in evolving domains

Abstract: We develop a unified theory for continuous-in-time finite element discretizations of partial differential equations posed in evolving domains, including the consideration of equations posed on evolving surfaces and bulk domains, as well as coupled surface bulk systems. We use an abstract variational setting with time-dependent function spaces and abstract time-dependent finite element spaces. Optimal a priori bounds are shown under usual assumptions on perturbations of bilinear forms and approximation properties of the abstract finite element spaces. The abstract theory is applied to evolving finite elements in both flat and curved spaces. Evolving bulk and surface isoparametric finite element spaces defined on evolving triangulations are defined and developed. These spaces are used to define approximations to parabolic equations in general domains for which the abstract theory is shown to apply. Numerical experiments are described, which confirm the rates of convergence.

On the Complexity of an Augmented Lagrangian Method for Nonconvex Optimization

Abstract: In this paper we study the worst-case complexity of an inexact Augmented Lagrangian method for nonconvex constrained problems. Assuming that the penalty parameters are bounded, we prove a complexity bound of $\mathcal{O}(|\log(\epsilon)|)$ outer iterations for the referred algorithm to generate an $\epsilon$-approximate KKT point, for $\epsilon\in (0,1)$. When the penalty parameters are unbounded, we prove an outer iteration complexity bound of $\mathcal{O}\left(\epsilon^{-2/(\alpha-1)}\right)$, where $\alpha>1$ controls the rate of increase of the penalty parameters. For linearly constrained problems, these bounds yield to evaluation complexity bounds of $\mathcal{O}(|\log(\epsilon)|^{2}\epsilon^{-2})$ and $\mathcal{O}\left(\epsilon^{-\left(\frac{2(2+\alpha)}{\alpha-1}+2\right)}\right)$, respectively, when appropriate first-order methods ($p=1$) are used to approximately solve the unconstrained subproblems at each iteration. In the case of problems having only linear equality constraints, the latter bounds are improved to $\mathcal{O}(|\log(\epsilon)|^{2}\epsilon^{-(p+1)/p})$ and $\mathcal{O}\left(\epsilon^{-\left(\frac{4}{\alpha-1}+\frac{p+1}{p}\right)}\right)$, respectively, when appropriate $p$-order methods ($p\geq 2$) are used as inner solvers.

Uniform error estimates for artificial neural network approximations for heat equations

Abstract: Recently, artificial neural networks (ANNs) in conjunction with stochastic gradient descent optimization methods have been employed to approximately compute solutions of possibly rather high-dimensional partial differential equations (PDEs). Very recently, there have also been a number of rigorous mathematical results in the scientific literature which examine the approximation capabilities of such deep learning based approximation algorithms for PDEs. These mathematical results from the scientific literature prove in part that algorithms based on ANNs are capable of overcoming the curse of dimensionality in the numerical approximation of high-dimensional PDEs. In these mathematical results from the scientific literature usually the error between the solution of the PDE and the approximating ANN is measured in the $L^p$-sense with respect to some $p \in [1,\infty)$ and some probability measure. In many applications it is, however, also important to control the error in a uniform $L^\infty$-sense. The key contribution of the main result of this article is to develop the techniques to obtain error estimates between solutions of PDEs and approximating ANNs in the uniform $L^\infty$-sense. In particular, we prove that the number of parameters of an ANN to uniformly approximate the classical solution of the heat equation in a region $ [a,b]^d $ for a fixed time point $ T \in (0,\infty) $ grows at most polynomially in the dimension $ d \in \mathbb{N} $ and the reciprocal of the approximation precision $ \varepsilon > 0 $. This shows that ANNs can overcome the curse of dimensionality in the numerical approximation of the heat equation when the error is measured in the uniform $L^\infty$-norm.

Deep neural network approximation for high-dimensional elliptic PDEs with boundary conditions

Abstract: In recent work it has been established that deep neural networks are capable of approximating solutions to a large class of parabolic partial differential equations without incurring the curse of dimension. However, all this work has been restricted to problems formulated on the whole Euclidean domain. On the other hand, most problems in engineering and the sciences are formulated on finite domains and subjected to boundary conditions. The present paper considers an important such model problem, namely the Poisson equation on a domain $D\subset \mathbb{R}^d$ subject to Dirichlet boundary conditions. It is shown that deep neural networks are capable of representing solutions of that problem without incurring the curse of dimension. The proofs are based on a probabilistic representation of the solution to the Poisson equation as well as a suitable sampling method.

On strong stability of explicit Runge–Kutta methods for nonlinear semibounded operators

Abstract: Explicit Runge-Kutta methods are classical and widespread techniques in the numerical solution of ordinary differential equations (ODEs). Considering partial differential equations, spatial semidiscretisations can be used to obtain systems of ODEs that are solved subsequently, resulting in fully discrete schemes. However, certain stability investigations of high-order methods for hyperbolic conservation laws are often conducted only for the semidiscrete versions. Here, strong stability (also known as monotonicity) of explicit Runge-Kutta methods for ODEs with nonlinear and semibounded (also known as dissipative) operators is investigated. Contrary to the linear case, it is proven that many strong stability preserving (SSP) schemes of order two or greater are not strongly stable for general smooth and semibounded nonlinear operators. Additionally, it is shown that there are first order accurate explicit SSP Runge-Kutta methods that are strongly stable (monotone) for semibounded (dissipative) and Lipschitz continuous operators.

Conforming and nonconforming virtual element methods for a Kirchhoff plate contact problem

Abstract: We establish a general framework to study the conforming and nonconforming virtual element methods (VEMs) for solving a Kirchhoff plate contact problem with friction, which is a fourth-order elliptic variational inequality (VI) of the second kind. This VI contains a non-differentiable term due to the frictional contact. This theoretical framework applies to the existing virtual elements such as the conforming element, the $C^0$-continuous nonconforming element and the fully nonconforming Morley-type element. In the unified framework we derive a priori error estimates for these virtual elements and show that they achieve optimal convergence order for the lowest-order case. For demonstrating the performance of the VEMs we present some numerical results that confirm the theoretical prediction of the convergence order.

Equivalence of local-and global-best approximations, a simple stable local commuting projector, and optimal hp approximation estimates in H(div)

Abstract: Given an arbitrary function in H(div), we show that the error attained by the global-best approximation by H(div)-conforming piecewise polynomial Raviart-Thomas-Nedelec elements under additional constraints on the divergence and normal flux on the boundary, is, up to a generic constant, equivalent to the sum of independent local-best approximation errors over individual mesh elements, without constraints on the divergence or normal fluxes. The generic constant only depends on the shape-regularity of the underlying simplicial mesh, the space dimension, and the polynomial degree of the approximations. The analysis also gives rise to a stable, local, commuting projector in H(div), delivering an approximation error that is equivalent to the local-best approximation. We next present a variant of the equivalence result, where robustness of the constant with respect to the polynomial degree is attained for unbalanced approximations. These two results together further enable us to derive rates of convergence of global-best approximations that are fully optimal in both the mesh size h and the polynomial degree p, for vector fields that only feature elementwise the minimal necessary Sobolev regularity. We finally show how to apply our findings to derive optimal a priori hp-error estimates for mixed and least-squares finite element methods applied to a model diffusion problem.

Conforming and nonconforming VEMs for the fourth-order reaction–subdiffusion equation: a unified framework

Abstract: We establish a unified framework to study the conforming and nonconforming virtual element methods (VEMs) for a class of time dependent fourth-order reaction–subdiffusion equations with the Caputo derivative. To resolve the intrinsic initial singularity we adopt the nonuniform Alikhanov formula in the temporal direction. In the spatial direction three types of VEMs, including conforming virtual element, $C^0$ nonconforming virtual element and fully nonconforming Morley-type virtual element, are constructed and analysed. In order to obtain the desired convergence results, the classical Ritz projection operator for the conforming virtual element space and two types of new Ritz projection operators for the nonconforming virtual element spaces are defined, respectively, and the projection errors are proved to be optimal. In the unified framework we derive a prior error estimate with optimal convergence order for the constructed fully discrete schemes. To reduce the computational cost and storage requirements, the sum-of-exponentials (SOE) technique combined with conforming and nonconforming VEMs (SOE-VEMs) are built. Finally, we present some numerical experiments to confirm the theoretical correctness and the effectiveness of the discrete schemes. The results in this work are fundamental and can be extended into more relevant models.

Energy-preserving methods for nonlinear Schrödinger equations

Abstract: This paper is concerned with the numerical integration in time of nonlinear Schrodinger equations using different methods preserving the energy or a discrete analog of it. The Crank-Nicolson method is a well known method of order 2 but is fully implicit and one may prefer a linearly implicit method like the relaxation method introduced in [10] for the cubic nonlinear Schrodinger equation. This method is also an energy preserving method and numerical simulations have shown that its order is 2. In this paper we give a rigorous proof of the order of this relaxation method and propose a generalized version that allows to deal with general power law nonlinearites. Numerical simulations for different physical models show the efficiency of these methods.

An unfitted Eulerian finite element method for the time-dependent Stokes problem on moving domains

Abstract: We analyse a Eulerian Finite Element method, combining a Eulerian time-stepping scheme applied to the time-dependent Stokes equations using the CutFEM approach with inf-sup stable Taylor-Hood elements for the spatial discretisation. This is based on the method introduced by Lehrenfeld \& Olshanskii [ESAIM: M2AN 53(2):585--614] in the context of a scalar convection-diffusion problems on moving domains, and extended to the non-stationary Stokes problem on moving domains by Burman, Frei \& Massing [arXiv:1910.03054 [math.NA]] using stabilised equal-order elements. The analysis includes the geometrical error made by integrating over approximated levelset domains in the discrete CutFEM setting. The method is implemented and the theoretical results are illustrated using numerical examples.

hp-FEM for the fractional heat equation

Abstract: We consider a time dependent problem generated by a nonlocal operator in space. Applying a discretization scheme based on $hp$-Finite Elements and a Caffarelli-Silvestre extension we obtain a semidiscrete semigroup. The discretization in time is carried out by using $hp$-Discontinuous Galerkin based timestepping. We prove exponential convergence for such a method in an abstract framework for the discretization in the original domain $\Omega$.

The Scharfetter–Gummel scheme for aggregation–diffusion equations

Abstract: In this paper, we propose a finite-volume scheme for aggregation-diffusion equations that is based on a Scharfetter--Gummel approximation of the nonlinear, nonlocal flux term. This scheme is analyzed concerning well-posedness and convergence towards solutions to the continuous problem. Also, it is proven that the numerical scheme has several structure-preserving features. More specifically, it is shown that the discrete solutions satisfy a free-energy dissipation relation analogous to the continuous model, and, as a consequence, the numerical solutions converge in the large time limit to stationary solutions, for which we provide a thermodynamic characterization.

The time-fractional Cahn–Hilliard equation: analysis and approximation

Abstract: We consider a time-fractional Cahn–Hilliard equation where the conventional first-order time derivative is replaced by a Caputo fractional derivative of order $\\alpha \\in (0,1)$. Based on an a priori bound of the exact solution, global existence of solutions is proved and detailed regularity results are included. A finite element method is then analyzed in a spatially discrete case and in a completely discrete case based on a convolution quadrature in time generated by the backward Euler method. Error bounds of optimal order are obtained for solutions with smooth and nonsmooth initial data, thereby extending earlier studies on the classical Cahn–Hilliard equation. Further, by proving a new result concerning the positivity of a discrete time-fractional integral operator, it is shown that the proposed numerical scheme inherits a discrete energy dissipation law at the discrete level. Numerical examples are presented to illustrate the theoretical results.

A numerical analysis focused comparison of several Finite Volume schemes for an Unipolar Degenerated Drift-Diffusion Model

Abstract: In this paper, we consider an unipolar degenerated drift-diffusion system where the relation between the concentration of the charged species c and the chemical potential h is h(c) = log c 1−c. We design four different finite volume schemes based on four different formulations of the fluxes. We provide a stability analysis and existence results for the four schemes. The convergence proof with respect to the discretization parameters is established for two of them. Numerical experiments illustrate the behaviour of the different schemes.

Space–time discontinuous Galerkin approximation of acoustic waves with point singularities

Abstract: We develop a convergence theory of space-time discretizations for the linear, 2nd-order wave equation in polygonal domains $\Omega\subset\mathbb{R}^2$, possibly occupied by piecewise homogeneous media with different propagation speeds. Building on an unconditionally stable space-time DG formulation developed in [Moiola, Perugia 2018], we (a) prove optimal convergence rates for the space-time scheme with local isotropic corner mesh refinement on the spatial domain, and (b) demonstrate numerically optimal convergence rates of a suitable \emph{sparse} space-time version of the DG scheme. The latter scheme is based on the so-called \emph{combination formula}, in conjunction with a family of anisotropic space-time DG-discretizations. It results in optimal-order convergent schemes, also in domains with corners, with a number of degrees of freedom that scales essentially like the DG solution of one stationary elliptic problem in $\Omega$ on the finest spatial grid. Numerical experiments for both smooth and singular solutions support convergence rate optimality on spatially refined meshes of the full and sparse space-time DG schemes.

Multivariate approximation of functions on irregular domains by weighted least-squares methods

Abstract: We propose and analyse numerical algorithms based on weighted least squares for the approximation of a bounded real-valued function on a general bounded domain Ω Ă R d. Given any n-dimensional approximation space Vn Ă L 2 pΩq, the analysis in [6] shows the existence of stable and optimally converging weighted least-squares estimators, using a number of function evaluations m of the order n ln n. When an L 2 pΩqorthonormal basis of Vn is available in analytic form, such estimators can be constructed using the algorithms described in [6, Section 5]. If the basis also has product form, then these algorithms have computational complexity linear in d and m. In this paper we show that, when Ω is an irregular domain such that the analytic form of an L 2 pΩq-orthonormal basis is not available, stable and quasi-optimally weighted leastsquares estimators can still be constructed from Vn, again with m of the order n ln n, but using a suitable surrogate basis of Vn orthonormal in a discrete sense. The computational cost for the calculation of the surrogate basis depends on the Christoffel function of Ω and Vn. Numerical results validating our analysis are presented.

A Banach space mixed formulation for the unsteady Brinkman–Forchheimer equations

Abstract: We propose and analyze a mixed formulation for the Brinkman-Forchheimer equations for unsteady flows. Our approach is based on the introduction of a pseudostress tensor related to the velocity gradient, leading to a mixed formulation where the pseudostress tensor and the velocity are the main unknowns of the system. We establish existence and uniqueness of a solution to the weak formulation in a Banach space setting, employing classical results on nonlinear monotone operators and a regularization technique. We then present well-posedness and error analysis for semidiscrete continuous-in-time and fully discrete finite element approximations on simplicial grids with spatial discretization based on the Raviart-Thomas spaces of degree $k$ for the pseudostress tensor and discontinuous piecewise polynomial elements of degree $k$ for the velocity and backward Euler time discretization. We provide several numerical results to confirm the theoretical rates of convergence and illustrate the performance and flexibility of the method for a range of model parameters.