Convergence in the Calculation of the Handoff Arrival Rate: A Log-Time Iterative Algorithm

Modeling to study the performance of wireless networks in recent years has produced sets of nonlinear equations with interrelated parameters. Because these nonlinear equations have no closed-form solution, the numerical values of the parameters are calculated by iterative algorithms. In a Markov chain model of a wireless cellular network, one commonly used expression for calculating the hando ﬀ arrival rate can lead to a sequence of oscillating iterative values that fail to converge. We present an algorithm that generates a monotonic sequence, and we prove that the monotonic sequence always converges. Lastly, we give a further algorithm that converges logarithmically, thereby permitting the hando ﬀ arrival rate to be calculated very quickly to any desired degree of accuracy.


INTRODUCTION
Wireless cellular networks provide service to mobile terminals, which can move from a given cell to any adjacent cell multiple times during the lifetime of a particular call. Therefore, a wireless network must take into account the rate at which ongoing calls arrive from neighboring cells, in addition to the arrival rate of new calls. When a user crosses the boundary from one cell to another, the network must react by handing off the call. However, there must be a channel available in the new cell for that call, or else the handoff fails and the service is abruptly terminated.
One approach for improving the likelihood that a free channel is available when a handoff call arrives is the dedication of a certain number of channels in each cell purely for handoff calls. These dedicated channels are called guard channels, and earlier works have focused on the benefit of determining the number of guard channels dynamically.
Models of cellular networks are very important for design as well as operation of the network. During the operation of a network, performance parameters can be estimated empirically by collecting data while the network is in operation. In fact, most of the current networks collect performance data and use it for decision making. However, if a network's performance is outside the desired range, some of the control parameters will need adjustment. The amount of adjustment to be made is determined from a model.
Simulation as well as analytical models are used for designing networks. Simulation models require a long computation time. However, in the absence of analytical models, simulation is the only available tool. Also, simulation models are necessary for the final evaluation of networks designed using approximate analytical models. For instance, even though call holding times and cell dwell times do not follow exponential distributions (see [1,2]), analytical models assume that they are exponentially distributed. Therefore, a cellular network designed using such a model can be finetuned using simulation.
On the other hand, analytical models are computationally efficient. One can estimate performance parameters very quickly. For instance, the (fuzzy associative memory) FAMbased call admission controller, reported in [3,4], used a simulation model for development of the FAM. It took about two months of simulation time on a Pentium IV PC to develop the FAM. However, the FAMs for the call admission controllers reported in [5,6] were developed using the algorithm presented in this paper requiring about a day on the same Pentium IV PC.
Since the late eighties, the modeling of wireless cellular networks for analysis of their performance has produced sets 2 EURASIP Journal on Wireless Communications and Networking of nonlinear equations with interrelated parameters. These nonlinear equations have no closed-form solution, so the numerical values of the parameters are calculated by iterative algorithms. When these iterations fail to converge, however, the precise values of the parameters cannot be determined.
The foregoing applies to wireless cellular networks, for which many studies have used Markov chains as models [7][8][9][10]. Some of these models treat all calls identically, while others create a priority status for handoff calls. With respect to the prioritization of handoff calls, there are two basic approaches: (a) the early reservation of channels and (b) the use of guard channels that are dedicated exclusively to handoff calls [8][9][10][11].
The number of guard channels can be established in advance (statically) or as an ongoing process (dynamically) (see [12,Chapter 2]). In the former case, bandwidth may be underutilized or handoff call failure rate may be too high. In the latter case, there must be a continual computation of the optimal number of guard channels [7,11]. This in turn creates a need for the computation of the handoff arrival rate. Note that although current handoff call arrival rate can be estimated from some "time averaging" of recent handoff call arrival records, the handoff call arrival rate that would result from the change of the number of guard channels must be determined from simulation or analytical model. Since simulations require a long time, analytical models are more desirable.
The absence of a closed-form expression for the handoff arrival rate requires an alternate method, which commonly involves the use of iterative algorithms [7]. One standard formula for the calculation of the handoff arrival rate generates a sequence of approximations that may oscillate around the actual value. When this sequence converges, a result is obtained within any desired degree of accuracy. Convergence is not guaranteed, however, and in that instance the sequence develops a bifurcation and oscillates repeatedly between two values above and below the actual handoff rate value.
In [13], fixed-point iteration for calculating the handoff arrival rate is proposed and used to overcome numerical overflow problems when a cell has a large number of channels. The paper also presents an algorithm for computing the optimal number of guard channels, but the optimization algorithm uses the proposed fixed-point iterative algorithm. The authors of that paper indicate that a proof of convergence of their algorithm is an open problem (last sentence of [13,Section VI]). The iterative algorithm in [7] attempts to avoid any potential nonconvergence by partitioning and bounding the solution interval. However, the process is rather slow-linear with the inverse of the desired accuracy.
In this paper, we present a novel iterative algorithm that always converges and which is logarithmic in nature (thereby assuring a relatively fast convergence). We also present proof of convergence of the algorithm. One can find further details of the work reported here in [14].

Definitions and notation
It is assumed that each cell in a network has a fixed number of channels, and at any given moment somewhere between none and all of them will be in use. Moreover, the cells are assumed to be identical, that is, the system is homogeneous. Calls arriving into a cell can be from one of two sources: (a) a call that was previously accepted by the network and that is now being handed off from an adjacent cell (a handoff ) and (b) a brand-new call that has just been received by the cell (a new call). Two time frames are relevant. The average length of time that a given call remains active from inception to uninterrupted completion is referred to as the holding time, whereas the average amount of time that a call remains in any given cell before departing is the dwell time. 1 Calls depart from a cell for one of two reasons: (a) the mobile terminal moves to an adjacent cell or (b) the customer completes the call and terminates the connection. These departures are distinct from calls that never enter the cell (although there is an attempt to enter). For a new call, if there is no available channel, then the call is simply not accepted. For a handoff, if similarly there is not an available channel, then the handoff fails and the existing call is forced to terminate.

Organization of the remaining sections
In Section 2, we first give a set of nonlinear equations for the parameters derived from a Markov model of a wireless cellular network. We then present one commonly used expression for iterative calculation of the handoff arrival rate and include an algorithm. Next, we use a straightforward example that shows that the iterations converge with one set of values, but fail to converge when a very slight change is made to one of the parameters. We finish the section by explaining the source of the oscillating nonconvergence and by proposing to use an alternative expression and an accompanying novel algorithm for calculating the handoff arrival rate that always converges. In Section 3, we give a rigorous proof that our novel approach always converges, both for the nonpriority case (no guard channels) and the priority case (a network with guard channels). In Section 4, we take the earlier results and give an algorithm that not only converges, but does so logarithmically. Section 6 contains our concluding remarks.

MARKOV MODEL AND CALCULATION OF HANDOFF ARRIVAL RATE
We first refine the definitions of two items and then express the steady state probabilities for a homogeneous cellular wireless network with C channels per cell, of which n are nonguard channels (see [8,9,13] for the derivation of the following equations). The offered load and the handoff 3 load are more precisely The steady state probabilities where one or more channels are in use can be split into two portions: (a) states in which any arriving call, whether a new call or a handoff one, will be accepted and (b) states in which only handoff calls will be accepted. These are P j = ρ j j! P 0 , for 0 < j ≤ n, The probability for state 0 (the state in which no channels are in use) is a normalization obtained from the fact that C j=0 P j = 1, The blocking probability of a new call (P b ) and the handoff failure probability of an ongoing call (P h f ) are given by

Existence of actual value for λ h
Before giving an expression for the calculation of the handoff arrival rate (λ h ), we note the possible range of values. Clearly the value cannot be negative, so zero is a lower bound. The quantity ηC is an upper bound, since the rate cannot exceed the number of channels C in the cell divided by the average dwell time 1/η. Also, since a finite, irreducible, positive recurrent Markov chain models a cell, it has a unique stationary distribution, and hence P b and P h f are uniquely determined. Since a standard expression for the handoff arrival rate is (see [8,9,13] for the details), for given values of λ 0 , μ, η, C, and n, the handoff rate is determined uniquely by (5). We will denote this unique value by λ * h . Iterative algorithms are typically used for the calculation of the handoff arrival rate, where the value from one iteration is then fed into the equation, thereby producing successive values (see (6)). The hope is that these iterations will converge, but some approaches do not always converge.

Standard approach
The iterative form is some small value > 0 λ h := new λ h := 0 do { the following steps } Step 1: λ h := new λ h Step 2: update values for the offered load ρ and handoff load ρ h per (1) Step 3: update the value of P 0 per (3) Step 4: update the values of state probabilities P 1 through P C per (2) Step 5: update the blocking probability P b and handoff failure probability P h f per (4) Step 6: compute the new value for λ h , that is, new λ h , per (6) while where P b (k) and P h f (k) are the values derived from using λ h (k) in the kth iteration. An algorithm that incorporates this approach is in Algorithm 1.
We will show that this approach works in some situations, but can also lead to oscillations that do not converge.
In our examples, assume that there are twenty channels in each cell, of which four are guard channels, and that for any given call the average duration (holding time) is 120 seconds and the average time in any given cell before departing (dwell time) is twelve seconds. Hence, we have the following values:

Works sometimes
Without loss of generality, we will choose zero as the initial value for the handoff arrival rate (i.e., λ h (0) = 0). If the new call arrival rate (λ 0 ) is a relatively low figure, such as 0.1, then using (6) will result in λ h = 0.899 027 7. Figure 1 shows the plot of the sequence of calculated values for the handoff arrival rate beginning with the initial value of λ h (0) = 0. The convergence occurs fairly quickly.

Can oscillate and not converge
On the other hand, increasing the value for λ 0 (the new call arrival rate) very slightly to 0.12 is sufficient to produce oscillations that do not converge. Once again using the initial value of λ h (0) = 0, we obtain from (6) the alternating pair of 1.170 851 55 and 0.712 858 as the calculated values for λ h . Figure 2 illustrates the oscillations.

(1) Why oscillation occurs
Referring to (6) the handoff failure probability P h f (k). The blocking probability is the sum of the steady state probabilities for those states where only handoff calls will be accepted (i.e., states n through C). When P b (k) is very low, the numerator of (6) becomes larger and results in a higher calculated value for λ h (k + 1). The handoff failure probability is the steady state probability for the final state (i.e., state C), and if P b (k) is low, then so will be P h f (k). The combination of low calculated values for P b (k) and P h f (k) produces a higher value for λ h (k + 1). When that higher value is then fed into (6), the system's general load (ρ(k + 1)) and handoff load (ρ h (k + 1)) are correspondingly higher. This shifts the weighted average of the state probabilities to the right, with the result that the guard states (states n through C) have higher probabilities. Thus, for this iteration, both P b (k + 1) and P h f (k + 1) increase. These increases result in a smaller numerator and larger denominator in (6), thereby producing a smaller calculated value for λ h (k + 2) for the next iteration.
This alternation between higher and lower values for the sequence λ h (k) can prevent convergence. The problem oc- States (i.e., number of channels in use)

Avoiding nonconverging oscillations
Rather than using (5) (or its iterative form, (6)) for the calculation of the handoff arrival rate (λ h ), we instead use the basic expression from which (5) is derived (see [8,13] for details). In general, the value for λ h is the expected number of channels in use (call it E(N)) divided by the average dwell time, that is, The value for E(N) is simply the weighted average of the number of channels in use: Just as there exists a unique steady state value for λ h given a set of values for the other system components (number of channels, number of guard channels, holding time, dwell time, and new call arrival rate), there is similarly a unique steady state value for E(N).
Combining these ideas, we obtain the following expression for the handoff arrival rate: The iterative form of the equation is The algorithm is identical to the one for the standard approach with the exception that new λ h is calculated using (11) (instead of (6)) and there is now no need to calculate the blocking probability (P b ) or handoff failure probability (P h f ).
As with the iterative algorithm using (6), Algorithm 2 is an iterative algorithm where the calculated value for λ h from one iteration is plugged into the next iteration. In contrast to the use of (6), however, the use of (11) always converges and does not experience the oscillations that plagued the first algorithm.
By way of illustration, Figures 4 and 5 show the results from using the set of values for C, n, μ, and η from (7) and the values 0.1 and 0.12, respectively, for λ 0 . In both cases we ob-  The reason that the iterative algorithm using (11) always converges is that two things occur simultaneously, and both are monotonic. If we begin with an initial value for λ h (0) that is less than the actual value λ * h , (1) each successive iteration produces values for E(N) and λ h that are larger than their immediate predecessor values; (2) no matter how many iterations are done, the calculated values for E(N) and λ h always remain less than the respective actual values.
If we start with an initial value for λ h (0) that is greater than the actual value, then the reverse holds true (i.e., the 6 EURASIP Journal on Wireless Communications and Networking iterations produce successively smaller calculated values for E(N) and λ h , but always greater than the actual values).

CONVERGENCE OF THE MONOTONIC APPROACH
We will denote by {P j (i)}, 0 ≤ j ≤ C, the steady state probability distribution of the one-dimensional finite, irreducible, positive recurrent Markov chain on {0, 1, . . . , C} determined by the parameters at the ith iteration of the algorithm. Our approach would be to show that these probability vectors satisfy a likelihood ratio ordering, which, as is well known (Lehmann [15, Section 3.3] and Shanthikumar [16]), implies strong stochastic ordering. The special case of this result that we use is stated for ease of reference and completeness.

Lemma 1. Suppose for nonnegative integers i 1 and i 2 ,
for all j, 0 ≤ j < C, with the convention that 0/0 = 0. Then for alll, 0≤l≤C, with strict inequality for at least one positive l.
Proof. Let j 0 be the least integer in {0, . . . , C} such that P j0 (i 1 ) > P j0 (i 2 ). Such an integer must exist because C j=0 P j (i 1 ) = C j=0 P j (i 2 ) = 1 and because of the ratio inequality. The remainder of the conditions imply that We cannot have j 0 = 0, because C j=0 P j (i 1 ) = C j=0 P j (i 2 ) = 1. Therefore we have For l < j 0 , the assertion follows by summing both sides of the inequality (15) over 0 ≤ j < l and subtracting from one. For l ≥ j 0 , the assertion follows by summing both sides of the inequality (16) over l ≤ j ≤ C. Moreover, we have strict inequality for all l ≥ j 0 .

Lemma 2. Suppose for nonnegative integers i 1 and i 2 ,
for all j, 0 ≤ j < C, with the convention that 0/0 = 0. Then Proof. By Lemma 1, we have for all l, 1 ≤ l ≤ C, with strict inequality for at least one l. Summing over all l, 1 ≤ l ≤ C, we obtain Interchanging the order of summation yields equivalently, This proves the result.

Lemma 3 (ratio lemma). For any nonnegative integers i 1 and
where the offered loads are Proof. The proof breaks down into two cases: (a) no guard channels and (b) the presence of guard channels. Where there are no guard channels, the proof follows essentially from the fact that in general the ratio of successive states in the same iteration results in If we have ρ(i 1 ) > ρ(i 2 ), then we can move from Similarly, if we start with the inequality between ratios of successive states in the same iteration, then we can end up with the inequality between loads ρ(i 1 ) and ρ(i 2 ). Hence the lemma holds in both directions where there are no guard channels.
In the presence of guard channels, there is an extra step involved in computing some of the ratios. Assume there are n nonguard channels. For P j+1 (k)/P j (k), where 0 ≤ j < n, the situation is identical to the one where there are no guard channels. For j = n, however, the numerator represents a guard channel state, whereas the denominator is a nonguard channel state. For n < j < C, both the numerator and denominator are guard channel states. We show that these ratios in fact lead to the same expression, which in turn verifies the lemma.
Dilip Sarkar et al. 7 Where j = n, we get Similarly, for n < j < C, we get If we have ρ h (i 1 ) > ρ h (i 2 ), then λ h (i 1 ) > λ h (i 2 ). We can add the new call arrival rate (λ 0 ) to each side and then divide by μ + η, giving us ρ(i 1 ) > ρ(i 2 ). We can then move from Similarly, if we start with the inequality between ratios of successive states in the same iteration, then we can end up with the inequality between loads ρ h (i 1 ) and ρ h (i 2 ). Hence the lemma also holds in both directions in the presence of guard channels.
We use these lemmas for showing convergence in our approach. In the next subsection we state and prove these theorems.

Convergence theorems
The technique of showing that the successive iterations of λ h (k) produce calculated values for the handoff arrival rate that monotonically approach the actual value λ * h demonstrates that (11) always converges.
If the initial value λ h (0) equals λ * h , we must have λ h (1) = λ * h , and the computation terminates. This is because in this case P j (0), 0 ≤ j ≤ C, are the steady state probabilities of the Markov chain. Since the steady state distribution is unique, any initial value λ h (0) not equal to λ * h would yield a λ h (1) that is not equal to λ h (0). Hence we must have the following: Here are the theorems on monotonic convergence of the proposed algorithm.
Theorem 1. Assume the use of (11) for the calculation of the successive values of λ h (k). If the initial value chosen for λ h (0) is not equal to λ * h , the sequence λ h (k), k = 1, 2, . . . , is monotonic.
If λ h (m+1) < λ h (m), by definition (see (1)) we have ρ(m+ 1) < ρ(m). Therefore, by Lemma 3 we have Now, by Lemma 2, we have The two cases considered above in the proof of Theorem 1 immediately, leads to the following two corollaries. Corollary 1. Assume the use of (11) for the calculation of the successive values of λ h (k). If the initial value chosen for λ h (0) is less than the actual value λ * h , the sequence λ h (k), k = 1, 2, . . . , is monotonically increasing. Corollary 2. Assume the use of (11) for the calculation of the successive values of λ h (k). If the initial value chosen for λ h (0) is greater than the actual value λ * h , the sequence λ h (k), k = 1, 2, . . . , is monotonically decreasing.
The following theorem asserts the convergence of the computation, in all cases, to the desired value.
Theorem 2. Assume the use of (11) for the calculation of the successive values of λ h (k). For any initial value of λ h (0), where {P j (k)}, 0 ≤ j ≤ C, is the steady state probability distribution of the one-dimensional finite, irreducible, positive recurrent Markov chain on {0, 1, . . . , C} determined by the parameters at the kth iteration of the algorithm.

EURASIP Journal on Wireless Communications and Networking
C j=l P j (k) is a monotone sequence in k for each l. Therefore all these sequences have a limit as k → ∞, and consequently P j (k) has a limit for all j. Therefore, taking limits in (11), we get Since the limits satisfy the balance equations, by uniqueness of the steady state distribution we must have lim k→∞ λ h (k) = λ * h .

FASTER CONVERGENCE BY A BISECTION ALGORITHM
Although the monotonic algorithm given in Section 2.5 does converge, the rate is much slower than necessary for practical applications in cellular networks. A faster approach makes use of the fact that the successive values of λ h (k) are monotonic. The basic idea is to take two values, lowλ h and hiλ h , that are known to be lower and higher, respectively, than λ * h . These two values are averaged, and the result is deemed to be the testValue for λ h . The testValue is then fed into the iterative process (11), which produces a resultValue.
By virtue of monotonicity, if the resultValue is less than the testValue, then we know that λ * h is less than the result-Value. In other words, In that case, we keep the same lowerValue and we make the resultValue the new higherValue. In the same manner, if a re-sultValue is greater than the testValue that produced it, then we know that λ * h is greater than the resultValue. Now the relationships are Here, the higherValue would remain the same, and the re-sultValue becomes the new lowerValue. The lower and higher values are averaged, which produces a new testValue. This continues until the difference between the lower and higher values is within the desired accuracy of the user. For original lower and higher values, we use the lower and higher bounds for λ * h , namely, 0 and ηC. The foregoing is captured in Algorithm 3.
This approach is an improvement over the monotonic algorithm, which merely used the result from one iteration as the initial value for the next iteration. By taking advantage of the knowledge given to us by Corollaries 1 and 2, we know from the relationship between testValue and result-Value whether the actual value λ * h is greater than or less than the resultValue, and we can adjust the lower or higher bound accordingly as we hone in on the actual value. In fact, our approach is even stronger than a pure bisection, because we are able to use resultValue (and not just the testValue) as the new lower or higher value for the following iteration. Hence, the some small value > 0 lowλ h := 0 hiλ h := ηC while (hiλ h − lowλ h > ) Step 1: testValue := (low λ h + hiλ h )/2 Step 2: update values for the offered load ρ and handoff load ρ h per (1) Step 3: update the value of P 0 per (3) Step 4: update the values of state probabilities P 1 through P C per (2) Step 5: compute the new value for λ h , that is, resultV alue, per (11) Step 6: if (resultValue < testValue) then range [lowerValue, higherValue] shrinks by more than onehalf with each iteration. We illustrate with two charts the speed with which the proposed bisection algorithm can achieve a very accurate approximation of λ * h quickly. In Figure 6, a value of 0.1 was used for λ 0 , and the result from the bisection algorithm is combined with results from Figures 1 and 4. Figure 7 is similar, using a value of 0.12 for λ 0 and combining with the results from Figures 2 (which did not converge) and 5 (which did converge, albeit somewhat slowly).
The convergence properties of the bisection algorithm can be expressed in a theorem. Theorem 3. The bisection algorithm converges. Moreover, for a given degree of accuracy > 0, the number of iterations required to achieve that level of accuracy is on the order of Proof. We begin with the maximum possible range of values for λ * h , which is [0, ηC]. Because of Corollaries 1 and 2, with each iteration one end of the range is adjusted in the direction of the actual value of λ * h , always keeping the actual value within the range, and hence the range continues to shrink as the number of iterations increases. We note that the initial gap is simply ηC − 0 = ηC. The gap is actually divided by more than a factor of 2 with each iteration. That can be observed from the fact that testValue is the average of the current range endpoints, but resultValue replaces one of the endpoints for the next iteration (leaving the other endpoint Dilip Sarkar et al.
Now the logarithmic convergence can be established from the following classical argument. For a given value of > 0, we need to keep dividing the gap until the range is within the desired degree of accuracy. The number of steps, call it m, needed to accomplish this can be expressed as which means The smallest integer n that satisfies this inequality is the maximum number of required iterations. This completes the proof of the theorem.

MODEL VALIDATION
For validation of the accuracy of handoff call arrival rates obtained from the algorithms presented in previous sections, we developed a simulation model. The cell layout for our simulation model is shown in Figure 8. The 49 white cells are part of the model and the shaded ones show the wraparound neighbors. The wraparound topology is used, since it eliminates the boundary effect keeping exactly six neighbors for each cell [9]. We assume a static channel allocation scheme for cells, that is, the number of channels allocated to a cell does not  change during the simulation. For the results reported here, all cells were assigned 20 channels. Mobility of terminals is modeled using a simple random walk, that is, a terminal moves to any of the current cell's neighbors with equal probability-1/6 for the hexagonal layout. New call arrivals into the network follow the Poisson distribution with mean λ calls/s. The call holding time and the cell dwell time follow exponential distributions with respective means 1/μ and 1/η seconds. For obtaining good estimates of the parameters, each simulation study was run for 1 000 000 new calls. Note that the assumptions for the simulations are identical to those to our analysis. Our extensive studies have shown a close match between the theoretical and simulation results.
Some typical results are shown in Table 1. We varied call arrival rate from 0.06 to 0.2 calls per second. The average call holding time and cell dwell time were kept constant at 120 seconds and 12 seconds, respectively. Out of the 20 channels, 4 were used as guard channels. As can be seen from the table, handoff call arrival rates calculated by the algorithm and obtained from simulations agree up through the hundredth place. Therefore, handoff call arrival rates calculated from the presented algorithm are very accurate.

CONCLUSION
Since the late eighties, the modeling and analysis of the performance of wireless networks have produced sets of nonlinear equations with interrelated parameters. These nonlinear equations have no closed-form solution, so the numerical values of the parameters are calculated by iterative algorithms. When these iterations fail to converge, however, the precise values of the parameters cannot be determined.
Using a Markov chain to model a wireless cellular network, we discussed a common expression for calculating the handoff arrival rate iteratively. We then provided for illustration an instance where the sequence of iterative values fails  to converge. After explaining the reason for the nonconverging oscillations, we gave an alternate simple iterative algorithm that generates a monotonic sequence and proved that the monotonic sequence always converges. Lastly, we refined this algorithm and, drawing upon the earlier results of this paper, set forth another algorithm that converges logarithmically.
The proposed algorithm can be used in existing cellular network optimization and call control algorithms [10,13].