697 research outputs found
Efficient Ehrlich–Aberth iteration for finding intersections of interpolating polynomials and rational functions
We analyze the problem of carrying out an efficient iteration to approximate the eigenvalues of some rank structured pencils obtained as linearization of sums of polynomials and rational functions expressed in (possibly different) interpolation bases. The class of linearizations that we consider has been introduced by Robol, Vandebril and Van Dooren in [17]. We show that a traditional QZ iteration on the pencil is both asymptotically slow (since it is a cubic algorithm in the size of the matrices) and sometimes not accurate (since in some cases the deflation of artificially introduced infinite eigenvalues is numerically difficult). To solve these issues we propose to use a specifically designed Ehrlich–Aberth iteration that can approximate the eigenvalues in O(kn2) flops, where k is the average number of iterations per eigenvalue, and n the degree of the linearized polynomial. We suggest possible strategies for the choice of the initial starting points that make k asymptotically smaller than O(n), thus making this method less expensive than the QZ iteration. Moreover, we show in the numerical experiments that this approach does not suffer of numerical issues, and accurate results are obtained
Geometric means of structured matrices
The geometric mean of positive definite matrices is usually identified
with the Karcher mean, which possesses all properties—generalized from the
scalar case— a geometric mean is expected to satisfy. Unfortunately,
the Karcher mean is typically not structure preserving, and destroys,
e.g., Toeplitz and band structures, which emerge in many
applications. For this reason, the Karcher mean is not always
recommended for modeling averages of structured matrices. In this
article a new definition of a geometric mean for structured matrices
is introduced, its properties are outlined, algorithms for its
computation, and numerical experiments are provided. In the Toeplitz case an existing mean based on the Kahler metric is analyzed for comparison.sponsorship: This work was
partially supported by MIUR grant number 2002014121; by the Research
Council KU Leuven, projects OT/11/055 (Spectral Properties of
Perturbed Normal Matrices and their Applications), CoE EF/05/006
Optimization in Engineering (OPTEC); by the Fund for Scientific
Research—Flanders (Belgium) project G034212N (Reestablishing
Smoothness for Matrix Manifold Optimization via Resolution of
Singularities); and by the Interuniversity Attraction Poles
Programme, initiated by the Belgian State, Science Policy Office,
Belgian Network DYSCO (Dynamical Systems, Control, and
Optimization).status: Publishe
Fast and backward stable computation of roots of polynomials, Part II: Backward error analysis; companion matrix and companion pencil
This work is a continuation of work by [J. L. Aurentz, T. Mach, R. Vandebril, and D. S. Watkins, J. Matrix Anal. Appl., 36 (2015), pp. 942–973]. In that paper we introduced a companion QR algorithm that finds the roots of a polynomial by computing the eigenvalues of the companion matrix in O(n2) time using O(n) memory. We proved that the method is backward stable. Here we introduce, as an alternative, a companion QZ algorithm that solves a generalized eigenvalue problem for a companion pencil. More importantly, we provide an improved backward error analysis that takes advantage of the special structure of the problem. The improvement is also due, in part, to an improvement in the accuracy (in both theory and practice) of the turnover operation, which is the key component of our algorithms. We prove that for the companion QR algorithm, the backward error on the polynomial coefficients varies linearly with the norm of the polynomial’s vector of coefficients. Thus, the companion QR algorithm has a smaller backward error than the unstructured QR algorithm (used by MATLAB’s roots command, for example), for which the backward error on the polynomial coefficients grows quadratically with the norm of the coefficient vector. The companion QZ algorithm has the same favorable backward error as companion QR, provided that the polynomial coefficients are properly scaled
A framework for structured linearizations of matrix polynomials in various bases
We present a framework for the construction of linearizations for scalar and matrix polynomials based on dual bases which, in the case of orthogonal polynomials, can be described by the associated recurrence relations. The framework provides an extension of the classical linearization theory for polynomials expressed in nonmonomial bases and allows us to represent polynomials expressed in product families, that is, as a linear combination of elements of the form φi(λ) ψj (λ), where φi(λ) andψj (λ) can either be polynomial bases or polynomial families which satisfy some mild assumptions. We show that this general construction can be used for many different purposes. Among them, we show how to linearize sums of polynomials and rational functions expressed in different bases. As an example, this allows us to look for intersections of functions interpolated on different nodes without converting them to the same basis. We then provide some constructions for structured linearizations for ∗-even and ∗-palindromic matrix polynomials. The extension of these constructions to ∗-odd and ∗-antipalindromic of odd degree is discussed and follows immediately from the previous results
An implicit multishift QR-algorithm for Hermitian plus low rank matrices
Hermitian plus possibly unhermitian low rank matrices can be efficiently reduced into Hessenberg form. The resulting Hessenberg matrix can still be written as the sum of a Hermitian plus low rank matrix. In this paper we develop a new implicit multishift QR-algorithm for Hessenberg matrices, which are the sum of a Hermitian plus a possibly non-Hermitian low rank correction.The proposed algorithm exploits both the symmetry and low rank structure to obtain a QR-step involving only O(n) floating point operations instead of the standard O(n^2) operations needed for performing a QR-step on a Hessenberg matrix. The algorithm is based on a suitable O(n) representation of the Hessenberg matrix. The low rank parts present in both the Hermitian and low rank part of the sum are compactly stored by a sequence of Givens transformations and few vectors. Due to the new representation, we cannot apply classical deflation techniques for Hessenberg matrices. A new, efficient technique is developed to overcome this problem. Some numerical experiments based on matrices arising in applications are performed.The experiments illustrate effectiveness and accuracy of both the -algorithm and the newly developed deflation technique
When is a matrix unitary or Hermitian plus low rank?
Hermitian and unitary matrices are two representatives of the class of normal matrices whose full eigenvalue decomposition can be stably computed in quadratic computing complexity once the matrix has been reduced, for instance, to tridiagonal or Hessenberg form. Recently, fast and reliable eigensolvers dealing with low-rank perturbations of unitary and Hermitian matrices have been proposed. These structured eigenvalue problems appear naturally when computing roots, via confederate linearizations, of polynomials expressed in, for example, the monomial or Chebyshev basis. Often, however, it is not known beforehand whether or not a matrix can be written as the sum of a Hermitian or unitary matrix plus a low-rank perturbation. In this paper, we give necessary and sufficient conditions characterizing the class of Hermitian or unitary plus low-rank matrices. The number of singular values deviating from 1 determines the rank of a perturbation to bring a matrix to unitary form. A similar condition holds for Hermitian matrices; the eigenvalues of the skew-Hermitian part differing from 0 dictate the rank of the perturbation. We prove that these relations are linked via the Cayley transform. Then, based on these conditions, we identify the closest Hermitian or unitary plus rank k matrix to a given matrix A, in Frobenius and spectral norm, and give a formula for their distance from A. Finally, we present a practical iteration to detect the low-rank perturbation. Numerical tests prove that this straightforward algorithm is effective
The seriation problem in the presence of a double Fiedler value
Seriation is a problem consisting of seeking the best enumeration order of a set of units whose interrelationship is described by a bipartite graph, that is, a graph whose nodes are partitioned in two sets and arcs only connect nodes in different groups. An algorithm for spectral seriation based on the use of the Fiedler vector of the Laplacian matrix associated to the problem was developed by Atkins et al., under the assumption that the Fiedler value is simple. In this paper, we analyze the case in which the Fiedler value of the Laplacian is not simple, discuss its effect on the set of the admissible solutions, and study possible approaches to actually perform the computation. Examples and numerical experiments illustrate the effectiveness of the proposed methods.sponsorship: Anna Concas, Caterina Fenu, and Giuseppe Rodriguez were partially supported by Regione Autonoma della Sardegna research project "Algorithms and Models for Imaging Science (AMIS)"(RASSR57257, intervento finanziato con risorse FSC 2014-2020 - Patto per lo Sviluppo della RegioneSardegna) and by the INdAM-GNCS research project "Tecniche numeriche per l'analisi delle reti com-plesse e lo studio dei problemi inversi." Caterina Fenu also gratefully acknowledges Regione Autonomadella Sardegna for the financial support provided under the Operational Programme P.O.R. Sardegna F.S.E.(European Social Fund 2014-2020 - Axis III Education and Formation, Objective 10.5, Line of Activity10.5.12). The research of Raf Vandebril was partially supported by the Research Council KU Leuven,project C16/21/002 (Manifactor: Factor Analysis for Maps into Manifolds). (Regione Autonoma della Sardegna research project "Algorithms and Models for Imaging Science (AMIS)", INdAM-GNCS research project "Tecniche numeriche per l'analisi delle reti com-plesse e lo studio dei problemi inversi."|RASSR57257, Regione Autonomadella Sardegna, Research Council KU Leuven, C16/21/002)status: Published onlin
Structured backward errors in linearizations
A standard approach to compute the roots of a univariate polynomial is to compute the eigenvalues of an associated confederate matrix instead, such as, for instance the companion or comrade matrix. The eigenvalues of the confederate matrix can be computed by Francis’s QR algorithm. Unfortunately, even though the QR algorithm is provably backward stable, mapping the errors back to the original polynomial coefficients can still lead to huge errors. However, the latter statement assumes the use of a non-structure exploiting QR algorithm. In [J. Aurentz et al., Fast and backward stable computation of roots of polynomials, SIAM J. Matrix Anal. Appl., 36(3), 2015] it was shown that a structure exploiting QR algorithm for companion matrices leads to a structured backward error on the companion matrix. The proof relied on decomposing the error into two parts: a part related to the recurrence coefficients of the basis (monomial basis in that case) and a part linked to the coefficients of the original polynomial. In this article we prove that the analysis can be extended to other classes of comrade matrices. We first provide an alternative backward stability proof in the monomial basis using structured QR algorithms; our new point of view shows more explicitly how a structured, decoupled error on the confederate matrix gets mapped to the associated polynomial coefficients. This insight reveals which properties must be preserved by a structure exploiting QR algorithm to end up with a backward stable algorithm. We will show that the previously formulated companion analysis fits in this framework and we will analyze in more detail Jacobi polynomials (Comrade matrices) and Chebyshev polynomials (Colleague matrices).sponsorship: The work of Vanni Noferini was supported by an Academy of Finland grant (Suomen Akatemian paatos 331240) ; the work of Leonardo Robol was supported by an INdAM/GNCS research grant "Metodi low-rank per problemi di algebra lineare con struttura data-sparse". (Academy of Finland|331240, INdAM/GNCS research grant "Metodi low-rank per problemi di algebra lineare con struttura data-sparse", Academy of Finland (AKA)|331240)status: Publishe
A unification of unitary similarity transforms to compressed representations
In this paper a new framework for transforming arbitrary matrices to compressed representations is presented. The framework provides a generic way of transforming a matrix via unitary similarity transformations to e.g.\ Hessenberg, Hessenberg\-/like and combinations of both. The new algorithms are deduced, based on the -factorization of the original matrix. Based on manipulations with Givens transformations, all the algorithms consists of eliminating the correct set of Givens transformations, resulting in a matrix obeying the desired structural constraints. Based on this new reduction procedure we investigate further correspondences such as irreducibility, unicity of the reduction procedure and the link with (rational) Krylov methods. The unitary similarity transform to Hessenberg\-/like form as presented here, differssignificantly from the one presented in earlier work. Not only does it use less Givens transformations to obtain the desired structure, also the convergence to rational Ritz values is not observed in the standard way
A unitary similarity transform of a normal matrix to complex symmetric form
In this work a new unitary similarity transformation of a normal matrix to complex symmetric form will be discussed. A constructive proof as well as some properties and examples will be given.sponsorship: The author has a grant as “Postdoctoraal Onderzoeker” from the Fund for Scientific Research–Flanders (Belgium). The research of the author was partially supported by the Research Council K.U. Leuven, CoE EF/05/006 Optimization in Engineering (OPTEC) and by the Interuniversity Attraction Poles Programme, initiated by the Belgian State, Science Policy Office, Belgian Network DYSCO (Dynamical Systems, Control, and Optimization).status: Publishe
- …
