Showing papers in "Mathematics of Computation in 2010"
TL;DR: In this paper, the authors introduce a new technique for the error analysis of hybridizable discontinuous Galerkin (HDG) methods, which relies on the use of a new projection whose design is inspired by the form of the numerical traces of the methods.
Abstract: We introduce a new technique for the error analysis of hybridizable discontinuous Galerkin (HDG) methods. The technique relies on the use of a new projection whose design is inspired by the form of the numerical traces of the methods. This renders the analysis of the projections of the discretization errors simple and concise. By showing that these projections of the errors are bounded in terms of the distance between the solution and its projection, our studies of influence of the stabilization parameter are reduced to local analyses of approximation by the projection. We illustrate the technique on a specific HDG method applied to a model second-order elliptic problem.
254 citations
TL;DR: A Jacobi-collocation spectral method is developed for Volterra integral equations of second kind with a weakly singular kernel so that the solution of the new equation possesses better regularity and the Jacobi orthogonal polynomial theory can be applied conveniently.
Abstract: In this paper, a Jacobi-collocation spectral method is developed for Volterra integral equations of second kind with a weakly singular kernel. We use some function transformation and variable transformations to change the equation into a new Volterra integral equation deflned on the standard interval (i1;1), so that the solution of the new equation possesses better regularity and the Jacobi orthogonal polynomial theory can be applied conveniently. In order to obtain high order accuracy for the approx- imation, the integral term in the resulting equation is approximated by using Jacobi spectral quadrature rules. The convergence analysis of this novel method is based on the Lebesgue constants corresponding to the Lagrange interpolation polynomials, poly- nomials approximation theory for orthogonal polynomials and the operator theory. The spectral rate of convergence for the proposed method is established in the L 1 -norm and weighted L 2 -norm. Numerical results are presented to demonstrate the efiectiveness of the proposed method.
231 citations
TL;DR: A new multiscale finite element method which is able to accurately capture solutions of elliptic interface problems with high contrast coefficients by using only coarse quasiuniform meshes, and without resolving the interfaces is introduced.
Abstract: We introduce a new multiscale finite element method which is
able to accurately capture solutions of elliptic interface problems with high
contrast coefficients by using only coarse quasiuniform meshes, and without
resolving the interfaces. A typical application would be the modelling of flow
in a porous medium containing a number of inclusions of low (or high) permeability
embedded in a matrix of high (respectively low) permeability. Our
method is H^1- conforming, with degrees of freedom at the nodes of a triangular
mesh and requiring the solution of subgrid problems for the basis functions on
elements which straddle the coefficient interface but which use standard linear
approximation otherwise. A key point is the introduction of novel coefficientdependent
boundary conditions for the subgrid problems. Under moderate
assumptions, we prove that our methods have (optimal) convergence rate of
O(h) in the energy norm and O(h^2) in the L_2 norm where h is the (coarse)
mesh diameter and the hidden constants in these estimates are independent
of the “contrast” (i.e. ratio of largest to smallest value) of the PDE coefficient.
For standard elements the best estimate in the energy norm would be
O(h^(1/2−e)) with a hidden constant which in general depends on the contrast.
The new interior boundary conditions depend not only on the contrast of the
coefficients, but also on the angles of intersection of the interface with the
element edges.
207 citations
TL;DR: A new error analysis of discontinuous Galerkin methods is developed using only the H k weak formulation of a boundary value problem of order 2k using a discrete energy norm that is well defined for functions in H k .
Abstract: The standard a priori error analysis of discontinuous Galerkin methods requires additional regularity on the solution of the elliptic boundary value problem in order to justify the Galerkin orthogonality and to handle the normal derivative on element interfaces that appear in the discrete energy norm. In this paper, a new error analysis of discontinuous Galerkin methods is developed using only the H k weak formulation of a boundary value problem of order 2k. This is accomplished by replacing the Galerkin orthogonality with estimates borrowed from a posteriori error analysis and by using a discrete energy norm that is well defined for functions in H k .
204 citations
TL;DR: It is shown that quasi-optimality is obtained under the conditions that $kh/ p$ is sufficiently small and the polynomial degree $p$ is at least $O(\log k).$
Abstract: A rigorous convergence theory for Galerkin methods for a model Helmholtz problem in $\Bbb R^d, d \in \{1,2,3\}$ is presented. General conditions on the approximation properties of the approximation space are stated that ensure quasi-optimality of the method. As an application of the general theory, a full error analysis of the classical $hp$-version of the finite element method is presented for the model problem where the dependence on the mesh width $h,$ the approximation order $p,$ and the wave number $k$ is given explicitly. In particular, it is shown that quasi-optimality is obtained under the conditions that $kh/p$ is sufficiently small and the polynomial degree $p$ is at least $O(\log k).$
202 citations
TL;DR: This paper presents weak second and third order schemes for the Cox-Ingersoll-Ross process, without any restriction on its parameters, and gives a general recursive construction method to get weak second-order schemes that extends the one introduced by Ninomiya and Victoir~\cite{NV.
Abstract: This paper presents weak second and third order schemes for the Cox-Ingersoll-Ross (CIR) process, without any restriction on its parameters. At the same time, it gives a general recursive construction method to get weak second-order schemes that extends the one introduced by Ninomiya and Victoir~\cite{NV}. Combining these both results, this allows to propose a second-order scheme for more general affine diffusions. Simulation examples are given to illustrate the convergence of these schemes on CIR and Heston models. Algorithms are stated in a pseudocode language.
194 citations
TL;DR: Two discrete functional analysis tools are established for spaces of piecewise polynomial functions on general meshes and the convergence of Discontinuous Galerkin approximations of the steady incompressible Navier--Stokes equations are proved.
Abstract: Two discrete functional analysis tools are established for spaces of piecewise polynomial functions on general meshes: (i) a discrete counterpart of the continuous Sobolev embeddings, in both Hilbertian and non-Hilbertian settings; (ii) a compactness result for bounded sequences in a suitable Discontinuous Galerkin norm, together with a weak convergence property for some discrete gradients. The proofs rely on techniques inspired by the Finite Volume literature, which differ from those commonly used in Finite Element analysis. The discrete functional analysis tools are used to prove the convergence of Discontinuous Galerkin approximations of the steady incompressible Navier--Stokes equations. Two discrete convective trilinear forms are proposed, a non-conservative one relying on Temam's device to control the kinetic energy balance and a conservative one based on a nonstandard modification of the pressure.
187 citations
TL;DR: This paper aims to demonstrate the efforts towards in-situ applicability of EMMARM, which aims to provide real-time information about the physical properties of E-modulus in the response of the immune system to shocks.
Abstract: Fujian NSF [S0750017]; National NSF of China [10531080]; Ministry of Education of China; 973 High Performance Scientific Computation Research Program [2005CB321703]
174 citations
TL;DR: It is proved that the approximation error is of order k + 1 in both the displacement and the stress, and that a postprocessed displacement approximation converging at order k - 2 can be computed element by element.
Abstract: We introduce a new mixed method for linear elasticity. The novelty is a simplicial element for the approximate stress. For every positive integer k, the row-wise divergence of the element space spans the set of polynomials of total degree k. The degrees of freedom are suited to achieve continuity of the normal stresses. What makes the element distinctive is that its dimension is the smallest required for enforcing a weak symmetry condition on the approximate stress. This is achieved using certain "bubble matrices", which are special divergence-free matrix-valued polynomials. We prove that the approximation error is of order k + 1 in both the displacement and the stress, and that a postprocessed displacement approximation converging at order k + 2 can be computed element by element. We also show that the globally coupled degrees of freedom can be reduced by hybridization to those of a displacement approximation on the element boundaries.
148 citations
TL;DR: In this article, the error of the Euler scheme applied to a stochastic partial dierential equation was studied and it was shown that the weak order of convergence is twice the strong order.
Abstract: We study the error of the Euler scheme applied to a stochastic partial dierential equation We prove that as it is often the case, the weak order of convergence is twice the strong order A key ingredient in our proof is Malliavin calculus which enables us to get rid of the irregular terms of the error We apply our method to the case a semilinear stochastic heat equation driven by a space-time white noise
123 citations
TL;DR: A general estimate is given for the condition number and least singular value of the matrix M + N n, generalizing an earlier result of Spielman and Teng for the case when x is gaussian and involves the norm ∥M∥.
Abstract: Let x be a complex random variable with mean zero and bounded variance. Let N n be the random matrix of size n whose entries are iid copies of x and let M be a fixed matrix of the same size. The goal of this paper is to give a general estimate for the condition number and least singular value of the matrix M + N n , generalizing an earlier result of Spielman and Teng for the case when x is gaussian. Our investigation reveals an interesting fact that the "core" matrix M does play a role on tail bounds for the least singular value of M + N n . This does not occur in Spielman-Teng studies when x is gaussian. Consequently, our general estimate involves the norm ∥M∥. In the special case when ∥M∥ is relatively small, this estimate is nearly optimal and extends or refines existing results.
TL;DR: The theory underlying the approach guarantees that homotopy paths lead to all isolated solutions, and this capability can be used to generate witness supersets for solution components at any dimension, the first step in computing an irreducible decomposition of the solution set of a system of polynomial equations.
Abstract: We present a new technique, based on polynomial continuation, for solving systems of n polynomials in N complex variables. The method allows equations to be introduced one-by-one or in groups, obtaining at each stage a representation of the solution set that can be extended to the next stage until finally obtaining the solution set for the entire system. At any stage where positive dimensional solution components must be found, they are sliced down to isolated points by the introduction of hyperplanes. By moving these hyperplanes, one may build up the solution set to an intermediate system in which a union of hyperplanes "regenerates" the intersection of the component with the variety of the polynomial (or system of polynomials) brought in at the next stage. The theory underlying the approach guarantees that homotopy paths lead to all isolated solutions, and this capability can be used to generate witness supersets for solution components at any dimension, the first step in computing an irreducible decomposition of the solution set of a system of polynomial equations. The method is illustrated on several challenging problems, where it proves advantageous over both the polyhedral homotopy method and the diagonal equation-by-equation method, formerly the two leading approaches to solving sparse polynomial systems by numerical continuation.
TL;DR: In this article, the authors presented a space-efficient algorithm to compute the Hilbert class polynomial H D (X) modulo a positive integer P, based on an explicit form of the Chinese Remainder Theorem.
Abstract: We present a space-efficient algorithm to compute the Hilbert class polynomial H D (X) modulo a positive integer P, based on an explicit form of the Chinese Remainder Theorem. Under the Generalized Riemann Hypothesis, the algorithm uses O(|D| 1/2+e log P) space and has an expected running time of O(|D| 1+e ). We describe practical optimizations that allow us to handle larger discriminants than other methods, with |D| as large as 10 13 and h(D) up to 10 6 . We apply these results to construct pairing-friendly elliptic curves of prime order, using the CM method.
TL;DR: This work proposes a new matrix geometric mean satisfying the ten properties given by Ando, Li and Mathias and provides a geometric interpretation and a generalization which includes as special cases the authors' mean and the Ando-Li-Mathias mean.
Abstract: We propose a new matrix geometric mean satisfying the ten properties given by Ando, Li and Mathias (Linear Alg. Appl. 2004]. This mean is the limit of a sequence which converges superlinearly with convergence of order 3 whereas the mean introduced by Ando, Li and Mathias is the limit of a sequence having order of convergence 1. This makes this new mean very easily computable. We provide a geometric interpretation and a generalization which includes as special cases our mean and the Ando-Li-Mathias mean.
TL;DR: A semilocal convergence analysis for directional Newton methods in n-variables with the following advantages: weaker convergence conditions; larger convergence domain; finer error estimates on the distances involved; and an at least as precise information on the location of the zero of the function.
Abstract: A semilocal convergence analysis for directional Newton methods in n-variables is provided in this study. Using weaker hypotheses than in the elegant related work by Y. Levin and A. Ben-Israel and introducing the center-Lipschitz condition we provide under the same computational cost as in Levin and Ben-Israel a semilocal convergence analysis with the following advantages: weaker convergence conditions; larger convergence domain; finer error estimates on the distances involved, and an at least as precise information on the location of the zero of the function. A numerical example where our results apply to solve an equation but not the ones in Levin and Ben-Israel is also provided in this study.
TL;DR: In this article, it was shown that the problem of integrating a polynomial function f over a rational simplex is NP-hard for arbitrary polynomials via a generalization of a theorem of Motzkin and Straus.
Abstract: This paper starts by settling the computational complexity of the problem of integrating a polynomial function f over a rational simplex. We prove that the problem is NP-hard for arbitrary polynomials via a generalization of a theorem of Motzkin and Straus. On the other hand, if the polynomial depends only on a fixed number of variables, while its degree and the dimension of the simplex are allowed to vary, we prove that integration can be done in polynomial time. As a consequence, for polynomials of fixed total degree, there is a polynomial time algorithm as well. We explore our algorithms with some experiments. We conclude the article with extensions to other polytopes and discussion of other available methods. 1.
TL;DR: An edge element adaptive strategy with error control is developed for wave scattering by biperiodic structures by truncating the unbounded computational domain by a perfectly matched layer technique.
Abstract: An edge element adaptive strategy with error control is developed for wave scattering by biperiodic structures. The unbounded computational domain is truncated to a bounded one by a perfectly matched layer (PML) technique. The PML parameters, such as the thickness of the layer and the medium properties, are determined through sharp a posteriori error estimates. Numerical experiments are presented to illustrate the competitive behavior of the proposed adaptive method.
TL;DR: A new method for rigorously computing smooth branches of zeros of nonlinear operators f: ℝ l 1 × B 1 → ℜ l 2 ×B 2 , where B 1 and B 2 are Banach spaces is presented.
Abstract: In this paper, we present a new method for rigorously computing smooth branches of zeros of nonlinear operators f: ℝ l 1 × B 1 → ℝ l 2 ×B 2 , where B 1 and B 2 are Banach spaces. The method is first introduced for parameter continuation and then generalized to pseudo-arclength continuation. Examples in the context of ordinary, partial and delay differential equations are given.
TL;DR: The resulting evolution equation is discretized by a Fourier spectral method with respect to the spatial variables and a modified Crank—Nicolson scheme in time.
Abstract: This paper is concerned with the analysis of a numerical algorithm for the approximate solution of a class of nonlinear evolution problems that arise as L 2 gradient flow for the Modica-Mortola regularization of the functional v ∈ BV(T d ; {-1,1}) ↦ E(v) := γ/2 ∫ Td ∇v + 1/2 ∑ k∈ℤd σ(k)v(k) 2 . Here γ is the interfacial energy per unit length or unit area, T d is the flat torus in ℝ d , and σ is a nonnegative Fourier multiplier, that is continuous on ℝ d , symmetric in the sense that σ(ξ) = σ(-ξ) for all ξ ∈ ℝ d and that decays to zero at infinity. Such functionals feature in mathematical models of pattern-formation in micromagnetics and models of diblock copolymers. The resulting evolution equation is discretized by a Fourier spectral method with respect to the spatial variables and a modified Crank—Nicolson scheme in time. Optimal-order a priori bounds are derived on the global error in the l ∞ (0, T; L 2 (T d )) norm.
TL;DR: It is proved that in parallel and simultaneously in converse to Galerkin finite element methods, under some circumstances, the weak solution is the orthogonal projection of the postprocessed mixed finite element approximation onto the H 1 0 (Ω) space and also establish several links between mixed finiteelement approximations and some generalized weak solutions.
Abstract: We derive in this paper a unified framework for a priori and a posteriori error analysis of mixed finite element discretizations of second-order elliptic problems. It is based on the classical primal weak formulation, the postprocessing of the potential proposed in [T. Arbogast and Z. Chen, On the implementation of mixed methods as nonconforming methods for second-order elliptic problems, Math. Comp. 64 (1995), 943-972], and the discrete Friedrichs inequality. Our analysis in particular avoids any explicit use of the uniform discrete inf-sup condition and in a straightforward manner and under minimal necessary assumptions, known convergence and superconvergence results are recovered. The same framework then turns out to lead to optimal a posteriori energy error bounds. In particular, estimators for all families and orders of mixed finite element methods on grids consisting of simplices or rectangular parallelepipeds are derived. They give a guaranteed and fully computable upper bound on the energy error, represent error local lower bounds, and are robust under some conditions on the diffusion-dispersion tensor. They are thus suitable for both overall error control and adaptive mesh refinement. Moreover, the developed abstract framework and a posteriori error estimates are quite general and apply to any locally conservative method. We finally prove that in parallel and simultaneously in converse to Galerkin finite element methods, under some circumstances, the weak solution is the orthogonal projection of the postprocessed mixed finite element approximation onto the H 1 0 (Ω) space and also establish several links between mixed finite element approximations and some generalized weak solutions.
TL;DR: A fractal scalar conservation law modified by a fractional power of the Laplace operator is considered, and a numerical method to approximate its solutions is proposed, showing the efficiency of the scheme and illustrating qualitative properties of the solution to the fractal conservation law.
Abstract: We consider a fractal scalar conservation law, that is to say, a conservation law modified by a fractional power of the Laplace operator, and we propose a numerical method to approximate its solutions. We make a theoretical study of the method, proving in the case of an initial data belonging to L ∞ ∩ BV that the approximate solutions converge in L ∞ weak-* and in L P strong for p < ∞, and we give numerical results showing the efficiency of the scheme and illustrating qualitative properties of the solution to the fractal conservation law.
TL;DR: In this paper, a priori energy estimates for the finite element method on shape regular grids are presented, an assumption which allows for highly graded meshes and which much more closely matches the typical practical situation.
Abstract: Local energy error estimates for the nite element method for el- liptic problems were originally proved in 1974 by Nitsche and Schatz. These estimates show that the local energy error may be bounded by a local ap- proximation term, plus a global \pollution" term that measures the inuence of solution quality from outside the domain of interest and is heuristically of higher order. However, the original analysis of Nitsche and Schatz is restricted to quasi-uniform grids. We present local a priori energy estimates that are valid on shape regular grids, an assumption which allows for highly graded meshes and which much more closely matches the typical practical situation. Our chief technical innovation is an improved superapproximation result.
TL;DR: The theoretical aspect of the Local Discontinuous Galerkin (LDG) method converges uniformly under the Shishkin mesh for singularly perturbed two-point boundary problems of the convection-diffusion type under a simplified ODE model is investigated.
Abstract: It has been observed from the authors' numerical experiments (2007) that the Local Discontinuous Galerkin (LDG) method converges uniformly under the Shishkin mesh for singularly perturbed two-point boundary problems of the convection-diffusion type. Especially when using a piecewise polynomial space of degree k, the LDG solution achieves the optimal convergence rate k+1 1 under the L 2 -norm, and a superconvergence rate 2k+1 1 for the one-sided flux uniformly with respect to the singular perturbation parameter ∈. In this paper, we investigate the theoretical aspect of this phenomenon under a simplified ODE model. In particular, we establish uniform convergence rates √∈ (1nN/N) k+1 for the L 2 -norm and ( ln N / N ) 2k+1 for the one-sided flux inside the boundary layer region. Here N (even) is the number of elements.
TL;DR: In this paper, an elementary way to distinguish between the twists of an ordinary elliptic curve E over F p in order to identify the one with p+1-2U points, when p = U 2 + dV 2 with 2U, 2V ∈ ℤ and E is constructed using the CM method, is given.
Abstract: We give an elementary way to distinguish between the twists of an ordinary elliptic curve E over F p in order to identify the one with p+1—2U points, when p = U 2 + dV 2 with 2U, 2V ∈ ℤ and E is constructed using the CM method for finding elliptic curves with a prescribed number of points. Our algorithms consist in most cases of reading off simple congruence conditions on U and V modulo 4.
TL;DR: For certain selfadjoint operators, it is proved that the rate of convergence of Hill's method to the least eigenvalue is faster than any polynomial power.
Abstract: Hill's method is a means to numerically approximate spectra of linear differential operators with periodic coefficients. In this paper, we address different issues related to the convergence of Hill's method. We show the method does not produce any spurious approximations, and that for selfadjoint operators, the method converges in a restricted sense. Furthermore, assuming convergence of an eigenvalue, we prove convergence of the associated eigenfunction approximation in the L 2 -norm. These results are not restricted to selfadjoint operators. Finally, for certain selfadjoint operators, we prove that the rate of convergence of Hill's method to the least eigenvalue is faster than any polynomial power.
TL;DR: This paper presents an efficient way of constructing infinite families of elliptic curves with given torsion group structures over cubic number fields with explicit examples of the theoretical result recently developed by the first two authors and A. Schweizer.
Abstract: In this paper we construct infinite families of elliptic curves with given torsion group structures over cubic number fields. This result provides explicit examples of the theoretical result recently developed by the first two authors and A. Schweizer; they determined all the group structures which occur infinitely often as the torsion of elliptic curves over cubic number fields. In fact, this paper presents an efficient way of constructing such families of elliptic curves with prescribed torsion group structures over cubic number fields.
TL;DR: It is shown that intermediate regularity is allowed for the Raviart-Thomas interpolation of arbitrary order when the vector field being approximated has components in W j+1,p, for triangles or tetrahedra, where 0 < j < k and 1 < p < ∞.
Abstract: We prove optimal order error estimates for the Raviart-Thomas interpolation of arbitrary order under the maximum angle condition for triangles and under two generalizations of this condition, namely, the so-called three-dimensional maximum angle condition and the regular vertex property, for tetrahedra. Our techniques are different from those used in previous papers on the subject, and the results obtained are more general in several aspects. First, intermediate regularity is allowed; that is, for the Raviart-Thomas interpolation of degree k > 0, we prove error estimates of order j + 1 when the vector field being approximated has components in W j+1,p , for triangles or tetrahedra, where 0 < j < k and 1 < p < ∞. These results are new even in the two-dimensional case. Indeed, the estimate was known only in the case j = k. On the other hand, in the three-dimensional case, results under the maximum angle condition were known only for k = 0.
TL;DR: This work considers the application of a perfectly matched layer (PML) technique to approximate solutions to the elastic wave scattering problem in the frequency domain and analyzes a Galerkin numerical approximation to the truncated PML problem to prove that it is well posed.
Abstract: We consider the application of a perfectly matched layer (PML) technique to approximate solutions to the elastic wave scattering problem in the frequency domain. The PML is viewed as a complex coordinate shift in spherical coordinates which leads to a variable complex coefficient equation for the displacement vector posed on an infinite domain (the complement of the scatterer). The rapid decay of the PML solution suggests truncation to a bounded domain with a convenient outer boundary condition and subsequent finite element approximation (for the truncated problem). We prove existence and uniqueness of the solutions to the infinite domain and truncated domain PML equations (provided that the truncated domain is sufficiently large). We also show exponential convergence of the solution of the truncated PML problem to the solution of the original scattering problem in the region of interest. We then analyze a Galerkin numerical approximation to the truncated PML problem and prove that it is well posed provided that the PML damping parameter and mesh size are small enough. Finally, computational results illustrating the efficiency of the finite element PML approximation are presented.
TL;DR: A class of hp-DG methods that is closely related to other DG schemes, however, combines both p-optimal jump penalty as well as lifting stabilization and it is proved that the resulting error estimates are optimal with respect to both the local element sizes and polynomial degrees.
Abstract: The aim of this paper is to overcome the well-known lack of p-optimality in hp-version discontinuous Galerkin (DG) discretizations for the numerical approximation of linear elliptic problems. For this purpose, we shall present and analyze a class of hp-DG methods that is closely related to other DG schemes, however, combines both p-optimal jump penalty as well as lifting stabilization. We will prove that the resulting error estimates are optimal with respect to both the local element sizes and polynomial degrees.
TL;DR: It is demonstrated that stationary statistical properties of the time discrete approximations, i.e., numerical scheme, converge to those of the underlying continuous dissipative infinite-dimensional dynamical system under three very natural assumptions as the time step approaches zero.
Abstract: We consider temporal approximation of stationary statistical properties of dissipative infinite-dimensional dynamical systems. We demonstrate that stationary statistical properties of the time discrete approximations, i.e., numerical scheme, converge to those of the underlying continuous dissipative infinite-dimensional dynamical system under three very natural assumptions as the time step approaches zero. The three conditions that are sufficient for the convergence of the stationary statistical properties are: (1) uniform dissipativity of the scheme in the sense that the union of the global attractors for the numerical approximations is pre-compact in the phase space; (2) convergence of the solutions of the numerical scheme to the solution of the continuous system on the unit time interval [0,1] uniformly with respect to initial data from the union of the global attractors; (3) uniform continuity of the solutions to the continuous dynamical system on the unit time interval [0,1] uniformly for initial data from the union of the global attractors. The convergence of the global attractors is established under weaker assumptions. An application to the infinite Prandtl number model for convection is discussed.