Subspace Recycling 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.

An accelerated deflation preconditioner for acoustic BEM systems based on subspace recycling

Following a classic Krylov subspace recycling [1] technique in the context of frequency sweep analyses performed with the acoustic BEM [2] might inflict unnecessary overhead costs. To speed up such a procedure, this work employs an accelerated recycling technique that is based on the recently proposed model order reduction framework for BEM systems [3]. Specifically, employing an affine approximation of a BEM system in combination with a Galerkin projection on the deflation basis, the operation of constructing the corresponding deflation based preconditioner is significantly accelerated. The efficiency of the deflation basis for the entire frequency range is guaranteed through an Automatic Krylov subspace Recycling (AKR) [4] algorithm, thus leading to a drastic decrease of the number of iterations required by an iterative solver. The proposed approach is tested on an academic BEM example and achieves significant acceleration of the computed frequency sweep by ensuring fast convergence of the employed iterative solver, avoiding overhead projection costs.

An accelerated subspaces recycling strategy for the deflation of parametric linear systems based on model order reduction

Krylov subspace recycling has been extensively used to facilitate the solution of sequences of linear systems by constructing a deflation subspace and accelerating the convergence of a corresponding iterative solver. However, most existing techniques update the recycled subspace sequentially for each system, thus inducing a potentially high computational cost. In that context, this work proposes a method to accelerate the above procedure for the case of multiresolution analyses of affine parametric systems, by decoupling the solution of the system from the construction of the recycled subspace. The proposed method follows the projection based Model Order Reduction (MOR) logic that splits operations into an offline and an online phase and therefore proves particularly beneficial in case of affine parametrizations that involve multiple affine coefficients. In the offline phase of the method an especially tailored version of the Automatic Krylov subspaces Recycling (AKR) algorithm (Panagiotopoulos et al., 2021), proposed within this work and denoted as AKR-D, is employed. In brief, AKR-D constructs a high quality recycling basis W for a predefined parameter interval Ψ by targeting a desirable convergence rate at an automatically generated set of parameter values Ω⊂Ψ. The upfront construction of W enables an a-priori Galerkin projection of the affinely described full parametric system to yield a reduced order model (ROM). This ROM can be leveraged in the online phase of the proposed recycling method to accelerate the construction of the corresponding projector P onto W=span{W}, whose implicit assembly is required in most recycling based deflation schemes. The effectiveness of the proposed strategy is investigated in terms of algorithmic complexity and is demonstrated on a densely parametrized system arising within the boundary integral solver of the Helmholtz equation and a random sparse parametrization example.

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 Flexible Global GCRO-DR Method for Shifted Linear Systems and General Coupled Matrix Equations

In this paper, we mainly focus on the development and study of a new global GCRO-DR method that allows both the flexible preconditioning and the subspace recycling for sequences of shifted linear systems. The novel method presented here has two main advantages: firstly, it does not require the right-hand sides to be related, and, secondly, it can also be compatible with the general preconditioning. Meanwhile, we apply the new algorithm to solve the general coupled matrix equations. Moreover, by performing an error analysis, we deduce that a much looser tolerance can be applied to save computation by limiting the flexible preconditioned work without sacrificing the closeness of the computed and the true residuals. Finally, numerical experiments demonstrate that the proposed method illustrated can be more competitive than some other global GMRES-type methods.

A survey of subspace recycling iterative methods

AbstractThis survey concerns subspace recycling methods, a popular class of iterative methods that enable effective reuse of subspace information in order to speed up convergence and find good initial vectors over a sequence of linear systems with slowly changing coefficient matrices, multiple right‐hand sides, or both. The subspace information that is recycled is usually generated during the run of an iterative method (usually a Krylov subspace method) on one or more of the systems. Following introduction of definitions and notation, we examine the history of early augmentation schemes along with deflation preconditioning schemes and their influence on the development of recycling methods. We then discuss a general residual constraint framework through which many augmented Krylov and recycling methods can both be viewed. We review several augmented and recycling methods within this framework. We then discuss some known effective strategies for choosing subspaces to recycle before taking the reader through more recent developments that have generalized recycling for (sequences of) shifted linear systems, some of them with multiple right‐hand sides in mind. We round out our survey with a brief review of application areas that have seen benefit from subspace recycling methods.

Accelerating linear system solvers for time-domain component separation of cosmic microwave background data

Component separation is one of the key stages of any modern cosmic microwave background data analysis pipeline. It is an inherently nonlinear procedure and typically involves a series of sequential solutions of linear systems with similar but not identical system matrices, derived for different data models of the same data set. Sequences of this type arise, for instance, in the maximization of the data likelihood with respect to foreground parameters or sampling of their posterior distribution. However, they are also common in many other contexts. In this work we consider solving the component separation problem directly in the measurement (time-) domain. This can have a number of important benefits over the more standard pixel-based methods, in particular if non-negligible time-domain noise correlations are present, as is commonly the case. The approach based on the time-domain, however, implies significant computational effort because the full volume of the time-domain data set needs to be manipulated. To address this challenge, we propose and study efficient solvers adapted to solving time-domain-based component separation systems and their sequences, and which are capable of capitalizing on information derived from the previous solutions. This is achieved either by adapting the initial guess of the subsequent system or through a so-called subspace recycling, which allows constructing progressively more efficient two-level preconditioners. We report an overall speed-up over solving the systems independently of a factor of nearly 7, or 5, in our numerical experiments, which are inspired by the likelihood maximization and likelihood sampling procedures, respectively.

Mstab: Stabilized Induced Dimension Reduction for Krylov Subspace Recycling

We introduce $\mathcal{M}$stab, a Krylov subspace recycling method for the iterative solution of sequences of linear systems, where the system matrix is fixed and is large, sparse, and nonsymmetric, and the right-hand-side vectors are available in sequence. $\mathcal{M}$stab utilizes the short-recurrence principle of induced dimension reduction-type methods, adapted to solve sequences of linear systems. Using IDRstab for solving the linear system with the first right-hand side, the proposed method then recycles the Petrov space constructed throughout the solution of that system, generating a larger initial space for subsequent linear systems. The richer space potentially produces a rapidly convergent scheme. Numerical experiments demonstrate that $\mathcal{M}$stab often enters the superlinear convergence regime faster than other Krylov-type recycling methods.

Robust and efficient adjoint solver for complex flow conditions

Abstract A key step in gradient-based aerodynamic shape optimisation using the Reynolds-averaged Navier–Stokes equations is to compute the adjoint solution. Adjoint equations inherit the linear stability and the stiffness of the nonlinear flow equations. Therefore for industrial cases with complex geometries at off-design flow conditions, solving the resulting stiff adjoint equation can be challenging. In this paper, Krylov subspace solvers enhanced by subspace recycling and preconditioned with incomplete lower-upper factorisation are used to solve the stiff adjoint equations arising from typical design and off-design conditions. Compared to the baseline matrix-forming adjoint solver based on the generalized minimal residual method, the proposed algorithm achieved memory reduction of up to a factor of two and convergence speedup of up to a factor of three, on industry-relevant cases. These test cases include the DLR-F6 and DLR-F11 configurations, a wing-body configuration in pre-shock buffet and a large civil aircraft with mesh sizes ranging from 3 to 30 million. The proposed method seems to be particularly effective for the more difficult flow conditions.

Enabling off-design linearised aerodynamics analysis using Krylov subspace recycling technique

The major computational challenge, when using frequency domain linearised computational fluid dynamics in the analysis of aeroelastic problems such as aircraft flutter, gust response or shock buffet, are the excessive memory and CPU time requirements to solve the large sparse linear systems of equations. To address these issues found with the generalised minimal residual linear equation solver, the generalised conjugate residual solver with deflated restarting is adopted here in which an invariant Krylov subspace is recycled both between restarts when solving a single linear frequency domain problem and for a sequence of equations when varying the system matrix and forcing terms. The proposed method is implemented in an industrial code and applied to three test cases including the forced excitation and buffet onset of a pitch-plunge aerofoil, a realistic passenger aircraft in inviscid transonic flow and a generic half wing-body model at a pre-buffet condition. The memory requirements for the problems investigated are reduced by up to an order of magnitude, while the CPU times are reduced by up to a factor of three.

A combination of the fast multipole boundary element method and Krylov subspace recycling solvers

The solution of the Helmholtz equation by the Boundary Element Method leads to a sequence of frequency dependent linear systems of equations, where each is typically solved independently. The Krylov Subspace Recycling algorithms, like the GCRO-DR and the GCROT, are based on the idea that the solutions of consecutive systems have similarities and the information of the previous cycle can be reused to accelerate the convergence. These solvers showed very good results for sparse matrices arising in the FEM and are now applied to the fully populated BEM matrices. Additionally, the solution of a single system of equations is accelerated by the Fast Multipole Method, which shows a mostly linear correlation between iterations and calculation time. Hence the newly proposed combination has a high potential of achieving a faster solution process. The 3D Fast Multipole Boundary Element Method additionally incorporates a Burton–Miller formulation and a halfspace formulation to be applicable to a wider range of engineering problems. The method is illustrated and discussed by two different numerical examples. The advantages and critical aspects of the combination are presented.

Block Krylov Subspace Recycling for Shifted Systems with Unrelated Right-Hand Sides

Many Krylov subspace methods for shifted linear systems take advantage of the invariance of the Krylov subspace under a shift of the matrix. However, exploiting this fact in the non-Hermitian case introduces restrictions; e.g., initial residuals must be collinear and this collinearity must be maintained at restart. Thus we cannot simultaneously solve shifted systems with unrelated right-hand sides using this strategy, and all shifted residuals cannot be simultaneously minimized over a Krylov subspace such that collinearity is maintained. It has been shown that this renders them generally incompatible with techniques of subspace recycling [K. M. Soodhalter, D. B. Szyld, and F. Xue, Appl. Numer. Math., 81C (2014), pp. 105--118]. This problem, however, can be overcome. By interpreting a family of shifted systems as one Sylvester equation, we can take advantage of the known “shift invariance” of the Krylov subspace generated by the Sylvester operator. Thus we can simultaneously solve all systems over one block Krylov subspace using the full orthogonalization method or GMRES-type methods, even when they have unrelated right-hand sides. Because residual collinearity is no longer a requirement at restart, these methods are fully compatible with subspace recycling techniques. Furthermore, we realize the benefits of block sparse matrix operations which arise in the context of high-performance computing applications. In this paper, we discuss exploiting this Sylvester equation point of view which has yielded methods for shifted systems which are compatible with unrelated right-hand sides. From this, we propose a recycled GMRES method for simultaneous solution of shifted systems. Numerical experiments demonstrate the effectiveness of the methods.

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.

Krylov subspace recycling for sequences of shifted linear systems

We study the use of Krylov subspace recycling for the solution of a sequence of slowly-changing families of linear systems, where each family consists of shifted linear systems that differ in the coefficient matrix only by multiples of the identity. Our aim is to explore the simultaneous solution of each family of shifted systems within the framework of subspace recycling, using one augmented subspace to extract candidate solutions for all the shifted systems. The ideal method would use the same augmented subspace for all systems and have fixed storage requirements, independent of the number of shifted systems per family. We show that a method satisfying both requirements cannot exist in this framework.As an alternative, we introduce two schemes. One constructs a separate deflation space for each shifted system but solves each family of shifted systems simultaneously. The other builds only one recycled subspace and constructs approximate corrections to the solutions of the shifted systems at each cycle of the iterative linear solver while only minimizing the base system residual. At convergence of the base system solution, we apply the method recursively to the remaining unconverged systems. We present numerical examples involving systems arising in lattice quantum chromodynamics.

Subspace recycling accelerates the parametric macro‐modeling of MEMS

SUMMARYA fast computational technique that speeds up the process of parametric macro‐model extraction is proposed. An efficient starting point is the technique of parametric model order reduction (PMOR). The key step in PMOR is the computation of a projection matrix V, which requires the computation of multiple moment matrices of the underlying system. In turn, for each moment matrix, a linear system with multiple right‐hand sides has to be solved. Usually, a considerable number of linear systems must be solved when the system includes more than two free parameters. If the original system is of very large size, the linear solution step is computationally expensive. In this paper, the subspace recycling algorithm outer generalized conjugate residual method combined with generalized minimal residual method with deflated restarting (GCRO‐DR), is considered as a basis to solve the sequence of linear systems. In particular, two more efficient recycling algorithms, G‐DRvar1 and G‐DRvar2, are proposed. Theoretical analysis and simulation results show that both the GCRO‐DR method and its variants G‐DRvar1 and G‐DRvar2 are very efficient when compared with the standard solvers. Furthermore, the presented algorithms overcome the bottleneck of a recently proposed subspace recycling method the modified Krylov recycling generalized minimal residual method. From these subspace recycling algorithms, a PMOR process for macro‐model extraction can be significantly accelerated. Copyright © 2013 John Wiley & Sons, Ltd.

Parallel finite element density functional computations exploiting grid refinement and subspace recycling

In this communication computational methods that facilitate finite element analysis of density functional computations are developed. They are: (i) h–adaptive grid refinement techniques that reduce the total number of degrees of freedom in the real space grid while improving on the approximate resolution of the wanted solution; and (ii) subspace recycling of the approximate solution in self-consistent cycles with the aim of improving the performance of the generalized eigenproblem solver. These techniques are shown to give a convincing speed-up in the computation process by alleviating the overhead normally associated with computing systems with many degrees-of-freedom.

Multi-dimensional gyrokinetic parameter studies based on eigenvalue computations

Plasma microinstabilities, which can be described in the framework of the linear gyrokinetic equations, are routinely computed in the context of stability analyses and transport predictions for magnetic confinement fusion experiments. The GENE code, which solves the gyrokinetic equations, has been coupled to the SLEPc package for an efficient iterative, matrix-free, and parallel computation of rightmost eigenvalues. This setup is presented, including the preconditioner which is necessary for the newly implemented Jacobi–Davidson solver. The fast computation of instabilities at a single parameter set is exploited to make parameter scans viable, that is to compute the solution at many points in the parameter space. Several issues related to parameter scans are discussed, such as an efficient parallelization over parameter sets and subspace recycling.

Multisplitting for regularized least squares with Krylov subspace recycling

SUMMARYThe method of multisplitting (MS), implemented as a restricted additive Schwarz type algorithm, is extended for the solution of regularized least squares problems. The presented non‐stationary version of the algorithm uses dynamic updating of the weights applied to the subdomains in reconstituting the global solution. Standard convergence results follow from extensive prior literature on linear MS schemes. Additional convergence results on nonstationary iterations yield convergence conditions for the presented nonstationary MS algorithm. The global iteration uses repeated solves of local problems with changing right hand sides but a fixed system matrix. These problems are solved inexactly using a conjugate gradient least squares algorithm which provides a seed Krylov subspace. Recycling of the seed system Krylov subspace to obtain the solutions of subsequent nearby systems of equations improves the overall efficiency of the MS algorithm, and is apparently novel in this context. The obtained projected solution is not always of sufficient accuracy to satisfy a reasonable inner convergence condition on the local solution. Improvements to accuracy may be achieved by reseeding the solution space either every few steps, or when the successive right hand sides are sufficiently close as measured by a provided tolerance. Restarting and augmenting the solution space are also discussed. Any time a new space is generated it is used for subsequent steps. Numerical simulations validate the use of the recycling algorithm. These numerical experiments use the standard reconstruction of the two dimensional Shepp–Logan phantom, as well as a two dimensional problem from seismic tomography. Copyright © 2011 John Wiley & Sons, Ltd.

A Flexible Generalized Conjugate Residual Method with Inner Orthogonalization and Deflated Restarting

This work is concerned with the development and study of a minimum residual norm subspace method based on the generalized conjugate residual method with inner orthogonalization (GCRO) method that allows flexible preconditioning and deflated restarting for the solution of nonsymmetric or non-Hermitian linear systems. First we recall the main features of flexible generalized minimum residual with deflated restarting (FGMRES-DR), a recently proposed algorithm of the same family but based on the GMRES method. Next we introduce the new inner-outer subspace method named FGCRO-DR. A theoretical comparison of both algorithms is then made in the case of flexible preconditioning. It is proved that FGCRO-DR and FGMRES-DR are algebraically equivalent if a collinearity condition is satisfied. While being nearly as expensive as FGMRES-DR in terms of computational operations per cycle, FGCRO-DR offers the additional advantage to be suitable for the solution of sequences of slowly changing linear systems (where both the matrix and right-hand side can change) through subspace recycling. Numerical experiments on the solution of multidimensional elliptic partial differential equations show the efficiency of FGCRO-DR when solving sequences of linear systems.

Fast inexact subspace iteration for generalized eigenvalue problems with spectral transformation

We study inexact subspace iteration for solving generalized non-Hermitian eigenvalue problems with spectral transformation, with focus on a few strategies that help accelerate preconditioned iterative solution of the linear systems of equations arising in this context. We provide new insights into a special type of preconditioner with “tuning” that has been studied for this algorithm applied to standard eigenvalue problems. Specifically, we propose an alternative way to use the tuned preconditioner to achieve similar performance for generalized problems, and we show that these performance improvements can also be obtained by solving an inexpensive least squares problem. In addition, we show that the cost of iterative solution of the linear systems can be further reduced by using deflation of converged Schur vectors, special starting vectors constructed from previously solved linear systems, and iterative linear solvers with subspace recycling. The effectiveness of these techniques is demonstrated by numerical experiments.