Table of Contents Author Guidelines Submit a Manuscript
Advances in Operations Research
Volume 2011 (2011), Article ID 216790, 26 pages
Research Article

Two Coupled Queues with Vastly Different Arrival Rates: Critical Loading Case

1Department of Mathematics, Statistics and Computer Science, University of Illinois at Chicago, 851 South Morgan Street Chicago, IL 60607-7045, USA
2Consultant Bell Laboratories, Alcatel-Lucent, 600 Mountain Avenue, Murray Hill, NJ 07974, USA

Received 8 September 2010; Accepted 18 January 2011

Academic Editor: A. Gómez-Corral

Copyright © 2011 Charles Knessl and John A. Morrison. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.


We consider two coupled queues with a generalized processor sharing service discipline. The second queue has a much smaller Poisson arrival rate than the first queue, while the customer service times are of comparable magnitude. The processor sharing server devotes most of its resources to the first queue, except when it is empty. The fraction of resources devoted to the second queue is small, of the same order as the ratio of the arrival rates. We assume that the primary queue is heavily loaded and that the secondary queue is critically loaded. If we let the small arrival rate to the secondary queue be 𝑂(𝜀), where 0𝜀1, then in this asymptotic limit the number of customers in the first queue will be large, of order 𝑂(𝜀1), while that in the second queue will be somewhat smaller, of order 𝑂(𝜀1/2). We obtain a two-dimensional diffusion approximation for this model and explicitly solve for the joint steady state probability distribution of the numbers of customers in the two queues. This work complements that in (Morrison, 2010), which the second queue was assumed to be heavily or lightly loaded, leading to mean queue lengths that were 𝑂(𝜀1) or 𝑂(1), respectively.

1. Introduction

The study of two coupled queues is a fundamental problem in queueing theory and applied probability. Classic examples include the shortest queue problem [13], the longer queue problem [4], fork-join models [57] and two coupled queues with generalized processor sharing [810], which is the subject of the present investigation. Computing the joint probability distribution for these models typically leads to functional equations that may sometimes be recast as boundary value problems [11], such as Dirichlet and Riemann-Hilbert problems.

Generalized processor sharing (GPS) models have become quite popular in recent years, as they provide scheduling algorithms that yield both service differentiation among different customer classes and also gains from statistical multiplexing. Some recent investigations and applications of such models appear in [1215], where they are used, for example, for flow control in integrated service networks.

We consider here two parallel queues with respective Poisson arrival rates 𝜆 and 𝜀𝜎, where 0<𝜀1. Thus, the arrivals to the second queue are much less frequent than those to the first, and we immediately scale the second arrival rate by 𝜀, thus introducing 𝜎. The service times are assumed to be exponentially distributed in both queues, with respective means 1/𝜇 and 1. Thus, we are taking the unit of time as the mean service time in the second queue. There is a single processor sharing server that works at unit rate and devotes 1𝜀𝜅=1𝑂(𝜀) of its capacity to the first queue and the remaining 𝜀𝜅 to the second queue, provided both queues are nonempty. If one queue is empty, the processor devotes all of its capacity to the other queue. The total load is given by 𝜆/𝜇+𝜀𝜎, and we assume that the system is in heavy traffic so that this quantity will be close to 1. Hence, we define 𝜔 from the relation 𝜆/𝜇+𝜀𝜎=1𝜀𝜔 and assume that 𝜔>0, so that the system is stable. This also means that the first queue is heavily loaded. The second queue has traffic intensity 𝜀𝜎/(𝜀𝜅)=𝜎/𝜅, and we say it is underloaded if 𝜎/𝜅<1, overloaded if 𝜎/𝜅>1 and critically loaded if 𝜎/𝜅1 (more precisely 𝜎/𝜅=1+𝑂(𝜀)). The underloaded and overloaded cases were analyzed in [16]. We denote by 𝑁1 (resp., 𝑁2) the number of customers in the first (resp., second) queue and the joint steady state probability distribution will be denoted by 𝑝(𝑚,𝑛)=Prob[𝑁1=𝑚,𝑁2=𝑛].

For the underloaded case most of the mass occurs on the scale 𝑚=𝑂(𝜀1) and 𝑛=𝑂(1), so there will tend to be only a few customers in the second queue. Asymptotically, 𝑝(𝑚,𝑛) has a product form behavior, with an exponential distribution in 𝜀𝑁1 and a geometric distribution in 𝑁2 (see [16]). This analysis was recently extended to an arbitrary number of parallel queues by Morrison and Borst [17], as long as one queue is heavily loaded and all the others are underloaded (with similar assumptions about arrival rates and processor-sharing factors as above). For the overloaded case most of the probability mass occurs for both 𝑚,𝑛=𝑂(𝜀1), and in [16] a diffusion limit of the form 𝑝(𝑚,𝑛)𝜀2𝒫(𝜁,𝜏)=𝜀2𝒫(𝜀𝑚,𝜀𝑛) is obtained. Here 𝒫 may be characterized as the solution to a parabolic PDE, in the variables 𝜁 and 𝜏. Here, we will analyze the critically loaded case, which also leads to a diffusion limit with now 𝑝(𝑚,𝑛)𝜀3/2𝜙0(𝜁,𝑤)=𝜀3/2𝜙0(𝜀𝑚,𝜀𝑛), where 𝜙0 will satisfy an elliptic PDE. Thus the critically loaded case has 𝑁2=𝑂(𝜀1/2) and leads to a somewhat more difficult problem than either the underloaded or overloaded cases.

We will obtain the PDE for 𝜙0(𝜁,𝑤) as a limiting case of the difference equation(s) satisfied by 𝑝(𝑚,𝑛), and explicitly solve the PDE by transform methods. We will obtain detailed results for the marginal distributions Prob[𝑁1=𝑚]=𝑛=0𝑝(𝑚,𝑛) and Prob[𝑁2=𝑛]=𝑚=0𝑝(𝑚,𝑛), as well as the mean queue lengths. We will also obtain other approximations to 𝑝(𝑚,𝑛) that are valid on scales where 𝑚=𝑜(𝜀1) and/or 𝑛=𝑜(𝜀1/2). In particular, we shall show that 𝑝(𝑚,𝑛) is 𝑂(𝜀) on the scale 𝑚=𝑂(𝜀1/2) and 𝑛=𝑂(1).

Previous work on this model includes Fayolle and Iasnogorodski [8] (see also [10]) and the more recent study of Guillemin and Pinchon [9]. There the authors consider the double generating function 𝐹(𝑥,𝑦)=𝑚,𝑛𝑥𝑚𝑦𝑛𝑝(𝑚,𝑛) and obtain a functional equation for the boundary values 𝐹(𝑥,0) and 𝐹(0,𝑦). This is ultimately converted to a Dirichlet problem, which is solved to yield the boundary values of 𝐹 in terms of elliptic integrals. One of the authors (J. A. Morrison) has verified that by analyzing the results in [8] and [9], in the asymptotic limit we consider and with the scaling 𝑥=1𝑂(𝜀) and 𝑦=1𝑂(𝜀), and then inverting the double transform, we also obtain our main approximation 𝑝(𝑚,𝑛)𝜀3/2𝜙0(𝜁,𝑤). However, this involves a lengthy calculation that takes far more work than our direct approach, which consists of deriving a limiting PDE and solving it, along with appropriate boundary conditions. Also, this direct approach should work for other models of this type, including ones with finite capacities of customers, and with 3 or more coupled queues.

Other recent work on diffusion approximations for generalized processor sharing models includes Ramanan and Reiman [18] (see also references therein). This work, however, is more concerned with theoretical aspects of the diffusion approximations, such as convergence of the discrete problem to a certain diffusion process. Here, our focus is on obtaining the explicit solution to the limiting diffusion equation that arises from the balance equations. It is highly likely that this equation can be interpreted as the Kolmogorov forward equation for some appropriate diffusion process, but we do not consider such “process level” aspects here, as our approach is largely analytical.

Our approach has the merit that it can be used to compute correction terms to the diffusion approximations, and we do this in some cases here. Also, we treat scales other than the basic diffusion scale, where, for example, some of the variables remain discrete, which we then relate to the diffusion scale by asymptotic matching. This type of analysis is needed, for example, to accurately compute boundary probabilities.

From a mathematical viewpoint, the diffusion approximation we obtain (i.e., 𝜙0(𝜁,𝑤)) is somewhat nonstandard in that the density vanishes as 𝜁0 and the approximation breaks down for 𝜁=𝑂(𝜀). Also, the corner behavior of the problem is much different than what is typical. In [19] we analyzed a more general version of this model in another heavy traffic limit, assuming that the arrival rates and processor-sharing factors were of comparable magnitude. That analysis led to an elliptic PDE that was more complicated than the one obtained here, but probably more representative of typical diffusion approximations to two coupled queues, such as those considered in [2022].

The present scaling limit leads to a separable, elliptic PDE in the variables 𝜀𝑚 and 𝜀𝑛. Since the boundary conditions are somewhat simpler than those for the diffusion model in [19], we are able to obtain a more explicit solution to this equation, using classical transform theory [23]. We then evaluate the solution in various limiting cases, to obtain even simpler results that yield more insight into model behavior.

Yet another analysis of the model considered here is done in [24], but there it was assumed that both the arrival and service rates of the secondary customers are small, while the server devotes comparable resources to each queue.

The paper is organized as follows. In Section 2, we state the problem more precisely and give the balance equations satisfied by 𝑝(𝑚,𝑛). In Section 3, we summarize all of our main results. The derivations are given in Section 4 for the scale 𝑚=𝑂(𝜀1),𝑛=𝑂(𝜀1/2) and in Section 5 for the other ranges of 𝑚,𝑛.

Throughout the paper we will use the notation 𝑓(𝑥)𝑔(𝑥) to mean lim𝑥𝑥0[𝑓(𝑥)/𝑔(𝑥)]=1, 𝑓(𝑥)=𝑜(𝑔(𝑥)) to mean lim𝑥𝑥0[𝑓(𝑥)/𝑔(𝑥)]=0, and 𝑓(𝑥)=𝑂(𝑔(𝑥)) to mean that |𝑓(𝑥)/𝑔(𝑥)| is bounded for 𝑥 sufficiently close to 𝑥0.

2. Formulation

We consider two parallel infinite capacity queues for different traffic classes. The jobs arrive as Poisson processes with rates 𝜆 for the primary class, and 𝜀𝜎 for the secondary class, where 0<𝜀1. Hence, the secondary jobs arrive much less frequently than the primary ones. Moreover, it is assumed that the primary and secondary jobs have exponentially distributed service requirements with mean service times 1/𝜇 and 1, respectively, where 𝜇=𝑂(1). The server works at unit rate, and if neither queue is empty devotes fractions 1𝜀𝜅 and 𝜀𝜅 of its effort to the primary and secondary queues, respectively, where 𝜅=𝑂(1). The corresponding service rates are (1𝜀𝜅)𝜇 and 𝜀𝜅. If one queue is empty the server works at unit rate on the other queue, so that this model is work conserving. It is assumed that the primary queue is heavily loaded, with𝜆𝜇+𝜀𝜎=1𝜀𝜔,0<𝜔=𝑂(1).(2.1) Moreover, we assume that the secondary queue is critically loaded, with𝜎=𝜅+𝛿𝜀,𝛿=𝑂(1),(2.2) where 𝛿 may have either sign. Thus, the asymptotic limit we consider has 𝜀0+ with 𝜔,𝜇,𝜅 and 𝛿 fixed, and then 𝜎 and 𝜆 vary with 𝜀 so that (2.1) and (2.2) hold.

Since 𝜔>0 the system is stable. Let 𝑝(𝑚,𝑛) denote the stationary probability that there are 𝑚 jobs in the primary queue and 𝑛 jobs in the secondary queue. Then, the balance equations satisfied by 𝑝(𝑚,𝑛) are[]+(𝜆+𝜀𝜎+(1𝜀𝜅)𝜇+𝜀𝜅𝑝(𝑚,𝑛)=𝜆𝑝(𝑚1,𝑛)+𝜀𝜎𝑝(𝑚,𝑛1)(1𝜀𝜅)𝜇𝑝(𝑚+1,𝑛)+𝜀𝜅𝑝(𝑚,𝑛+1),𝑚1,𝑛1,(2.3)(𝜆+𝜀𝜎+𝜇)𝑝(𝑚,0)=𝜆𝑝(𝑚1,0)+𝜇𝑝(𝑚+1,0)+𝜀𝜅𝑝(𝑚,1),𝑚1,(2.4)𝜆+𝜀𝜎+1)𝑝(0,𝑛)=𝜀𝜎𝑝(0,𝑛1)+(1𝜀𝜅)𝜇𝑝(1,𝑛)+𝑝(0,𝑛+1),𝑛1,(2.5)(𝜆+𝜀𝜎)𝑝(0,0)=𝜇𝑝(1,0)+𝑝(0,1).(2.6) The normalization condition is𝑚=0𝑛=0𝑝(𝑚,𝑛)=1.(2.7)

The mean number of jobs in the primary and secondary queues are𝐄𝑁1=𝑚=0𝑛=0𝐄𝑁𝑚𝑝(𝑚,𝑛),2=𝑚=0𝑛=0𝑛𝑝(𝑚,𝑛),(2.8) From Little’s Law, the corresponding mean waiting times are 𝐄(𝑁1)/𝜆 and 𝐄(𝑁2)/(𝜀𝜎). From [8], the conservation law of Kleinrock implies, using (2.1) and (2.2), that𝐄𝑁1𝜇𝑁+𝐄2=1𝜆𝜀𝜔𝜇2=1+𝜀𝜎𝜀𝜇𝜔1+𝜀(𝜇𝜅𝜔𝜅)+𝜀3/2.(𝜇1)𝛿(2.9)

3. Summary of Results

We consider 𝑚=𝜁/𝜀=𝑂(1/𝜀) and 𝑛=𝑤/𝜀=𝑂(1/𝜀), with𝑝𝜁𝜀,𝑤𝜀=𝜀3/2𝜙0(𝜁,𝑤)+𝜀𝜙1,(𝜁,𝑤)+𝑂(𝜀)0<𝜁=𝑂(1),0<𝑤=𝑂(1).(3.1) We then have

Proposition 3.1. 𝜙0(𝜁,𝑤)+𝐼(𝛿<0)𝛿𝜔𝜅exp𝜔𝜁+𝛿𝑤𝜅=2𝜔𝜋exp𝜔𝜁2+𝛿𝜔2𝜅0𝛽[]𝛽cos(𝑤𝛽)+(𝛿/2𝑘)sin(𝑤𝛽)𝛽2+𝛿2/4𝜅2exp𝜔2+𝛿2+𝜅𝜇4𝜅𝛽2𝜇𝜁2𝑑𝛽.(3.2)

This gives the limiting density of (𝑁1,𝑁2), which applies for fixed values of 𝜁 and 𝑤. We next evaluate this density in various limiting cases, such as 𝑤 and 𝛿±, to gain more insight into its structure, and to verify consistency with results in [16].

Corollary 3.2. If 𝛿=𝑂(1), 𝑤1 and 0<𝜁/𝑤1, then 𝜙0𝜔𝛿(𝜁,𝑤)2+𝜅𝜇𝜔23/4𝜁𝜋𝜇𝛿+𝛿2+𝜅𝜇𝜔2𝑤3/2×exp𝜔𝜁2𝑤exp2𝜅𝛿2+𝜅𝜇𝜔2.𝛿(3.3)

Corollary 3.3. If 𝛿1,𝑤/|𝛿|>0 and 0<𝜁=𝑂(1), then 𝜙0𝜔||𝛿||(𝜁,𝑤)2𝜅exp𝜔𝜁2𝑤exp2𝜅𝛿2+𝜅𝜇𝜔2+||𝛿||exp𝜇𝜔24||𝛿||𝑤×exp𝜔𝜁2||𝛿||𝜁Erfc𝜇𝜔𝑤2||𝛿||𝜇𝑤exp𝜔𝜁2||𝛿||𝜁Erfc𝜇𝜔𝑤+2||𝛿||,𝜇𝑤(3.4) where the complementary error function is given by [25] 2Erfc(𝑧)=𝜋𝑧𝑒𝑢2𝑑𝑢.(3.5)

Corollary 3.4. If 𝑤𝛿1 and 0<𝜁=𝑂(1), then 𝜙02||𝛿||(𝜁,𝑤)5/2𝜁𝜋𝜇𝜅𝜇𝜔𝑤3/2exp𝜔𝜁2𝑤exp2𝜅𝛿2+𝜅𝜇𝜔2+||𝛿||.(3.6)

Corollary 3.5. If 𝛿1,0<𝑤=𝑂(1) and 0<𝜁=𝑂(1), then 𝜙0𝜔||𝛿||(𝜁,𝑤)𝜅||𝛿||𝑤exp𝜔𝜁𝜅.(3.7)

Remark 3.6. This matches with [16, Result 7], with 𝑛=𝑤/𝜀,0<𝜀1 and |𝛿|𝜀1.

Corollary 3.7. If 𝛿1,0<𝑤/𝛿=𝑂(1) and 0<𝜁=𝑂(1), then 𝜙0𝜔(𝜁,𝑤)𝛿𝜁2𝜋𝜇𝑤3/2𝛿exp4𝜇𝑤𝜁+𝜇𝜔𝑤𝛿2.(3.8)

Remark 3.8. This matches with [16, Result 4], with 𝜏=𝑤𝜀,0<𝜀1 and 𝛿𝜀1.

In Propositions 3.9 and 3.16 below, we give the limiting marginals of 𝑁2 and 𝑁1, which apply for 𝑤 and 𝜁 fixed, respectively. Then, we simplify these marginal densities in various limiting cases.

Proposition 3.9. The scaled lowest order asymptotic approximation to the stationary distribution of the number of jobs in the secondary queue is 𝑚=0𝑝𝑤𝑚,𝜀𝜀0𝜙0(𝜁,𝑤)𝑑𝜁,(3.9) where 0𝜙0𝛿(𝜁,𝑤)𝑑𝜁+𝐼(𝛿<0)𝜅exp𝛿𝑤𝜅=4𝜔𝜋exp𝛿𝑤×2𝜅0𝛽[]𝛽cos(𝑤𝛽)+(𝛿/2𝜅)sin(𝑤𝛽)𝛽2+𝛿2/4𝜅2𝜔+𝜔2+𝛿2/(𝜅𝜇)+4𝜅𝛽2/𝜇𝑑𝛽.(3.10)

Corollary 3.10. If 𝛿=𝑂(1) and 𝑤1, then 0𝜙04𝛿(𝜁,𝑤)𝑑𝜁2+𝜅𝜇𝜔23/4𝜋𝜇𝜔𝛿+𝛿2+𝜅𝜇𝜔2𝑤3/2𝑤exp2𝜅𝛿2+𝜅𝜇𝜔2,𝛿(3.11) which is consistent with (3.3).

Corollary 3.11. If 𝛿1 and 0<𝑤/|𝛿|=𝑂(1), then 0𝜙0||𝛿||(𝜁,𝑤)𝑑𝜁𝜅𝑤exp2𝜅𝛿2+𝜅𝜇𝜔2+||𝛿||×1+𝜇𝜔2𝑤2||𝛿||exp𝜇𝜔2𝑤4||𝛿||𝜔Erfc2𝜇𝑤||𝛿||𝜔𝜇𝑤𝜋||𝛿||.(3.12)

Corollary 3.12. If 𝑤𝛿1, then 0𝜙08||𝛿||(𝜁,𝑤)𝑑𝜁5/2𝜋𝜇𝜅𝜇𝜔3𝑤3/2𝑤exp2𝜅𝛿2+𝜅𝜇𝜔2+||𝛿||,(3.13) which is consistent with (3.6).

Corollary 3.13. If 𝛿1 and 0<𝑤=𝑂(1), then 0𝜙0||𝛿||(𝜁,𝑤)𝑑𝜁𝜅||𝛿||𝑤exp𝜅,(3.14) which is consistent with (3.7).

Corollary 3.14. If 𝛿1 and 0<𝑤/𝛿=𝑂(1), then 0𝜙0(𝜁,𝑤)𝑑𝜁𝜔𝜇𝛿1𝜋𝑤exp𝜇𝜔2𝑤𝜔4𝛿2𝜇𝛿𝜔Erfc2𝜇𝑤𝛿,(3.15) which is consistent with (3.8).

Remark 3.15. This matches with [16, Result 3], with 𝜏=𝑤𝜀,0<𝜀1 and 𝛿𝜀1.

Proposition 3.16. The scaled asymptotic approximation to the stationary distribution of the number of jobs in the primary queue is 𝑛=0𝑝𝜁𝜀,𝑛=𝜀𝜔exp(𝜔𝜁)+𝜀3/2𝑄1(𝜀𝜁)+𝑂2,(3.16) where 𝑄1=(𝜁)𝐼(𝛿>0)𝛿(1𝜔𝜁)exp(𝜔𝜁)𝜔𝜇𝜋exp𝜔𝜁20𝛽2𝛽2+𝛿2/4𝜅22×𝜔+𝜔2+𝛿2+𝜅𝜇4𝜅𝛽2𝜇exp𝜔2+𝛿2+𝜅𝜇4𝜅𝛽2𝜇𝜁22𝜔exp𝜔𝜁2𝑑𝛽.(3.17)

Corollary 3.17. If |𝛿|1, then 𝑛=0𝑝𝜁𝜀,𝑛=𝜀exp(𝜔𝜁)𝜔+𝜀𝜀𝐼(𝛿>0)𝛿(1𝜔𝜁)+𝑂2.(3.18)

Remark 3.18. This matches with [16, Result 7] for 𝛿<0, and with [16, Result 2] for 𝛿>0 since, for 𝛿𝜀1, [](𝜔+𝜎𝜅)exp(𝜔+𝜎𝜅)𝜁=exp(𝜔𝜁)𝜔+𝜀𝛿(1𝜔𝜁)+𝑂(𝜀).(3.19)

We next give expansions for the mean queue lengths, for 𝜀0 (with 𝜔,𝜅,𝜇,𝛿 fixed).

Proposition 3.19. The lowest order asymptotic approximation to the stationary mean secondary queue length is 𝐄𝑁2𝑠𝜀,𝛿𝑠=𝜇𝜔211+2𝜋𝑐41𝑐4cos1(𝑐)1+𝑐4cos1𝑐2𝑐1𝑐212𝑐2,(3.20) where 𝛿1<𝑐=𝛿2+𝜅𝜇𝜔2<1,0<cos1(±𝑐)<𝜋.(3.21)

Remark 3.20. We have verified from (3.17) that 0𝜁𝑄1(𝜁)𝑑𝜁=𝜇𝑠 so that from (3.16), 𝐄𝑁1=1𝜀𝜔𝜇𝑠𝜀+𝑂(1),(3.22) which is consistent with (2.9).

Corollary 3.21. If |𝛿|𝜀1, then 𝛿𝑠𝜇𝜔2𝜅,if𝛿1;𝑠𝛿,if𝛿1.(3.23) These match with [16, Result 1 and Corollary 10], respectively.

We next give some asymptotic results for 𝜙0(𝜁,𝑤) that apply for 𝛿 fixed and 𝜁 and/or 𝑤. We also give the “corner” behavior as (𝜁,𝑤)(0,0).

Proposition 3.22. (i)   𝜁,𝑤 with 𝛿<0 and 0𝑤/𝜁<|𝛿|/(𝜇𝜔), 𝜙0𝜔||𝛿||(𝜁,𝑤)𝜅||𝛿||exp𝜔𝜁𝜅𝑤,(3.24)
(ii)𝜁,𝑤 with 0<𝑤/𝜁< for 𝛿>0, or |𝛿|/(𝜇𝜔)<𝑤/𝜁< for 𝛿0,𝜙0𝜔(𝜁,𝑤)𝐾(𝜁,𝑤)exp2𝛿𝜁+12𝜅𝑤2𝜔2+𝛿2𝜇𝜅𝜇𝜅𝑤2+𝜁2,(3.25)𝐾(𝜁,𝑤)=𝜇𝜋𝜅𝜔𝑏𝑠𝑏𝑠𝜔+𝛿/(2𝜅)2+𝛿2𝜅𝜇1/4𝜁𝜁2+𝜇𝜅𝑤23/4,𝑏𝑠=𝑏𝑠𝑤𝜁=𝜇𝜔2𝜅2+𝛿2𝜅𝜇1/2𝑤𝜁2+(𝜇/𝜅)𝑤2,(3.26)
(iii)𝜁,𝑤=𝑂(1),𝛿>0, 𝜙0𝜔(𝜁,𝑤)𝛿𝜇3/2𝜅𝜔2+𝛿2𝜅𝜇3/41𝜋𝜁3/2𝑤+2𝜅𝛿exp𝛿𝑤12𝜅exp2𝜔+𝜔2+𝛿2𝜁𝜅𝜇,(3.27)
(iv)𝑤,𝜁,𝜁=𝑂(𝑤), 𝜙0(𝜁,𝑤)2𝑏+𝜅𝜇𝑏+𝜔𝑏+𝜁+𝛿/(2𝜅)𝜋𝑤3/2𝛿exp2𝜅𝑏+𝜔𝑤2𝑏𝜁+𝜅𝜁2𝜇2𝑤,𝑏(3.28)+=𝜇2𝜅𝜔2+𝛿2,𝜅𝜇(3.29)
(v)𝜁,𝑤0 with 0<𝑤/𝜁<,𝜙0(𝜁,𝑤)2𝜔𝜋𝜇𝜅𝜁𝜁2+(𝜇/𝜅)𝑤2,(3.30)
(vi)𝜁0 with 0<𝑤<,𝜙0(𝜁,𝑤)𝜅𝜇𝜔𝜋𝜁𝑎(𝑤)exp𝛿𝑤𝛿2𝜅,(3.31)𝑎(𝑤)=12𝜅𝑤𝑒𝑖𝑤𝛽𝜔(𝜇/4𝜅)2+𝛿2/𝜇𝜅𝑖(𝛽𝛿/2𝜅)(𝛽+𝑖𝛿/2𝜅)2𝛽2+𝜔(𝜇/4𝜅)2+𝛿2/𝜅𝜇1/2+𝜇𝑑𝛽𝜔4𝜅2+𝛿21𝜇𝜅𝑤2𝑒𝑖𝑤𝛽𝛽2+𝜇𝜔4𝜅2+𝛿2𝜅𝜇3/2𝑑𝛽.(3.32)

Remark 3.23. The results show that for 𝛿<0 the density 𝜙0 is asymptotically of product form in a sector, and distinctly nonproduct form in the complimentary sector. For 𝛿>0 the product form behavior is absent. Item (v) shows that 𝜙0 has an integrable singularity near the corner, while item (vi) shows that 𝜙0=𝑂(𝜁) as 𝜁0. The second integral in (3.32) may be expressed in terms of a modified Bessel function, using the identity 𝑒𝑖𝑍𝑈𝑍2+13/2𝑑𝑍=2𝑈𝐾1(𝑈).(3.33)

The approximation 𝑝(𝑚,𝑛)𝜀3/2𝜙0(𝜁,𝑤) is only valid for 𝑚=𝑂(𝜀1) and 𝑛=𝑂(𝜀1/2). For other ranges of (𝑚,𝑛), other expansions must be constructed, and these we summarize below.

Proposition 3.24. (i)   𝑚=𝑆/𝜀=𝑂(𝜀1/2), 𝑛=𝑤/𝜀=𝑂(𝜀1/2), 𝑝(𝑚,𝑛)𝜀2𝑆𝒫+(𝑤)+𝒫,𝒫(𝑤)+𝜔(𝑤)=𝜋𝜅𝜇exp𝛿𝑤𝒫2𝜅𝑎(𝑤),(𝑤)=𝜇𝑤𝒫+(𝑢)𝑑𝑢,(3.34)

Remark 3.25. We comment that for 𝑛=𝑂(1) and 𝜁>0(𝑚=𝑂(𝜀1))𝑝(𝑚,𝑛)𝜀3/2𝜙0(𝜁,0) so that the diffusion approximation still applies. For 𝑚=𝑂(1) and 𝑤>0 item (i) still applies, and then 𝑝(𝑚,𝑛)𝜀2𝒫(𝑤), which is independent of 𝑚 and can be used to estimate the boundary probabilities 𝑝(0,𝑛), which are 𝑂(𝜀2) for 𝑛=𝑂(𝜀1/2). Note that in item (ii) and for 𝑛=0 in item (iii), 𝑝(𝑚,𝑛) is 𝑂(𝜀), which is larger than the order of magnitude (𝑂(𝜀3/2)) of the diffusion approximation on the (𝜁,𝑤) scale. But, the total mass in the range in (ii) is 𝑂(𝜀), while that in ranges (i) and (iii) is 𝑂(𝜀).

4. Analysis of the Main Diffusion Approximation

If we let 𝑚=𝜁/𝜀,𝑛=𝑤/𝜀 and use (3.1), (2.1), and (2.2), in (2.3) and (2.4), we obtain to lowest order𝜇𝜕2𝜙0𝜕𝜁2+𝜇𝜔𝜕𝜙0𝜕𝜕𝜁+𝜅2𝜙0𝜕𝑤2𝛿𝜕𝜙0𝜕𝑤=0,𝜁>0,𝑤>0,(4.1) and the boundary condition𝜅𝜕𝜙0𝜕𝑤(𝜁,0)=𝛿𝜙0(𝜁,0),𝜁>0.(4.2) We will discuss the second boundary condition, along 𝜁=0, after (4.16).

We let𝜙0(𝜁,𝑤)=exp𝜔𝜁2+𝛿𝑤Φ2𝜅0(𝜁,𝑤).(4.3) It follows that𝜇𝜕2Φ0𝜕𝜁2𝜕+𝜅2Φ0𝜕𝑤214𝜇𝜔2+𝛿2𝜅Φ0=0,𝜁>0,𝑤>0,(4.4)𝜕Φ0𝛿𝜕𝑤(𝜁,0)=Φ2𝜅0(𝜁,0).(4.5)

To solve (4.4) and (4.5) we apply a transform in the 𝑤 variable. Using the theory of distributions and Green’s functions for ordinary differential equations (see [23, p. 294, exercise 4.24]), we have the following transform pair:𝐺(𝐵)=0[]𝐵cos(𝐵𝑥)+𝐴sin(𝐵𝑥)𝐹(𝑥)𝑑𝑥,(4.6)𝐹(𝑥)=𝐼(𝐴<0)𝐶𝐴𝑒𝐴𝑥+2𝜋0𝐵cos(𝐵𝑥)+𝐴sin(𝐵𝑥)𝐴2+𝐵2𝐺(𝐵)𝑑𝐵.(4.7) Here, the constant 𝐶 appears only when 𝐴<0 and the term 𝐴𝑒𝐴𝑥 corresponds to a single discrete eigenvalue in the spectral theory. By multiplying (4.7) by 𝑒𝐴𝑥 and integrating from 𝑥=0 and 𝑥=, we find that 𝐶=20𝐹(𝑥)𝑒𝐴𝑥𝑑𝑥. Applying (4.6) with 𝐵=𝛽 and 𝑥=𝑤, we letΩ0(𝜁,𝛽)=0𝛿𝛽cos(𝑤𝛽)+Φ2𝜅sin(𝛽𝑤)0(𝜁,𝑤)𝑑𝑤,(4.8) then integration by parts and the boundary condition (4.5) leads to0𝛿𝛽cos(𝛽𝑤)+𝜕2𝜅sin(𝛽𝑤)2Φ0𝜕𝑤2(𝜁,𝑤)𝑑𝑤=𝛽2Ω0(𝜁,𝛽).(4.9) Hence, from (4.4),𝜕2Ω0𝜕𝜁214𝜔2+𝛿2+𝜅𝜇4𝜅𝛽2𝜇Ω0=0,(4.10) so thatΩ0(𝜁,𝛽)=Ω0(0,𝛽)exp𝜔2+𝛿2+𝜅𝜇4𝜅𝛽2𝜇𝜁2.(4.11) Applying the inversion formula (4.7)Φ0𝛿(𝜁,𝑢)+𝐼(𝛿<0)𝜅exp𝛿𝑢2𝜅0exp𝛿𝑤Φ2𝜅0=2(𝜁,𝑤)𝑑𝑤𝜋0[]𝛽cos(𝑢𝛽)+(𝛿/2𝜅)sin(𝑢𝛽)𝛽2+𝛿2/4𝜅2Ω0(𝜁,𝛽)𝑑𝛽.(4.12) Let𝑃0(𝜁)=0𝜙0(𝜁,𝑤)𝑑𝑤.(4.13) Then, from (4.1) and (4.2),𝑑2𝑃0𝑑𝜁2+𝜔𝑑𝑃0𝑑𝜁=0,𝜁>0.(4.14) But, from the normalization condition (2.7), (3.1) and the Euler-Maclaurin summation formula [26], 0𝑃0(𝜁)𝑑𝜁=1, so that𝑃0(𝜁)=𝜔𝑒𝜔𝜁.(4.15) Hence, from (4.3) and (4.13),0exp𝛿𝑤Φ2𝜅0(𝜁,𝑤)𝑑𝑤=𝜔𝑒𝜔𝜁/2.(4.16) Also, if we use (3.1), (2.1), and (2.2), in (2.5), we obtain the lowest order boundary condition 𝜕𝜙0/𝜕𝑤(0,𝑤)=0, for 𝑤>0. We conclude that 𝜙0(0,𝑤) and Φ0(0,𝑤) are proportional to a delta function at 𝑤=0+ and hence, from (4.8) and (4.16), that Ω0(0,𝛽)=𝜔𝛽. Proposition 3.1 follows from (4.3), (4.11), (4.12), and (4.15).

We may rewrite the integral in (3.2) as12𝛽𝑒𝑖𝑤𝛽(𝛽+𝑖𝛿/(2𝜅))exp𝜔2+𝛿2+𝜅𝜇4𝜅𝛽2𝜇𝜁2𝑑𝛽.(4.17) We then deform the contour of integration to one around a cut in the 𝛽-plane from (𝑖/2𝜅)𝛿2+𝜅𝜇𝜔2 to 𝑖, and let1𝛽=𝑖2𝜅𝛿2+𝜅𝜇𝜔2+𝑦.(4.18) For 𝛿<0 there is a contribution from the pole at 𝛽=𝑖𝛿/(2𝜅), and we obtain𝜙0(𝜁,𝑤)=2𝜔𝜋exp𝜔𝜁2+𝛿𝑤𝑤2𝜅exp2𝜅𝛿2+𝜅𝜇𝜔2×0𝛿2+𝜅𝜇𝜔2𝑒+2𝜅𝑦𝑤𝑦𝛿+𝛿2+𝜅𝜇𝜔2+2𝜅𝑦sin𝑦𝜅𝑦+𝛿2+𝜅𝜇𝜔2𝜁𝜇𝑑𝑦.(4.19) For 𝛿=𝑂(1),𝑤1 and 0<𝜁/𝑤1, the main contribution to the integral comes from 𝑦=𝑂(1/𝑤), and Corollary 3.2 follows.

If 𝛿1,𝑤>0 and 0<𝜁=𝑂(1), then the integral in (4.19) is approximated by||𝛿||2𝜅0𝑒𝑤𝑦𝑦+𝜇𝜔2||𝛿||𝜁/4sin||𝛿||𝑦𝜇||𝛿||𝑑𝑦=𝜅0||𝛿||𝑥𝑥exp𝜇𝑤/2𝑥2+𝜔2/4sin(𝜁𝑥)𝑑𝑥,(4.20) which leads [27] to Corollary 3.3. Corollaries 3.4 and 3.5 follow from the asymptotic approximation [25] Erfc(𝑧)exp(𝑧2)/(𝜋𝑧),𝑧1, and the limiting value Erfc()=2. If 𝛿1,0<𝑤/𝛿=𝑂(1) and 0<𝜁=𝑂(1), then from (4.19), we obtain the approximation𝜙0(𝜔𝜁,𝑤)𝜋exp𝜔𝜁2𝜇𝜔2𝑤4𝛿0𝑒𝑤𝑦sin𝛿𝑦𝜇𝜁𝑑𝑦,(4.21) and hence Corollary 3.7.

Proposition 3.9 follows from Proposition 3.1, (3.1) and the Euler-Maclaurin summation formula [26]. Next, from (4.19)0𝜙0(𝜁,𝑤)𝑑𝜁=8𝜔𝜋𝜇exp𝛿𝑤𝑤2𝑘exp2𝜅𝛿2+𝜅𝜇𝜔2×0𝛿2+𝜅𝜇𝜔2+2𝜅𝑦𝑦𝜅𝑦+𝛿2+𝜅𝜇𝜔2𝑒𝑤𝑦𝛿+𝛿2+𝜅𝜇𝜔2+2𝜅𝑦𝜇𝜔2+4𝑦𝜅𝑦+𝛿2+𝜅𝜇𝜔2𝑑𝑦.(4.22) For 𝛿=𝑂(1) and 𝑤1, the main contribution to the integral comes from 𝑦=𝑂(1/𝑤), and Corollary 3.10 follows. If 𝛿1 and 0<𝑤/|𝛿|=𝑂(1), then0𝜙0𝜔(𝜁,𝑤)𝑑𝜁𝜇||𝛿||𝑤𝜋𝜅exp2𝜅𝛿2+𝜅𝜇𝜔2+||𝛿||0𝑦𝑒𝑤𝑦𝑦+𝜇𝜔2/4||𝛿||2𝑑𝑦,(4.23) which leads [27] to Corollary 3.11. Corollaries 3.12, 3.13, and 3.14 follow directly from Corollaries 3.4, 3.5, and 3.7, respectively.

We now consider the first order correction term 𝜙1(𝜁,𝑤) in (3.1). If we use (3.1), (2.1), and (2.2), in (2.3) and (2.4), we obtain𝜇𝜕2𝜙1𝜕𝜁2+𝜇𝜔𝜕𝜙1𝜕𝜕𝜁+𝜅2𝜙1𝜕𝑤2𝛿𝜕𝜙1𝜇𝜕𝑤+𝛿𝜕𝜙0+1𝜕𝜁2𝜕2𝜙0𝜕𝑤2=0,(4.24) for 𝜁>0 and 𝑤>0, and the boundary condition𝜅𝜕𝜙1𝜕𝑤(𝜁,0)𝛿𝜙1𝜕(𝜁,0)+𝜇2𝜙0𝜕𝜁2(𝜁,0)+𝜇(𝜔+𝜅)𝜕𝜙0𝜅𝜕𝜁(𝜁,0)+2𝜕2𝜙0𝜕𝑤2(𝜁,0)=0.(4.25) We let𝑃1(𝜁)=0𝜙1(𝜁,𝑤)𝑑𝑤.(4.26) It follows, from (4.15) and (4.24)–(4.26), that𝜇𝑑2𝑃1𝑑𝜁2+𝜇𝜔𝑑𝑃1𝜕𝑑𝜁+𝜇2𝜙0𝜕𝜁2(𝜁,0)+𝜇(𝜔+𝜅)𝜕𝜙0𝜅𝜕𝜁(𝜁,0)+2𝜕2𝜙0𝜕𝑤2𝛿(𝜁,0)2𝜕𝜙0𝜕𝑤(𝜁,0)=𝛿𝜇𝜔2𝑒𝜔𝜁.(4.27) Hence, from (4.1) at 𝑤=0+,𝑃1(𝜁)=𝑄11(𝜁)2𝜙0(𝜁,0),(4.28) where𝑑2𝑄1𝑑𝜁2+𝜔𝑑𝑄1𝑑𝜁=𝛿𝜔2𝑒𝜔𝜁𝜅𝜕𝜙0𝜕𝜁(𝜁,0).(4.29)

From Proposition 3.1,𝜙0(𝜁,0)+𝐼(𝛿<0)𝛿𝜔𝜅𝑒𝜔𝜁=2𝜔𝜋exp𝜔𝜁20𝛽2𝛽2+𝛿2/4𝜅2exp𝜔2+𝛿2+𝜅𝜇4𝜅𝛽2𝜇𝜁2𝑑𝛽.(4.30) It follows that𝑑2𝑄1𝑑𝜁2+𝜔𝑑𝑄1𝑑𝜁𝐼(𝛿>0)𝛿𝜔2𝑒𝜔𝜁=𝜔𝜅𝜋𝑒𝜔𝜁/20𝛽2𝜔+𝜔2+𝛿2/(𝜅𝜇)+4𝜅𝛽2/𝜇𝛽2+𝛿2/4𝜅2exp𝜔2+𝛿2+𝜅𝜇4𝜅𝛽2𝜇𝜁2𝑑𝛽.(4.31) From (3.1), (4.13), (4.26), (4.28), and the Euler-Maclaurin summation formula [26],𝑛=0𝑝𝜁𝜀,𝑛=𝜀𝑃0(𝜁)+𝜀3/2𝑄1(𝜀𝜁)+𝑂2.(4.32) Hence,0𝑄1(𝜁)𝑑𝜁=0.(4.33) Proposition 3.16 follows in an elementary manner from (4.31) and (4.33).

If we let 2𝜅𝑦=𝛿2+𝜅𝜇𝜔2(cosh𝑢1) in (4.22), we obtain0𝜙0(𝜁,𝑤)𝑑𝜁=2𝜔𝜋𝜇𝜅0sinh2𝑢cosh𝑢exp(𝑤/2𝜅)𝛿2+𝜅𝜇𝜔2(cosh𝑢𝑐)(cosh𝑢+𝑐)cosh2𝑢𝑐2𝑑𝑢,(4.34) where 𝑐 is given by (3.21). Hence,00𝑤𝜙0(𝜁,𝑤)𝑑𝜁𝑑𝑤=8𝛿1𝑐23/2𝜋𝜇𝜔2𝑐0sinh2𝑢cosh𝑢(cosh𝑢+𝑐)2(cosh𝑢𝑐)3𝑑𝑢.(4.35) The evaluation of the integral in (4.35) is routine, but tedious, and Proposition 3.19 follows from (3.1), since𝐄𝑁21𝜀00𝑤𝜙0(𝜁,𝑤)𝑑𝜁𝑑𝑤.(4.36)

We next establish the various asymptotic formulas in Proposition 3.22. We first note that the integrand in Proposition 3.1 is an even function of 𝛽, and if 𝛽 is viewed as complex it has a simple pole at 𝛽=𝑖𝛿/(2𝜅) and branch points at𝛽=±𝑖𝜇2𝜅𝜔2+𝛿2𝜅𝜇𝑖𝑏±.(4.37) Then, we can represent, for any 𝛿, 𝜙0 as the contour integral𝜙0(𝜔𝜁,𝑤)=exp2𝛿𝜁+𝑤12𝜅𝜋𝒞𝛽𝜔1𝛽+𝑖𝛿/2𝜅×exp𝑖𝑤𝛽2𝜁𝜔2+𝛿2+𝜅𝜇4𝜅𝜇𝛽2𝑑𝛽.(4.38) Here, 𝒞 is a horizontal contour in the 𝛽-plane, on which𝛿max2𝜅,0<Im(𝛽)<𝑏+.(4.39) The condition in (4.39) insures that if 𝛿<0 the pole at 𝛽=𝑖|𝛿|/(2𝜅) lies below the contour 𝒞.

If 𝛿>0 we can shift the contour to the real 𝛽-axis, and then (4.38) becomes the same as (3.2). If 𝛿<0 the pole must be taken into account in making this shift, and the residue from this pole yields the exponential terms in the left hand side of (3.2).

To evaluate (4.38) for 𝜁,𝑤 we employ the saddle point method. There is a saddle point where𝑑1𝑑𝛽𝑖𝑤𝛽2𝜁𝜔2+𝛿2+𝜅𝜇4𝜅𝜇𝛽2=0,(4.40) so that𝛽=𝑖𝑏𝑠=𝑖𝑏𝑠𝑤𝜁𝜇𝑖2𝜅𝜔2+𝛿2𝑤𝜅𝜇(𝜇/𝜅)𝑤2+𝜁2.(4.41) The saddle is on the imaginary axis and the directions of steepest descent are arg(𝛽𝑖𝑏𝑠)=0,𝜋. Then shifting the contour 𝒞 into another horizontal contour through 𝑖𝑏𝑠 leads to (3.24). Such a shift is always permissible if 𝛿>0, but if 𝛿<0 we need the saddle to lie above the pole, that is, 𝑏𝑠>|𝛿|/(2𝜅), and this occurs precisely when 𝑤/𝜁>|𝛿|/(𝜇𝜔). We thus obtain the condition in item (i) of Proposition 3.22. If 𝛿<0 and 𝑤/𝜁<|𝛿|/(𝜇𝜔) the pole dominates the saddle point contribution, and we obtain (3.23).

For 𝜁 and 𝑤=𝑂(1) a different analysis is needed, as now 𝑏𝑠0 so the saddle approaches the real axis, where the integrand in (4.38) has as simple zero. In this case (which applies only if 𝛿>0), we shift 𝒞 back to the real axis and expand the integrand for 𝛽0. Using𝛽𝜔𝑒𝛽+𝑖𝛿/(2𝜅)𝑖𝑤𝛽=2𝛽𝜅𝜔𝑖𝛿1+𝑖𝛽𝑤+2𝜅𝛿𝛽+𝑂2,𝛽0,(4.42) we thus obtain𝜙0𝜔(𝜁,𝑤)exp2+12𝜔2+𝛿2𝜁𝛿𝜅𝜇exp𝑤×2𝜅2𝜅𝜔𝜋𝛿𝛽𝑖+𝛽2𝑤+2𝜅𝛿𝛽+𝑂3exp𝜅𝜁𝛽𝜇Δ2𝑑𝛽,(4.43) where Δ=𝜔2+𝛿2/(𝜅𝜇). Evaluating the integral(s) in (4.43) leads to (3.25). If we consider the opposite limit, where 𝜁=𝑂(1) and 𝑤, then the saddle point approaches the upper branch point at 𝑖𝑏+. But by deforming 𝒞 to an integral about the branch cut we can show that the final result coincides with the expansion of (3.25) for 𝑤/𝜁1, which is given by (3.28).

Next, we consider the corner behavior of 𝜙0 as 𝜁,𝑤0. Now, the main contribution to the integral will come from where |𝛽| is large. From (3.2) for 𝛿>0, we then obtain 𝜙0𝜔(𝜁,𝑤)𝜋𝑒𝑖𝑤𝛽𝑒(𝜅/𝜇)𝜁|𝛽|=𝜔𝑑𝛽𝜋𝜇𝜅2𝜁𝜇1+𝜅𝑤2𝜁21,(4.44) which yields (3.30), and this can be shown to remain valid for 𝛿<0.

Finally, we fix 𝑤 and let 𝜁0. Simply, setting 𝜁=0 in (3.2) leads to a divergent integral. However, an integration by parts leads to, for 𝛿>0, 𝜙0𝜔(𝜁,𝑤)=𝜋𝛿exp𝑤2𝜅exp𝑖𝑤𝛽𝜅𝜇𝛽2+𝜇Δ4𝜅2𝜁𝑔0(𝛽)𝑑𝛽,(4.45) where Δ=𝜔2+𝛿2/𝜅𝜇 and𝑔0𝛽(𝛽)=𝛽𝛽+𝑖𝛿/2𝜅𝑖𝑤𝛽2+(𝜇/4𝜅)Δ2𝜅𝜇𝜁1.(4.46) Expanding (4.46) for 𝜁0 and noting that, by contour integration (if 𝛿>0),𝑒𝑖𝑤𝛽𝑑𝛽𝑑𝛽𝛽+𝑖𝛿/(2𝜅)𝑑𝛽=0,(4.47) we write the integrand as𝑒𝑖𝑤𝛽1𝜅𝜇𝛽2+𝜇Δ4𝜅2𝜁𝜁+𝑂2×𝑖𝑤𝑑𝛽𝑑𝛽+𝜁𝛽+𝑖𝛿/(2𝜅)𝑤2𝑑𝛽𝑑𝛽𝛽𝛽+𝑖𝛿/(2𝜅)𝛽2+(𝜇/(4𝜅))Δ2𝜁+𝑂2.(4.48) By using (4.47), identifying the 𝑂(𝜁) terms in (4.48), and explicitly performing the differentiation with respect to 𝛽, we ultimately obtain (3.31) and (3.32). This completes the (sketched) derivation of Proposition 3.22.

5. Analysis of Boundary and Corner Regions

We analyze cases where 𝑚=𝑜(𝜀1) and/or 𝑛=𝑜(𝜀1/2). While these carry mass that is asymptotically small, they must be considered to insure that 𝑝(𝑚,𝑛) is properly normalized to higher orders in 𝜀, and to compute higher order approximations to the moments. Also, to determine 𝜙0(𝜁,𝑤) we used the boundary condition 𝜙0(0,𝑤)=𝜔𝛿(𝑤), and analysis of cases where 𝜁 and 𝑤 are small will allow us to examine this condition more carefully.

First we observe from (3.31) that 𝜙0(𝜁,𝑤) vanishes linearly as 𝜁0, which indicates a nonuniformity in the asymptotics. We first consider the scale 𝑚=𝑂(1) with 𝑤>0 and set𝑝(𝑚,𝑛)=𝜀𝜈1𝒫(𝑚,𝑤;𝜀)𝜀𝜈1𝒫(𝑚,𝑤).(5.1) Here, 𝜈1 is a constant that will be determined by asymptotic matching. From (2.3), in terms of the variables 𝑚 and 𝑤, we have[]𝒫𝜆+𝜇+𝜀(𝜎+𝜅𝜅𝜇)𝒫(𝑚,𝑤;𝜀)=𝜆𝒫(𝑚1,𝑤;𝜀)+𝜀𝜎𝑚,𝑤𝒫𝜀;𝜀+𝜇(1𝜀𝜅)𝒫(𝑚+1,𝑤;𝜀)+𝜀𝜅𝑚,𝑤+.𝜀;𝜀(5.2) Thus, the leading term must satisfy, since 𝜆=𝜇+𝑂(𝜀),2𝜇𝒫(𝑚,𝑤)=𝜇𝒫(𝑚1,𝑤)+𝜇𝒫(𝑚+1,𝑤), and we write 𝒫 as𝒫𝒫(𝑚,𝑤)=𝑚+(𝒫𝑤)+(𝑤),(5.3) which is a linear function of 𝑚. The expansion for 𝑚=𝑂(1) must satisfy the boundary condition (2.5), which is𝒫(𝜆+𝜀𝜎+1)𝒫(0,𝑤;𝜀)=𝜀𝜎0,𝑤𝒫𝜀;𝜀+𝜇(1𝜀𝜅)𝒫(1,𝑤;𝜀)+0,𝑤+𝜀;𝜀.(5.4) To leading order, (5.4) implies that 𝒫(0,𝑤)=𝒫(1,𝑤) so that (5.3) becomes 𝒫𝒫(𝑚,𝑤)=(𝑤). But, as 𝑚,𝜀𝜈1𝒫(𝑤) cannot match to 𝜀3/2𝜙0(𝜁,𝑤), since 𝜙0 vanishes as 𝜁=𝑚𝜀0. This indicates that we must analyze another scale, which has 𝑚1 and 𝑚1/𝜀 (with 𝑤 fixed).

We thus set 𝑚=𝑆/𝜀 and𝑝(𝑚,𝑛)=𝜀𝜈2𝒫(𝑆,𝑤;𝜀)𝜀𝜈2𝒫(𝑆,𝑤).(5.5) By rewriting (2.3) on the (𝑆,𝑤) scale, we obtain𝜇1𝜀(𝜔+𝜅)𝜀3/2𝛿𝒫(𝑆,𝑤;𝜀)𝒫𝑆𝒫𝜀,𝑤;𝜀=𝜀𝜎𝑆,𝑤𝒫𝜀;𝜀𝒫(𝑆,𝑤;𝜀)+𝜀𝜅𝑆+𝒫𝜀,𝑤;𝜀𝒫(𝑆,𝑤;𝜀)+𝜇(1𝜀𝜅)𝑆+.𝜀,𝑤;𝜀𝒫(𝑆,𝑤;𝜀)(5.6) The limiting form of (5.6) as 𝜀0 is 𝒫𝑆𝑆=0, so we write𝒫(𝑆,𝑤)=𝑆𝒫+(𝑤)+𝒫(𝑤).(5.7) On the (𝑆,𝑤) scale the boundary condition (2.5), using also 𝜆=𝜇[1𝜀(𝜔+𝜎)] and 𝜎=𝜅+𝜀𝛿, becomes[]𝜇+1+𝜀(𝜎𝜎𝜇𝜔𝜇)𝒫(0,𝑤;𝜀)=𝜀𝜎𝒫0,𝑤𝜀;𝜀+𝒫0,𝑤+𝜀;𝜀+𝜇(1𝜀𝜅)𝒫.𝜀,𝑤;𝜀(5.8) To leading order as 𝜀0, (5.8) implies that𝒫𝑤(0,𝑤)+𝜇𝒫𝑆(0,𝑤)=0.(5.9) Combining (5.7) and (5.9), we thus write the approximation on the (𝑆,𝑤) scale as𝑝(𝑚,𝑛)𝜀𝜈2𝑆𝒫+(𝑤)+𝜇𝑤𝒫+(𝑢)𝑑𝑢.(5.10) We then determine 𝜈2 and 𝒫+(𝑤) by asymptotically matching (5.10), as 𝑆, to 𝜀3/2𝜙0(𝜁,𝑤) as 𝜁0. Noting that 𝑆=𝜁/𝜀, we find, in view of (3.31), that𝜈2=2,𝒫+𝜔(𝑤)=𝜋𝜅𝜇exp𝛿𝑤2𝜅𝑎(𝑤).(5.11) This establishes (3.34) in Proposition 3.24. Expression (5.10) remains valid for 𝑚=𝑂(1) as then 𝑆=𝜀𝑚0 and we obtain 𝑝(𝑚,𝑛)𝜀2𝜇𝑤𝒫+(𝑢)𝑑𝑢, which is independent of 𝑚 and consistent with our previous analysis for 𝑚=𝑂(1).

We consider some limiting cases of the (𝑆,𝑤) scale result in (3.34). As 𝛿 from (3.28) we obtain 𝑏+𝛿/(2𝜅) and then (3.32) yields, after some computations,𝑎(𝑤)exp𝑤𝛿2𝜅2𝜋𝑤3/212𝛿2𝜅;𝛿+,𝑤=𝑂(1).(5.12) It follows from (3.34) that𝒫+𝜔(𝑤)2𝜋𝛿𝜅1𝑤3/2,𝒫𝜔(𝑤)𝜋1𝛿𝜇𝑤;𝛿+,𝑤=𝑂(1).(5.13) The above results are consistent with those of Morrison in [16] for the case 𝜎>𝜅, where it was shown that𝑚=0𝑝(𝑚,𝑛)𝜀𝜔𝜇𝜋(𝜎𝜅)𝜀𝑛𝑝(𝑚,𝑛)𝜀(𝜎𝜅);𝜎>𝜅,𝑛1.(5.14) Hence, for 𝛿=(𝜎𝜅)/𝜀 our approximation 𝑝(𝑚,𝑛)𝜀2𝒫(𝑤) agrees with (5.14).

We will next consider the scale 𝑚,𝑛=𝑂(1) and matching this to the (𝑆,𝑤) (or (𝑚,𝑤)) scale(s) will require the behavior of (3.34) as 𝑤0. For 𝑤0 the second integral in the definition of 𝑎(𝑤) in (3.32) dominates, and we obtain2𝑎(𝑤)𝑤2,𝑤0+,(5.15) so that𝒫+𝜔(𝑤)𝜋𝜅𝜇2𝑤2,𝑤0+.(5.16) Also, 𝒫(𝑤) will be less singular (𝑂(𝑤1)) than 𝒫+(𝑤), and, since 𝑤=𝑛/𝜀, for 𝑤0+ the approximation 𝜀2[𝜀𝑚𝒫+(𝑤)+𝒫(𝑤)] becomes of order 𝑂(𝜀3/2) on the (𝑚,𝑛) scale.

On the scale 𝑚,𝑛=𝑂(1) we will show that 𝑝(𝑚,0) is asymptotically larger than 𝑝(𝑚,𝑛) for 𝑛1, which we just inferred to be 𝑂(𝜀3/2). We thus set𝜀𝑝(𝑚,𝑛)=3/2𝑄𝜀(𝑚,𝑛)+𝑂2,𝑛1,𝜀𝑄(𝑚)+𝜀3/2𝑄(1)𝜀(𝑚)+𝑂2,𝑛=0.(5.17) On the (𝑚,𝑛) scale, the limiting form of (2.3), with (5.17), becomes2𝜇𝑄(𝑚,𝑛)=𝜇𝑄(𝑚1,𝑛)+𝜇𝑄(𝑚+1,𝑛),(5.18) while the boundary condition (2.5) yields(𝜇+1)𝑄(0,𝑛)=𝜇𝑄(1,𝑛)+𝑄(0,𝑛+1),𝑛2.(5.19) Hence, 𝑄(𝑚,𝑛) is linear in 𝑚, and writing𝑄(𝑚,𝑛)=𝑄+(𝑛)𝑚+𝑄(𝑛),(5.20) we then obtain from (2.5), 𝑄(𝑛)𝑄(𝑛+1)=𝜇𝑄+(𝑛) which we sum to get𝑄(𝑛)=𝜇𝑗=𝑛𝑄+(𝑗).(5.21) Next, we consider the boundary condition (2.4) and the corner condition (2.6). Since 𝜆=𝜇+𝑂(𝜀), we find that2𝜇𝑄(𝑚)=𝜇𝑄(𝑚1)+𝜇𝑄(𝑚+1),(5.22) and 𝑄(1)(𝑚) also satisfies (5.22). Hence, we have𝑄(𝑚)=𝛾+𝛾𝑚,(5.23) and 𝑄(1)(𝑚)=𝛾(1)+𝛾(1)𝑚. The corner condition (2.6), with (5.17), yields to leading order (𝑂(𝜀))𝜇𝑄(0)=𝜇𝑄(1),(5.24) and at the next order (𝑂(𝜀3/2))𝜇𝑄(1)(0)=𝜇𝑄(1)(1)+𝑄(0,1).(5.25) It follows from (5.24) that 𝛾=0. Also, (2.5) with 𝑛=1 becomes (𝜇+1)𝑝(0,1)=𝜇𝑝(1,1)+𝑝(0,2)+𝑂(𝜀𝑝(0,0)) so that (5.19) holds to leading order even if 𝑛=1. To summarize, we have obtained𝜀𝑝(𝑚,𝑛)3/2𝑚𝑄+(𝑛)+𝜇𝑗=𝑛𝑄+𝜀(𝑗),𝑛1,𝜀𝛾+𝑂3/2,𝑛=0.(5.26) It remains to determine 𝑄+(𝑛) and the constant 𝛾. By letting 𝑚=𝑆/𝜀 and 𝑛 and matching to the (𝑆,𝑤) scale, we conclude that as 𝑛𝑄+𝜔(𝑛)𝜋𝜅𝜇2𝑛2,𝑄𝜔(𝑛)𝜋2𝜇𝜅𝑛.(5.27) While these relations are consistent with (5.21), they do not determine 𝑄+(𝑛) for 𝑛=𝑂(1).

We thus examine another scale, that has 𝑚=𝑆/𝜀=𝑂(𝜀1/2) and 𝑛=𝑂(1). By matching to (5.26) as 𝑚, we expect that 𝑝(𝑚,𝑛) will be 𝑂(𝜀) on this scale, so we set𝑝(𝑚,𝑛)𝜀𝑄(𝑆,𝑛).(5.28) Then, from (2.3), we obtain to leading order𝜇𝑄𝑆𝑆(𝑆,𝑛)+𝜅𝑄(𝑆,𝑛+1)+𝑄(𝑆,𝑛1)2𝑄(𝑆,𝑛)=0,𝑛1,(5.29) while the boundary condition (2.4) leads to𝜇𝑄𝑆𝑆(𝑆,0)+𝜅𝑄(𝑆,1)𝑄(𝑆,0)=0.(5.30)

Thus (5.29) corresponds to driftless diffusion in the 𝑆 variable and a driftless random walk in the discrete 𝑛 variable. We may replace (5.30) by the “artificial” boundary condition 𝑄(𝑆,1)=𝑄(𝑆,0), by extending (5.29) to hold also at 𝑛=0.

To analyze (5.29) consider the discrete Green’s function problem𝐺(𝑛+1)+𝐺(𝑛1)+(𝜈2)𝐺(𝑛)=𝛿𝑛,𝑛0,𝐺(1)=𝐺(0),(5.31) where 𝜈 is a “spectral” parameter. This is easily solved to yield𝐺=𝐺𝑛;𝑛0=1