On the numerical analysis of infinite multi-dimensional Markov chains

Abstract

A system with multiple interacting subsystems that exhibits the Markov property can be represented as a multi-dimensional Markov chain (MC). Usually the reachable state space of this MC is a proper subset of its product state space, that is, Cartesian product of its subsystem state spaces. Compact storage of the infinitesimal generator matrix underlying such a MC and efficient implementation of analysis methods using Kronecker operations require the set of reachable states to be represented as a union of Cartesian products of subsets of subsystem state spaces. We first show that the problem of partitioning the reachable state space of a three or higher dimensional system with a minimum number of partitions into Cartesian products of subsets of subsystem state spaces is NP-complete. Two algorithms, one merge based the other refinement based, that yield possibly nonoptimal partitionings are presented. Results of experiments on a set of problems from the literature and those that are randomly generated indicate that, although it may be more time and memory consuming, the refinement based algorithm almost always computes partitionings with a smaller number of partitions than the merge based algorithm. When the infinitesimal generator matrix underlying the MC is represented compactly using Kronecker products, analysis methods based on vector– Kronecker product multiplication need to be employed. When the factors in the Kronecker product terms are relatively dense, vector–Kronecker product multiplication can be performed efficiently by the shuffle algorithm. When the factors are relatively sparse, it may be more efficient to obtain nonzero elements of the generator matrix in Kronecker form on-the-fly and multiply them with corresponding elements of the vector. We propose a modification to the shuffle algorithm that multiplies relevant elements of the vector with submatrices of factors in which zero rows and columns are omitted. This approach avoids unnecessary floating-point operations that evaluate to zero during the course of the multiplication. Numerical experiments on a large number of models indicate that, in many cases the modified shuffle algorithm performs a smaller number of floating-point operations than the shuffle algorithm and the algorithm that generates nonzeros on-the-fly, sometimes with minimum number of floating-point operations and amount of memory possible. Although the generator matrix is stored compactly using Kronecker products, solution vectors used in the analysis still require memory proportional to the size of the reachable state space. This becomes a bigger problem as the number of dimensions increases. We show that it is possible to use the hierarchical Tucker decomposition (HTD) to store the solution vectors during Kroneckerbased Markovian analysis relatively compactly and still carry out the basic operation of vector–matrix multiplication in Kronecker form relatively efficiently. The time evolution of a stochastic chemical system modelled as a continuoustime MC (CTMC) can be described as a system of ordinary differential equations (ODEs) known as the chemical master equation (CME). The CME can be analyzed by discretizing time and solving a linear system obtained by truncating the countably infinite state space at each time step. However, it is not trivial to choose a truncated state space that includes few states with negligible probabilities and leaves out only a small probability mass. We show that it is possible to decrease the memory requirement of the ODE solver using HTD with adaptive truncation strategies and we propose a novel approach to truncate the countably infinite state space using prediction vectors. Numerical experiments indicate that adaptive truncation strategies improve time and memory efficiency significantly when fixed truncation strategies are inefficient. Finally, we consider a multi-class multi-server retrial queueing system in which customers arrive according to a class-dependent Markovian arrival process (MAP). Service and retrial times follow class-dependent phase-type (PH) distributions with the further assumption that PH distributions of retrial times are acyclic. Here, we obtain a necessary and sufficient condition for ergodicity from criteria based on drifts. The countably infinite state space of the model is truncated with an appropriately chosen Lyapunov function. The truncated model is described as a multi-dimensional MC and a Kronecker representation of its infinitesimal generator matrix is numerically analyzed.