Federico Poloni

Vanni Noferini, Lauri Nyman, Federico Poloni

Given a square complex matrix $A$, we tackle the problem of finding the nearest matrix with multiple eigenvalues or, equivalently when $A$ had distinct eigenvalues, the nearest defective matrix. To this goal, we extend the general framework described in [M. Gnazzo, V. Noferini, L. Nyman, F. Poloni, \emph{Riemann-Oracle: A general-purpose Riemannian optimizer to solve nearness problems in matrix theory}, Found. Comput. Math., To appear] and based on variable projection and Riemannian optimization, allowing the ambient manifold to simultaneously track left and right eigenvectors. Our method also allows us to impose arbitrary complex-linear constraints on either the perturbation or the perturbed matrix; this can be useful to study structured eigenvalue condition numbers. We present numerical experiments, comparing with preexisting algorithms.

Dario A. Bini, Beatrice Meini, Federico Poloni

We present some advances, both from a theoretical and from a computational point of view, on a quadratic vector equation (QVE) arising in Markovian Binary Trees. Concerning the theoretical advances, some irreducibility assumptions are relaxed, and the minimality of the solution of the QVE is expressed in terms of properties of the Jacobian of a suitable function. From the computational point of view, we elaborate on the Perron vector-based iteration proposed in [http://arxiv.org/abs/1006.0577]. In particular we provide a condition which ensures that the Perron iteration converges to the sought solution of the QVE. Moreover we introduce a variant of the algorithm which consists in applying the Newton method instead of a fixed-point iteration. This method has the same convergence behaviour as the Perron iteration, since its convergence speed increases for close-to-critical problems. Moreover, unlike the Perron iteration, the method has a quadratic convergence. Finally, we show that it is possible to alter the bilinear form defining the QVE in several ways without changing the solution. This modification has an impact on convergence speed of the algorithms.

Beatrice Meini, Federico Poloni

We consider the multilinear pagerank problem studied in [Gleich, Lim and Yu, Multilinear Pagerank, 2015], which is a system of quadratic equations with stochasticity and nonnegativity constraints. We use the theory of quadratic vector equations to prove several properties of its solutions and suggest new numerical algorithms. In particular, we prove the existence of a certain minimal solution, which does not always coincide with the stochastic one that is required by the problem. We use an interpretation of the solution as a Perron eigenvector to devise new fixed-point algorithms for its computation, and pair them with a homotopy continuation strategy. The resulting numerical method is more reliable than the existing alternatives, being able to solve a larger number of problems.

Beatrice Meini, Federico Poloni

We propose a novel numerical method for solving a quadratic vector equation arising in Markovian Binary Trees. The numerical method consists in a fixed point iteration, expressed by means of the Perron vectors of a sequence of nonnegative matrices. A theoretical convergence analysis is performed. The proposed method outperforms the existing methods for close-to-critical problems.

Federico Poloni

We study in an unified fashion several quadratic vector and matrix equations with nonnegativity hypotheses. Specific cases of such problems (QBD equations, nonsymmetric algebraic Riccati equations, Lu's simple equation, Markovian binary trees equations) have been studied extensively in the past by several authors. Many of the results appearing here have already been proved for one or more of the single instances of the problems, resorting to specific characteristics of the problem. In some cases the proofs we present here are mere rewriting of the original proofs with a little change of notation to adapt them to our framework, but in some cases we are effectively able to remove some hypotheses and generalize the results by abstracting the specific aspects of each problem.

Federico Greco, Bruno Iannazzo, Federico Poloni

It is proved that among the rational iterations locally converging with order s>1 to the sign function, the Padé iterations and their reciprocals are the unique rationals with the lowest sum of the degrees of numerator and denominator.

Federico Poloni, Timo Reis

We introduce a numerical method for the numerical solution of the so-called Lur'e matrix equations that arise in balancing-related model reduction and linear-quadratic infinite time horizon optimal control. Based on the fact that the set of solutions can be characterized in terms of deflating subspaces of even matrix pencils, an iterative scheme is derived that converges linearly to the maximal solution.

Tayyebe Haqiri, Federico Poloni

We describe a procedure based on the Krawczyk method to compute a verified enclosure for the stabilizing solution of a continuous-time algebraic Riccati equation $A^*X+XA+Q=XGX$ building on the work of [B.~Hashemi, \emph{SCAN} 2012] and adding several modifications to the Krawczyk procedure. We show that after these improvements the Krawczyk method reaches results comparable with the current state-of-the-art algorithm [Miyajima, \emph{Jpn. J. Ind. Appl. Math} 2015], and surpasses it in some examples. Moreover, we introduce a new direct method for verification which has a cubic complexity in term of the dimension of $X$, employing a fixed-point formulation of the equation inspired by the ADI procedure. The resulting methods are tested on a number of standard benchmark examples.

Bruno Iannazzo, Federico Poloni

The worst situation in computing the minimal nonnegative solution of a nonsymmetric algebraic Riccati equation associated with an M-matrix occurs when the corresponding linearizing matrix has two very small eigenvalues, one with positive and one with negative real part. When both these eigenvalues are exactly zero, the problem is called critical or null recurrent. While in this case the problem is ill-conditioned and the convergence of the algorithms based on matrix iterations is slow, there exist some techniques to remove the singularity and transform the problem to a well-behaved one. Ill-conditioning and slow convergence appear also in close-to-critical problems, but when none of the eigenvalues is exactly zero the techniques used for the critical case cannot be applied. In this paper, we introduce a new method to accelerate the convergence properties of the iterations also in close-to-critical cases, by working on the invariant subspace associated with the problematic eigenvalues as a whole. We present a theoretical analysis and several numerical experiments which confirm the efficiency of the new method.

Nigel Bean, Giang T. Nguyen, Federico Poloni

Fluid queues are mathematical models frequently used in stochastic modelling. Their stationary distributions involve a key matrix recording the conditional probabilities of returning to an initial level from above, often known in the literature as the matrix $Ψ$. Here, we present a probabilistic interpretation of the family of algorithms known as \emph{doubling}, which are currently the most effective algorithms for computing the return probability matrix $Ψ$. To this end, we first revisit the links described in \cite{ram99, soares02} between fluid queues and Quasi-Birth-Death processes; in particular, we give new probabilistic interpretations for these connections. We generalize this framework to give a probabilistic meaning for the initial step of doubling algorithms, and include also an interpretation for the iterative step of these algorithms. Our work is the first probabilistic interpretation available for doubling algorithms.

Fernando De Terán, Bruno Iannazzo, Federico Poloni et al.

We provide necessary and sufficient conditions for the generalized $\star$-Sylvester matrix equation, $AXB + CX^\star D = E$, to have exactly one solution for any right-hand side E. These conditions are given for arbitrary coefficient matrices $A, B, C, D$ (either square or rectangular) and generalize existing results for the same equation with square coefficients. We also review the known results regarding the existence and uniqueness of solution for generalized Sylvester and $\star$-Sylvester equations.

Federico Poloni

The aim of this paper is twofold. First, we introduce a new class of linearizations, based on the generalization of a construction used in polynomial algebra to find the zeros of a system of (scalar) polynomial equations. We show that one specific linearization in this class, which is constructed naturally from the QR factorization of the matrix obtained by stacking the coefficients of $A(x)$, has good conditioning and stability properties. Moreover, while analyzing this class, we introduce a general technique to derive new linearizations from existing ones. This technique generalizes some ad-hoc arguments used in dealing with the existing linearization classes, and can hopefully be used to derive a simpler and more general theory of linearizations. This technique relates linearizations to \emph{pencil arithmetic}, a technique used in solving matrix equations that allows to extend some algebraic operations from matrix to matrix pencils.

Elias Jarlebring, Federico Poloni

The delay Lyapunov equation is an important matrix boundary-value problem which arises as an analogue of the Lyapunov equation in the study of time-delay systems $\dot{x}(t) = A_0x(t)+A_1x(t-τ)+B_0u(t)$. We propose a new algorithm for the solution of the delay Lyapunov equation. Our method is based on the fact that the delay Lyapunov equation can be expressed as a linear system of equations, whose unknown is the value $U(τ/2)\in\mathbb{R}^{n\times n}$, i.e., the delay Lyapunov matrix at time $τ/2$. This linear matrix equation with $n^2$ unknowns is solved by adapting a preconditioned iterative method such as GMRES. The action of the $n^2\times n^2$ matrix associated to this linear system can be computed by solving a coupled matrix initial-value problem. A preconditioner for the iterative method is proposed based on solving a T-Sylvester equation $MX+X^TN=C$, for which there are methods available in the literature. We prove that the preconditioner is effective under certain assumptions. The efficiency of the approach is illustrated by applying it to a time-delay system stemmingfrom the discretization of a partial differential equation with delay. Approximate solutions to this problem can be obtained for problems of size up to $n\approx 1000$, i.e., a linear system with $n^2\approx 10^6$ unknowns, a dimension which is outside of the capabilities of the other existing methods for the delay Lyapunov equation.

Giang T. Nguyen, Federico Poloni

Markov-modulated fluid queues are popular stochastic processes frequently used for modelling real-life applications. An important performance measure to evaluate in these applications is their steady-state behaviour, which is determined by the stationary density. Computing it requires solving a (nonsymmetric) M-matrix algebraic Riccati equation, and indeed computing the stationary density is the most important application of this class of equations. Xue, Xu and Li [Numer. Math., 2012] provided a componentwise first-order perturbation analysis of this equation, proving that the solution can be computed to high relative accuracy even in the smallest entries, and suggested several algorithms for computing it. An important step in all proposed algorithms is using so-called triplet representations, which are special representations for M-matrices that allow for a high-accuracy variant of Gaussian elimination, the GTH-like algorithm. However, triplet representations for all the M-matrices needed in the algorithm were not found explicitly. This can lead to an accuracy loss that prevents the algorithms to converge in the componentwise sense. In this paper, we focus on the structured doubling algorithm, the most efficient among the proposed methods in Xue et al., and build upon their results, providing (i) explicit and cancellation-free expressions for the needed triplet representations, allowing the algorithm to be performed in a really cancellation-free fashion; (ii) an algorithm to evaluate the final part of the computation to obtain the stationary density; and (iii) a componentwise error analysis for the resulting algorithm, the first explicit one for this class of algorithms. We also present numerical results to illustrate the accuracy advantage of our method over standard (normwise-accurate) algorithms.

Federico Poloni

In this paper, we analyze the process of "assembling" new matrix geometric means from existing ones, through function composition or limit processes. We show that for n=4 a new matrix mean exists which is simpler to compute than the existing ones. Moreover, we show that for n>4 the existing proving strategies cannot provide a mean computationally simpler than the existing ones.

Federico Poloni

We propose a new O(n)-space implementation of the GKO-Cauchy algorithm for the solution of linear systems with Cauchy-like matrix. Despite its slightly higher computational cost, this new algorithm makes a more efficient use of the processor cache memory. Thus, for matrices of size larger than about 500-1000, it outperforms the existing algorithms. We present an applicative case of Cauchy-like matrices with non-reconstructible main diagonal. In this special instance, the O(n) space algorithms can be adapted nicely to provide an efficient implementation of basic linear algebra operations in terms of the low displacement-rank generators.

Dario A. Bini, Federico Poloni

We review some known inclusion results for the roots of a polynomial, and adapt them to a conjecture recently presented by S. A. Vavasis. In particular, we provide strict upper and lower bounds to the distance of the closest root of a polynomial p(z) from a given root of p'(z).