scispace - formally typeset
Search or ask a question

Showing papers in "SIAM Journal on Numerical Analysis in 1983"


Journal ArticleDOI
TL;DR: A class of iterative algorithms for solving systems of linear equations where the coefficient matrix is nonsymmetric with positive-definite symmetric part, modelled after the conjugate gradient method, are considered.
Abstract: We consider a class of iterative algorithms for solving systems of linear equations where the coefficient matrix is nonsymmetric with positive-definite symmetric part. The algorithms are modelled after the conjugate gradient method, and are well suited for large sparse systems. They do not make use of any associated symmetric problems. Convergence results and error bounds are presented.

919 citations


Journal ArticleDOI
TL;DR: It is shown in this paper that an approximate solution of the trust region problem may be found by the preconditioned conjugate gradient method, and it is shown that the method has the same convergence properties as existing methods based on the dogleg strategy using an approximate Hessian.
Abstract: Algorithms based on trust regions have been shown to be robust methods for unconstrained optimization problems. All existing methods, either based on the dogleg strategy or Hebden-More iterations, require solution of system of linear equations. In large scale optimization this may be prohibitively expensive. It is shown in this paper that an approximate solution of the trust region problem may be found by the preconditioned conjugate gradient method. This may be regarded as a generalized dogleg technique where we asymptotically take the inexact quasi-Newton step. We also show that we have the same convergence properties as existing methods based on the dogleg strategy using an approximate Hessian.

829 citations


Journal ArticleDOI
TL;DR: Statically and kinematically admissible fields are explicitly derived from the finite element solution of the primal form of linear models as mentioned in this paper, and the contribution of each element to this error allows to implement an automatic mesh refinement procedure leading to a uniform distribution of a given accuracy.
Abstract: Statically and kinematically admissible fields are explicitly derived from the finite element solution of the primal form of linear models The error on constitutive law for these fields yields an expression of the finite element error Moreover, the contribution of each element to this error allows to implement an automatic mesh refinement procedure leading to a uniform distribution of a given accuracy

670 citations


Journal ArticleDOI
TL;DR: It is proved that two of the algorithms investigated are optimal for band graphs, andumerical evidence indicates that these two algorithms are nearly optimal on practical problems.
Abstract: Given a mapping with a sparse Jacobian matrix, the problem of minimizing the number of function evaluations needed to estimate the Jacobian matrix by differences is investigated. This problem can be attacked as a graph coloring problem and this approach leads to very efficient algorithms. The behavior of these algorithms is studied and, in particular, it is proved that two of the algorithms are optimal for band graphs. Numerical evidence is presented which indicates that these two algorithms are nearly optimal on practical problems.

467 citations


Journal ArticleDOI
TL;DR: The notion of a generalized finite element method is introduced and this class of methods is analyzed and their relation to mixed methods is discussed.
Abstract: The notion of a generalized finite element method is introduced. This class of methods is analyzed and their relation to mixed methods is discussed. The class of generalized finite element methods offers a wide variety of computational procedures from which particular procedures can be selected for particular problems. A particular generalized finite element method which is very effective for problems with rough coefficients is discussed in detail.

426 citations


Journal ArticleDOI
TL;DR: The design of algorithms for interpolating discrete data using $C^1 $-quadratic splines in such a way that the monotonicity and/or convexity of the data is preserved is discussed.
Abstract: In this paper we discuss the design of algorithms for interpolating discrete data using $C^1 $-quadratic splines in such a way that the monotonicity and/or convexity of the data is preserved. The analysis culminates in an interactive algorithm which takes full advantage of the flexibility which quadratic splines permit.

255 citations


Journal ArticleDOI
TL;DR: The presented proof applies to procedures with any number of smoothing iterations and to the V-cycle and proves convergence under natural assumptions on the discretization and the elliptic problem.
Abstract: For a positive definite finite element equation we describe a multigrid iteration and prove convergence under natural assumptions on the discretization and the elliptic problem. Hitherto existing convergence proofs require a sufficiently large number of smoothing iterations and exclude the “V-cycle”. The presented proof applies to procedures with any number of smoothing iterations and to the V-cycle.

225 citations


Journal ArticleDOI
TL;DR: In this article, the incomplete inverse is used to form a preconditioning matrix whose inverse is a polynomial in G with unit coefficients, and then the parameters are chosen to minimize the condition number of the product of the inverse and the pseudo-Neumann series.
Abstract: Dubois, Greenbaum and Rodrigue proposed using a truncated Neumann series as an approximation to the inverse of a matrix A for the purpose of preconditioning conjugate gradient iterative approximations to $Ax = b$. If we assume that A has been symmetrically scaled to have unit diagonal and is thus of the form $(I - G)$, then the Neumann series is a power series in G with unit coefficients. The incomplete inverse was thought of as a replacement of the incomplete Cholesky decomposition suggested by Meijerink and van der Vorst in the family of methods ICCG $(n)$. The motivation for the replacement was the desire to have a preconditioned conjugate gradient method which only involved vector operations and which utilized long vectors.We here suggest parameterizing the incomplete inverse to form a preconditioning matrix whose inverse is a polynomial in G. We then show how to select the parameters to minimize the condition number of the product of the polynomial and $(I - G)$. Theoretically the resulting algorithm...

205 citations


Journal ArticleDOI
TL;DR: In this article, the authors developed a general framework for calculating the eigenvalues of a symmetric matrix using ordinary differential equations and suggested new algorithms and old algorithms, including $QR$, are interpreted.
Abstract: In this paper the authors develop a general framework for calculating the eigenvalues of a symmetric matrix using ordinary differential equations. New algorithms are suggested and old algorithms, including $QR$, are interpreted.

183 citations


Journal ArticleDOI
TL;DR: In this article, a finite difference procedure that reflects the dominance of convection in incompressible flow in porous media is developed. But this method is not suitable for the case of two-phase, incompressibly flow.
Abstract: Two-phase, incompressible flow in porous media is governed by a system of nonlinear partial differential equations. Convection physically dominates diffusion, and the object of this paper is to develop a finite difference procedure that reflects this dominance. The pressure equation, which is elliptic in appearance, is discretized by a standard five-point difference method. The concentration equation is treated by an implicit finite difference method that applies a form of the method of characteristics to the transport terms. A convergence analysis is given for the method.

163 citations


Journal ArticleDOI
TL;DR: In extensive computational tests, a tensor algorithm is significantly more efficient than a similar algorithm based on the standard linear model, both on standard nonsingular test problems and on problems where the Jacobian at the solution is singular.
Abstract: A new class of methods for solving systems of nonlinear equations, called tensor methods, is introduced. Tensor methods are general purpose methods intended especially for problems where the Jacobian matrix at the solution is singular or ill-conditioned. They base each iteration on a quadratic model of the nonlinear function, the standard linear model augmented by a simple second order term. The second order term is selected so that the model interpolates function values from several previous iterations, as well as the current function value and Jacobian. The tensor method requires no more function and derivative information per iteration and hardly more storage or arithmetic per iteration, than a standard method based on Newton’s method. In extensive computational tests, a tensor algorithm is significantly more efficient than a similar algorithm based on the standard linear model, both on standard nonsingular test problems and on problems where the Jacobian at the solution is singular.

Journal ArticleDOI
TL;DR: In this paper, a finite dimensional stability test for checking velocity/pressure finite element trial spaces is presented and applications are made to a new class of element pairs proposed in this paper as well as to existing spaces.
Abstract: A finite dimensional stability test for checking velocity/pressure finite element trial spaces is presented. Applications are made to a new class of element pairs proposed in this paper as well as to existing spaces.

Journal ArticleDOI
TL;DR: It is pointed out that the kernel of the error functional, at any complex point outside the interval of integration, can be evaluated accurately and efficiently by a recursive procedure.
Abstract: For Gaussian quadrature rules over a finite interval, applied to analytic or meromorphic functions, we develop error bounds from contour integral representations of the remainder term. As in previous work on the subject, we consider both circular and elliptic contours. In contrast with earlier work, however, we attempt to determine exactly where on the contour the kernel of the error functional attains its maximum modulus. We succeed in answering this question for a large class of weight distributions (including all Jacobi weights) when the contour is a circle. In the more difficult case of elliptic contours, we can settle the question for certain special Jacobi weight distributions with parameters $ \pm \frac{1} {2}$, and we provide empirical results for more general Jacobi weights. We further point out that the kernel of the error functional, at any complex point outside the interval of integration, can be evaluated accurately and efficiently by a recursive procedure. The same procedure is useful also t...

Journal ArticleDOI
TL;DR: In this paper, the condition number of a nonsingular matrix in the 2-norm has been estimated with probability at least 1 − ε - ε ≥ 0 using probabilistic algorithms.
Abstract: The paper describes “probabilistic” algorithms which may be used to make rough estimates of the largest and smallest eigenvalues of a positive definite matrix and the condition number of a nonsingular matrix in the 2-norm. Given $\varepsilon > 0$ and a prescribed relative error, the algorithms compute estimates which, with probability at least $1 - \varepsilon $, have relative errors less than that prescribed. In particular, the method gives a reliable way to estimate the condition number of a matrix of large degree.

Journal ArticleDOI
TL;DR: In this article, it was shown that the pseudospectral Chebyshev second derivative operator with separated boundary conditions is real, negative and distinct, and an application to the full potential equation is also presented.
Abstract: It is shown that the eigenvalues of the pseudospectral Chebyshev second derivative operator with separated boundary conditions are real, negative and distinct. An application to the full potential equation is also presented.

Journal ArticleDOI
TL;DR: In this article, the collocation method for integral equations of the second kind is surveyed and analyzed for the case in which the approximate solutions are only piecewise continuous, and a satisfactory sense of point evaluation is given for elements of the integral equation.
Abstract: The collocation method for integral equations of the second kind is surveyed and analyzed for the case in which the approximate solutions are only piecewise continuous. Difficulties with the usual function space setting of $L_\infty (D)$ are discussed, and a satisfactory sense of point evaluation is given for elements of $L_\infty (D)$. Other approaches which are discussed include (i) reformulation as a degenerate kernel method, (ii) the prolongation-restriction framework of Noble, (iii) other function space settings, and (iv) reformulation as a continuous approximation problem by iterating the piecewise continuous approximate solution in the original integral equation.

Journal ArticleDOI
TL;DR: The method is related to inverse iteration and to Newton's method applied to the eigenvalue problem and is analogous to iterative improvement for the solution of linear systems.
Abstract: This paper describes a computational method for improving the accuracy of a given eigenvalue and its associated eigenvector. The method is analogous to iterative improvement for the solution of linear systems. An iterative algorithm using working precision arithmetic is applied to increase the accuracy of the eigenpair. The only extended precision computation is the residual calculation. The method is related to inverse iteration and to Newton's method applied to the eigenvalue problem.

Journal ArticleDOI
TL;DR: In this paper, a special nonpolynomial spline space is introduced to model the structure of these solutions near the point of nonsmooth behavior, where collocation in these spaces will once more lead to high-order methods.
Abstract: Volterra integral equations of the second kind with weakly singular kernels possess, in general, solutions which are not smooth near the left endpoint of the interval of integration. Since ordinary polynomial spline collocation cannot lead to high-order convergence we introduce special nonpolynomial spline spaces which are modelled after the structure of these solutions near the point of nonsmooth behavior; collocation in these spaces will once more lead to high-order methods. Analogous results are derived for Volterra integro-differential equations with weakly singular kernels.

Journal ArticleDOI
TL;DR: The conjugate gradient method is shown to have a superlinear rate of convergence when applied to this formulation of the problem and is an order of magnitude more efficient than previously known methods with the possible exception of multi-grid.
Abstract: A new method for the numerical solution of the first biharmonic Dirichlet problem in a rectangular domain is presented. For an $N \times N$ mesh the complexity of this algorithm is on the order of ...

Journal ArticleDOI
TL;DR: In this paper, the authors give a detailed analysis of the linear convergence rates for several types of singular problems and describe modifications of Newton's method which will restore quadratic convergence.
Abstract: If Newton’s method is employed to find a root of a map from a Banach space into itself and the derivative is singular at that root, the convergence of the Newton iterates to the root is linear rather than quadratic. In this paper we give a detailed analysis of the linear convergence rates for several types of singular problems. For some of these problems we describe modifications of Newton’s method which will restore quadratic convergence.

Journal ArticleDOI
TL;DR: New Runge–Kutta algorithms are developed which determine the solution of a system of ordinary differential equations at any point within a given integration step, as well as at the end of each step.
Abstract: New Runge–Kutta algorithms are developed which determine the solution of a system of ordinary differential equations at any point within a given integration step, as well as at the end of each step...

Journal ArticleDOI
TL;DR: In this article, Newton's method is analyzed in the neighborhood of irregular singularities, which include all minimizers at which the Hessian has a one-dimensional null space, and it is shown that depending on a parameter specific to the underlying optimization problem or system of simultaneous equations, the method is either to converge with a limiting ratio of about $\frac{2} {3}$, or to diverge from arbitrarily close starting points, or to behave in a certain sense chaotically.
Abstract: Recent work on Newton’s method at singularities of the Jacobian has established linear convergence under certain regularity assumptions. Here, Newton’s method is analyzed in the neighborhood of irregular singularities which include all minimizers at which the Hessian has a one-dimensional null space. Depending on a parameter specific to the underlying optimization problem or system of simultaneous equations, Newton’s method is found either to converge with a limiting ratio of about $\frac{2} {3}$, or to diverge from arbitrarily close starting points, or to behave in a certain sense chaotically.

Journal ArticleDOI
TL;DR: In this paper, an imbedded family of fully symmetric rules for numerical integration over an n-dimensional hypercube was studied in theory and shown to exist in practice for arbitrary degree.
Abstract: We consider an imbedded family of $m + 1$ generator $2m + 1$ degree fully symmetric rules for numerical integration over an n-dimensional hypercube. The rules are shown to exist in theory for arbitrary degree. In practice, the rules are most economical and useful for $3 \leqq m \leqq 6$, $3 \leqq n \leqq 3m$.

Journal ArticleDOI
TL;DR: A “Caratheodory–Fejer method” is presented for near-best real rational approximation on intervals, based on the eigenvalue (or singular value) analysis of a Hankel matrix of Chebyshev coefficients.
Abstract: A “Caratheodory–Fejer method” is presented for near-best real rational approximation on intervals, based on the eigenvalue (or singular value) analysis of a Hankel matrix of Chebyshev coefficients. In approximation of a smooth function F, the CF approximant $R^{cf} $ frequently differs from the best approximation $R^ * $ by only one part in millions or billions. To account for this we show here under weak assumptions that if F is approximated on $[ - \varepsilon ,\varepsilon ]$, then as $\varepsilon \to 0$, $||F - R^ * || = O(\varepsilon ^{m + n + 1} )$ while $||R^{cf} - R^ * || = O(\varepsilon ^{3m + 2n + 3} )$. In contrast, the latter figure would be $O(\varepsilon ^{m + n + 2} )$ for the Chebyshev economization approximant of Maehly or the Chebyshev–Pade approximant of Gragg. It follows that as $\varepsilon \to 0$, best approximation error curves approach the real parts of $m + n + 1$-winding rational functions of constant modulus to within $O(\varepsilon ^{3m + 2n + 3} )$. Numerical examples are given...

Journal ArticleDOI
TL;DR: In this paper, it was shown that the best condition matrix SR has a condition number approximately equal to the cosecant of the smallest angle between right subspaces belonging to different diagonal blocks of Θ.
Abstract: How ill-conditioned must a matrix S be if its columns are constrained to span certain subspaces? We answer this question in order to find nearly best conditioned matrices SR and SL that block diagonalize a given matrix pencil T=A+λB, i.e. S L −1 TSR=Θ is bloc diagonal. We show that the best conditioned SR has a condition number approximately equal to the cosecant of the smallest angle between right subspaces belonging to different diagonal blocks of Θ. Thus, the more nearly the right subspaces overlap the more ill-conditioned SR must be. The same is true of SL and the left subspaces. For the standard eigenproblem (T=A−λI), SL = SR and the cosecant of the angle between subspaces turns out equal to an earlier estimate of the smallest condition number, namely the norm of the projection matrix associated with one of the subspaces. We apply this result to bound the error in an algorithm to compute analytic functions of matrices, for instance exp(T).

Journal ArticleDOI
TL;DR: This work proves an error estimate for a fully discrete method for the numerical solution of a two-dimensional model problem in neutron transport theory based on using the discrete ordinates method forThe angular variable and the discontinuous Galerkin finite element method with piecewise linear trial functions for the space variable.
Abstract: We prove an error estimate for a fully discrete method for the numerical solution of a two-dimensional model problem in neutron transport theory based on using the discrete ordinates method for the...

Journal ArticleDOI
TL;DR: In this article, the authors present a collection of stability results for finite difference approximations to the advection-diffusion equation, which are derived from a uniform framework based on the Schur-Cohn theory of Simple von Neumann Polynomials and are necessary and sufficient for the stability of Cauchy problem.
Abstract: We present a collection of stability results for finite difference approximations to the advection-diffusion equation $u_t\ = a u_x\ + b u_{xx}$. The results are for centered difference schemes in space and include explicit and implicit schemes in time up to fourth order and schemes that use different space and time discretizations for the advective and diffusive terms. The results are derived from a uniform framework based on the Schur-Cohn theory of Simple von Neumann Polynomials and are necessary and sufficient for the stability of the Cauchy problem. Some of the results are believed to be new.

Journal ArticleDOI
TL;DR: The suitability of B-splines as a basis for piecewise polynomial solution representation for solving differential equations is challenged in this article, where two alternative local solution representations are considered in the context of collocating ordinary differential equations: "Hermite-type" and "Lmonomial" representations.
Abstract: The suitability of B-splines as a basis for piecewise polynomial solution representation for solving differential equations is challenged. Two alternative local solution representations are considered in the context of collocating ordinary differential equations: “Hermite-type” and “Lmonomial”. Both are much easier and shorter to implement and somewhat more efficient than B-splines.A new condition number estimate for the B-splines and Hermite-type representations is presented. One choice of the Hermite-type representation is experimentally determined to produce roundoff errors at most as large as those for B-splines. The monomial representation is shown to have a much smaller condition number than the other ones, and correspondingly produces smaller roundoff errors, especially for extremely nonuniform meshes. The operation counts for the two local representations considered are about the same, the Hermite-type representation being slightly cheaper. It is concluded that both representations are preferable,...

Journal ArticleDOI
TL;DR: In this paper, the convergence of the discrete ordinates method in the numerical solution of the transport equation for slab geometry is studied, and the combined effect of spatial and angular approximations is analyzed.
Abstract: We study the convergence of the discrete ordinates method in the numerical solution of the transport equation for slab geometry. In particular, the combined effect of spatial and angular approximations is analyzed. General conditions are given which guarantee the stability of the combined approximations. Error estimates are derived for both the scalar flux in the $L^p $-norm, $1 \leqq p \leqq \infty $, and for the approximate critical eigenvalues. A comparison with numerical results indicates that the eigenvalue error estimates are optimal.

Journal ArticleDOI
TL;DR: In this article, the authors obtained new bounds for the error in polynomial approximation in Sobolev spaces in an open bounded set which is star-shaped with respect to every point in a set of positive measure $B \subset \Omega $.
Abstract: We obtain new bounds for the error in the polynomial approximation in Sobolev spaces in an open bounded set which is star-shaped with respect to every point in a set of positive measure $B \subset \Omega $. This estimate follows by applying the Hardy–Littlewood maximal function.