Search

We present some accelerated variants of fixed point iterations for computing the minimal non-negative solution of the unilateral matrix equation associated with an M/G/1-type Markov chain. These variants derive from certain staircase regular splittings of the block Hessenberg M-matrix associated with the Markov chain. By exploiting the staircase profile, we introduce a two-step fixed point iteration. The iteration can be further accelerated by computing a weighted average between the approximations obtained at two consecutive steps. The convergence of the basic two-step fixed point iteration and of its relaxed modification is proved. Our theoretical analysis, along with several numerical experiments, shows that the proposed variants generally outperform the classical iterations.

We consider the problem of computing the minimal non-negative solution $G$ of the nonlinear matrix equation $X=\sum _{i=-1}^\infty A_iX^{i+1}$ where $A_i$, for $i\geqslant -1$, are non-negative square matrices such that $\sum _{i=-1}^\infty A_i$ is stochastic. This equation is fundamental in the analysis of M/G/1-type Markov chains, since the matrix $G$ provides probabilistic measures of interest. A new family of fixed point iterations for the numerical computation of $G$, which includes the classical iterations, is introduced. A detailed convergence analysis proves that the iterations in the new class converge faster than the classical iterations. Numerical experiments confirm the effectiveness of our extension.

In discrete time, $\ell$-blocks of red lights are separated by $\ell$-blocks of green lights. Cars arrive at random. \ We seek the distribution of maximum line length of idle cars, and justify conjectured probabilistic asymptotics algebraically for $2\leq\ell\leq3$ and numerically for $\ell\geq4$.

A new measure $c(e)$ of the centrality of an edge $e$ in an undirected graph $G$ is introduced. It is based on the variation of the Kemeny constant of the graph after removing the edge $e$. The new measure is designed in such a way that the Braess paradox is avoided. A numerical method for computing $c(e)$ is introduced and a regularization technique is designed in order to deal with cut-edges and disconnected graphs. Numerical experiments performed on synthetic tests and on real road networks show that this measure is particularly effective in revealing bottleneck roads whose removal would greatly reduce the connectivity of the network.

A Markov modulated Brownian motion (MMBM) is a substantial generalization of the classical Brownian motion and is obtained by allowing the Brownian parameters to be modulated by an underlying Marko...

We identify a relationship between the solutions of a nonsymmetric algebraic$$T$$T-Riccati equation ($$T$$T-NARE) and the deflating subspaces of a palindromic matrix pencil, obtained by arranging the coefficients of the$$T$$T-NARE. The interplay between$$T$$T-NAREs and palindromic pencils allows one to derive both theoretical properties of the solutions of the equation, and new methods for its numerical solution. In particular, we propose methods based on the (palindromic) QZ algorithm and the doubling algorithm, whose effectiveness is demonstrated by several numerical tests.

Let $a(z)=\sum_{i\in\mathbb Z}a_iz^i$ be a complex valued function defined for $|z|=1$, such that $\sum_{i\in\mathbb Z}|ia_i

We consider the problem of solving a rational matrix equation arising in the solution of G-networks. We propose and analyze two numerical methods: a fixed point iteration and the Newton–Raphson method. The fixed point iteration is shown to be globally convergent with linear convergence rate, while the Newton method is shown to have a local convergence, with quadratic convergence rate. Numerical experiments show the effectiveness of the proposed methods.

Quadratic matrix equations of the kind $A_1X^2+A_0X+A_{-1}=X$ are encountered in the analysis of Quasi--Birth-Death stochastic processes where the solution of interest is the minimal nonnegative solution $G$. In many queueing models, described by random walks in the quarter plane, the coefficients $A_1,A_0,A_{-1}$ are infinite tridiagonal matrices with an almost Toeplitz structure. Here, we analyze some fixed point iterations, including Newton's iteration, for the computation of $G$ and introduce effective algorithms and acceleration strategies which fully exploit the Toeplitz structure of the matrix coefficients and of the current approximation. Moreover, we provide a structured perturbation analysis for the solution $G$. The results of some numerical experiments which demonstrate the effectiveness of our approach are reported.

The treatment of two-dimensional random walks in the quarter plane leads to Markov processes which involve semi-infinite matrices having Toeplitz or block Toeplitz structure plus a low-rank correction. Finding the steady state probability distribution of the process requires to perform operations involving these structured matrices. We propose an extension of the framework of [5] which allows to deal with more general situations such as processes involving restart events. This is motivated by the need for modeling processes that can incur in unexpected failures like computer system reboots. Algebraically, this gives rise to corrections with infinite support that cannot be treated using the tools currently available in the literature. We present a theoretical analysis of an enriched Banach algebra that, combined with appropriate algorithms, enables the numerical treatment of these problems. The results are applied to the solution of bidimensional Quasi-Birth-Death processes with infinitely many phases which model random walks in the quarter plane, relying on the matrix analytic approach. This methodology reduces the problem to solving a quadratic matrix equation with coefficients of infinite size. We provide conditions on the transition probabilities which ensure that the solution of interest of the matrix equation belongs to the enriched algebra. The reliability of our approach is confirmed by extensive numerical experimentation on some case studies.

Denote by $\mathcal{W}_1$ the set of complex valued functions of the form $a(z)=\sum_{i=-\infty}^{+\infty}a_iz^i$ which are continuous on the unit circle, and such that $\sum_{i=-\infty}^{+\infty}|ia_i

Given a square matrix A, Brauer’s theorem [Duke Math. J. 19 (1952), 75–91] shows how to modify one single eigenvalue of A via a rank-one perturbation, without changing any of the remaining eigenvalues. We reformulate Brauer’s theorem in functional form and provide extensions to matrix polynomials and to matrix Laurent series A(z) together with generalizations to shifting a set of eigenvalues. We provide conditions under which the modified function \( \tilde{A}(z) \) has a canonical factorization \( \tilde{A}(z)=\tilde{U}(z)\tilde{L}(z^{-1}) \) and we provide explicit expressions of the factors \( \tilde{U}(z) \) and \( \tilde{L}(z) \). Similar conditions and expressions are given for the factorization of \( \tilde{A}(z^{-1}) \). Some applications are discussed.

In their 1960 book on finite Markov chains, Kemeny and Snell established that a certain sum is invariant. The value of this sum has become known as {\it Kemeny's constant}. Various proofs have been given over time, some more technical than others. We give here a very simple physical justification, which extends without a hitch to continuous-time Markov chains on a finite state space. For Markov chains with denumerably infinite state space, the constant may be infinite and even if it is finite, there is no guarantee that the physical argument will hold. We show that the physical interpretation does go through for the special case of a birth-and-death process with a finite value of Kemeny's constant. Keywords: Kemeny's constant; discrete-time Markov chains; continuous-time Markov chains; passage times; deviation matrix.

We revisit the shift technique applied to Quasi-Birth and Death (QBD) processes He et al. (2001) [13] in functional form by bringing the attention to the existence and properties of canonical factorizations. To this regard, we prove new results concerning the solutions of the quadratic matrix equations associated with the QBD. These results find applications to the solution of the Poisson equation for QBDs.

Let $a(z)=\sum_{i\in\mathbb Z}a_iz^i$ be a complex valued continuous function, defined for $|z|=1$, such that $\sum_{i=-\infty}^{+\infty}|ia_i

A class of nonsymmetric algebraic Riccati equation, where one of the two linear coefficients is block diagonal, is studied. These equations arise in the modeling of an adaptive MMAP[K]/PH[K]/1 queue. Some theoretical results are proved, and two new algorithms are introduced that exploit the diagonal structure of the linear coefficient.

We consider the Poisson equation $(I-P)\boldsymbol{u}=\boldsymbol{g}$, where $P$ is the transition matrix of a Quasi-Birth-and-Death (QBD) process with infinitely many levels, $\bm g$ is a given infinite dimensional vector and $\bm u$ is the unknown. Our main result is to provide the general solution of this equation. To this purpose we use the block tridiagonal and block Toeplitz structure of the matrix $P$ to obtain a set of matrix difference equations, which are solved by constructing suitable resolvent triples.

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.

We study the properties of palindromic quadratic matrix polynomials φ ( z ) = P + Qz + Pz 2 , i.e., quadratic polynomials where the coefficients P and Q are square matrices, and where the constant and the leading coefficients are equal. We show that, for suitable choices of the matrix coefficients P and Q, it is possible to characterize by means of φ ( z ) well known matrix functions, namely the matrix square root, the matrix polar factor, the matrix sign and the geometric mean of two matrices. Finally we provide some integral representations of these matrix functions.

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.

Cyclic reduction is an algorithm invented by G. H. Golub and R. W. Hockney in the mid 1960s for solving linear systems related to the finite differences discretization of the Poisson equation over a rectangle. Among the algorithms of Gene Golub, it is one of the most versatile and powerful ever created. Recently, it has been applied to solve different problems from different applicative areas. In this paper we survey the main features of cyclic reduction, relate it to properties of analytic functions, recall its extension to solving more general finite and infinite linear systems, and different kinds of nonlinear matrix equations, including algebraic Riccati equations, with applications to Markov chains, queueing models and transport theory. Some new results concerning the convergence properties of cyclic reduction and its applicability are proved under very weak assumptions. New formulae for overcoming breakdown are provided.

We provide a simple convergence proof for the cyclic reduction algorithm for M/G/1 type Markov chains together with a probabilistic interpretation which helps to understand better the relationships between logarithmic reduction and cyclic reduction.

The cyclic reduction algorithm is specialized to palindromic matrix polynomials and a complete analysis of applicability and convergence is provided. The resulting iteration is then related to other algorithms as the evaluation/interpolation at the roots of unity of a certain Laurent matrix polynomial, the trapezoidal rule for a certain integral and an algorithm based on the finite sections of a tridiagonal block Toeplitz matrix.

The Erlangian approximation of Markovian fluid queues leads to the problem of computing the matrix exponential of a subgenerator having a block-triangular, block-Toeplitz structure. To this end, we propose some algorithms which exploit the Toeplitz structure and the properties of generators. Such algorithms allow us to compute the exponential of very large matrices, which would otherwise be untreatable with standard methods. We also prove interesting decay properties of the exponential of a generator having a block-triangular, block-Toeplitz structure.

Different kinds of nonlinear matrix equations are analyzed and their interplay with certain structured matrices is pointed out. By using the point of view and the tools of structured linear algebra, both theoretical results and numerical methods are derived for solving different classes of nonlinear matrix equations.

Non-skip-free Markov chains of the M/G/1 type are revisited in functional form. The problem of the computation of the steady state vector is reduced to inverting a Laurent matrix power series A ( z ) which is singular for z =1. This problem is related to the Wiener–Hopf factorization and to solving matrix equations. A way for removing the singularity is presented and some algorithms for inverting a Laurent matrix power series are shown. A generalization of Ramaswami's formula is derived from the Wiener–Hopf factorization of A ( z ).

We give a new characterization of the matrix square root and a new algorithm for its computation. We show how the matrix square root is related to the constant block coefficient of the inverse of a suitable matrix Laurent polynomial. This fact, besides giving a new interpretation of the matrix square root, allows one to design an efficient algorithm for its computation. The algorithm, which is mathematically equivalent to Newton's method, is quadratically convergent and numerically insensitive to the ill-conditioning of the original matrix and works also in the special case where the original matrix is singular and has a square root.

Let p(z) be a polynomial of degree n having zeros |ξ1|≤⋅⋅⋅≤|ξ m

In this paper, based on matrix structure analysis, we derive and analyze efficient algorithms to solve nonlinear matrix equations of the form X +∑ 1⩽ i ⩽ d A i X −1 D i = C . This class of equations is encountered in the solution of Tree-Like stochastic processes which are a generalization of Quasi-Birth-and-Death (QBD) processes.

The problem of the computation of the stochastic matrix G associated with discrete-time quasi-birth-death (QBD) Markov chains is analyzed. We present a shifted cyclic reduction algorithm and show that the speed of convergence of the latter modified algorithm is always faster than that of the original cyclic reduction.

Summary Matrix equations of the kind A1X2+A0X+A−1=X, where both the matrix coefficients and the unknown are semi-infinite matrices belonging to a Banach algebra, are considered. These equations, where coefficients are quasi-Toeplitz matrices, are encountered in certain quasi-birth–death processes as the tandem Jackson queue or in any other processes that can be modeled as a reflecting random walk in the quarter plane. We provide a numerical framework for approximating the minimal nonnegative solution of these equations that relies on semi-infinite quasi-Toeplitz matrix arithmetic. In particular, we show that the algorithm of cyclic reduction can be effectively applied and can approximate the infinite-dimensional solutions with quadratic convergence at a cost that is comparable to that of the finite case. This way, we may compute a finite approximation of the sought solution and of the invariant probability measure of the associated quasi-birth–death process, within a given accuracy. Numerical experiments, performed on a collection of benchmarks, confirm the theoretical analysis.

We analyze the problem of the computation of the solution of the nonlinear matrix equation $X=\sum_{i=0}^{+\infty}X^iA_i$, arising in queueing models. We propose a technique based on regular splittings, that one hand leads to a new method for computing the solution, on the other hand it may be used to construct nonlinear matrix equations equivalent to starting one, that can be possibly solved by applying different algorithms. TEL:: 050593448 EMAIL:: favati@imc.pi.cnr.it

We propose new algorithms for solving n x n banded Toeplitz systems with bandwidth m. If the function associated with the Toeplitz matrix has no zero in the unit circle, then $O(n\log m + m\log ^2 m\log\log \epsilon^{-1})$ arithmetic operations (ops) are sufficient to approximate the solution of the system up to within the error $\epsilon$; otherwise the cost becomes $O(n\log m +m\log^2 m\log {n\over m})$ ops. Here $m=o(n)$ and $n>\log \epsilon^{-1}$. Some applications are presented. The methods can be applied to infinite and bi-infinite systems and to block matrices.

Several multiple access control (MAC) protocols of industrial interest can be modeled by discrete-state discrete-time Markov chains, with finite lower block Hessenberg probability transition matrices P/sub N/, where the diagonal blocks have a size which decreases while passing from the top left block (A/sub 0,0/) down to the right bottom block (A/sub N,N/). Such matrices have been identified in the literature as funnel matrices. In general, for meaningful systems of engineering interest, the size of P/sub N/ is so large that the computation of the P/sub N/ stationary probabilities (/spl pi/) with common numerical methods cannot be performed due to the computational cost involved. To calculate the stationary probabilities of P/sub N/ we have developed an innovative computational method which fully exploits the block Hessenberg structure of the matrix I-P/sub N/. In this way, we drastically reduce the overall computational cost with respect to the customarily used LU factorization, while still keeping the strong numerical stability of Gaussian elimination with diagonal adjustment. The potential of the new method is exploited to assess the performance of an access protocol for third-generation mobile systems, called PRMA++ (packet reservation multiple access protocol). This was designed within the European project RACE and has so far been studied, to the best of our knowledge, via simulative analysis. PRMA++ has received a lot of attention from manufacturers in the ongoing fifth framework program of the European Community.

A survey on numerical methods, based on Toeplitz matrix computation, for the solution of M/G/l type Markov chains is presented together with some new advances. We show how the block Toeplitz structure, shared by the transition matrices associated with M/G/l type Markov chains, can be exploited as the basis for devising fast and numerically stable algorithms, for the computation of the probability invariant vector and for the solution of the nonlinear matrix equation .

The cyclic reduction technique (Buzbee et al., 1970), rephrased in functional form (Bini and Meini, 1996), provides a numerically stable, quadratically convergent method for solving the matrix equation X = ∑+ ∞ i=0 Xi Ai, where the Ai's are nonnegative k × k matrices such that ∑+ ∞ i=0 Ai is column stochastic. In this paper we propose a further improvement of the above method, based on a point-wise evaluation/interpolation at a suitable set of Fourier points, of the functional relations defining each step of cyclic reduction (Bini and Meini,1996). This new technique allows us to devise an algorithm based on FFT having a lower computational cost and a higher numerical stability. Numerical results and comparisons are provided.

In this chapter we consider quasi-birth and death processes with low rank downward and upward transitions. We show how such structure can be exploited to reduce the computational cost of the cyclic reduction iteration. The proposed algorithm saves computation by performing ultiplications and inversions of matrices of small size (equal to the rank instead of to the phase space dimension) and inherits the stability property of the customary cyclic reduction. Numerical experiments show the gain of the new algorithm in terms of computational cost.

In this chapter we consider quasi-birth and death processes with low rank downward and upward transitions. We show how such structure can be exploited to reduce the computational cost of the cyclic reduction iteration. The proposed algorithm saves computation by performing multiplications and inversions of matrices of small size (equal to the rank instead of to the phase space dimension) and inherits the stability property of the customary cyclic reduction. Numerical experiments show the gain of the new algorithm in terms of computational cost. Quasi-birth and death process, Low rank matrix, Cyclic reduction

A technique for deflating two eigenvalues λ 1 , λ 2 of an n × n quadratic matrix polynomial A ( z ) is proposed. The method requires the knowledge of a right eigenvector corresponding to λ 1 , and of a left eigenvector corresponding to λ 2 . This technique consists of two stages: the shift of λ 1 and λ 2 to 0 and ∞ , respectively, and the deflation of 0 and ∞ . The final result is the construction of an ( n - 1 ) × ( n - 1 ) quadratic matrix polynomial, which shares with A ( z ) all the eigenvalues, except for λ 1 and λ 2 . The particular case of ∗ -palindromic quadratic matrix polynomials is treated.

Matrix-analytic methods have advanced considerably since the pioneering work of Marcel Neuts [6, 5] on Quasi-Birth-Death (QBD), GI/M/1- and M/G/1- type Markov chains (MCs). Especially the algorithms involved to (iteratively) solve these structured Markov chains have matured a lot, which has resulted in more efficient, but also more complex algorithms [4, 1]. While the first algorithms were straightforward to implement---as they were based on simple functional iterations---more advanced algorithms/features like cyclic-reduction, the Newton iteration or the shift technique (to accelerate convergence), require more effort; in particular for GI/M/1- and M/G/1-type Markov chains. This has motivated us to develop the Structured Markov Chain Solver (SMCSolver) tool [2], which implements a large number of basic and more advanced algorithms for solving QBD, GI/M/1- and M/G/1-type MCs1 (as well as the more general Non-Skip-Free M/G/1-type MCs). The MATLAB version of the tool consists of a collection of MATLAB functions, while the Fortran version is accompanied by a graphical user-interface (GUI). Apart from making these more advanced algorithms accessible to non-specialists, the tool is also useful as a platform for the development and study of new algorithms and acceleration techniques. Since its initial release in 2006, various extensions have been made. In [3] different transformation techniques and shift strategies are incorporated in order to speed up and optimize the algorithms, while even more recently an efficient Newton iteration for GI/M/1- and M/G/1-type Markov chains was included [8]. Matrix-analytic methods have also been very effective in the analysis of many queueing systems in both discrete- and continuous-time. The Q-MAM tool [7] is a collection of MATLAB functions that allows one to compute the queue length, waiting time and delay distribution of various queueing systems of infinite size. It includes amongst others implementations of the PH/PH/1, MAP/MAP/1, MAP/M/c, MAP/D/c, RAP/RAP/1, MMAP[K]/PH[K]/1, MMAP[K]/SM[K]/1, SM[K]/PH[K]/1 (many in both discrete- and continuous-time), where state-of-the-art solution techniques are used to solve these models efficiently. The Matlab version of the SMCSolver and Q-MAM tool is available at http://win.ua.ac.be/?vanhoudt/ while the Fortran 90 version of the SMCSolver tool with the GUI can be downloaded from http://bezout.dm.unipi.it/SMCSolver.

Consider a quasi-birth-and-death (QBD) Markov chain [6], having probability transition matrix where Bi, Ai, i = ?1, 0, 1, are m x m matrices. In the numerical solution of QBD Markov chains a crucial step is the efficient computation of the minimal nonnegative solution R of the quadratic matrix equation X = X 2 A?1 + XA 0 + A 1 . (1) To this purpose, many numerical methods, with different properties, have been designed in the last years (see for instance [1, 2, 3, 4]). However, many of these numericalmethods are defined for general block coefficients A?1, A0 and A1, and do not exploit the possible structure of these blocks. Recently, some attention has been addressed to the case where A ?1 has only few non-null columns, or A1 has only few non-null rows. These properties are satisfied when the QBD has restricted transitions to higher (or lower) levels. In particular, in [7] the authors exploit these properties of the matrix A ?1 , or A 1 , to formulate the QBD in terms of an M/G/1 type Markov chain, where the block matrices have size smaller than m; in particular, when both A?1 and A1 have the desired property, the latter M/G/1 type Markov chain reduces to a QBD. In [5] the structure of A ?1 is used in order to reduce the computational cost of some algorithms for computing R. Here we assume that both the matrices A ?1 and A 1 have small rank with respect to their size m. In particular, if A?1 and A1 have only few non-null columns and rows, respectively, they have small rank. We show that, under this assumption, the matrix R can be computed by using the cyclic reduction algorithm, where the matrices A(k) i , i = ?1, 0, 1, generated at the kth step of the algorithm, can be represented by small rank matrices. In particular, if r ?1 is the rank of A ?1 , and if r 1 is the rank of A 1 , then each step of cyclic reduction can be performed by means of O((r ?1+r1 ) 3 ) arithmetic operations. This cost estimate must be compared with the cost of O(m3) arithmetic operations, needed without exploiting the structure of A ?1 and A 1 . Therefore, if r 1 and r 1 /are much smaller than m, the advantage is evident. It remains an open issue to understand how the structure can be exploited in the case where only one between A ?1 and A 1 has low rank.

This article reports the main properties of M/G/1-type Markov chains together with the classical and the most advanced algorithms for their analysis. Keywords: Markov chains; M/G/1-type; non-skip-free chain; QBD process