The computation of isotropic vectors

We describe algorithms to compute isotropic vectors for matrices with real or complex entries. These are unit vectors b satisfying b ^* Ab = 0. For real matrices the algorithm uses only the eigenvectors of the symmetric part corresponding to the extreme eigenvalues. For complex matrices, we first use the eigenvalues and eigenvectors of the Hermitian matrix K = (A − A ^*)/2i. This works in many cases. In case of failure we use the Hermitian part H or a combination of eigenvectors of H and K. We give some numerical experiments comparing our algorithms with those proposed by R. Carden and C. Chorianopoulos, P. Psarrakos and F. Uhlig.     

Micro-local analysis for the Metropolis algorithm

We prove sharp rates of convergence to stationarity for a simple case of the Metropolis algorithm: the placement of a single disc of radius h randomly into the interval [ − 1 − h, 1 + h], with h > 0 small. We find good approximations for the top eigenvalues and eigenvectors. The analysis gives rigorous proof for the careful numerical work (in Exp. Math. 13, 207–213). The micro-local techniques employed offer promise for the analysis of more realistic problems.     

Homogenization of eigenvalue problems in perforated domains

In this paper, we treat some eigenvalue problems in periodically perforated domains and study the asymptotic behaviour of the eigenvalues and the eigenvectors when the number of holes in the domain increases to infinity Using the method of asymptotic expansion, we give explicit formula for the homogenized coefficients and expansion for eigenvalues and eigenvectors. If we denote by ε the size of each hole in the domain, then we obtain the following aysmptotic expansion for the eigenvalues: Dirichlet: λ_ε = ε⁻² λ + λ₀ +O (ε), Stekloff: λ_ε = ελ₁ +O (ε²), Neumann: λ_ε = λ₀ + ελ₁ +O (ε²). Using the method of energy, we prove a theorem of convergence in each case considered here. We briefly study correctors in the case of Neumann eigenvalue problem.     

Parameter-robust numerical method for a system of singularly perturbed initial value problems

In this work we study a system of M( ≥ 2) first-order singularly perturbed ordinary differential equations with given initial conditions. The leading term of each equation is multiplied by a distinct small positive parameter, which induces overlapping layers. A maximum principle does not, in general, hold for this system. It is discretized using backward Euler difference scheme for which a general convergence result is derived that allows to establish nodal convergence of O(N ^− 1ln N) on the Shishkin mesh and O(N ^− 1) on the Bakhvalov mesh, where N is the number of mesh intervals and the convergence is robust in all of the parameters. Numerical experiments are performed to support the theoretical results.     

Anharmonic Oscillators in the Complex Plane, $\boldsymbol{\mathcal{PT}}$-symmetry,and Real Eigenvalues

For integers m ≥ 3 and 1 ≤ ℓ ≤ m − 1, we study the eigenvalue problems − u ^″(z) + [( − 1)^ℓ(iz) ^m  − P(iz)]u(z) = λu(z) with the boundary conditions that u(z) decays to zero as z tends to infinity along the rays argz=-\fracp2±\frac(l+1)pm+2\arg z=-\frac{\pi}{2}\pm \frac{(\ell+1)\pi}{m+2} in the complex plane, where P is a polynomial of degree at most m − 1. We provide asymptotic expansions of the eigenvalues λ _n . Then we show that if the eigenvalue problem is PT\mathcal{PT}-symmetric, then the eigenvalues are all real and positive with at most finitely many exceptions. Moreover, we show that when gcd(m,l)=1\gcd(m,\ell)=1, the eigenvalue problem has infinitely many real eigenvalues if and only if one of its translations or itself is PT\mathcal{PT}-symmetric. Also, we will prove some other interesting direct and inverse spectral results.     

An algorithm for second order initial and boundary value problems with an automatic error estimate based on a third derivative method

A third derivative method (TDM) with continuous coefficients is derived and used to obtain a main and additional methods, which are simultaneously applied to provide all approximations on the entire interval for initial and boundary value problems of the form y′′ = f(x, y, y′). The convergence analysis of the method is discussed. An algorithm involving the TDMs is developed and equipped with an automatic error estimate based on the double mesh principle. Numerical experiments are performed to show efficiency and accuracy advantages.     

A particle approximation of the solution of the Kushner–Stratonovitch equation

We construct a sequence of branching particle systems α _n convergent in measure to the solution of the Kushner–Stratonovitch equation. The algorithm based on this result can be used to solve numerically the filtering problem. We prove that the rate of convergence of the algorithm is of order n ^?. This paper is the third in a sequence, and represents the most efficient algorithm we have identified so far. Received: 4 February 1997 / Revised version: 26 October 1998     

Spectral estimates for Schrödinger operators with sparse potentials on graphs

A construction of “sparse potentials,” suggested by the authors for the lattice \mathbbZ^d {\mathbb{Z}^d} , d > 2, is extended to a large class of combinatorial and metric graphs whose global dimension is a number D > 2. For the Schr?dinger operator − Δ − αV on such graphs, with a sparse potential V, we study the behavior (as α → ∞) of the number N_(−Δ − αV) of negative eigenvalues of − Δ − αV. We show that by means of sparse potentials one can realize any prescribed asymptotic behavior of N_(−Δ − αV) under very mild regularity assumptions. A similar construction works also for the lattice \mathbbZ² {\mathbb{Z}^2} , where D = 2. Bibliography: 13 titles.     

A Low-Complexity Divide-and-Conquer Method for Computing Eigenvalues and Eigenvectors of Symmetric Band Matrices

A framework for an efficient low-complexity divide-and-conquer algorithm for computing eigenvalues and eigenvectors of an n × n symmetric band matrix with semibandwidth b n is proposed and its arithmetic complexity analyzed. The distinctive feature of the algorithm—after subdivision of the original problem into p subproblems and their solution—is a separation of the eigenvalue and eigenvector computations in the central synthesis problem. The eigenvalues are computed recursively by representing the corresponding symmetric rank b(p–1) modification of a diagonal matrix as a series of rank-one modifications. Each rank-one modifications problem can be solved using techniques developed for the tridiagonal divide-and-conquer algorithm. Once the eigenvalues are known, the corresponding eigenvectors can be computed efficiently using modified QR factorizations with restricted column pivoting. It is shown that the complexity of the resulting divide-and-conquer algorithm is O (n ² b ²) (in exact arithmetic).This revised version was published online in October 2005 with corrections to the Cover Date.     

A spectral method for elliptic equations: the Neumann problem

Let Ω be an open, simply connected, and bounded region in ℝ ^d , d ≥ 2, and assume its boundary ∂Ω is smooth. Consider solving the elliptic partial differential equation − Δu + γu = f over Ω with a Neumann boundary condition. The problem is converted to an equivalent elliptic problem over the unit ball B, and then a spectral method is given that uses a special polynomial basis. In the case the Neumann problem is uniquely solvable, and with sufficiently smooth problem parameters, the method is shown to have very rapid convergence. Numerical examples illustrate exponential convergence.     

A new algorithm for computing the Geronimus transformation with large shifts

A monic Jacobi matrix is a tridiagonal matrix which contains the parameters of the three-term recurrence relation satisfied by the sequence of monic polynomials orthogonal with respect to a measure. The basic Geronimus transformation with shift α transforms the monic Jacobi matrix associated with a measure dμ into the monic Jacobi matrix associated with dμ/(x − α) + Cδ(x − α), for some constant C. In this paper we examine the algorithms available to compute this transformation and we propose a more accurate algorithm, estimate its forward errors, and prove that it is forward stable. In particular, we show that for C = 0 the problem is very ill-conditioned, and we present a new algorithm that uses extended precision.     

On the local convergence of a family of Euler-Halley type iterations with a parameter

This paper aims to study the local convergence of a family of Euler-Halley type methods with a parameter α for solving nonlinear operator equations under the second-order generalized Lipschitz assumption. The radius r _α of the optimal convergence ball and the error estimation of the method corresponding to α are estimated for each α ∈ ( − ∞ , + ∞ ). For each α > 0, we get r _α  ≥ r _− α and the upper bound of the error estimation of the method with α > 0 is not larger than the one with α < 0. For each α ≤ 0, we get the precise value of r _α , which is closely linked to the dynamical property of the method applied to a real or a complex function, and the optimal error estimation, which decreases when α→0⁻. Results show that the method corresponding to α is better than the one corresponding to − α for each α > 0 and the Chebyshev-Euler method is the best among all methods in the family with α ∈ ( − ∞ , 0] from the view of both safe choice of the initial point and error estimation.     

The two-level local projection stabilization as an enriched one-level approach

The two-level local projection stabilization is considered as a one-level approach in which the enrichments on each element are piecewise polynomial functions. The dimension of the enrichment space can be significantly reduced without losing the convergence order. On triangular meshes, for example, using continuous piecewise polynomials of degree r ≥ 1, only 2r − 1 functions per macro-cell are needed for the enrichment compared to r ² in the two-level approach. In case of the Stokes problem r − 1 functions per macro-cell are already sufficient to guarantee stability and to preserve convergence order. On quadrilateral meshes the corresponding reduction rates are even higher. We give examples of “reduced” two-level approaches and study how the constant in the local inf-sup condition for the one-level and different two-level approaches, respectively, depends on the polynomial degree r.     

Presentation geometrique des methodes de calcul des valeurs propres

Summary We present a survey of recent work on the convergence of methods for computing eigenvalues and eigenvectors of matrices. We try to maintain a geometric point of view and give pride of place to the R algorithm.

---

---

Cet article a été ecrit pendant le sejour de l'auteur en laboratoire d'Analyse Numérique de l'Université de Paris 6     

Piecewise-linear,discontinuous Galerkin method for a fractional diffusion equation

We use a piecewise-linear, discontinuous Galerkin method for the time discretization of a fractional diffusion equation involving a parameter in the range − 1 < α < 0. Our analysis shows that, for a time interval (0,T) and a spatial domain Ω, the error in L_￥((0,T);L₂(W))L_\infty\bigr((0,T);L_2(\Omega)\bigr) is of order k ^2 + α , where k denotes the maximum time step. Since derivatives of the solution may be singular at t = 0, our result requires the use of non-uniform time steps. In the limiting case α = 0 we recover the known O(k ²) convergence for the classical diffusion (heat) equation. We also consider a fully-discrete scheme that employs standard (continuous) piecewise-linear finite elements in space, and show that the additional error is of order h ²log(1/k). Numerical experiments indicate that our O(k ^2 + α ) error bound is pessimistic. In practice, we observe O(k ²) convergence even for α close to − 1.     

Weighted compound integration rules with higher order convergence for all N

Quasi-Monte Carlo integration rules, which are equal-weight sample averages of function values, have been popular for approximating multivariate integrals due to their superior convergence rate of order close to 1/N or better, compared to the order 1/?N1/\sqrt{N} of simple Monte Carlo algorithms. For practical applications, it is desirable to be able to increase the total number of sampling points N one or several at a time until a desired accuracy is met, while keeping all existing evaluations. We show that although a convergence rate of order close to 1/N can be achieved for all values of N (e.g., by using a good lattice sequence), it is impossible to get better than order 1/N convergence for all values of N by adding equally-weighted sampling points in this manner. We then prove that a convergence of order N ^− α for α > 1 can be achieved by weighting the sampling points, that is, by using a weighted compound integration rule. We apply our theory to lattice sequences and present some numerical results. The same theory also applies to digital sequences.     

Projection–difference method for controlled Fourier filtering

A projection–difference method is developed for approximating controlled Fourier filtering for quasilinear parabolic functional-differential equations. The method relies on a projection–difference scheme (PDS) for the approximation of the differential problem and derives a O(τ^1/2 + h) bound on the rate of convergence of PDS in the weighted energy norm without prior assumptions of additional smoothness of the generalized solutions. The PDS leads to a natural approximation of the objective functional in the optimal Fourier filtering problem. A bound of the same order is obtained for the rate of convergence in the functional of the problems approximating the Fourier filter control problem.     

Extending the eigCG algorithm to nonsymmetric Lanczos for linear systems with multiple right‐hand sides

The technique that was used to build the eigCG algorithm for sparse symmetric linear systems is extended to the nonsymmetric case using the BiCG algorithm. We show that, similar to the symmetric case, we can build an algorithm that is capable of computing a few smallest magnitude eigenvalues and their corresponding left and right eigenvectors of a nonsymmetric matrix using only a small window of the BiCG residuals while simultaneously solving a linear system with that matrix. For a system with multiple right‐hand sides, we give an algorithm that computes incrementally more eigenvalues while solving the first few systems and then uses the computed eigenvectors to deflate BiCGStab for the remaining systems. Our experiments on various test problems, including Lattice QCD, show the remarkable ability of eigBiCG to compute spectral approximations with accuracy comparable with that of the unrestarted, nonsymmetric Lanczos. Furthermore, our incremental eigBiCG followed by appropriately restarted and deflated BiCGStab provides a competitive method for systems with multiple right‐hand sides. Copyright © 2013 John Wiley & Sons, Ltd.     

The <Emphasis Type="Italic">k</Emphasis>-Centrum Straight-line Location Problem

A line is sought in the plane which minimizes the sum of the k largest (Euclidean) weighted distances from n given points. This problem generalizes the known straight-line center and median problems and, as far as the authors are aware, has not been tackled up to now. By way of geometric duality it is shown that such a line may always be found which either passes through two of the given points or lying at equal weighted distance from three of these. This allows construction of an algorithm to find all t-centrum lines for 1 ≤ t ≤ k in O((k + logn)n ³). Finally it is shown that both, the characterization of an optimal line and the algorithm, can be extended to any smooth norm.     

Estimates for weighted eigenvalues of a fourth-order elliptic operator with variable coefficients

We investigate the Dirichlet weighted eigenvalue problem for a fourth-order elliptic operator with variable coefficients in a bounded domain in \mathbbRⁿ {\mathbb{R}^n} . We establish a sharp inequality for its eigenvalues. It yields an estimate for the upper bound of the (k + 1)th eigenvalue in terms of the first k eigenvalues. Moreover, we also obtain estimates for some special cases of this problem. In particular, our results generalize the Wang–Xia inequality (J. Funct. Anal., 245, No. 1, 334–352 (2007)) for the clamped-plate problem to a fourth-order elliptic operator with variable coefficients.