A Levin-type algorithm for accelerating the convergence of Fourier series

A nonlinear sequence transformation is presented which is able to accelerate the convergence of Fourier series. It is tailored to be exact for a certain model sequence. As in the case of the Levin transformation and other transformations of Levin-type, in this model sequence the partial sum of the series is written as the sum of the limit (or antilimit) and a certain remainder, i.e., it is of Levin-type. The remainder is assumed to be the product of a remainder estimate and the sum of the first terms oftwo Poincaré-type expansions which are premultiplied by two different phase factors. This occurrence of two phase factors is the essential difference to the Levin transformation. The model sequence for the new transformation may also be regarded as a special case of a model sequence based on several remainder estimates leading to the generalized Richardson extrapolation process introduced by Sidi. An algorithm for the recursive computation of the new transformation is presented. This algorithm can be implemented using only two one-dimensional arrays. It is proved that the sequence transformation is exact for Fourier series of geometric type which have coefficients proportional to the powers of a numberq, |q|<1. It is shown that under certain conditions the algorithm indeed accelerates convergence, and the order of the convergence is estimated. Finally, numerical test data are presented which show that in many cases the new sequence transformation is more powerful than Wynn's epsilon algorithm if the remainder estimates are properly chosen. However, it should be noted that in the vicinity of singularities of the Fourier series the new sequence transformation shows a larger tendency to numerical instability than the epsilon algorithm.     

On nonsymmetric saddle point matrices that allow conjugate gradient iterations

Linear systems in saddle point form are usually highly indefinite,which often slows down iterative solvers such as Krylov subspace methods. It has been noted by several authors that negating the second block row of a symmetric indefinite saddle point matrix leads to a nonsymmetric matrix ${{\mathcal A}}Linear systems in saddle point form are usually highly indefinite,which often slows down iterative solvers such as Krylov subspace methods. It has been noted by several authors that negating the second block row of a symmetric indefinite saddle point matrix leads to a nonsymmetric matrix whose spectrum is entirely contained in the right half plane. In this paper we study conditions so that is diagonalizable with a real and positive spectrum. These conditions are based on necessary and sufficient conditions for positive definiteness of a certain bilinear form,with respect to which is symmetric. In case the latter conditions are satisfied, there exists a well defined conjugate gradient (CG) method for solving linear systems with . We give an efficient implementation of this method, discuss practical issues such as error bounds, and present numerical experiments. In memory of Gene Golub (1932–2007), our wonderful friend and colleague, who had a great interest in the conjugate gradient method and the numerical solution of saddle point problems. The work of J?rg Liesen was supported by the Emmy Noether-Program and the Heisenberg-Program of the Deutsche Forschungsgemeinschaft.     

Preconditioning Landweber iteration in Hilbert scales

In this paper we investigate convergence of Landweber iteration in Hilbert scales for linear and nonlinear inverse problems. As opposed to the usual application of Hilbert scales in the framework of regularization methods, we focus here on the case s≤0, which (for Tikhonov regularization) corresponds to regularization in a weaker norm. In this case, the Hilbert scale operator L^−2s appearing in the iteration acts as a preconditioner, which significantly reduces the number of iterations needed to match an appropriate stopping criterion. Additionally, we carry out our analysis under significantly relaxed conditions, i.e., we only require instead of which is the usual condition for regularization in Hilbert scales. The assumptions needed for our analysis are verified for several examples and numerical results are presented illustrating the theoretical ones. supported by the Austrian Science Foundation (FWF) under grant SFB/F013     

A modified harmonic block Arnoldi algorithm with adaptive shifts for large interior eigenproblems

The harmonic block Arnoldi method can be used to find interior eigenpairs of large matrices. Given a target point or shift τ$\tau$ to which the needed interior eigenvalues are close, the desired interior eigenpairs are the eigenvalues nearest τ$\tau$ and the associated eigenvectors. However, it has been shown that the harmonic Ritz vectors may converge erratically and even may fail to do so. To do a better job, a modified harmonic block Arnoldi method is coined that replaces the harmonic Ritz vectors by some modified harmonic Ritz vectors. The relationships between the modified harmonic block Arnoldi method and the original one are analyzed. Moreover, how to adaptively adjust shifts during iterations so as to improve convergence is also discussed. Numerical results on the efficiency of the new algorithm are reported.     

A <Emphasis Type="Italic">V</Emphasis>-cycle Multigrid for multilevel matrix algebras: proof of optimality

We analyze the convergence rate of a multigrid method for multilevel linear systems whose coefficient matrices are generated by a real and nonnegative multivariate polynomial f and belong to multilevel matrix algebras like circulant, tau, Hartley, or are of Toeplitz type. In the case of matrix algebra linear systems, we prove that the convergence rate is independent of the system dimension even in presence of asymptotical ill-conditioning (this happens iff f takes the zero value). More precisely, if the d-level coefficient matrix has partial dimension n _r at level r, with , then the size of the system is , , and O(N(n)) operations are required by the considered V-cycle Multigrid in order to compute the solution within a fixed accuracy. Since the total arithmetic cost is asymptotically equivalent to the one of a matrix-vector product, the proposed method is optimal. Some numerical experiments concerning linear systems arising in 2D and 3D applications are considered and discussed.     

The block grade of a block Krylov space

The aim of the paper is to compile and compare basic theoretical facts on Krylov subspaces and block Krylov subspaces. Many Krylov (sub)space methods for solving a linear system Ax=b have the property that in exact computer arithmetic the true solution is found after ν iterations, where ν is the dimension of the largest Krylov subspace generated by A from r₀, the residual of the initial approximation x₀. This dimension is called the grade of r₀ with respect to A. Though the structure of block Krylov subspaces is more complicated than that of ordinary Krylov subspaces, we introduce here a block grade for which an analogous statement holds when block Krylov space methods are applied to linear systems with multiple, say s, right-hand sides. In this case, the s initial residuals are bundled into a matrix R₀ with s columns. The possibility of linear dependence among columns of the block Krylov matrix , which in practical algorithms calls for the deletion (or, deflation) of some columns, requires extra care. Relations between grade and block grade are also established, as well as relations to the corresponding notions of a minimal polynomial and its companion matrix.     

Incomplete block-matrix factorization iterative methods for convection-diffusion problems

Standard Galerkin finite element methods or finite difference methods for singular perturbation problems lead to strongly unsymmetric matrices, which furthermore are in general notM-matrices. Accordingly, preconditioned iterative methods such as preconditioned (generalized) conjugate gradient methods, which have turned out to be very successful for symmetric and positive definite problems, can fail to converge or require an excessive number of iterations for singular perturbation problems.This is not so much due to the asymmetry, as it is to the fact that the spectrum can have both eigenvalues with positive and negative real parts, or eigenvalues with arbitrary small positive real parts and nonnegligible imaginary parts. This will be the case for a standard Galerkin method, unless the meshparameterh is chosen excessively small. There exist other discretization methods, however, for which the corresponding bilinear form is coercive, whence its finite element matrix has only eigenvalues with positive real parts; in fact, the real parts are positive uniformly in the singular perturbation parameter.In the present paper we examine the streamline diffusion finite element method in this respect. It is found that incomplete block-matrix factorization methods, both on classical form and on an inverse-free (vectorizable) form, coupled with a general least squares conjugate gradient method, can work exceptionally well on this type of problem. The number of iterations is sometimes significantly smaller than for the corresponding almost symmetric problem where the velocity field is close to zero or the singular perturbation parameter =1.The 2 ^nd author's research was sponsored by Control Data Corporation through its PACER fellowship program.The 3 ^rd author's research was supported by the Netherlands organization for scientific research (NWO).On leave from the Institute of Mathematics, Academy of Science, 1090 Sofia, P.O. Box 373, Bulgaria.     

Algebraic multilevel preconditioning methods. I

Summary A recursive way of constructing preconditioning matrices for the stiffness matrix in the discretization of selfadjoint second order elliptic boundary value problems is proposed. It is based on a sequence of nested finite element spaces with the usual nodal basis functions. Using a nodeordering corresponding to the nested meshes, the finite element stiffness matrix is recursively split up into two-level block structures and is factored approximately in such a way that any successive Schur complement is replaced (approximated) by a matrix defined recursively and thereform only implicitely given. To solve a system with this matrix we need to perform a fixed number (v) of iterations on the preceding level using as an iteration matrix the preconditioning matrix already defined on that level. It is shown that by a proper choice of iteration parameters it suffices to use \left( {1 - \gamma ^2 } \right)^{ - \tfrac{1}{2}} $$ " align="middle" border="0"> iterations for the so constructedv-foldV-cycle (wherev=2 corresponds to aW-cycle) preconditioning matrices to be spectrally equivalent to the stiffness matrix. The conditions involve only the constant in the strengthened C.-B.-S. inequality for the corresponding two-level hierarchical basis function spaces and are therefore independent of the regularity of the solution for instance. If we use successive uniform refinements of the meshes the method is of optimal order of computational complexity, if . Under reasonable assumptions of the finite element mesh, the condition numbers turn out to be so small that there are in practice few reasons to use an accelerated iterative method like the conjugate gradient method, for instance.Dedicated to the memory of Peter HenriciThe research of the second author reported here was supported in part by the Committee of Science, Bulgaria, under Grant No. 55/26.03.87     

A proximal point algorithm for minimax problems

This paper presents a proximal point algorithm for solving discretel _∞ approximation problems of the form minimize ∥Ax−b∥_∞. Let ε_∞ be a preassigned positive constant and let ε _l ,l = 0,1,2,... be a sequence of positive real numbers such that 0 < ε _l < ε_∞. Then, starting from an arbitrary pointz ₀, the proposed method generates a sequence of points z _l ,l= 0,1,2,..., via the rule . One feature that characterizes this algorithm is its finite termination property. That is, a solution is reached within a finite number of iterations. The smaller are the numbers ε _l the smaller is the number of iterations. In fact, if ε ₀ is sufficiently small then z₁ solves the original minimax problem. The practical value of the proposed iteration depends on the availability of an efficient code for solving a regularized minimax problem of the form minimize where ∈ is a given positive constant. It is shown that the dual of this problem has the form maximize , and ify solves the dual thenx=A ^T y solves the primal. The simple structure of the dual enables us to apply a wide range of methods. In this paper we design and analyze a row relaxation method which is suitable for solving large sparse problems. Numerical experiments illustrate the feasibility of our ideas.     

QMR: a quasi-minimal residual method for non-Hermitian linear systems

Summary The biconjugate gradient (BCG) method is the natural generalization of the classical conjugate gradient algorithm for Hermitian positive definite matrices to general non-Hermitian linear systems. Unfortunately, the original BCG algorithm is susceptible to possible breakdowns and numerical instabilities. In this paper, we present a novel BCG-like approach, the quasi-minimal residual (QMR) method, which overcomes the problems of BCG. An implementation of QMR based on a look-ahead version of the nonsymmetric Lanczos algorithm is proposed. It is shown how BCG iterates can be recovered stably from the QMR process. Some further properties of the QMR approach are given and an error bound is presented. Finally, numerical experiments are reported.This work was supported in part by DARPA via Cooperative Agreement NCC 2-387 between NASA and the Universities Space Research Association (USRA).     

A Jacobi eigenreduction algorithm for definite matrix pairs

Summary We propose a Jacobi eigenreduction algorithm for symmetric definite matrix pairsA, J of small to medium-size real symmetric matrices withJ ²=I,J diagonal (neitherJ norA itself need be definite). Our Jacobi reduction works only on one matrix and usesJ-orthogonal elementary congruences which include both trigonometric and hyperbolic rotations and preserve the symmetry throughout the process. For the rotation parameters only the pivotal elements of the current matrix are needed which facilitates parallelization. We prove the global convergence of the method; the quadratic convergence was observed in all experiments. We apply our method in two situations: (i) eigenreducing a single real symmetric matrix and (ii) eigenreducing an overdamped quadratic matrix pencil. In both cases our method is preceded by a symmetric indefinite decomposition and performed in its one-sided variant on the thus obtained factors. Our method outdoes the standard methods like standard Jacobi orqr/ql in accuracy in spite of the use of hyperbolic transformations which are not orthogonal (a theoretical justification of this behaviour is made elsewhere). The accuracy advantage of our method can be particularly drastic if the eigenvalues are of different order. In addition, in working with quadratic pencils our method is shown to either converge or to detect non-overdampedness.     

On acceleration methods for coupled nonlinear elliptic systems

Summary We compare both numerically and theoretically three techniques for accelerating the convergence of a nonlinear fixed point iterationuT(u), arising from a coupled elliptic system: Chebyshev acceleration, a second order stationary method, and a nonlinear version of the Generalized Minimal Residual Algorithm (GMRES) which we call NLGMR. All three approaches are implemented in Jacobian-free mode, i.e., only a subroutine which returnsT(u) as a function ofu is required.We present a set of numerical comparisons for the drift-diffusion semiconductor model. For the mappingT which corresponds to the nonlinear block Gauß-Seidel algorithm for the solution of this nonlinear elliptic system, NLGMR is found to be superior to the second order stationary method and the Chebychev acceleration. We analyze the local convergence of the nonlinear iterations in terms of the spectrum [T _u (u ^(*))] of the derivativeT _u at the solutionu ^(*). The convergence of the original iteration is governed by the spectral radius [T _U (u ^(*))]. In contrast, the convergence of the two second order accelerations are related to the convex hull of [T _u (u ^(*))], while the convergence of the GMRES-based approach is related to the local clustering in [I–T _u (u ^(*))]. The spectrum [I–T _u (u ^(*))] clusters only at 1 due to the successive inversions of elliptic partial differential equations inT. We explain the observed superiority of GMRES over the second order acceleration by its ability to take advantage of this clustering feature, which is shared by similar coupled nonlinear elliptic systems.     

Finding a Nash equilibrium in noncooperativeN-person games by solving a sequence of linear stationary point problems

In this paper we present an algorithm for finding a Nash equilibrium in a noncooperative normal formN-person game. More generally, the algorithm can be applied for solving a nonlinear stationary point problem on a simplotope, being the Cartesian product of several simplices. The algorithm solves the problem by solving a sequence of linear stationary point problems. Each problem in the sequence is solved in a finite number of iterations. Although the overall convergence cannot be proved, the method performs rather well. Computational results suggest that this algorithm performs at least as good as simplicial algorithms do.For the special case of a bi-matrix game (N=2), the algorithm has an appealing game-theoretic interpretation. In that case, the problem is linear and the algorithm always finds a solution. Furthermore, the equilibrium found in a bi-matrix game is perfect whenever the algorithm starts from a strategy vector at which all actions are played with positive probability.This research is part of the VF-program Co-operation and Competition, which has been approved by the Netherlands Ministery of Education and Sciences.     

Nonlinear hybrid procedures and fixed point iterations

Let (x_n) and (?_n) be two vector sequences converging to a common limit. First, we shall define nonlinear hybrid procedures which consist of constructing a new vector sequence (y_n) with better convergence properties than (x_n) and (?_n). Then, this procedure is used for accelerating the convergence of a given sequence and applied to the construction of fixed point methods. New methods are derived. Finally, the connection between fixed point iterations and methods for the numerical integration of differential equations is also exploited. Numerical results are given.     

A spectral Galerkin approximation of the Orr-Sommerfeld eigenvalue problem in a semi-infinite domain

Summary The convergence of a Galerkin approximation of the Orr-Sommerfeld eigenvalue problem, which is defined in a semi-infinite domain, is studied theoretically. In case the system of trial functions is based on a composite of Jacobi polynomials and an exponential transform of the semi-infinite domain, the error of the Galerkin approximation is estimated in terms of the transformation parametera and the numberN of trial functions. Finite or infinite-order convergence of the spectral Galerkin method is obtained depending on how the transformation parameter is chosen. If the transformation parameter is fixed, then convergence is of finite order only. However, ifa is varied proportional to 1/N with an exponent 0<<1, then the approximate eigenvalue converges faster than any finite power of 1/N asN. Some numerical examles are given.     

On an iterative algorithm with superquadratic convergence for solving nonlinear operator equations

We study an iterative method with order for solving nonlinear operator equations in Banach spaces. Algorithms for specific operator equations are built up. We present the received new results of the local and semilocal convergence, in case when the first-order divided differences of a nonlinear operator are Hölder continuous. Moreover a quadratic nonlinear majorant for a nonlinear operator, according to the conditions laid upon it, is built. A priori and a posteriori estimations of the method’s error are received. The method needs almost the same number of computations as the classical Secant method, but has a higher order of convergence. We apply our results to the numerical solving of a nonlinear boundary value problem of second-order and to the systems of nonlinear equations of large dimension.     

Nonlinear coordinate representations of smooth optimization problems

Nonlinear coordinate representations of smooth optimization problems are investigated from the point of view of variable metric algorithms. In other words, nonlinear coordinate systems, in the sense of differential geometry, are studied by taking into consideration the structure of smooth optimization problems and variable metric methods.Both the unconstrained and constrained cases are discussed. The present approach is based on the fact that the nonlinear coordinate transformation of an optimization problem can be replaced by a suitable Riemannian metric belonging to the Euclidean metric class. In the case of equality and inequality constraints, these questions are related closely to the right inverses of full-rank matrices; therefore, their characterization is a starting point of the present analysis. The main results concern a new subclass of nonlinear transformations in connection with the common supply of coordinates to two Riemannian manifolds, one immersed in the other one. This situation corresponds to the differentiable manifold structure of nonlinear optimization problems and improves the insight into the theoretical background of variable metric algorithms. For a wide class of variable metric methods, a convergence theorem in invariant form (not depending on coordinate representations) is proved. Finally, a problem of convexification by nonlinear coordinate transformations and image representations is studied.This research was supported by the Hungarian National Research Foundation, Grant Nos. OTKA-2568 and OTKA-2116, and by the Project Trasporti of the Italian National Research Council (CNR).     

On the Vector ε-Algorithm for Solving Linear Systems of Equations

The four vector extrapolation methods, minimal polynomial extrapolation, reduced rank extrapolation, modified minimal polynomial extrapolation and the topological epsilon algorithm, when applied to linearly generated vector sequences are Krylov subspace methods and it is known that they are equivalent to some well-known conjugate gradient type methods. However, the vector -algorithm is an extrapolation method, older than the four extrapolation methods above, and no similar results are known for it. In this paper, a determinantal formula for the vector -algorithm is given. Then it is shown that, when applied to a linearly generated vector sequence, the algorithm is also a Krylov subspace method and for a class of matrices the method is equivalent to a preconditioned Lanczos method. A new determinantal formula for the CGS is given, and an algebraic comparison between the vector -algorithm for linear systems and CGS is also given.     

Error estimates for approximating non-hyperbolic heteroclinic orbits of maps

Summary. In this paper we consider heteroclinic orbits in discrete time dynamical systems that connect a hyperbolic fixed point to a non-hyperbolic fixed point with a one-dimensional center direction. A numerical method for approximating the heteroclinic orbit by a finite orbit sequence is introduced and a detailed error analysis is presented. The loss of hyperbolicity requires special tools for proving the error estimate – the polynomial dichotomy of linear difference equations and a (partial) normal form transformation near the non-hyperbolic fixed point. This situation appears, for example, when one fixed point undergoes a flip bifurcation. For this case, the approximation method and the validity of the error estimate is illustrated by an example.Supported by the DFG Research Group Spectral analysis, asymptotic distributions and stochastic dynamicsMathematics Subject Classification (2000): 37C29, 65P30Acknowledgement The authors thank the referees for several helpful suggestions which improved the first version of this paper.