text:
# Computing stationary probability distributions and large deviation rates for constrained random walks. The undecidability results.

## 1 Introduction

The main model considered in this paper is a constrained homogeneous random walk in a $`d`$-dimensional nonnegative orthant $`𝒵_+^d`$, where $`𝒵_+^d`$ is the space of $`d`$-dimensional vectors with integral nonnegative components. Specifically, the transitions with positive probabilities can occur only to neighboring states and the transition probabilities depend only on the face that the current state of the random walk belongs to, but not on the size of the components of the state.

Ever since the appearance of the papers by Malyshev , ,  and Menshikov , random walks in $`𝒵_+^d`$ have assumed a prominent role in modelling and analysis of queueing networks of certain type, for example Markovian queueing networks. Specifically, the question of positive recurrence or stability was analyzed. One of the main techniques used for the stability analysis of these type of random walks is a Lyapunov function technique also known as Foster’s criteria. A comprehensive study of constrained random walks in $`𝒵_+^d`$ was conducted by Fayolle, Malyshev and Menshikov in  and many additional results appeared after the book was published. Specifically, a very interesting connection between constrained random walks and general dynamical systems on compact manifolds was established by Malyshev . Exact conditions for positive recurrence for the case $`d4`$ were obtained by Ignatyuk and Malyshev in . The large deviation principle for special cases and modifications of random walks in $`𝒵_+^d`$ was established by Ignatyuk, Malyshev and Scherbakov in . This followed by efforts to actually compute the large deviation rates, which turned out to be a very complicated problem. See for example Kurkova and Suhov , where large deviations limits are computed for a random walk in $`𝒵_+^2`$ arising from joint-the-shortest queueing system. The analysis uses a fairly complicated complex-analytic techniques. The goal of the current paper is to explain the difficulty in obtaining such results for general dimensions.

Analysis of random walks arising from special types of multiclass queueing networks became a subject of particularly aggressive research efforts during the previous decade. Many interesting and deep results were established which connect stability of such queueing networks with stability of corresponding fluid models, obtained by Law of Large Numbers type of rescaling. This research direction was initiated in pioneering works by Rybko and Stolyar  and Dai , where it was shown that stability of a fluid model implies stability of the underlying queueing system. The converse of this result is not true in general, see Dai, Hasenbein and Vande Vate , Bramson , but is true under some stronger conditions, Dai , Meyn . Despite these results, to the day no full characterization of stable queueing networks is available. Stability was characterized only for special types of queueing networks or special scheduling polices. For example feedforward networks are known to be stable for all work-conserving policies, Down and Meyn , Dai . Stability of fluid networks with two processing stations operating under arbitrary work-conserving policies is fully characterized in Bertsimas, Gamarnik and Tsitsiklis  by means of a linear programming and in Dai and Vande Vate  by direct methods. The question of computing stationary distributions comes naturally after the question of stability. Several results are available again in the context of multiclass queueing networks. Some of the results were obtained using quadratic Lyapunov functions, Bertsimas, Paschalidis and Tsitsiklis , Kumar and Kumar , Kumar and Meyn , Kumar and Morrison , using piece-wise linear Lyapunov functions in Bertsimas, Gamarnik and Tsitsiklis , and using more direct methods, Bertsimas and Nino-Mora ,  All of the results obtain only bounds on the stationary probabilities. Computing exactly the stationary probabilities seems beyond the existing techniques.

## 2 Our results

It was established by the author in  that positive recurrence of a constrained homogeneous random walk in $`𝒵_+^d`$ is an undecidable property. Meaning, no algorithm can possibly exists which given the description of the random walks (given the dimension and the transition matrix) will be able to check whether the walk is positive recurrent. This result was also established for queueing systems operating under the class of so called generalized priority policies. This result explains the difficulty in stability analysis by stating that these problems are simply insolvable. It was conjectured in the same paper that the stability of multiclass queueing networks operating under the class of much studied priority or First-In-First-Out policies, is undecidable as well. The conjecture remains unproven.

In the current paper we continue the decidability analysis of constrained random walks by asking the following question: given a constrained homogeneous random walk, can we compute its stationary distribution, provided that the existence of a stationary distribution can be checked, for example, by constructing a Lyapunov function? To put this question into a proper computation theoretic framework, we ask the following question. Given a constrained homogeneous random walk, which possesses a unique stationary distribution $`\pi `$, given a state $`q𝒵_+^d`$, for example $`q=0`$, and given a rational value $`r>0`$, is it true that the stationary probability $`\pi (q)`$ of this state satisfies $`\pi (q)r`$? In this paper we prove that this problem is undecidable,  *even if* a Lyapunov function witnessing positive recurrence is available. Thus, no algorithm can exist which given a positive recurrent constrained homogeneous random walk computes its stationary distribution. Specifically, the stationary distribution cannot be written down in any constructive way using some formulas. Contrast this with random walks corresponding to product form type networks, for example Jackson networks, for which a very simple formula is available. We then prove that computing large deviations limits for the same model is an undecidable problem as well. In particular, we show that given a random walk in $`𝒵_+^d`$ with a unique stationary distribution (witnessed, for example, by a Lyapunov function), and given a vector $`v\mathrm{}^d`$, the problem of deciding whether $`lim_n\mathrm{}\mathrm{log}(\pi (vn))/n`$ is finite, is undecidable. We extend these results to queueing systems operating under a class of generalized priority policies. Finally, we observe that, nevertheless, estimating stationary distribution of a constrained random walk is a decidable problem, if one is willing to tolerate a two-sided error and a Lyapunov function exists. Specifically, given such a random walk with the unique stationary distribution $`\pi `$, given a Lyapunov function, given a state $`q𝒵_+^d`$ and any value $`ϵ>0`$, an interval $`(r,r+ϵ)`$ can be constructed which contains $`\pi (q)`$. This result is an easy consequence of a powerful result obtained by Meyn and Tweedie , which obtains exponential bounds on the mixing rate of Markov chains, using Lyapunov function methods. We note that such approximation result cannot be obtained for large deviations rates since, as we mentioned above, even determining whether a given large deviation rate is finite is an undecidable problem.

The remainder of the paper is organized as follows. In the following section we describe our model – constrained homogeneous random walk in $`𝒵_+^d`$ and introduce Lyapunov functions. In Section 4 we introduce a counter machine – a modification of a Turing machine which for us is the main tool for establishing the undecidability results. In Section 5 we prove that computing a stationary distribution of a positive recurrent random walk in $`𝒵_+^d`$ is an undecidable problem. In Section 6 we prove that computing large deviations rates for positive recurrent random walks in $`𝒵_+^d`$ is an undecidable problem as well. Extension of these results to queueing systems is established in Section 7. In Section 8 we show how stationary distribution can be computed with a two-sided error using a Lyapunov function technique and Meyn and Tweedie results from . Conclusions and open problems are discussed in Section 9.

## 3 Constrained homogeneous random walk in $`𝒵_+^d`$. Lyapunov function and stationary distribution

Let $`𝒵_+^d`$ denote the space of $`d`$-dimensional vectors with nonnegative integer components. Our model is a random walk $`Q(t),t=0,1,2,\mathrm{}`$ which has $`𝒵_+^d`$ as a state-space. For each $`\mathrm{\Lambda }\{1,2,\mathrm{},d\}`$ let $`𝒵_\mathrm{\Lambda }`$ denote the corresponding face:

$$𝒵_\mathrm{\Lambda }=\{(z_1,z_2,\mathrm{},z_d)𝒵_+^d:z_i>0\mathrm{for}i\mathrm{\Lambda },z_i=0\mathrm{for}i\mathrm{\Lambda }\}.$$

The transition probabilities are face-homogeneous – they depend entirely on the face the random walk is currently on. In addition the transition vectors have at most unit length in $`\mathrm{max}`$ norm. In other words, for each $`\mathrm{\Lambda }\{1,2,\mathrm{},d\}`$ and each $`\mathrm{\Delta }\{1,0,1\}^d`$ a certain value $`p(\mathrm{\Lambda },\mathrm{\Delta })`$ (the transition probability) is defined. These values satisfy

$$\underset{\mathrm{\Delta }\{1,0,1\}^d}{}p(\mathrm{\Lambda },\mathrm{\Delta })=1$$

for each $`\mathrm{\Lambda }`$ and $`p(\mathrm{\Lambda },\mathrm{\Delta })=0`$ if for some $`i\mathrm{\Lambda },\mathrm{\Delta }_i=1`$. The latter condition is simply a consistency condition which prevents transitions into states with negative components. Given a current state $`Q(t)𝒵_+^d`$ of the random walk, the next state is chosen to be $`Q(t)+\mathrm{\Delta }`$ with probability $`p(\mathrm{\Lambda },\mathrm{\Delta })`$, if the state $`Q(t)`$ belongs to the face $`𝒵_\mathrm{\Lambda }`$. We will also write $`p(q,q^{})`$ instead of $`p(\mathrm{\Lambda },\mathrm{\Delta })`$ if the state $`q𝒵_\mathrm{\Lambda }`$ and $`q^{}q=\mathrm{\Delta }`$. We denote by $`p^{(t)}(q,q^{})`$ the $`t`$-step transition probabilities: $`\mathrm{Prob}\{Q(t)=q^{}|Q(0)=q\}`$. The model above will be referred to as constrained homogeneous random walk in $`𝒵_+^d`$. We will say that our walk is deterministic if $`p(\mathrm{\Lambda },\mathrm{\Delta })\{0,1\}`$ for all $`\mathrm{\Lambda }`$ and $`\mathrm{\Delta }`$. In other words, the transition vector $`\mathrm{\Delta }=\mathrm{\Delta }(\mathrm{\Lambda })`$ deterministically depends on the face. The set of parameters $`p(\mathrm{\Lambda },\mathrm{\Delta })`$ is finite, and, in particular, it contains $`6^d`$ elements corresponding to $`2^d`$ faces $`𝒵_\mathrm{\Lambda }`$ and $`3^d`$ transition vectors $`\mathrm{\Delta }`$ per face (with some transitions occurring with zero probability). Let $`Q(t)`$ denote $`L_1`$ norm. That is $`Q(t)=_{id}Q_i(t)`$.

For any state $`q𝒵_+^d`$ and subset $`X𝒵_+^d`$, let $`T=T(q,X)`$ denote the first hitting time for the set $`X`$ when the initial state of the walk is $`q`$, including the possibility $`T=\mathrm{}`$. That is

$$T=\mathrm{min}\{t:Q(t)X|Q(0)=q\}.$$
(1)
The following definition is standard in the theory or infinite Markov chains.

###### Definition 1

A homogeneous random walk is defined to be positive recurrent or stable if there exist some $`C>0`$ such that the random walk visits the set $`X_C\{z𝒵_+^d:_{i=1}^dz_iC\}`$ infinitely often with probability one, and $`E[T(q,X_C)]`$ is finite for all $`q𝒵_+^d`$.

Stability of a constrained homogeneous random walk $`Q(t)`$ can be checked, for example, by constructing a suitable Lyapunov function.

###### Definition 2

A function $`\mathrm{\Phi }:𝒵_+^d\mathrm{}_+`$ is defined to be a Lyapunov function with drift $`\gamma <0`$ and exception set $`𝒵_+^d`$ if $`||<\mathrm{}`$ and for every state $`q`$

$$E[\mathrm{\Phi }(Q(t+1)|Q(t)=q]\mathrm{\Phi }(q)=\underset{q𝒵_+^d}{}\mathrm{\Phi }(q^{})p(q,q^{})\mathrm{\Phi }(q)\gamma .$$
(2)
In other words, the expected value of the Lyapunov function should decrease at each time step, whenever the random walk is outside of the exception set. Existence of a Lyapunov function under some additional assumptions implies stability. For a comprehensive survey of Lyapunov function methods see Meyn and Tweedie . Various forms of Lyapunov functions, specifically piece-wise linear and quadratic Lyapunov functions were used to prove stability of random walks corresponding to Markovian queueing networks , , , , , . In some simple cases even linear Lyapunov function of the form $`\mathrm{\Phi }(q)=_{i=1}^dw_iq_i`$ with $`w_i0`$ can prove stability of a constrained homogeneous random walk. It is easy to see that a linear function $`\mathrm{\Phi }(q)=w^Tq`$ is a Lyapunov function if and only if for some $`\gamma >0`$ and every nonempty set $`\mathrm{\Lambda }\{1,2,\mathrm{},d\}`$ the following inequality holds

$$E[w^TQ(t+1)w^TQ(t)|Q(t)𝒵_\mathrm{\Lambda }]=\underset{\mathrm{\Delta }\{1,0,1\}^d}{}(w^T\mathrm{\Delta })p(\mathrm{\Lambda },\mathrm{\Delta })\gamma .$$
(3)
The existence of a linear Lyapunov function is only sufficient but not necessary for stability of the constrained random walk. It is also useful sometimes to consider a geometric Lyapunov function, defined as follows.

###### Definition 3

A function $`\mathrm{\Phi }_g:𝒵_+^d[1,+\mathrm{})`$ is defined to be a geometric Lyapunov function with drift $`0<\gamma _g<1`$ and exception set $`𝒵_+^d`$ if $`||<\mathrm{}`$ and for every state $`q`$

$$\frac{E[\mathrm{\Phi }_g(Q(t+1)|Q(t)=q]}{\mathrm{\Phi }_g(q)}=\underset{q𝒵_+^d}{}\frac{\mathrm{\Phi }_g(q^{})}{\mathrm{\Phi }_g(q)}p(q,q^{})\gamma _g<1.$$
(4)
A geometric Lyapunov function is used, for example in Meyn and Tweedie , to prove exponentially fast mixing of a Markov chain which admits a geometric Lyapunov function. The precise statement of this result will be given below in Section 8. If the condition (3) is met for some $`w`$ and $`\gamma `$, then a function of the form $`\mathrm{\Phi }_g(q)=\mathrm{exp}(\delta w^Tq)`$ is a geometric Lyapunov function for a suitable value of $`\delta >0`$.

Throughout the paper we will assume all the states $`q`$ communicate with the state $`0`$, that is $`p^{(t)}(q,0)>0`$ for some $`t0`$. As a consequence, the random walk is irreducible. If it is in addition positive recurrent, then it possesses a unique stationary distribution $`\pi :𝒵_+^d[0,1]`$. Namely $`_{q𝒵_+^d}\pi (q)=1`$ and for any state $`q`$

$$\underset{q^{}𝒵_+^d}{}\pi (q^{})p(q^{},q)=\pi (q)$$
(5)
This stationary distribution is defined completely by the set of transition parameters $`p(\mathrm{\Lambda },\mathrm{\Delta })`$. Computing the stationary probability distribution for these walks is the main focus of this paper. It was established by the author in  that checking positive recurrence of a constrained homogeneous random walk is an undecidable problem - no algorithm can exist to achieve this task. However, if one is lucky to construct a Lyapunov function, for example by checking condition 3 for some nonnegative vector $`w\mathrm{}_+^d`$, then the random walk is in fact positive recurrent. One might be tempted to believe that in this case the analysis of the random walk is simplified significantly. In Section 5 we show even if a linear Lyapunov function exists, computing the stationary probability distribution is still an undecidable problem. As in the case of stability analysis, our main tool for establishing this undecidability result is a counter machine and the halting problem defined in Section 4.

## 4 Counter Machines, Halting Problem and Undecidability

A counter machine (see , ) is a deterministic computing machine which is a simplified version of a Turing machine – a general description of an algorithm working on a particular input. In his classical work on the Halting Problem, Turing showed that certain decision problems simply cannot have a corresponding solving algorithm, and thus are undecidable. For a definition of a Turing machine and the Turing Halting Problem see . Ever since many quite natural problems, in mathematics and computer science were found to be undecidable. Some of the undecidability results in control theory were obtained by reduction from a counter machine, see Blondel et al. . For a survey of decidability results in control theory area see Blondel and Tsitsiklis .

A counter machine is described by 2 counters $`R_1,R_2`$ and a finite collection of states $`S`$. Each counter contains some nonnegative integer in its register. Depending on the current state $`sS`$ and depending on whether the content of the registers is positive or zero, the counter machine is updated as follows: the current state $`s`$ is updated to a new state $`s^{}S`$ and one of the counters has its number in the register incremented by one, decremented by one or no change in the counters occurs.

Formally, a counter machine is a pair $`(S,\mathrm{\Gamma })`$. $`S=\{s_0,s_1,\mathrm{},s_{m1}\}`$ is a finite set of states and $`\mathrm{\Gamma }`$ is configuration update function $`\mathrm{\Gamma }:S\times \{0,1\}^2S\times \{2,1,0,1,2\}`$. A configuration of a counter machine is an arbitrary triplet $`(s,z_1,z_2)S\times 𝒵_+^2`$. A configuration $`(s,z_1,z_2)`$ is updated to a configuration $`(s^{},z_1^{},z_2^{})`$ as follows. First a binary vector $`b=(b_1,b_2)`$ is computed were $`b_i=1`$ if $`z_i>0`$ and $`b_i=0`$ if $`z_i=0`$, $`i=1,2`$. If $`\mathrm{\Gamma }(s,b)=(s^{},1)`$, then the current state is changed from $`s`$ to $`s^{}`$, the content of the first counter is incremented by one and the second counter does not change: $`z_1^{}=z_1+1,z_2^{}=z_2`$. We will also write $`\mathrm{\Gamma }:(s,z_1,z_2)(s^{},z_1+1,z_2)`$ and $`\mathrm{\Gamma }:ss^{},\mathrm{\Gamma }:z_1z_1+1,\mathrm{\Gamma }:z_2z_2`$. If $`\mathrm{\Gamma }(s,b)=(s^{},1)`$, then the current state becomes $`s^{}`$, $`z_1^{}=z_11,z_2^{}=z_2`$. Similarly, if $`\mathrm{\Gamma }(s,b)=(s^{},2)`$ or $`\mathrm{\Gamma }(s,b)=(s^{},2)`$, the new configuration becomes $`(s^{},z_1,z_2+1)`$ or $`(s^{},z_1,z_21)`$, respectively. If $`\mathrm{\Gamma }(s,b)=(s^{},0)`$ then the state is updated to $`s^{}`$, but the contents of the counters do not change. This definition can be extended to the one which incorporates more than two counters, but, in most cases, such an extension is not necessary for our purposes.

Given an initial configuration $`(s^0,z_1^0,z_2^0)`$ the counter machine uniquely determines subsequent configurations $`(s^1,z_1^1,z_2^1),(s^2,z_1^2,z_2^2),\mathrm{},(s^t,z_1^t,z_2^t),\mathrm{}.`$ We fix a certain configuration $`(s^{},z_1^{},z_2^{})`$ and call it a *halting* configuration. If this configuration is reached then the process halts and no additional updates are executed. The following theorem establishes the undecidability of the halting property.

###### Theorem 1

Given a counter machine $`(S,\mathrm{\Gamma })`$, initial configuration $`(s^0,z_1^0,z_2^0)`$ and the halting configuration $`(s^{},z_1^{},z_2^{})`$, the problem of determining whether the halting configuration is reached in finite time is undecidable. It remains undecidable even if the initial and the halting configurations are the same with both counters equal to zero: $`s^0=s^{},z_1^0=z_2^0=z_1^{}=z_2^{}=0`$.

The first part of this theorem is a classical result and can be founded in . The restricted case of $`s^0=s^{},z_i^0=z_i^{},i=1,2`$ can be proven similarly by extending the set of states and the set of transition rules. It is the restricted case of the theorem which will be used in the current paper.

## 5 Computing the stationary probability distribution. The undecidability result

Theorem 1 was used in  to prove that the stability of a constrained random walk in $`𝒵_+^d`$ is undecidable. Naturally, the problem of stability comes before the problem of computing the stationary distribution of a stable random walk. As we mentioned in Section 3, stability can be checked sometimes by constructing a Lyapunov function. In this section we prove our main result: even if such a Lyapunov function, witnessing stability, is available and is provided as a part of the data parameters, computing stationary distribution is an undecidable problem.

We now give an informal outline of the proof. The proof uses a reduction from a halting problem for a counter machine. We embed a counter machine with initial and halting configuration $`(s^{},0,0)`$ into a deterministic walk in $`𝒵_+^d`$. The state space and the transition rules of this walk are then extended in some way that incorporates an independent Bernoulli process with some fixed parameter $`p`$. We then show that

* If the original counter machine never returns to the initial configuration $`(s^{},0,0)`$, then the constructed random walk, when started from the origin, returns into the origin in $`2t+2`$ steps with probability $`(1p)p^t`$ for $`t=0,1,2,\mathrm{}.`$ In particular, the expected return time to the origin is $`2/(1p)`$.
* If the original counter machine returns to the initial configuration in $`T`$ steps, then the modified random walk returns into the origin in $`2t+2`$ steps with probability $`(1p)p^t`$ for $`tT1`$ and in $`2T+2`$ steps with the remaining probability $`1_{tT1}(1p)p^t=p^T`$. In particular, the expected return time to the origin is $`(22p^{T+1})/(1p)`$.

The stationary probability distribution of any state is exactly the expected return time to this state. Therefore, the stationary probability of the origin is $`(1p)/2`$ if the counter machine halts and is strictly greater, if the counter machine does not halt. Since the value $`p`$ is our control parameter, and since checking whether the counter machine halts is an undecidable problem, then computing stationary probability is undecidable as well. We now state and prove rigorously this result. As before, let $`\pi `$ denote the unique stationary distribution of an irreducible positive recurrent random walk. Let also $`0`$ denote the origin of the nonnegative lattice $`𝒵_+^d`$.

###### Theorem 2

Given an irreducible constrained random walk with transition probabilities $`p(\mathrm{\Lambda },\mathrm{\Delta })`$, given a linear vector $`w\mathrm{}_+^d`$ satisfying (3) and given a rational value $`0r1`$, the problem of checking whether $`\pi (0)r`$ is undecidable – no algorithm exists which achieves this task.

Remarks :

1. The stationary distribution can in principle take non-rational values. In order to put the problem into a framework suitable for algorithmic analysis we modified the question into the one of checking whether $`\pi ()r`$ for rational values $`r`$. This is a standard method in the theory of Turing decidable numbers, see .

2. A simple example where computing the stationary probability distribution is a decidable problem is Jackson networks, . For such a network with $`d`$ stations the stationary probability of the state $`m=(m_1,\mathrm{},m_d)`$ is given by $`_{j=1}^d(1\rho _j)\rho _j^{m_j}`$, where $`\rho _j`$ is the traffic intensity in station $`j`$. Specifically, the stationary probability of the state $`0`$ is $`_{j=1}^d(1\rho _j)`$. Given any rational value $`0r1`$, it is a trivial computation to check whether this product is at least $`r`$.

Proof of Theorem 2: we start with a construction used in . Namely, we embed a given counter machine with states $`s_0,s_1,\mathrm{},s_{m1}`$ into a deterministic walk in $`𝒵_+^{m+1}`$ as follows. Without the loss of generality, assume that $`s^{}=s_0`$. Let configuration $`(s_i,z_1,z_2),1im1`$ correspond to the state $`q=(e_i,z_1,z_2)𝒵_+^{m+1}`$, where $`e_i`$ is unit vector with $`1`$ in $`i`$-th coordinate and zero everywhere else. Also, let configurations $`(s_0,z_1,z_2)`$ correspond to $`(0,z_1,z_2)`$, with zeros in first $`m1`$ coordinates. Specifically, the initial and halting configuration $`(s_0,0,0)`$ corresponds to the origin $`0`$. We now describe the set of transition vectors $`\mathrm{\Delta }=\mathrm{\Delta }(\mathrm{\Lambda })`$. We describe it first for subsets $`\mathrm{\Lambda }\{1,2,\mathrm{},m+1\}`$ which correspond to an encoding of some configuration of a counter machine. Specifically, $`\mathrm{\Lambda }\{1,2,\mathrm{},m1\}=\mathrm{}`$ (corresponding to configurations with state $`s_0`$) or $`\mathrm{\Lambda }\{1,2,\mathrm{},m1\}=\{i\}`$ for some $`1im1`$, corresponding to configurations with state $`s_i`$. Fix any configuration $`(s_i,z_1,z_2)`$. Suppose the corresponding update rule is $`\mathrm{\Gamma }((s_i,z_1,z_2))=(s_j,+1)`$ for some $`1jm1`$. That is, the state is changed into $`s_j`$, the first counter is incremented by $`1`$ and second counter remains unchanged. We make the corresponding transition vector to be $`\mathrm{\Delta }=\mathrm{\Delta }(\mathrm{\Lambda })`$, where the $`i`$-th coordinate of $`\mathrm{\Delta }`$ is $`1`$, the $`j`$-th coordinate is $`+1`$, the $`m`$-th coordinate is $`+1`$ and all the other coordinates are zeros. It is easy to see that if at time $`t`$, the state $`Q(t)`$ corresponds to some configuration $`(s_i,z_1,z_2)`$, that is $`Q(t)=(e_i,z_1,z_2)`$, then $`Q(t+1)=Q(t)+\mathrm{\Delta }`$ corresponds to the configuration $`(s_j,z_1+1,z_2)`$ obtained by applying rule $`\mathrm{\Gamma }`$. We construct transition vectors similarly for other cases of configuration updates. In particular, if the state $`s_i`$ is changed to state $`s_0`$, then the corresponding $`\mathrm{\Delta }`$ has $`1`$ in the $`i`$-th coordinate and zeros in all the coordinates $`1jm1,ji`$. As we will see later, if $`Q(t)`$ corresponds to some configuration of a counter machine at time $`t`$, then it does so for all the later time $`t^{}t`$. Now if $`Q(t)`$ belongs to some face $`𝒵_\mathrm{\Lambda }`$ which does not correspond to some configuration, then we simply set $`\mathrm{\Delta }(\mathrm{\Lambda })=e_i`$ where $`i`$ is the smallest coordinate which belongs to $`\mathrm{\Lambda }`$. Then at some later time $`t^{}>t`$ the state $`Q(t^{})`$ will correspond to some configuration.

Construction above is exactly the one used in  to analyze stability. We now modify the construction by adding two additional coordinates. Our new state at time $`t`$ is thus denoted by $`\overline{Q}(t)=(Q(t),q_1(t),q_2(t))𝒵_+^{m+3}`$. Also a parameter $`0<p<1`$ is fixed. The transition rules are modified as follows.

1. When $`q_2(t)=1`$, the first part $`Q(t)`$ of the state is updated exactly as above. Also, if $`Q(t)>0`$, in other words, $`Q(t)`$ does not represent the halting configuration $`(s_0,0,0)`$, then the value of $`q_2(t)`$ stays $`1`$ with probability $`p`$ and switches to $`0`$ with probability $`1p`$. If, on the other hand $`Q(t)=0`$ then we set $`q_2(t+1)=0`$ with probability $`1`$. Finally, the value of $`q_1(t)\{2,1,0,1,2\}`$ is selected in such a way that $`(Q(t+1),q_1(t+1))=(Q(t),q_1(t))+1`$, where $`(Q(t),q_1(t))=_{i=1}^{m+1}Q_i(t)+q_1(t)`$. It is easy to see that such a value of $`q_1(t)`$ always exists. For example if $`Q(t)`$ encodes $`(s_i,z_1,z_2),i0`$ and the configuration is changed into $`(s_j,z_1,z_21),j0`$, then we put $`q_1(t)=2`$.

Remark : We stipulated before that the transition vectors $`\mathrm{\Delta }`$ must belong to $`\{1,0,1\}^{m+4}`$ for our constrained random walk, whereas above the value of $`q_1(t)`$ can change by $`2`$ and $`2`$. It is easy to satisfy this constraint for $`q_1(t)`$ by splitting it into two coordinates $`q_1(t),q_1^{}(t)`$ and making $`q_1(t)=q_1^{}(t)=1`$ in case $`q_1(t)`$ was assigned $`2`$ before, and $`q_1(t)=q_1^{}(t)=1`$ in case $`q_1(t)`$ was assigned $`2`$. We keep only one $`q_1(t)`$ for simplicity, allowing it to take values $`2,2`$.
2. When $`q_2(t)=0`$, we set $`\mathrm{\Delta }_k=1,\mathrm{\Delta }_i=0,ik,1im+3`$, where $`k`$ is the smallest coordinate such that $`\overline{Q}_k(t)>0`$. In particular, $`q_2(t)`$ stays equal to $`0`$. If $`Q(t)=q_1(t)=0`$, (in particular $`Q(t)`$ encodes the initial-terminal configuration $`(s_0,0,0)`$) then $`Q(t)`$ and $`q_1(t)`$ are updated as in the case $`q_2(t)=1`$ above. Also $`q_2(t)`$ in this case is switched to $`1`$ with probability $`p`$ and stays $`0`$ with probability $`1p`$.

Note, that the only stochastic part in our random walk is the last component $`q_2(t)`$.

###### Proposition 1

The constructed random walk $`\overline{Q}(t)`$ is irreducible and positive recurrent with the unique stationary distribution $`\pi `$. Moreover,

1. If the counter machine with the initial configuration $`(s_0,0,0)`$ does not halt, then the random walk $`\overline{Q}(t)`$ with the initial state $`\overline{Q}(0)=0`$ returns to the origin in $`2t+2`$ steps with probability $`(1p)p^t`$, for $`t=0,1,2,\mathrm{}.`$ As a result, the expected recurrence time of the state $`0`$ is $`1/\pi (0)=2/(1p)`$.
2. If the counter machine with the initial configuration $`(s_0,0,0)`$ halts in $`T1`$ steps, then the random walk $`\overline{Q}(t)`$ with the initial state $`\overline{Q}(0)=0`$ returns to the origin in $`2+2t`$ steps with probability $`(1p)p^t`$ for $`t<T`$, and in $`2+2T`$ steps with the remaining probability $`p^T`$. As a result, the expected recurrence time of the state $`0`$ is $`1/\pi (0)=(22p^{T+1})/(1p)`$.
3. For any $`C2/(1p)`$ the function $`_{i=1}^{m+1}Q_i+q_1(t)+Cq_2(t)`$ is a linear Lyapunov function with drift $`\gamma =1`$ and an exception set $`=\{0\}`$.

We first show that the proposition above implies the theorem. Suppose, we had an algorithm $`𝒜`$ which given an irreducible constrained random walk $`Q(t)`$, with a linear Lyapunov function $`w^TQ(t)`$ and given a rational value $`0r1`$ could determine whether the unique stationary distribution $`\pi `$ satisfies $`\pi (0)r`$. We take a counter machine and construct a random walk $`\overline{Q}(t)`$ as described above. Proposition 1 implies that this walk is a valid input for the algorithm $`𝒜`$. We use $`𝒜`$ to determine whether $`\pi (0)r(1p)/2`$. From Proposition 1, this is the case if and only if the underlying counter machine does not halt. In this fashion, we obtain an algorithm for checking halting property for counter machines. This is a contradiction to Theorem 1. $`\mathrm{}`$

Proof of Proposition 1: Suppose the underlying counter machine does not halt. Let us trace the dynamics of our random walk $`\overline{Q}(t)`$ starting from $`\overline{Q}(0)=0`$. Initially, by applying rule 2, it moves into some state $`(Q(1),0,1)`$ with probability $`p`$ or state $`(Q(1),0,0)`$ with probability $`1p`$. An independent Bernoulli process for $`q_2(t)`$ with parameter $`p`$ is continued in the first case. Suppose this process succeeds exactly $`t0`$ times (including the transition from initial state $`0`$), which occurs with probability $`(1p)p^t`$. Then, applying rule 1, at times $`t`$ and $`t+1`$ we have states $`(Q(t),q_1(t),1),(Q(t+1),q_1(t+1),0)`$ with $`Q(t)+q_1(t)=t,Q(t+1)+q_1(t+1)=t+1`$. At this moment rule 2 becomes applicable. Since at each step the norm $`\overline{Q}(t)`$ decreases exactly by one, the origin is reached at time $`(t+1)+(t+1)`$. We conclude that the return time is $`2+2t`$ with probability $`(1p)p^t,t=0,1,2,\mathrm{}`$. The expected return time is then $`2/(1p)`$ and the stationary probability of the state $`0`$ is $`(1p)/2`$.

Suppose, now, the underlying counter machine reaches the terminal state $`(s_0,0,0)`$ in exactly $`T1`$ steps. Suppose also the Bernoulli process for $`q_2(t)`$ succeeds exactly $`t0`$ times for $`t<T`$. Then, exactly as above, the origin is reached in $`2+2t`$ steps and this occurs with probability $`(1p)p^t`$. If, however, by the time $`T`$ the Bernoulli process does not fail, which occurs with probability $`p^T`$, then the state $`\overline{Q}(T)=(Q(T),q_1(t),1)=(0_{m+1},q_1(t),1)`$ is reached at time $`T`$, where $`0_{m+1}`$ denotes a $`m+1`$-dimensional zero vector. By the choice of $`q_1(t)`$ in rule 1, $`q_1(T)=T`$. At time $`T+1`$, by rule 1, we have a state $`\overline{Q}(T+1)=(0_{m+2},T+1,0)`$ and rule 2 applies. At time $`(T+1)+(T+1)`$ the origin is reached. We conclude that the random walk returns to the origin in $`2+2T`$ steps. Combining, the expected return time to the origin is then $`(22p^{T+1})/(1p)`$, if the counter machine halts in $`T`$ steps, and the stationary probability of the state $`0`$ is $`(1p)/(22p^{T+1})`$.

To complete the proof of the proposition, we analyze the expected change of the function $`_{i=1}^{m+1}Q_i+q_1(t)+Cq_2(t)`$. When $`q_2(t)=0`$ and $`\overline{Q}(t)0`$, the sum decreases deterministically by $`1`$. When $`q_2(t)=1`$, the value of $`_{i=1}^{m+1}Q_i+q_1(t)`$ increases deterministically by $`1`$, and the value of $`Cq_2(t)`$ stays the same with probability $`p`$ or decreases by $`C`$ with probability $`1p`$. Therefore, the expected change of the sum is $`1C(1p)`$. When $`C2/(1p)`$, the expected change is at most $`1`$. $`\mathrm{}`$

An important implication of Theorem 2 is that it is impossible to express the stationary distribution of a positive recurrent random walk $`Q(t)`$ as a function of the parameters $`p(\mathrm{\Lambda },\mathrm{\Delta })`$ via some computable function $`f()`$. For example, the stationary distribution cannot be expressed as roots of some polynomial equations with rational coefficient, as inequalities $`xr`$ can be checked for any root $`x`$ of such a polynomial and any rational value $`r`$. This is a startling contrast to a simple expression $`(1\rho _j)\rho _j^{n_j}`$ corresponding to a stationary distribution of a Jackson network.

## 6 Large Deviation Rates. The undecidability result

In this section we discuss the question of computing large deviation rates for our model. Specifically, we focus on computing large deviation rates for the stationary distribution $`\pi `$ of our random walk $`Q(t)`$ in $`𝒵_+^d`$. Let $`\mathrm{}`$ and $`\mathrm{}_+`$ denote the set of real values and the set of nonnegative real values, respectively. For any $`x\mathrm{}`$ let $`x`$ denote largest integer not bigger than $`x`$, and for any $`x\mathrm{}^d`$ let $`x=(x_1,\mathrm{},x_d)`$. We say that a function $`L:\mathrm{}_+^d\mathrm{}_+\{\mathrm{}\}`$ is a large deviation rate function for a given irreducible positive recurrent random walk $`Q(t)`$ in $`𝒵_+^d`$ if for any vector $`v\mathrm{}_+^d`$, the stationary distribution $`\pi `$ satisfies

$$\underset{n\mathrm{}}{lim}\frac{\mathrm{log}(\pi (vn))}{n}=L(v),$$
(6)
In other words, the stationary probability of being in state $`vn`$ is asymptotically $`\mathrm{exp}(L(v)n)`$ for large $`n`$. For results on large deviations for specific types of constrained random walks in $`𝒵_+^d`$ see . There are numerous works on large deviation in the context of queueing systems, see Shwartz and Weiss  for a survey. Specifically, Kurkova and Suhov  study large deviation rates for a two dimensional random walk corresponding to join-the-shortest queue. The analysis is quite intricate and uses complex-analytic techniques developed by Malyshev , ,  back in 70’s. To the best of our knowledge, the existence of the large deviations limits (6) is not fully proved for general constrained homogeneous random walks $`Q(t)`$ in $`𝒵_+^d`$. One can instead consider limits

$$L_{}(v)\underset{n\mathrm{}}{lim\; inf}\frac{\mathrm{log}(\pi (vn))}{n},L_+(v)\underset{n\mathrm{}}{lim\; sup}\frac{\mathrm{log}(\pi (vn))}{n}.$$
(7)
The goal of the present section is to prove that computing the large deviation rate function $`L(v)`$ is an undecidable problem, even if the walk is known to be a priori positive recurrent via, for example, existence of a linear Lyapunov function, and even if the large deviation limit function $`L(v)`$ is known to exist. The following is the main result of this section.

###### Theorem 3

Given an irreducible constrained random walk, given a linear vector $`w\mathrm{}_+^d`$ satisfying (3), given a rational value $`0r1`$ and a vector $`vZ_+^d`$, the problems of determining whether $`L_{}(v)r,L_+(v)r`$ are undecidable.

Remark : As we will see below, the large deviations limit function $`L(v)=L_{}(v)=L_+(v)`$ exists for the subclass of random walks we consider. As before, the reason for including a linear Lyapunov function into the condition of the theorem is to provide a simple way of insuring that the walk is positive recurrent.

Proof : The proof is again based on reduction from a halting problem for a counter machine. Given a counter machine with $`m`$ states consider the extended $`m+3`$-dimensional random walk $`\overline{Q}(t)`$ constructed in the proof of Theorem 2. We extend it even further by adding an additional coordinate $`(\overline{Q}(t),q_3(t))`$. Recall, that the rules for updating $`q_1(t)`$ were such that $`Q(t)+q_1(t)=t`$ as long as $`q_2(t^{})=1`$ for $`1t^{}t`$. Construct the rules for updating $`q_3(t)`$ as follows. If $`q_2(t)=1`$, then $`q_3(t+1)=q_3(t)+1`$. Also, if $`Q(t)=q_1(t)=q_2(t)=0`$ then again $`q_3(t+1)=q_3(t)+1`$. In other words, as long as the random walk starts from the origin and as long as $`q_2(t)`$ remains equal to $`1`$, $`q_3(t)=Q(t)+q_1(t)=t`$. Once $`q_2(t)`$ becomes zero, the value of $`q_3(t)`$ stays the same as long as $`\overline{Q}(t)0`$ and decreases by one when $`\overline{Q}(t)`$ becomes zero, and continues decreasing until it itself becomes zero.

Let $`v=(0,\mathrm{},0,1)`$ be an $`m+4`$ dimensional vector with the last coordinate equal to unity and all other coordinates equal to zero. We now analyze the large deviation rate $`L(v)`$ for this vector with respect to the unique stationary distribution $`\pi `$. Specifically, we show that the value of $`L(v)`$ depends on whether the counter machine halts. Indeed, if the counter machine halts in $`T`$ steps, then the value of $`Q(t)+q_1(t)+q_2(t)+q_3(t)`$ is bounded by $`2T+1`$ and, as a result, the stationary probability of the state $`nv`$, for large $`n`$ becomes zero. That is $`L(v)=+\mathrm{}`$. Now we show that if the counter machine does not halt, then $`L(v)=\mathrm{log}p`$. We compute $`\pi (nv)`$ by computing the expected return time to state $`nv`$, when the random walk is in this state at time $`0`$. Thus, we have $`q_3(0)=n,Q(0)=q_1(0)=q_2(0)=0`$. By the update rules of $`q_3`$, it decreases by one at each time step and at time $`t=n`$ it becomes zero. All the other components remain equal to zero. Beginning from this time $`t=n`$, the random walks keeps returning to the origin $`0`$ after some random time intervals. It is easy to see that the probability that the state $`nv`$ is visited in between any given two visits to the origin is exactly $`p^n`$ – the probability that the Bernoulli process survives at least $`n`$ steps. Let $`R_1,R_2,\mathrm{}`$ denote the random time intervals between successive visits to the origin. For a fixed $`m1`$ the probability that $`R_m`$ is the first interval during which $`vn`$ is visited is $`p^n[1p^n]^{m1},m=1,2,\mathrm{}`$, and the expected number of intervals $`R_m`$ before state $`nv`$ is visited for the first time is $`(p^n)^1`$. Let $`I_n`$ denote the indicator function for the event ”state $`vn`$ is visited between visits to the origin”. In particular, $`\mathrm{Prob}\{I_n\}=p^n`$. We now compute $`E[R_m|I_n]`$ and $`E[R_m|\overline{I}_n]`$. Note, that the time $`n`$ which takes to get from $`nv`$ to the origin plus the expected time it takes to get from the origin to $`nv`$, conditioned on $`I_n`$ is exactly $`E[R_m|I_n]`$. We then obtain that the expected recurrence time of the state $`nv`$ is $`(p^n)^1E[R_m|\overline{I}_n]+E[R_m|I_n]`$.

To compute $`E[R_m|]`$, recall from Proposition 1 that

$$E[R_m]=E[R_m|I_n]\mathrm{Prob}\{I_n\}+E[R_m|\overline{I}_n]\mathrm{Prob}\{\overline{I}_n\}=2/(1p).$$
(8)
If the Bernoulli process survives $`tn`$ steps then $`R_m=3t+1`$ steps. Then

$$E[R_m|I_n]=\frac{E[R_mI_n]}{\mathrm{Prob}\{I_n\}}=\frac{\underset{tn}{}(1p)p^t(3t+1)}{p^n}=\frac{(4p)p^n}{p^n}=4p$$
(9)
Then we obtain from (8)

$$E[R_m|\overline{I}_n]=\frac{2/(1p)(4p)p^n}{1p^n}$$

We conclude that the expected return time to the state $`nv`$ is

$$\frac{1}{\pi (nv)}=(p^n)^1\frac{2/(1p)(4p)p^n}{1p^n}+(4p),$$

and

$$\underset{n\mathrm{}}{lim}\frac{\mathrm{log}(\pi (vn))}{n}=\mathrm{log}p$$

as we claimed. We see that the value of $`L(v)`$ depends on whether the underlying counter machine halts or not. Specifically, by taking any rational value $`r>\mathrm{log}p`$, we conclude that the problem of checking whether $`L(v)r`$ is undecidable, by appealing again to Theorem 1. $`\mathrm{}`$

Remarks :

1. Note that we cannot determine the value of $`L(v)`$ even approximately, as we cannot distinguish between the cases $`L(v)<+\mathrm{}`$ and $`L(v)=+\mathrm{}`$. Contrast this with the results of Section 8.

2. We would not need extra coordinate $`q_3(t)`$ if we were interested in large deviation rate
$`lim_n\mathrm{}\pi (\overline{Q}(t))/n`$ of the stationary distribution of the norm of the state. The analysis would be identical to the one above.

## 7 Application to queueing systems

The results of the previous sections have implications to a certain type of queueing systems. A queueing system consisting of a single station processor and operating under a certain class of generalized priority policies was introduced in . It was shown that, similar to constrained random walks, determining stability for these queueing systems is an undecidable problem. In this section we consider the same class of system and show that computing stationary probabilities and large deviation rates are undecidable problems as well.

We start with the description of the system. Consider a single station queueing system $`𝒬`$ consisting of a single server and $`I`$ types of parts arriving externally. The parts corresponding to type $`i=1,2,\mathrm{},I`$ visit the station $`J_i`$ times. On each visit each part must receive a service before proceeding to the next visit. Only one part among all the types can receive service at a time. While waiting for service for the $`j`$-th time, the type $`i`$ part is stored in buffer $`B_{ij}`$. We denote by $`n`$ the total number of buffers $`n=_{i=1}^IJ_i`$. The service time for each part in each visit is assumed to be equal to unity. Each part can arrive into the system only in times which are multiples of some fixed integer value $`M`$. Specifically, certain values $`0p_i1`$ are fixed for each type $`i`$. For each type $`i`$ and each $`m=0,1,2,\mathrm{},`$ exactly one part arrives at time $`mM`$ with probability $`p_i`$ and no part arrives with probability $`1p_i`$, independently for all $`m`$ and all other types. In particular, interarrival times are geometrically distributed with expected interarrival time equal to $`1/\lambda _i=M/p_i`$, where, correspondingly, $`\lambda _i`$ is the arrival rate for type $`i`$.

A scheduling policy $`u`$ is defined to be a generalized priority policy if it operates in the following manner. A function $`u:\{0,1\}^n\{0,1,2,\mathrm{},n\}`$ is fixed. At each time $`t=0,1,2,\mathrm{}`$ the scheduler looks at the system and computes the binary vector $`b=(b_1,b_2,\mathrm{},b_n)\{0,1\}^n`$, where $`b_i=1`$ if there are parts in the $`i`$-th buffer and $`b_i=0`$, otherwise. Then the value $`k=u(b),0kn`$ is computed. If $`k>0`$ then the station processes a part in the $`k`$-th buffer. If $`k=0`$ the server idles. The map $`u`$ is assumed to satisfy the natural consistency condition: $`u(b)=k>0`$ only when $`b_k=1`$. That is, processing can be done in buffer $`k`$ only when there are jobs in buffer $`k`$. Note that the generalized priority scheduling policy is defined in finitely many terms and is completely state dependent - the scheduling decision at time $`t`$ does not depend on the state of the queueing system at times $`t^{}<t`$. A usual priority policy corresponds to the case when there is some permutation $`\theta `$ of the buffers $`\{1,2,\mathrm{},n\}`$ and $`u(b)=k`$ if and only if $`b_k=1`$ and $`b_i=0`$ for all $`i`$ such that $`\theta (i)<\theta (k)`$. In words, priority scheduling policy processes parts in buffers with lowest value (highest priority) $`\theta `$, which still has parts. Once we specify the queueing system $`𝒬`$ and some generalized priority policy $`u`$ we have specified some discrete time discrete space stochastic process. This process considered in times $`t=mM,m=0,1,2,\mathrm{}`$ is in fact a Markov chain.

Given a generalized priority policy $`u`$, a pair $`(𝒬,u)`$ is defined to be stable if there exists a finite number $`C>0`$ such that the total number of parts in the queueing system $`𝒬`$ at time $`t`$ does not exceed $`C`$ for infinitely many $`t`$ with probability $`1`$. In other words, the underlying Markov chain is positive recurrent. In this case there exists at least one stationary probability distribution. It is known that the necessary condition for stability is the following load condition

$$\rho \underset{i=1}{\overset{I}{}}\underset{j=1}{\overset{J_i}{}}\lambda _i<1.$$
(10)
This condition is also sufficient for stability if the policy is work conserving, which does not apply here, since we allow idling $`u(b)=0`$. We assume that the load condition above holds. We define a Lyapunov function and large deviations rates $`L(v)`$ for this queueing system in the same way we did for constrained homogeneous random walks in Section 6. As for constrained random walks, we now show that computing stationary probability distributions and computing large deviations rates for queueing systems operating under generalized priority policies is not possible. As for constrained random walks, we show that these problems are impossible to solve even if the underlying Markov chain is known to be irreducible and a linear Lyapunov function is available. Let $`\pi `$ denote the unique stationary distribution of a given irreducible positive recurrent queueing system $`(𝒬,u)`$. Let also $`0`$ denote the state of the system with all buffers empty.

###### Theorem 4

Given a queueing system $`𝒬`$ operating under some generalized priority policy $`u`$, given a linear Lyapunov function $`\mathrm{\Phi }`$ and given a rational value $`0r1`$, the problem of determining whether $`\pi (0)r`$ is undecidable. Likewise, given a vector $`v`$, the problem of determining whether $`L(v)r`$ is undecidable.

Proof : A reduction from a counter machine to a queueing system operating under some generalized priority policy was constructed in . This reduction had the following features. Given a counter machine with $`m`$ states, the corresponding queueing system had $`24`$ buffers and $`I=3m+7`$ arrival streams. There is a one-to-one correspondence between the configurations of the counter machine and states of the queueing system. In particular, if a counter machine has configuration $`(s_i,z_1,z_2)`$ at time $`t`$, then the queueing system at time $`(3m+26)t`$ is in a state which corresponds to this configuration in some well-defined way. We omit the details of this reduction and instead refer the reader to . We now modify the reduction to incorporate the extended random walk $`\overline{Q}(t)=(Q(t),q_1(t),q_2(t))`$ that was constructed in Section 5. Recall, that the part $`Q(t)`$ of this walk represented exactly $`m`$ states and the two counters of the underlying counter machine. We add two additional streams of arrivals which correspond to coordinates $`q_1(t)`$ and $`q_2(t)`$. We also construct additional buffers for $`q_1`$ and $`q_2`$ exactly in the way we did in  for counters $`z_1,z_2`$. The interarrival times for all the arrival streams, except for the stream corresponding to $`q_2`$, are deterministic and equal to some integer $`M`$ which is selected to be bigger than the number of buffers. For the stream corresponding to $`q_2`$, at most one part arrives at times $`Mt,t=0,1,2,\mathrm{}`$ independently for all $`t`$, and the probability that a part does arrive at time $`t`$ is equal to $`p`$, where $`p`$ is the parameter selected in construction of the random walk $`\overline{Q}(t)`$. Thus, $`p_i=p`$ for the arrival streams corresponding to $`q_2`$ and $`p_i=1`$ for all the other arrival streams. Finally, we modify the rules of the generalized priority policy to incorporate the rules by which the values of $`q_1(t),q_2(t)`$ are updated. This can be done in a way similar to the rules corresponding to $`z_1,z_2`$ in . We thus obtain a system which mimics the dynamics of $`\overline{Q}(t)`$ at times $`Mt,t=0,1,2,\mathrm{}`$. A linear Lyapunov function can be constructed again, provided that the parameter $`p`$ is sufficiently small. Arguing as in the proof of Theorem 2, we show that the problem of checking whether $`\pi (0)r`$ is undecidable. Similarly, we show that the problem of checking whether $`L(v)r`$ is undecidable, where $`v`$ is the unit vector with one in the coordinate corresponding to $`q_2(t)`$ and zero in all the other coordinates. For the latter case of computing large deviations rates, we add an additional arrival stream and buffers to represent the part $`q_3`$. $`\mathrm{}`$

## 8 Computing stationary probabilities approximately using a Lyapunov function

In this section we show that, despite the results of Section 5, computing the stationary probability is possible, if we are willing to tolerate some two-sided error and a computable geometric Lyapunov function $`\mathrm{\Phi }_g`$ exists. Our result is a simple consequence of the following result established by Meyn and Tweedie , which shows that infinite state Markov chain mixes exponentially fast when a geometric Lyapunov function can be constructed. The following is Theorem 2.3 proven in .

###### Theorem 5

Given an irreducible Markov chain $`Q(t)`$, suppose $`\mathrm{\Phi }_g`$ is a geometric Lyapunov function with a geometric drift $`\gamma _g<1`$ and the exception set $``$. Suppose also that $`\pi `$ is the unique stationary distribution. Then, there exist constants $`R>0,0<\rho <1`$ such that for any state $`x𝒳`$ and any function $`\varphi :𝒳\mathrm{}`$ satisfying $`\varphi (x)\mathrm{\Phi }(x),x𝒳`$, the following bound holds

$$\left|\underset{y𝒳}{}\varphi (y)\left(\mathrm{Prob}\{Q(t)=y|Q(0)=x\}\pi (y)\right)\right|\mathrm{\Phi }_g(x)R\rho ^t.$$
(11)
The constants $`R,\rho `$ are computable functions which depend on $`\gamma _g,\mathrm{max}_x\mathrm{\Phi }_g(x)`$ and

$$\nu _g^\mathrm{\Phi }=\underset{x,x^{}𝒵_+^d}{\mathrm{max}}\{\frac{\mathrm{\Phi }_g(x^{})}{\mathrm{\Phi }_g(x)}:p(x,x^{})>0\},p_{\mathrm{min}}^{}\underset{x,y}{\mathrm{min}}p(x,y).$$
(12)
Exact formulas for computing $`R,\rho `$ are provided in . They are quite lengthy and we do not repeat them here. These formulas give meaningful bounds only in case $`0<\gamma _g<1;\nu _g^\mathrm{\Phi }<\mathrm{};p_{\mathrm{min}}^{}>0`$.

Given a fixed state $`x_0𝒵_+^d`$, consider the function $`\varphi (x_0)=1/\mathrm{\Phi }_g(x_0),\varphi (x)=0,xx_0`$. This function satisfies the conditions of the theorem and one obtains a computable bound on the difference $`|\mathrm{Prob}\{Q(t)=x_0|Q(0)=x\}\pi (x_0)|`$, which decreases exponentially fast with $`t`$. This bound can be used for computing stationary probability distribution $`\pi `$.

###### Theorem 6

Given a constrained random walk $`Q(t)`$ in $`𝒵_+^d`$, given a state $`x_0𝒵_+^d`$ and an arbitrary value $`ϵ>0`$, under the conditions of Theorem 5, there exists an computable value $`\widehat{x}`$ which satisfies $`\pi (x)[\widehat{x}ϵ,\widehat{x}+ϵ]`$. In other words, the stationary probability of the state $`x_0`$ can be computed approximately with an arbitrary degree of accuracy.

Proof : The proof is a simple consequence of Theorem 5. We fix an arbitrary initial state $`Q(0)`$, say $`Q(0)=0`$. Compute the values $`R,\rho ,1/\mathrm{\Phi }_g(x_0)`$. Select $`t`$ large enough, so that $`\mathrm{\Phi }_g(Q(0))R\rho ^t<ϵ`$. Compute the transient probability $`Q(t)=x_0`$ conditioned on $`Q(0)=0`$. This can be done by direct calculation since $`t`$ is finite and from any state there are only finitely many neighboring states that can be entered with positive probability. The value $`\mathrm{Prob}\{Q(t)=x_0|Q(0)=0\}`$ can be taken as $`\widehat{x}`$, using inequality (11) and by the choice of $`t`$. $`\mathrm{}`$

As we mentioned above, a similar result cannot be established for large deviations rates $`L(v)`$, since the value of $`L(v)`$ changes between $`\mathrm{log}p`$ and $`+\mathrm{}`$ depending on whether the underlying counter machine halts or not. Therefore, computing the value of $`L(v)`$ even approximately still is an undecidable problem.

## 9 Conclusions

We considered in this paper the problems of computing stationary probability distributions and large deviations rates for constrained homogeneous random walks in $`𝒵_+^d`$. Both problems were shown to be undecidable – no algorithmic procedure for solving these problems can exist. An implication of these results is that no useful formulas for computing these quantities, for example along the lines of formulas for product form networks, can exist. For the problems of computing stationary probabilities, we showed that an approximate computation is possible with arbitrary degree of accuracy if a suitable geometric Lyapunov function can be constructed. Yet the problem of computing large deviation rates remains to be undecidable even in approximation sense as even checking whether a large deviation rate along a given vector is finite or not, is an undecidable problems.

We conjecture that these problems remain to be undecidable in more restrictive and interesting class of Markov chains corresponding to multiclass queueing networks operating under more conventional scheduling policies like First-In-First-Out or priority polices.


format: text

text:
# 1 Introduction

## 1 Introduction

Seven-dimensional manifolds of $`G_2`$ holonomy have long been known to exist. The construction of explicit non-compact $`G_2`$ metrics began ten years ago, when asymptotically conical metrics of cohomogeneity one were found . The physical interest of $`G_2`$ manifolds has increasingly significantly with the discovery of M-theory, because they are the most natural compactifying spaces from the eleven-dimensional point of view. It is expected that M-theory compactified on a $`G_2`$ manifold gives rise to an $`𝒩=1`$ super Yang-Mills theory in $`D=4`$ . The $`G_2`$ manifold with principal orbits $`S^3\times S^3`$ provides a geometrical demonstration of the classical flop of the type IIA superstring theory . In , M-theory dynamics on a $`G_2`$ manifold were discussed.

Recently, a large class of new metrics of $`G_2`$ holonomy have been obtained -, following the construction of the first examples of asymptotically locally conical spin(7) manifolds . These examples have non-abelian isometry groups. $`G_2`$ metrics with nilpotent isometry groups were also constructed in , which can be obtained by taking the Heisenberg or Euclidean limits of the non-abelian examples. Whilst it is of great interest to construct regular $`G_2`$ metrics, physically, it is essential to have an appropriate singularity structure to give rise to chiral fermions in $`D=4`$ .

In section 2, we construct a new non-compact cohomogeneity three metric with $`G_2`$ holonomy. The metric has three non-compact radial-type coordinates, with the principal orbits being a $`T^3`$ bundle over $`S^1`$. The isometry group of the metric is a four-dimensional nilpotent Lie group. The metric has either power-law singularities or delta-function singularities. The solution can be viewed as the intersection of three smeared Taub-NUTs. When one of the Taub-NUT charges is set to zero, the metric describes a product of $`S^1`$ with a six-dimensional non-compact Calabi-Yau manifold.

In section 3, we consider an M-theory vacuum (Minkowski)$`{}_{4}{}^{}\times _7`$, where $`_7`$ is the $`G_2`$ manifold. We show that by dimensionally reducing the solution on one of the circles in the $`T^3`$, we obtain a type IIA configuration with one D6-brane, one Taub-NUT and one 6-brane with an R-R 2-form flux. On the other hand, if we reduce the solution on the base space $`S^1`$, we obtain an intersection of three 6-branes all carrying R-R 2-form flux. These two configurations can be viewed as the classical flop in type IIA string theory on a non-compact six-dimensional Kähler manifold with a nilpotent isometry group. The origin of the flop is that the $`T^3`$ bundle over $`S^1`$ principal orbits can also be viewed as $`S^1`$ bundle over $`T^3`$.

In section 4, we perform a Kaluza-Klein reduction on the full principal orbits and the spatial world-volume. We obtain three perpendicularly intersecting membranes in $`D=4`$, describing a lattice universe. In this picture, the three non-compact coordinates of the $`G_2`$ manifold are identified with the spatial coordinates of our universe. We conclude the letter in section 5.

## 2 New $`G_2`$ metric

The metric ansatz is given by

$`ds_7^2`$ $`=`$ $`H_1dx_1^2+H_2dx_2^2+H_3dx_3^2+H_1H_2H_3dz_4^2+H_3H_1^1(dz_1+H_1^{}z_2dz_4)^2`$ (1)
$`+H_1H_2^1(dz_2+H_2^{}z_3dz_4)^2+H_2H_3^1(dz_3+H_3^{}z_1dz_4)^2,`$
where $`H_1`$, $`H_2`$ and $`H_3`$ are functions of $`x_1`$, $`x_2`$ and $`x_3`$ respectively. The prime on $`H_i^{}`$ denotes a derivative with respect to the argument of $`H_i`$;

$$H_1^{}=_{x_1}H_1,H_2^{}=_{x_2}H_2,H_3^{}=_{x_3}H_3.$$
(2)
The natural vielbein basis is

$`e^0=\sqrt{H_1H_2H_3}dz_4,e^1=\sqrt{H_1}dx_1,e^2=\sqrt{H_2}dx_2,e^3=\sqrt{H_3}dx_3,`$
$`e^4=\sqrt{H_3H_1^1}(dz_1+H_1^{}z_2dz_4),e^5=\sqrt{H_1H_2^1}(dz_2+H_2^{}z_3dz_4),`$
$`e^6=\sqrt{H_2H_3^1}(dz_1+H_3^{}z_1dz_4).`$ (3)
The associative 3-form in this basis is given by

$$\mathrm{\Phi }=e^{016}+e^{024}+e^{035}+e^{125}e^{134}+e^{236}e^{456},$$
(4)
where $`e^{ijk}=e^ie^je^k`$. The metric (1) has $`G_2`$ holonomy if and only if $`\mathrm{\Phi }`$ is closed and co-closed. We find that the closure and co-closure of $`\mathrm{\Phi }`$ implies that

$$H_i^{\prime \prime }=0,i=1,2,3,$$
(5)
implying that

$$H_1=1+m_1x_1,H_2=1+m_2x_2,H_3=1+m_3x_3.$$
(6)
Here the constant 1 is included so that $`H_i`$ does not vanish when $`m_i=0`$. Clearly the metric has a power-law singularity whenever any of the $`H_i`$ vanishes. The metric can also be recast in a “co-moving” frame,

$`ds_7^2`$ $`=`$ $`dr_1^2+dr_2^2+dr_3^2+\frac{9}{4}(m_1m_2m_3r_1r_2r_3)^{2/3}dz_4^2+\left({\displaystyle \frac{m_3r_3}{m_1r_1}}\right)^{2/3}(dz_1+m_1z_2dz_4)^2`$ (7)
$`+\left({\displaystyle \frac{m_1r_1}{m_2r_2}}\right)^{2/3}(dz_2+m_2z_3dz_4)^2+\left({\displaystyle \frac{m_2r_2}{m_3r_3}}\right)^{2/3}(dz_3+m_3z_1dz_4)^2.`$
The fibration in the $`z_1`$, $`z_2`$ and $`z_3`$ coordinates implies that the constants $`m_i`$ are quantised, namely

$$m_1=n_1\frac{L_1}{L_2L_4},m_2=n_2\frac{L_2}{L_3L_4},m_3=n_3\frac{L_3}{L_1L_4},$$
(8)
where $`n_i`$ are integers and $`L_i`$ are the periods of the $`z_i`$. For simplicity, we can set $`L_i=\mathrm{}_p`$ where $`\mathrm{}_p`$ is the Plank length, and then $`m_i=n_i/\mathrm{}_p`$.

The metric has a power-law singularity when any of the $`H_i`$ vanishes. This can be avoided by instead taking

$$H_i=1+\underset{\alpha }{}m_i^\alpha |x_ix_i^\alpha |.$$
(9)
such that $`H_i`$ is positive definite. However, in doing so, we have introduced delta function singularities.

When all three of the $`H_i^{}`$ are non-vanishing, the metric describes three intersecting Taub-NUTs with three independent non-vanishing smeared charges. The metric has three non-compact coordinates $`x_1`$, $`x_2`$ and $`x_3`$. The principal orbits are $`T^3`$ bundle over $`S^1`$; they are parameterised by the coordinates $`(z_1,z_2,z_3)`$ and $`z_4`$ respectively. The metric is of cohomogeneity three since it depends explicitly on the three non-compact coordinates $`x_i`$. Despite the dependence on the $`(z_1,z_2,z_3)`$ coordinates, they, together with $`z_4`$, parameterise a four-dimensional nilpotent Lie group $`𝒢`$, which is the isometry group of the metric, and thus the four-dimensional principal orbits are homogeneous.

When two of the $`H_i^{}`$ vanish, the metric describes a direct product of Euclidean 3-space and a smeared Taub-NUT. if instead only one of $`H_i^{}`$ vanishes, in which case the metric was obtained in , it describes a product of an $`S^1`$ with a Calabi-Yau 6-manifold. To see this in detail, let us set $`H_3=1`$. The metric of the Calabi-Yau manifold is then given by

$`ds_6^2`$ $`=`$ $`H_1dx_1^2+H_2dx_2^2+H_1H_2dz_4^2+H_2dz_3^2`$ (10)
$`+H_1^1(dz_1+H_1^{}z_2dz_4)^2+H_1H_2^1(dz_2+H_2^{}z_3dz_4)^2.`$
and the Kähler form is given by

$$J=e^0e^5e^1e^4+e^2e^6,$$
(11)
where the vielbein is given by (3) with $`H_3=1`$.

## 3 Intersecting D6-branes

Having obtained the new $`G_2`$ metric, one may consider an M-theory vacuum solution given by the direct product of Minkowski 4-spacetime and the $`G_2`$ manifold, namely

$$ds_{11}^2=dt^2+dw_1^2+dw_2^2+dw_3^2+ds_7^2.$$
(12)
The solution can be viewed as a triple intersection of smeared Taub-NUTs, with the metric represented by the diagram

$`t`$ $`w_1`$ $`w_2`$ $`w_3`$ $`x_1`$ $`x_2`$ $`x_3`$ $`z_1`$ $`z_2`$ $`z_3`$ $`z_4`$    $`H_1(x_1)`$ $`\times `$ $`\times `$ $`\times `$ $`\times `$ $``$ $`\times `$ $`\times `$ $``$ $``$ $`\times `$ $``$  $`H_2(x_2)`$ $`\times `$ $`\times `$ $`\times `$ $`\times `$ $`\times `$ $``$ $`\times `$ $`\times `$ $``$ $``$ $``$  $`H_3(x_3)`$ $`\times `$ $`\times `$ $`\times `$ $`\times `$ $`\times `$ $`\times `$ $``$ $``$ $`\times `$ $``$ $``$

Diagram 1. Triple intersections of Taub-NUTs. Here $`\times `$, $``$ and $``$ denote

the world-volume, transverse space, and fibre coordinates respectively.X

There is a $`U(1)`$ isometry for each of the $`z_i`$ coordinates, and so we can reduce the metric on any $`z_i`$ to obtain a solution in type IIA theory. The $`z_i`$ for $`i=1,2,3`$ are equivalent, and hence the reduction can be discussed using $`z_1`$ as a representative. The resulting type IIA solution is given by

$`e^{{\scriptscriptstyle \frac{1}{6}}\varphi }ds_{10}^2`$ $`=`$ $`dt^2+dw_1^2+dw_2^2+dw_3^2+H_1dx_1^2+H_2dx_2^2+H_3dx_3^2+H_1H_2H_3dz_4^2`$
$`+H_1H_2^1(dz_2+H_2^{}z_3dz_4)^2+H_2H_3^1dz_3^2+H_3H_1(H_1^{}z_2)^2dz_4^2`$
$`W^1(H_3H_1^1H_1^{}z_2dz_4H_2H_3^1H_3^{}z_4dz_3)^2,`$
$`e^\varphi `$ $`=`$ $`W^{{\scriptscriptstyle \frac{3}{4}}},W=H_3H_1^1+H_2H_3^1(H_3^{}z_4)^2,`$
$`A_{\left(1\right)}`$ $`=`$ $`W^1(H_3H_1^1H_1^{}z_2dz_4H_2H_3^1H_3^{}z_4dz_3).`$ (13)
Note that before performing the Kaluza-Klein reduction, we have made a coordinate transformation $`z_3z_3H_3^{}z_1z_4`$ in the metric (1). Clearly, the solution describes an intersection of three objects. The one parameterised by $`H_1`$ is a smeared D6-brane, and the one parameterised by $`H_2`$ is a Taub-NUT. The one associated with $`H_3`$ is a 6-brane carrying an R-R 2-form flux, but it differs from a standard D6-brane.

We can instead reduce the solution on the $`z_4`$ coordinate, giving the type IIA solution

$`e^{{\scriptscriptstyle \frac{1}{6}}\varphi }ds_{10}^2=dt^2+dw_1^2+dw_2^2+dw_3^2+H_1dx_1^2+H_2dx_2^2+H_3dx_3^2`$
$`+H_3H_1^1dz_1^2+H_1H_2^1dz_2^2+H_2H_3^1dz_3^2`$
$`W^1(H_3H_1^1H_1^{}z_2dz_1+H_1H_2^1H_2^{}z_3dz_2+H_2H_3^1H_3^{}z_1dz_3)^2,`$
$`e^\varphi =W^{{\scriptscriptstyle \frac{3}{4}}},W=H_1H_2H_3+H_3H_1^1(H_1^{}z_2)^2+H_1H_2^1(H_2^{}z_3)^2+H_2H_3^1(H_3^{}z_1)^2,`$
$`A_{\left(1\right)}=W^1(H_3H_1^1H_1^{}z_2dz_1+H_1H_2^1H_2^{}z_3dz_2+H_2H_3^1H_3^{}z_1dz_3).`$ (14)
The solution describes three intersecting 6-branes all carrying R-R 2-form flux. These 6-branes are different from the usual D6-brane coming from the reduction of the fibre coordinate of a Taub-NUT in $`D=11`$.

The two configurations arising from the reduction on $`z_1`$ or $`z_4`$ can be viewed as a classical flop in the type IIA string theory on the non-compact Kähler manifold. The flop in $`D=10`$ can be geometrically explained by the fact that the $`T^3`$ bundle over $`S^1`$ principal orbits of the four-dimensional nilpotent Lie group $`𝒢`$ can also be described as an $`S^1`$ bundle over $`T^3`$. However, the two descriptions are somewhat different. In the latter case, the fibre is the circle group in $`𝒢`$ generated by $`\frac{}{z_4}`$. In the former case, the fibre is not the orbit of a three-dimensional subgroup of $`𝒢`$ because $`\frac{}{z_i}`$ for $`i=1,2,3`$ are not themselves Killing vectors; we must add a multiple of $`\frac{}{z_4}`$. In fact the flop involves interchanging the fibre and base spaces of the $`U(1)`$ fibration. This is analogous to the flop discussed in .

## 4 Lattice universe

The new $`G_2`$ metric (1) that we have obtained is in fact inspired by the four-dimensional intersecting membrane solution that describes the lattice universe .111In , a triply quasi-periodic Gibbons-Hawking metric was obtained. There has been experimental evidence suggesting that the network of galaxy superclusters and voids seems to form a three-dimensional lattice with a spacing of about $`120h^1Mpc`$ (where $`h^1`$ is the Hubble constant in units of $`100kms^1Mpc^1`$) . In , an M-theory solution was constructed to describe such a lattice structure, which can be realised by considering non-standard brane intersections of two M5-branes and one Taub-NUT, or two Taub-NUTs and one M5-brane. In the latter case, turning off the M5-brane charge causes the solution to reduce to the product of 5-dimensional Minkowski spacetime and the non-compact Calabi-Yau manifold given in (10).

In section 2, we obtained the new $`G_2`$ metric (1) by adding an extra fibration on the seventh coordinate. This procedure follows the general prescription of obtaining $`G_2`$ manifolds from six-dimensional Kähler manifolds, described in detail in .

If we reduce the M-theory solution on the world-volume spatial coordinates $`w_i`$ and also on the $`T^3`$ bundle over $`S^1`$ principal orbits, we obtain three perpendicularly intersecting membranes in $`D=4`$, with the metric

$$ds_4^2=(H_1H_2H_3)^{{\scriptscriptstyle \frac{1}{2}}}(dt^2+H_1dx_1^2+H_2dx_2^2+H_3dx_3^2).$$
(15)
This metric was first obtained in , although it was supported by very different field strength. The functions of $`H_i`$ in this case are given by (9) describing periodic arrays of intersecting membranes.

In this static cosmological model the spatial world-volume and the $`T^3`$ bundle over $`S^1`$ principal orbits are viewed as an internal space, whilst the three non-compact coordinates $`x_1`$, $`x_2`$ and $`x_3`$ of the $`G_2`$ manifold are identified with the spatial coordinates of our universe. This is rather natural since the principal orbits are clearly compact, and the spatial world-volume can be wrapped on a compact space such as $`T^3`$. Although it is not likely that the metric (15) describes our actual universe, since it preserves $`𝒩=1`$ supersymmetry; it is nevertheless rather suggestive that the lattice structure should emerge from a metric with $`G_2`$ holonomy.

## 5 Conclusions

In this letter, we constructed a new cohomogeneity-three metric with $`G_2`$ holonomy. It has three radial-type coordinates, with the principle orbits being a $`T^3`$ bundle over $`S^1`$. The solution can be viewed as three intersecting Taub-NUTs. We performed Kaluza-Klein reduction on the $`S^1`$ and instead on a circle in the $`T^3`$. The two resulting type IIA configurations can be viewed as a classical flop in type IIA string theory on a non-compact Kähler six-manifold. Although the type IIA solutions do not describe the triply intersecting D6-branes advocated in  for the realisation of chiral fermions, it is nevertheless of interest to investigate further if chirality could arise from the singularities of our $`G_2`$ metric. The metric provides a concrete example for studying such an issue, since it can be viewed as the lifting of a D6-brane configurations in $`D=10`$.

If we perform a Kaluza-Klein reduction on the entire four-dimensional full principal orbits and the spatial world-volume, we obtain triply intersecting membranes in $`D=4`$, describing a lattice universe. The construction takes full advantage of the non-compact nature of the manifold, in that the three non-compact coordinates are precisely identified with the spatial coordinates of our universe. It is of great interest to investigate further the significance of such a configuration arising from a $`G_2`$ manifold.

## Acknowledgement

We are grateful to Mirjam Cvetič for extensive discussions on triply intersecting D6-branes and chiral fermions, to Gary Gibbons for pointing out the isometry group of the metric, and to Chris Pope and Justin Vázquez-Poritz for discussions.