This paper reviews the Fourier-series method for calculating cumulative distribution functions (cdf's) and probability mass functions (pmf's) by numerically inverting characteristic functions, Laplace transforms and generating functions. Some variants of the Fourier-series method are remarkably easy to use, requiring programs of less than fifty lines. The Fourier-series method can be interpreted as numerically integrating a standard inversion integral by means of the trapezoidal rule. The same formula is obtained by using the Fourier series of an associated periodic function constructed by aliasing; this explains the name of the method. This Fourier analysis applies to the inversion problem because the Fourier coefficients are just values of the transform. The mathematical centerpiece of the Fourier-series method is the Poisson summation formula, which identifies the discretization error associated with the trapezoidal rule and thus helps bound it. The greatest difficulty is approximately calculating the infinite series obtained from the inversion integral. Within this framework, lattice cdf's can be calculated from generating functions by finite sums without truncation. For other cdf's, an appropriate truncation of the infinite series can be determined from the transform based on estimates or bounds. For Laplace transforms, the numerical integration can be made to produce a nearly alternating series, so that the convergence can be accelerated by techniques such as Euler summation. Alternatively, the cdf can be perturbed slightly by convolution smoothing or windowing to produce a truncation error bound independent of the original cdf. Although error bounds can be determined, an effective approach is to use two different methods without elaborate error analysis. For this...

The Laguerre method for numerically inverting Laplace transforms is an old established method based on the 1935 Tricomi-Widder theorem, which shows (under suitable regularity conditions) that the desired function can be represented as a weighted sum of Laguerre functions, where the weights are coefficients of a generating function constructed from the Laplace transform using a bilinear transformation. We present a new variant of the Laguerre method based on: (1) using our previously developed variant of the Fourier-series method to calculate the coefficients of the Laguerre gener-ating function, (2) developing systematic methods for scaling, and (3) using Wynn’s ɛ-algorithm to accelerate convergence of the Laguerre series when the Laguerre coefficients do not converge to zero geometrically fast. These contributions significantly expand the class of transforms that can be effectively inverted by the Laguerre method. We provide insight into the slow convergence of the Laguerre coefficients as well as propose a remedy. Before acceleration, the rate of convergence can often be determined from the Laplace transform by applying Darboux’s theorem. Even when the Laguerre coefficients converge to zero geometrically fast, it can be difficult to calculate the desired functions for large arguments because of roundoff errors. We solve this problem by calculating very small Laguerre coefficients with low relative error through appropriate scaling. We also develop another acceleration technique for the case in which the Laguerre coefficients converge to zero ge-ometrically fast. We illustrate the effectiveness of our algorithm through numerical examples. Subject classifications: Mathematics, functions: Laplace transforms. Probability, distributions: calculation by transform inversion. Queues, algorithms: Laplace transform inversion.

In this paper we investigate operators that map one or more probability distributions on the positive real line into another via their Laplace-Stieltjes transforms. Our goal is to make it easier to construct new transforms by manipulating known transforms. We envision the results here assisting modelling in conjunction with numerical transform inversion software. We primarily focus on operators related to infinitely divisible distributions and Le vy ´ processes, drawing upon Feller (1971). We give many concrete examples of infinitely divisible distributions. We consider a cumulant-moment-transfer operator that allows us to relate the cumulants of one distribution to the moments of another. We consider a power-mixture operator corresponding to an independently stopped Lévy process. The special case of exponential power mixtures is a continuous analog of geometric random sums. We introduce a further special case which is remarkably tractable, exponential mixtures of inverse Gaussian distributions (EMIGs). EMIGs arise naturally as approximations for busy periods in queues. We show that the steady-state waiting time in an M/G/1 queue is the difference of two EMIGs when the service-time distribution is an EMIG. We consider several transforms related to first passage times, e.g., for the M/M/1 queue, reflected Brownian motion and Lévy processes. Some of the associated probability density functions involve Bessel functions and theta functions. We describe properties of the operators, including how they transform moments.

In this paper we describe the time-dependent moments of the workload process in the M/G/1 queue. The k th moment as a function of time can be characterized in terms of a differential equation involving lower moment functions and the time-dependent server-occupation probability. For general initial conditions, we show that the first two moment functions can be represented as the difference of two nondecreasing functions, one of which is the moment function starting at zero. The two nondecreasing components can be regarded as probability cumulative distribution functions (cdf's) after appropriate normalization. The normalized moment functions starting empty are called moment cdf's; the other normalized components are called moment-difference cdf's. We establish relations among these cdf's using stationary-excess relations. We apply these relations to calculate moments and derivatives at the origin of these cdf's. We also obtain results for the covariance function of the stationary workload process. It is interesting that these various time-dependent characteristics can be described directly in terms of the steady-state workload distribution. Subject classification: queues, transient results: M/G/1 workload process. queues, busy-period analysis: M/G/1 queue. In this paper, we derive some simple descriptions of the transient behavior of the classical M/G/1 queue. In particular, we focus on the workload process {W(t) : t 0} (also known as the unfinished work process and the virtual waiting time process), which is convenient to analyze because it is a Markov process. Our main results describe the time-dependent probability that the server is busy, P(W(t) > 0), the time-dependent moments of the workload process, E[W(t) k ], and the covariance function of the stationary ...

Choudhury and Lucantoni recently developed an algorithm for calculating moments of a probability distribution by numerically inverting its moment generating function. They also showed that high-order moments can be used to calculate asymptotic parameters of the complementary cumulative distribution function when an asymptotic form is assumed, such as F c (x) ~ ax b e -hx as x . Moment-based algorithms for computing asymptotic parameters are especially useful when the transforms are not available explicitly, as in models of busy periods or polling systems. Here we provide additional theoretical support for this moment-based algorithm for computing asymptotic parameters and new refined estimators for the case b 0. The new refined estimators converge much faster (as a function of moment order) than the previous estimators, which means that fewer moments are needed, thereby speeding up the algorithm. We also show how to compute all the parameters in a multi-term asymptote of the form F c (x) ~ k = 1 S m a k x b - k + 1 e -hx . We identify conditions under which the estimators converge to the asymptotic parameters and we determine rates of convergence, focusing especially on the case b 0. Even when b = 0, we show that it is necessary to assume the asymptotic form for the complementary distribution function; the asymptotic form is not implied by convergence of the moment-based estimators alone. In order to get good estimators of the asymptotic decay rate h and the asymptotic power b when b 0, a multiple-term asymptotic expansion is required. Such asymptotic expansions typically hold when b 0, corresponding to the dominant singularity of the transform being a multiple pole (b a positive integer) or an algebraic singularity (branch point, b non-integer)...

Motivated by extreme-value engineering in service systems, we develop and evaluate simple approximations for the distributions of maximum values of queueing processes over large time intervals. We provide approximations for several different processes, such as the waiting times of successive customers, the remaining workload at an arbitrary time, and the queue length at an arbitrary time, in a variety of models. All our approximations are based on extreme-value limit theorems. Our first approach is to approximate the queueing process by one-dimensional reflected Brownian motion (RBM). We then apply the extreme-value limit for RBM, which we derive here. Our second approach starts from exponential asymptotics for the tail of the steadystate distribution. We obtain an approximation by relating the given process to an associated sequence of i.i.d. random variables with the same asymptotic exponential tail. We use estimates of the asymptotic variance of the queueing process to determine an approximate number of variables in this associated i.i.d. sequence. Our third approach is to simplify GI/G/1 extremevalue limiting formulas in Iglehart (1972) by approximating the distribution of an idle period by the stationary-excess distribution of an interarrival time. We use simulation to evaluate the quality of these approximations for the maximum workload. From the simulations, we obtain a rough estimate of the time when the extreme value limit theorems begin to yield good approximations.

Motivated by interest in probability density functions (pdf’s) with nonexponential tails in queueing and related areas, we introduce and investigate two classes of beta mixtures of exponential pdf’s. These classes include distributions introduced by Boxma and Cohen (1997) and Gaver and Jacobs (1998) to study queues with long-tail service-time distributions. When the standard beta pdf is used as the mixing pdf, we obtain pdf’s with an exponentially damped power tail, i.e., f(t) ∼ αt −q e −ηt as t → ∞. This pdf decays exponentially, but analysis is complicated by the power term. When the beta pdf of the second kind is used as the mixing pdf, we obtain pdf’s with a power tail, i.e., f(t) ∼ αt −q as t → ∞. We obtain explicit representations for the cumulative distributions functions, Laplace transforms, moments and asymptotics by exploiting connections to the Tricomi function. Properties of the power-tail class can be deduced directly from properties of the other class, because the power-tail pdf’s are undamped versions of the other pdf’s. The power-tail class can also be represented as gamma mixtures of Pareto pdf’s. Both classes of pdf’s have simple explicit Laguerre-series expansions.

This paper is an extension of Abate and Whitt (1988b), in which we studied the M/M/1 busyperiod distribution and proposed approximations for busy-period distributions in more general single-server queues. Here we provide additional theoretical and empirical support for two approximations proposed in Abate and Whitt (1988b), the natural generalization of the asymptotic normal approximation in (4.3) there and the inverse Gaussian approximation in (6.6), (8.3) and (8.4) there. These approximations yield convenient closed-form expressions depending on only a few parameters, and they help reveal the general structure of the busy-period distribution. The busy-period distribution is known to be important for determining system behavior.

Abstract It is known that probability density functions and probability mass functions can usually be calculated quite easily by numerically inverting their transforms (Laplace transforms and generating functions, respectively) with the Fourier-series method, but other more general functions can be substantially more difficult to invert, because the aliasing and roundoff errors tend to be more difficult to control. In this paper we propose a simple new scaling procedure for non-probability functions that is based on transforming the given function into a probability density function or a probability mass function and transforming the point of inversion to the mean. This new scaling is even useful for probability functions, because it enables us to compute very small values at large arguments with controlled relative error. Subject classifications: Mathematics, functions: scaling for numerical transform inversion. Queues, algorithms: scaling for numerical transform inversion.

A mathematical model is developed to help analyze the benefit in contact-center perfor-mance obtained from increasing employee (agent) retention, by increasing agent job satisfac-tion. The contact-center “performance ” may be restricted to a traditional productivity measure such as the number of calls answered per hour or it may include a broader measure of the qual-ity of service, e.g., revenue earned per hour or the number of problems successfully resolved per hour. The analysis is based on an idealized model of a contact center, in which the number of employed agents is constant over time, assuming that a new agent is immediately hired to replace each departing agent. The agent employment periods are assumed to be independent and identically distributed random variables with a general agent-retention probability distri-bution, which depends upon management policy and actions. The steady-state staff-experience distribution is obtained from the agent-retention distribution by applying renewal theory. An increasing real-valued function specifies the average performance as a function of agent expe-rience. Convenient closed-form expressions for the overall performance as a function of model elements are derived when either the agent-retention distribution or the performance function has exponential structure. Management actions may cause the agent-retention distribution to change. The model describes the consequences of such changes upon the long-run average staff experience and the long-run average performance.