Approximate Invariant Subspaces Research Articles

Recycling Krylov subspaces for efficient partitioned solution of aerostructural adjoint systems

Robust and efficient solvers for coupled-adjoint linear systems are crucial to successful aerostructural optimization. Monolithic and partitioned strategies can be applied. The monolithic approach is expected to offer better robustness and efficiency for strong fluid-structure interactions. However, it requires a higher implementation cost and convergence may depend on appropriate scaling and initialization strategies. On the other hand, the modularity of the partitioned method enables a straightforward implementation while its convergence may require relaxation. In addition, a partitioned solver often leads to a higher number of iterations to get the same level of convergence as the monolithic one.The objective of this paper is to accelerate the fluid-structure coupled-adjoint partitioned solver by considering techniques borrowed from approximate invariant subspace recycling strategies adapted to sequences of linear systems with varying right-hand sides. Indeed, in a partitioned framework, the structural source term attached to the fluid block of equations affects the right-hand side with the nice property of quickly converging to a constant value. For the fluid block of equations, we also consider deflation of approximate eigenvectors in conjunction with advanced inner-outer Krylov solvers.To demonstrate the effectiveness of our proposed recycling strategy, we compute the coupled derivatives for two emblematic problems in transonic viscous flow. The first one is an aeroelastic configuration of the ONERA-M6 fixed wing. For this exercise, the fluid grid was coupled to a structural model specifically designed to exhibit a high flexibility. The second one is a wing-body aeroelastic configuration of the Common Research Model. All computations are performed using a fully linearized one-equation Spalart-Allmaras turbulence model. Numerical simulations show up to 39% reduction in matrix-vector products for GCRO-DR and up to 25% for the nested FGCRO-DR solver.

A Refinement of Approximate Invariant Subspaces of Matrices Based on SVD in High Dimensionality Reduction and Image Compression

A Refinement of Approximate Invariant Subspaces of Matrices Based on SVD in High Dimensionality Reduction and Image Compression

Krylov subspace recycling for evolving structures

Krylov subspace recycling is a powerful tool when solving a long series of large, sparse linear systems that change only slowly over time. In PDE constrained shape optimization, these series appear naturally, as typically hundreds or thousands of optimization steps are needed with only small changes in the geometry. In this setting, however, applying Krylov subspace recycling can be a difficult task. As the geometry evolves, in general, so does the finite element mesh defined on or representing this geometry, including the numbers of nodes and elements and element connectivity. This is especially the case if re-meshing techniques are used. As a result, the number of algebraic degrees of freedom in the system changes, and in general the linear system matrices resulting from the finite element discretization change size from one optimization step to the next. Changes in the mesh connectivity also lead to structural changes in the matrices. In the case of re-meshing, even if the geometry changes only a little, the corresponding mesh might differ substantially from the previous one. Obviously, this prevents any straightforward mapping of the approximate invariant subspace of the linear system matrix (the focus of recycling in this paper) from one optimization step to the next; similar problems arise for other selected subspaces. In this paper, we present an algorithm to map an approximate invariant subspace of the linear system matrix for the previous optimization step to an approximate invariant subspace of the linear system matrix for the current optimization step, for general meshes. This is achieved by exploiting the map from coefficient vectors to finite element functions on the mesh, combined with interpolation or approximation of functions on the finite element mesh. We demonstrate the effectiveness of our approach numerically with several proof of concept studies for a specific meshing technique.

A Posteriori Error Estimates for Elliptic Eigenvalue Problems Using Auxiliary Subspace Techniques

We propose an a posteriori error estimator for high-order p- or hp-finite element discretizations of selfadjoint linear elliptic eigenvalue problems that is appropriate for estimating the error in the approximation of an eigenvalue cluster and the corresponding invariant subspace. The estimator is based on the computation of approximate error functions in a space that complements the one in which the approximate eigenvectors were computed. These error functions are used to construct estimates of collective measures of error, such as the Hausdorff distance between the true and approximate clusters of eigenvalues, and the subspace gap between the corresponding true and approximate invariant subspaces. Numerical experiments demonstrate the practical effectivity of the approach.

On the Kahan-Parlett-Jiang theorem—A globally optimal backward perturbation error for two-sided invariant subspaces

A crucial question in large non-Hermitian eigenproblems is whether or not any residual bounds can be very useful in practical situations. For this problem, the famous Kahan-Parlett-Jiang theorem Kahan et al. (1982) [15] seems to be the best compromise between generality and sharpness. In essence, it presents an a posteriori error bound on approximate solutions to eigenproblems, which provides a powerful tool to evaluate quality of the computed two-sided invariant subspaces of large-scale non-Hermitian matrices. However, the perturbation error determined by this theorem is only locally optimal rather than globally optimal. In this work, we revisit this problem and derive a globally optimal backward perturbation error for given two-sided approximate invariant subspaces. Our theorem enhances the Kahan-Parlett-Jiang theorem, and numerical experiments demonstrate the effectiveness of our theoretical results.

A spectrally preconditioned and initially deflated variant of the restarted block GMRES method for solving multiple right-hand sides linear systems

The solution of large linear systems with multiple right-hand sides given simultaneously is required in many large-scale scientific and engineering applications modelled by either partial differential or boundary integral equations. Block Krylov subspace methods are attractive to use for this problem class as they can overcome the memory bottleneck of direct methods and they perform block matrix-vector products, achieving high computational efficiency on modern cache-based computer architectures. In this paper we introduce variants of the block GMRES method that combine initial deflation and eigenvalue recycling strategies to remedy some typical convergence problems of block Krylov solvers. The new class of block iterative solvers has the ability to handle the approximate linear dependence of the block of right-hand sides and exploits approximate invariant subspaces recycled over the iterations to mitigate the bad effects that small eigenvalues can have on the convergence, by adapting an existing preconditioner. We illustrate the numerical behavior of the spectrally updated and initially deflated block GMRES method on a set of linear systems arising from the discretization of the Dirac equation and of boundary integral equations in electromagnetics scattering.

A new preconditioner update strategy for the solution of sequences of linear systems in structural mechanics: application to saddle point problems in elasticity

Many applications in structural mechanics require the numerical solution of sequences of linear systems typically issued from a finite element discretization of the governing equations on fine meshes. The method of Lagrange multipliers is often used to take into account mechanical constraints. The resulting matrices then exhibit a saddle point structure and the iterative solution of such preconditioned linear systems is considered as challenging. A popular strategy is then to combine preconditioning and deflation to yield an efficient method. We propose an alternative that is applicable to the general case and not only to matrices with a saddle point structure. In this approach, we consider to update an existing algebraic or application-based preconditioner, using specific available information exploiting the knowledge of an approximate invariant subspace or of matrix-vector products. The resulting preconditioner has the form of a limited memory quasi-Newton matrix and requires a small number of linearly independent vectors. Numerical experiments performed on three large-scale applications in elasticity highlight the relevance of the new approach. We show that the proposed method outperforms the deflation method when considering sequences of linear systems with varying matrices.

GMRES Convergence Bounds for Eigenvalue Problems

Abstract The convergence of GMRES for solving linear systems can be influenced heavily by the structure of the right-hand side. Within the solution of eigenvalue problems via inverse iteration or subspace iteration, the right-hand side is generally related to an approximate invariant subspace of the linear system. We give detailed and new bounds on (block) GMRES that take the special behavior of the right-hand side into account and explain the initial sharp decrease of the GMRES residual. The bounds motivate the use of specific preconditioners for these eigenvalue problems, e.g., tuned and polynomial preconditioners, as we describe. The numerical results show that the new (block) GMRES bounds are much sharper than conventional bounds and that preconditioned subspace iteration with either a tuned or polynomial preconditioner should be used in practice.

Recycling BiCGSTAB with an Application to Parametric Model Order Reduction

Krylov subspace recycling is a process for accelerating the convergence of sequences of linear systems. Based on this technique, the recycling BiCG algorithm has been developed recently. Here, we now generalize and extend this recycling theory to BiCGSTAB. Recycling BiCG focuses on efficiently solving sequences of dual linear systems, while the focus here is on efficiently solving sequences of single linear systems (assuming non-symmetric matrices for both recycling BiCG and recycling BiCGSTAB). As compared with other methods for solving sequences of single linear systems with non-symmetric matrices (e.g., recycling variants of GMRES), BiCG based recycling algorithms, like recycling BiCGSTAB, have the advantage that they involve a short-term recurrence, and hence, do not suffer from storage issues and are also cheaper with respect to the orthogonalizations. We modify the BiCGSTAB algorithm to use a recycle space, which is built from left and right approximate invariant subspaces. Using our algorithm for a parametric model order reduction example gives good results. We show about 40% savings in the number of matrix-vector products and about 35% savings in runtime.

Rayleigh--Ritz Approximation For the Linear Response Eigenvalue Problem

Large scale eigenvalue computation is about approximating certain invariant subspaces associated with the interesting part of the spectrum, and the interesting eigenvalues are then extracted through projecting the problem through approximate invariant subspaces into a much smaller eigenvalue problem. In the case of the linear response eigenvalue problem (aka the random phase eigenvalue problem), it is the pair of deflating subspaces associated with the first few smallest positive eigenvalues that needs to be computed. This paper is concerned with approximation accuracy relationships between a pair of approximate deflating subspaces and approximate eigenvalues extracted by the pair. Lower and upper bounds on eigenvalue approximation errors are obtained in terms of canonical angles between an exact and computed pair of deflating subspaces. These bounds can also be interpreted as lower/upper bounds on the canonical angles in terms of eigenvalue approximation errors. They are useful in analyzing numerical sol...

Structured mapping problems for linearly structured matrices

Given an appropriate class of structured matrices S, we characterize matrices X and B for which there exists a matrix A∈S such that AX=B and determine all matrices in S mapping X to B. We also determine all matrices in S mapping X to B and having the smallest norm. We use these results to investigate structured backward errors of approximate eigenpairs and approximate invariant subspaces, and structured pseudospectra of structured matrices.

Recycling Krylov subspaces for efficient large-scale electrical impedance tomography

Electrical impedance tomography (EIT) captures images of internal features of a body. Electrodes are attached to the boundary of the body, low intensity alternating currents are applied, and the resulting electric potentials are measured. Then, based on the measurements, an estimation algorithm obtains the three-dimensional internal admittivity distribution that corresponds to the image. One of the main goals of medical EIT is to achieve high resolution and an accurate result at low computational cost. However, when the finite element method (FEM) is employed and the corresponding mesh is refined to increase resolution and accuracy, the computational cost increases substantially, especially in the estimation of absolute admittivity distributions. Therefore, we consider in this work a fast iterative solver for the forward problem, which was previously reported in the context of structural optimization. We propose several improvements to this solver to increase its performance in the EIT context. The solver is based on the recycling of approximate invariant subspaces, and it is applied to reduce the EIT computation time for a constant and high resolution finite element mesh. In addition, we consider a powerful preconditioner and provide a detailed pseudocode for the improved iterative solver. The numerical results show the effectiveness of our approach: the proposed algorithm is faster than the preconditioned conjugate gradient (CG) algorithm. The results also show that even on a standard PC without parallelization, a high mesh resolution (more than 150,000 degrees of freedom) can be used for image estimation at a relatively low computational cost.

Riccati-based preconditioner for computing invariant subspaces of large matrices

This paper introduces and analyzes the convergence properties of a method that computes an approximation to the invariant subspace associated with a group of eigenvalues of a large not necessarily diagonalizable matrix. The method belongs to the family of projection type methods. At each step, it refines the approximate invariant subspace using a linearized Riccati's equation which turns out to be the block analogue of the correction used in the Jacobi-Davidson method. The analysis conducted in this paper shows that the method converges at a rate quasi-quadratic provided that the approximate invariant subspace is close to the exact one. The implementation of the method based on multigrid techniques is also discussed and numerical experiments are reported.

On eigenvalue variations of rayleigh quotient matrix pencils of a definite pencil

Let A − λ B be a definite matrix pencil of order n, i.e., both A and B are n × n Hermitian and c(A,B)= def min xε C n,‖x‖ 2=1 |x H(A+iB)x|>0 . Suppose Y is an n × l matrix with full column rank whose column vectors span an approximate invariant subspace for A − λ B. This note investigates the relation between eigenvalues of A− λB and those of Y H AY− λY H BY. Our result for the spectral norm improves Sun's [ Linear Algebra Appl. 139:253–267 (1990)]. We also present bounds in Frobenius norms and in general unitarily invariant norms.

Bounds on the error of an approximate invariant subspace for non-self-adjoint matrices

Suppose one approximates an invariant subspace of an\(n\times n\) matrix in\({\Bbb C}^{n\times n}\) which in not necessarily self--adjoint. Suppose that one also has an approximation for the corresponding eigenvalues. We consider the question of how good the approximations are. Specifically, we develop bounds on the angle between the approximating subspace and the invariant subspace itself. These bounds are functions of the following three terms: (1) the residual of the approximations; (2) singular--value separation in an associated matrix; and (3) the goodness of the approximations to the eigenvalues.

Two Simple Residual Bounds for the Eigenvalues of a Hermitian Matrix

Let A be Hermitian and let the orthonormal columns of X span an approximate invariant subspace of X. Then the residual $R = AX - XM\, ( M = X^H AX )$ will be small. The theorems of this paper bound the distance of the spectrum of M from the spectrum of A in terms of appropriate norms of R.

Quantum Theory of Dynamical Collective Subspace for Large-Amplitude Collective Motion

By placing emphasis on conceptual correspondence to the theory which has been developed within the framework of the time· dependent Hartree·Fock theory, a full quantum theory appropriate for describing large-amplitude collective motion is proposed. A central problem of the quantum theory is how to determine an optimal representation called a dynamical representation; the representation is specific for the collective subspace where the large-amplitude collective motion is replicated as satisfactorily as possible. As an extension of the classical theory where the concept of an approximate integral surface plays an important role, the dynamical representation is properly characterized by introducing a concept of an approximate invariant subspace of the Hamiltonian.

Residual Bounds on Approximate Eigensystems of Nonnormal Matrices

For nonnormal matrices the norms of the residuals of approximate eigenvectors are not by themselves sufficient information to bound the error in the approximate eigenvalue. It is sufficient however to give a bound on the distance to the nearest matrix for which the given approximations are exact. This result is extended to cover approximate invariant subspaces and their residuals. The theorems are used to derive a useful set of criteria for terminating the two-sided Lanczos algorithm. The study begins with a list of error bounds for eigenvalues of nonnormal matrices.

Error Bounds for Approximate Invariant Subspaces of Closed Linear Operators

Let A be a closed linear operator on a separable Hilbert space $\mathcal{H}$ whose domain is dense in $\mathcal{H}$ Let $\mathcal{X}$ be a subspace of $\mathcal{H}$ contained in the domain of A and let $\mathcal{Y}$ be its orthogonal complement. Let B and C be the compressions of A to $\mathcal{Z}$ and $\mathcal{Y}$ respectively, let $G = Y^ * AX$, where X and Y are the injections of $\mathcal{X}$ and $\mathcal{Y}$ into $\mathcal{H}$. It is shown that if B and C have disjoint spectra and $\| G \|$ is sufficiently small, then there is an invariant subspace $\mathcal{X}'$ of A near $\mathcal{X}$. Bounds for the distance between $\mathcal{X}'$ and $\mathcal{X}$ are given, and the spectrum of A is related to the spectra of B and C. In the development a measure of the separation of the spectra of B and C which is insensitive to small perturbations in B and C is introduced and analyzed.