Melina A. Freitag

Shift-Invert Arnoldi's Method with Preconditioned Iterative Solves

We consider the computation of a few eigenvectors and corresponding eigenvalues of a large sparse nonsymmetric matrix using shift-invert Arnoldis method with and without implicit restarts. For the inner iterations we use preconditioned GMRES as the inexact iterative solver. The costs of the solves are measured by the number of inner iterations needed by the iterative solver at each outer step of the algorithm. We first extend the relaxation strategy developed by Simoncini [SIAM J. Numer. Anal., 43 (2005), pp. 1155-1174] to implicitly restarted Arnoldis method, which yields an improvement in the overall costs of the method. Secondly, we apply a new preconditioning strategy to the inner solver. We show that small rank changes to the preconditioner can produce significant savings in the total number of iterations. The combination of the new preconditioner with the relaxation strategy in implicitly restarted Arnoldi produces enhancement in the overall costs of around 50 percent in the examples considered here. Numerical experiments illustrate the theory throughout the paper.

The origin of power-law emergent scaling in large binary networks

We study the macroscopic conduction properties of large but finite binary networks with conducting bonds. By taking a combination of a spectral and an averaging based approach we derive asymptotic formulae for the conduction in terms of the component proportions p and the total number of components N. These formulae correctly identify both the percolation limits and also the emergent power-law behaviour between the percolation limits and show the interplay between the size of the network and the deviation of the proportion from the critical value of p=1/2. The results compare excellently with a large number of numerical simulations.

The calculation of the distance to a nearby defective matrix

SUMMARY The distance of a matrix to a nearby defective matrix is an important classical problem in numerical linear algebra, as it determines how sensitive or ill-conditioned an eigenvalue decomposition of a matrix is. The concept has been discussed throughout the history of numerical linear algebra, and the problem of computing the nearest defective matrix first appeared in Wilkinsons famous book on the algebraic eigenvalue problem. In this paper, a new fast algorithm for the computation of the distance of a matrix to a nearby defective matrix is presented. The problem is formulated following Alam and Bora introduced in (2005) and reduces to finding when a parameter-dependent matrix is singular subject to a constraint. The solution is achieved by an extension of the implicit determinant method introduced by Spence and Poulton in (2005). Numerical results for several examples illustrate the performance of the algorithm. Copyright

Calculating the

We propose a fast algorithm to calculate the

H_{\infty}

H_{\infty}

-norm Using the Implicit Determinant Method

-norm of a transfer matrix. The method builds on a well-known relationship between singular values of the transfer function and pure imaginary eigenvalues of a certain Hamiltonian matrix. Using this property we construct a two-parameter eigenvalue problem, where, in the generic case, the critical value corresponds to a two-dimensional Jordan block. We use the implicit determinant method which replaces the need for eigensolves by the solution of linear systems, a technique recently used in [M. A. Freitag and A. Spence, Linear Algebra Appl., 435 (2011), pp. 3189--3205] for finding the distance to instability. In this paper the method takes advantage of the structured linear systems that arise within the algorithm to obtain efficient solves using the staircase reduction. We give numerical examples and compare our method to the algorithm proposed in [N. Guglielmi, M. Gurbuzbalaban, and M. L. Overton, SIAM J. Matrix Anal. Appl., 34 (2013), pp. 709--737].

The effect of numerical model error on data assimilation

Strong constraint 4D-Variational data assimilation (4D-Var) is a method used to create an initialisation for a numerical model, that best replicates subsequent observations of the system it aims to recreate. The method does not take into account the presence of errors in the model, using the model equations as a strong constraint. This paper gives a rigorous and quantitative analysis of the errors introduced into the initialisation through the use of finite difference schemes to numerically solve the model equations. The 1D linear advection equation together with circulant boundary conditions, are chosen as the model equations of interest as they are representative of the advective processes relevant to numerical weather prediction, where 4D-Var is widely used. We consider the deterministic error introduced by finite difference approximations in the form of numerical dissipation and numerical dispersion and identify the relationship between these properties and the error in the 4D-Var initialisation. In particular, we find that a solely numerically dispersive scheme has the potential to introduce destructive interference resulting in the loss of some wavenumber components in the initialisation. Bounds for the error in the initialisation due to finite difference approximations are determined with and without observation errors. The bounds are found to depend on the smoothness of the true initial condition we wish to recover and the numerically dissipative and dispersive properties of the scheme. Numerical results are presented to demonstrate the effectiveness of the bounds. These lead to the conclusion that there exists a critical number of discretisation points when considering full sets of observations, where the effects of both the considered numerical model error and observational errors on the initialisation are minimised. The numerically dissipative and dispersive properties of the finite difference schemes also have the potential to alter the properties of the noise found in observations. Correlated noise structures may be introduced into the 4D-Var initialisation as a result. We determine when this occurs for observational errors in the form of additive white noise and find that the effect is reduced through the use of numerically non-dissipative finite difference schemes.

A low-rank approach to the solution of weak constraint variational data assimilation problems

Abstract Weak constraint four-dimensional variational data assimilation is an important method for incorporating data (typically observations) into a model. The linearised system arising within the minimisation process can be formulated as a saddle point problem. A disadvantage of this formulation is the large storage requirements involved in the linear system. In this paper, we present a low-rank approach which exploits the structure of the saddle point system using techniques and theory from solving large scale matrix equations. Numerical experiments with the linear advection–diffusion equation, and the non-linear Lorenz-95 model demonstrate the effectiveness of a low-rank Krylov subspace solver when compared to a traditional solver.

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.

Tuned preconditioners for inexact two-sided inverse and Rayleigh quotient iteration

SUMMARY Convergence results are provided for inexact two-sided inverse and Rayleigh quotient iteration, which extend the previously established results to the generalized non-Hermitian eigenproblem and inexact solves with a decreasing solve tolerance. Moreover, the simultaneous solution of the forward and adjoint problem arising in two-sided methods is considered, and the successful tuning strategy for preconditioners is extended to two-sided methods, creating a novel way of preconditioning two-sided algorithms. Furthermore, it is shown that inexact two-sided Rayleigh quotient iteration and the inexact two-sided Jacobi-Davidson method (without subspace expansion) applied to the generalized preconditioned eigenvalue problem are equivalent when a certain number of steps of a Petrov–Galerkin–Krylov method is used and when this specific tuning strategy is applied. Copyright