Showing papers in "Mathematics of Computation in 1982"

Scattered data interpolation: tests of some methods

Abstract: Absract. This paper is concerned with the evaluation of methods for scattered data interpolation and some of the results of the tests when applied to a number of methods. The process involves evaluation of the methods in terms of timing, storage, accuracy, visual pleasantness of the surface, and ease of implementation. To indicate the flavor of the type of results obtained, we give a summary table and representative perspective plots of several surfaces.

Upwind difference schemes for hyperbolic systems of conservation laws

Abstract: We derive a new upwind finite difference approximation to systems of nonlinear hyperbolic conservation laws. The scheme has desirable properties for shock calculations. Under fairly general hypotheses we prove that limit solutions satisfy the entropy condition and that discrete steady shocks exist which are unique and sharp. Numerical examples involving the Euler and Lagrange equations of compressible gas dynamics in one and two space dimensions are given.

Approximation results for orthogonal polynomials in Sobolev spaces

Abstract: We analyze the approximation properties of some interpolation operators and some L2-orthogonal projection operators related to systems of polynomials which are orthonormal with respect to a weight function o(x1, . . ., Xd), d > 1. The error estimates for the Legendre system and the Chebyshev system of the first kind are given in the norms of the Sobolev spaces H'. These results are useful in the numerical analysis of the approximation of partial differential equations by spectral methods. 0. Introduction. Spectral methods are a classical and largely used technique to solve differential equations, both theoretically and numerically. During the years they have gained new popularity in automatic computations for a wide class of physical problems (for instance in the fields of fluid and gas dynamics), due to the use of the Fast Fourier Transform algorithm. These methods appear to be competitive with finite difference and finite element methods and they must be decisively preferred to the last ones whenever the solution is highly regular and the geometric dimension of the domain becomes large. Moreover, by these methods it is possible to control easily the solution (filtering) of those numerical problems affected by oscillation and instability phenomena. The use of spectral and pseudo-spectral methods in computations in many fields of engineering has been matched by deeper theoretical studies; let us recall here the pioneering works by Orszag [25], [26], Kreiss and Oliger [14] and the monograph by Gottlieb and Orszag [13]. The theoretical results of such works are mainly concerned with the study of the stability of approximation of parabolic and hyperbolic equations; the solution is assumed to be infinitely differentiable, so that by an analysis of the Fourier coefficients an infinite order of convergence can be achieved. More recently (see Pasciak [27], Canuto and Quarteroni [10], [11], Maday and Quarteroni [20], [211, [22], Mercier [23]), the spectral methods have been studied by the variational techniques typical of functional analysis, to point out the dependence of the approximation error (for instance in the L2-norm, or in the energy norm) on the regularity of the solution of continuous problems and on the discretization parameter (the dimension of the space in which the approximate solution is sought). Indeed, often the solution is not infinitely differentiable; on the other hand, sometimes even if the solution is smooth, its derivatives may have very Received August 9, 1980; revised June 12, 1981. 1980 Mathematics Subject Classification. Primary 41A25; Secondary 41A 10, 41A05. ? 1982 American Mathematical Society 0025-571 8/82/0000-0470/$06.00 (67 This content downloaded from 207.46.13.111 on Tue, 09 Aug 2016 06:29:39 UTC All use subject to http://about.jstor.org/terms 68 C. CANUTO AND A. QUARTERONI large norms which affect negatively the rate of convergence (for instance in problems with boundary layers). Both spectral and pseudo-spectral methods are essentially Ritz-Galerkin methods (combined with some integration formulae in the pseudo-spectral case). It is well known that when Galerkin methods are used the distance between the exact and the discrete solution (approximation error) is bounded by the distance between the exact solution and its orthogonal projection upon the subspace (projection error), or by the distance between the exact solution and its interpolated polynomial at some suitable points (interpolation error). This upper bound is often realistic, in the sense that the asymptotic behavior of the approximation error is not better than the one of the projection (or even the interpolation) error. Even more, in some cases the approximate solution coincides with the projection of the true solution upon the subspace (for instance when linear problems with constant coefficients are approximated by spectral methods). This motivates the interest in evaluating the projection and the interpolation errors in differently weighted Sobolev norms. So we must face a situation different from the one of the classical approximation theory where the properties of approximation of orthogonal function systems, polynomial and trigonometric, are studied in the LP-norms, and mostly in the maximum norm (see, e.g., Butzer and Berens [6], Butzer and Nessel [7], Nikol'skiT [24], Sansone [291, Szego [30], Triebel [31], Zygmund [32]; see also Bube [5]). Approximation results in Sobolev norms for the trigonometric system have been obtained by Kreiss and Oliger [15]. In this paper we consider the systems of Legendre orthogonal polynomials, and of Chebyshev orthogonal polynomials of the first kind in dimension d > 1. The reason for this interest must be sought in the applications to spectral approximations of boundary value problems. Indeed, if the boundary conditions are not periodic, Legendre approximation seems to be the easiest to be investigated (the weight w is equal to 1). On the other hand, the Chebyshev approximation is the most effective for practical computations since it allows the use of the Fast Fourier Transform algorithm. The techniques used to obtain our results are based on the representation of a function in the terms of a series of orthogonal polynomials, on the use of the so-called inverse inequality, and finally on the operator interpolation theory in Banach spaces. For the theory of interpolation we refer for instance to Calderon [8], Lions [17], Lions and Peetre [19], Peetre [28]; a recent survey is given, e.g., by Bergh and Lofstrom [4]. An outline of the paper is as follows. In Section 1 some approximation results for the trigonometric system are recalled; the presentation of the results to the interpolation is made in the spirit of what will be its application to Chebyshev polynomials. In Section 2 we consider the La-projection operator upon the space of polynomials of degree at most N in any variable (w denotes the Chebyshev or Legendre weight). In Section 3 a general interpolation operator, built up starting by integration formulas which are not necessarily the same in different spatial dimensions, is considered, and its approximation properties are studied. In [22] Maday and Quarteroni use the results of Section 2 to study the approximation properties of some projection operators in higher order Sobolev norms. Recently, an interesting method which lies inbetween finite elements and This content downloaded from 207.46.13.111 on Tue, 09 Aug 2016 06:29:39 UTC All use subject to http://about.jstor.org/terms ORTHOGONAL POLYNOMIALS IN SOBOLEV SPACES 69 spectral methods has been investigated from the theoretical point of view by Babuska, Szabo and Katz [3]. In particular they obtain approximation properties of polynomials in the norms of the usual Sobolev spaces. Acknowledgements. Some of the results of this paper were announced in [9]; we thank Professor J. L. Lions for the presentation to the C. R. Acad. Sci. of Paris. We also wish to express our gratitude to Professors F. Brezzi and P. A. Raviart for helpful suggestions and continuous encouragement. Notations. Throughout this paper we shall use the following notations: I will be an open bounded interval c R, whose variable is denoted by x; Q the product Id C Rd (d integer > 1) whose variable is denoted by x = (x(.')I_ d; for a multi-integer k E Zd, we set ikV = jd X I'12 and IkloK = m x 1, Dj = a/ax@). The symbol X'J=p (q eventually + oo) will denote the summation over all integral k such that p 0 in U. Set L2(Q) = ({: Q -C I 0 is measurable and ( 0, set Hs ( () = C E L(Q) I 1111ksI, < +?}, where /d 2 11I412I= kENd f DI L/)4 D w dx.

Some Optimal Error Estimates for Piecewise Linear Finite Element Approximations

Abstract: It is shown that the Ritz projection onto spaces of piecewise linear finite elements is bounded in the Sobolev space, Wl, for 2 - p < oc. This implies that for functions in wf n W2 the error in approximation behaves like 0(h) in Wp, for 2 ? p < oo, and like 0(h2) in Lp, for 2 - p < oo. In all these cases the additional logarithmic factor previously included in error estimates for linear finite elements does not occur.

Vortex methods. II. Higher order accuracy in two and three dimensions

Abstract: In an earlier paper the authors introduced a new version of the vortex method for three-dimensional, incompressible flows and proved that it converges to arbitrarily high order accuracy, provided we assume the consistency of a discrete approximation to the Biot-Savart Law. We prove this consistency statement here, and also derive substantially sharper results for two-dimensional flows. A complete, simplified proof of convergence in two dimensions is included.

Procedures for computing one- and two-dimensional integrals of functions with rapid irregular oscillations

Abstract: A collocation procedure for efficient integration of rapidly oscillatory functions is presented. The integration problem is transformed into a certain O.D.E. problem, and this is solved by a collocation technique. The method is also extended to two-dimensional integration, and some numerical results are appended showing the efficiency of the method in handling difficult cases of rapid irregular oscillations. 1. The Procedure for One-Dimensional Integrals. We consider integrals of the form a (1.1)~~~~~~~~ 1 f (x)e iq(x) dx, where f is smooth and "nonoscillatory" and Iq'(x) I? (b a)-1. Two practical methods for evaluating rapidly oscillatory integrals are described in [1], the use of approximation as in Filon's method [3], [4] and the speedup method of Longman [5]. Formally both methods are applicable to any integral of the form (1.1), but their best performance is for the case of a constant frequency q' = W. In this note we present an efficient method which is applicable for cases of varying frequency q' using only a small number of values of f and q' in [a, b] and the values q(a) and q(b). The proposed method follows the spirit of Filon's method. It is based upon the fact that if f were of the form (1.2) f(x) = iq'(x)p(x) + p'(x) =_ L(l1p(x), a < x < b, then the integral could be evaluated directly as (1.3 J b (iq'(x)p(x) + pt(x))eiq(x) dx f | d (p(x)eiq(x)) dx = (b)e qbp(a)eiqa Equation (1.2) can be considered as a differential equation for p(x), and any solution of this equation can be used in (1.3) for evaluating I. The general solution of this equation is (1.4) p(x) e-iq(x)[ f(t)eiq(t) dt + cj Received March 24, 1981; revised July 29, 1981. 1980 Mathematics Subject Classification. Primary 65D30; Secondary 65D32. ?1982 American Mathematical Society 0025-5718/8 1 /OOO-0 1 23/$03.25 531 This content downloaded from 157.55.39.104 on Sun, 19 Jun 2016 06:36:01 UTC All use subject to http://about.jstor.org/terms

Vortex methods. I. Convergence in three dimensions

Abstract: Recently several different approaches have been developed for the simulation of three-dimensional incompressible fluid flows using vortex methods. Some versions use detailed tracking of vortex filament structures and often local curvatures of these filaments, while other methods require only crude information, such as the vortex blobs of the two-dimensional case. Can such "crude" algorithms accurately account for vortex stretching and converge? We answer this question affirmatively by constructing a new class of "crude" three-dimensional vortex methods and then proving that these methods are stable and convergent, and can even have arbitrarily high order accuracy without being more expensive than other "crude" versions of the vortex algorithm.

A finite element method for solving Helmholtz type equations in waveguides and other unbounded domains

Abstract: A finite element method is described for solving Helmholtz type boundary value problems in unbounded regions, including those with infinite boundaries. Typical examples include the propagation of acoustic or electromagnetic waves in waveguides. The radiation condition at infinity is based on separation of variables and differs from the classical Sommerfeld radiation condition. It is shown that the problem may be replaced by a boundary value problem on a fixed bounded domain. The behavior of the solution near infinity is incorporated in a nonlocal boundary condition. This problem is given a weak or variational formulation, and the finite element method is then applied. It is proved that optimal error estimates hold.

A $p+1$ method of factoring

Abstract: Let N have a prime divisor p such that p + 1 has only small prime divisors. A method is described which will allow for the determination of p, given N. This method is analogous to the p — 1 method of factoring which was described in 1974 by Pollard. The results of testing this method on a large number of composite numbers are also presented.

Numerical Comparisons of Nonlinear Convergence Accelerators

Abstract: As part of a continuing program of numerical tests of convergence accelerators, the iterated Aitken's Delta-squared method, Wynn's epsilon algorithm, Brezinski's theta algorithm, and Levin's u transform are compared on a broad range of test problems: linearly convergence alternating, monotone, and irregular-sign series, logarithmically convergent series, power method and Bernoulli method sequences, alternating and monotone asymptotic series, and some perturbation series arising in applications. In each category either the epsilon algorithm or the u transform gives the best results of the four methods tested. In some cases differences among methods are slight, and in others they are quite striking.

Analysis of some mixed finite element methods related to reduced integration

Abstract: We prove error estimates for the following two mixed finite element methods related to reduced integration: A method for Stokes' problem using rectangular elements with piecewise bilinear approximations for the velocities and piecewise constants for the pressure, and one method for a plate problem using bilinear approximations for transversal displacement and rotations and piecewise constants for the shear stress. The main idea of the proof in the case of Stokes' problem is to combine a weak Babuska-Brezzi type stability estimate for the pressure with a superapproximability property for the velocities. A similar technique is used in the case of the plate problem.

On a higher order accurate fully discrete Galerkin approximation to the Navier-Stokes equations

Abstract: We consider approximating the solution of the initial and boundary value problem for the Navier-Stokes equations in bounded twoand three-dimensional domains using a nonstandard Galerkin (finite element) method for the space discretization and the third order accurate, three-step backward differentiation method (coupled with extrapolation for the nonlinear terms) for the time stepping. The resulting scheme requires the solution of one linear system per time step plus the solution of five linear systems for the computation of the required initial conditions; all these linear systems have the same matrix. The resulting approximations of the velocity are shown to have optimal rate of convergence in L2 under suitable restrictions on the discretization parameters of the problem and the size of the solution in an appropriate function space.

Analysis of a multilevel iterative method for nonlinear finite element equations

Abstract: The multilevel iterative technique is a powerful technique for solving the systems of equations associated with discretized partial differential equations. We describe how this technique can be combined with a globally convergent approximate Newton method to solve nonlinear partial differential equations. We show that asymptotically only one Newton iteration per level is required; thus the complexity for linear and nonlinear problems is essentially equal.

On estimating the largest eigenvalue with the Lanczos algorithm

Abstract: The Lanczos algorithm applied to a positive definite matrix produces good approximations to the eigenvalues at the extreme ends of the spectrum after a few iterations. In this note we utilize this behavior and develop a simple algorithm which computes the largest eigenvalue. The algorithm is especially economical if the order of the matrix is large and the accuracy requirements are low. The phenomenon of misconvergence is discussed. Some simple extensions of the algorithm are also indicated. Finally, some numerical examples and a comparison with the power method are given.

The numerical evaluation of very oscillatory infinite integrals by extrapolation

Abstract: Recently the author has given two modifications of a nonlinear extrapolation method due to Levin and Sidi, which enable one to accurately and economically compute certain infinite integrals whose integrands have a simple oscillatory behavior at infinity. In this work these modifications are extended to cover the case of very oscillatory infinite integrals whose integrands have a complicated and increasingly rapid oscillatory behavior at infinity. The new method is applied to a number of complicated integrals, among them the solution to a problem in viscoelasticity. Some convergence results for this method are presented.

Mesh modification for evolution equations

Abstract: Finite element methods for which the underlying function spaces change with time are studied. The error estimates produced are all in norms that are very naturally associated with the problems. In some cases the Galerkin solution error can be seen to be quasi-optimal. K. Miller's moving finite element method is studied in one space dimension; convergence is proved for the case of smooth solutions of parabolic problems. Most, but not all, of the analysis is done on linear problems. Although second order parabolic equations are emphasized, there is also some work on first order hyperbolic and Sobolev equations.

On the number of Markoff numbers below a given bound

Abstract: According to a famous theorem of Markoff, the indefinite quadratic forms with exceptionally large minima (greater than f of the square root of the discriminant) are in 1 : 1 correspondence with the solutions of the Diophantine equation p2 + q2 + r1 = ~ipqr. By relating Markoffs algorithm for finding solutions of this equation to a problem of count- ing lattice points in triangles, it is shown that the number of solutions less than x equals Clog2 3x + 0(log x log log2 x) with an explicitly computable constant C = 0.18071704711507.... Numerical data up to 101300 is presented which suggests that the true error term is considerably smaller.

Class number computations of real abelian number fields

Abstract: In this paper we describe the calculation of the class numbers of most real abelian number fields of conductor ? 200. The technique is due to J. M. Masley and makes use of discriminant bounds of A. M. Odlyzko. In several cases we have to assume the generalized Riemann hypothesis. Introduction. It is well known that the class number h of an abelian number field can be written as h = h h -, where h ? is the class number of the maximal real subfield K? of K and h is an integer. We can determine the relative class number h in a straightforward way, using the complex analytic class number formula (see [7, Kap. III], or [9, Chapter 3, Section 3]). For the full cyclotomic fields Q(n), with 4O(n) < 256, and their subfields, one can deduce hfrom the tables of G. Schrutka von Rechtenstamm [15]; here Dn denotes a primitive nth root of unity, and 40 is the Euler function. For the class number factor h + the complex analytic class number formula is less useful, since it requires that the units of K? be known. Alternative techniques have been developed by J. M. Masley [13], who computed the class number of almost all real cyclic number fields of conductor < 100; here the conductor of K is the least f for which K C Q(tf ). In this paper we apply Masley's techniques, with a few additions, to determine the class numbers of a large collection of real abelian fields of conductor < 200; see Section 1 for a precise statement of our results, some of which assume the generalized Riemann hypothesis. An important ingredient of Masley's method is the use of discriminant lower bounds proved by A. M. Odlyzko [14]. These lead to an upper bound for the class number of a real abelian number field, provided that its conductor, or more precisely its root discriminant (see [13, Section 1]), is sufficiently small. It follows that this method can only be used for a finite number of real abelian number fields. The existence of infinite class field towers shows that this remains true after any future improvement of Odlyzko's bounds. In fact, examples of J. Martinet [12] show that the method will never apply to fields whose root discriminant is larger than five times the present bound, under assumption of the generalized Riemann hypothesis. The structure of this paper is as follows. Section 1 contains our results and Section 2 lists the theorems used in the proofs. The proofs themselves are largely suppressed. Received October 24, 1980; revised December 24, 1981. 1980 Mathematics Subject Classification. Primary 12-04, 12A35, 12A55.

Error estimates for the multidimensional two-phase Stefan problem

Abstract: In this paper we derive rates of convergence for regularizations of the multidimen- sional two-phase Stefan problem and use the regularized problems to define backward-dif- ference in time and C° piecewise-linear in space Galerkin approximations. We find an L2 rate of convergence of order \^ in the e-regularization and an L rate of convergence of order (h2/e + Ai/ \6F) in the Galerkin estimates which leads to the natural choices £ ~ A4/3, A; ~ A4/3, and a resulting 0(/i2//3) L2 rate of convergence of the numerical scheme to the solution of the differential equation. An essentially 0(h) rate is demonstrated when e = 0 and A/ ~ h2 in our Galerkin scheme under a boundedness hypothesis on the Galerkin approxima- tions. The latter result is consistent with computational experience. 1. Introduction. Given a smoothly bounded domain s C RN, we consider the equation, in distribution form, (l.H) ^ "Au +fiu) = 0

The sequence of radii of the Apollonian packing

Abstract: We consider the distribution function N(x) of the curvatures of the disks in the Apollonian packing of a curvilinear triangle. That is, N(x) counts the number of disks in the packing whose curvatures do not exceed x. We show that log N( x )/log x approaches the limit S as x tends to infinity, where S is the exponent of the packing. A numerical fit of a curve of the form y = Ans to the values of N-(1000n) for n 1, 2,..., 6400 produces the estimate S 1.305636 which is consistent with the known bounds 1.300197 < S < 1.314534.

Absorbing boundary conditions for the discretization schemes of the one-dimensional wave equation

Abstract: When computing a partial differential equation, it is often necessary to introduce artificial boundaries. Here we explain a systematic method to obtain boundary conditions for the wave equation in one dimension, fitting to the discretization scheme and stable. Moreover, we give error estimates on the reflected part.

Runge-Kutta theory for Volterra integral equations of the second kind

Abstract: The present paper develops the theory of general Runge-Kutta methods for Volterra integral equations of the second kind. The order conditions are derived by using the theory of P-series, which for our problem reduces to the theory of V-series. These results are then applied to two special classes of Runge-Kutta methods introduced by Pouzet and by Bel'tyukov.