ON THE OPTIMALITY OF AN ALGORITHM

OF REINGOLD AND SUPOWIT

Peter J. Grabner and Helmut Prodinger

Department of Mathematics Technical University of Graz, Austria and Department of Algebra and Discrete Mathematics Technical University of Vienna, Austria

Abstract. In [5], Reingold and Supowit have analyzed a divide–and–conquer heuristic to obtain a matching with a small cost. In this paper we show that – among a large class of similar strategies – the version of Reingold and Supowit produces the minimal expected costs.

Let n be an even integer, and P a set of n points in the plane. A matching is a set of n/2 edges such that each point of P is an endpoint of exactly one edge. The sum of the lengths of the edges is called the cost of the matching. In [5], Reingold and Supowit have analyzed a divide–and–conquer heuristic to obtain a matching with a small cost. In two previous papers [6,7] they have performed a worst-case analysis of several minimal matching algorithms.

Here, we want to demonstrate, that – among a large class of similar strategies – the version of Reingold and Supowit produces the minimal expected costs.

For the sake of shortness, we assume a certain knowledge of [5].

Reingold and Supowit start with a rectangle with sides \sqrt{2} and 1. Then, it is bisected vertically to form two rectangles of sides 1 and 1/\sqrt{2}, respectively. If both rectangles contain an even number of points, the matching is constructed in both rectangles separately. If both rectangles contain an odd number of points, the recursive strategy leaves one point in each rectangle (the “stranded point”), and the two lonely points are connected.

Now we consider – more generally – a rectangle with sides a and b. The first vertical cut produces two rectangles with respective sides \frac{a}{2} and b. They are only similar to the original one in the Reingold/Supowit setting. However, the next cut will be horizontal, and the produced rectangles are similar to the original one.

Let us now consider the expected costs according to n random points. They are of course dependent on the size of the rectangle from which we start. If we multiply both sides by a factor λ, the expected cost will also multiply by the factor λ. And if

Typeset by AMS-TEX

This research was supported by the Austrian-Hungarian cooperation 10U3.

The first author is supported by the Schrödinger scholarship J00936-PHY.

we want to compare the expected costs according to different rectangles, we should normalize them by dividing by \sqrt{ab}.

Instead of one recursion as in the original case, we must consider F_{n}, the expected cost, based on a rectangle with sides a and b, and G_{n}, the expected cost, based on a rectangle with sides b and \frac{a}{2}. The fact that one recursion was sufficient to compute the expected cost might have been – apart from some ingenious intuition – the motivation of the original setting of Reingold and Supowit.

Otherwise, the argument to obtain the recursion is pretty much the same as in [5].

We have F_{0} = G_{0} = 0, and for n ≥ 1

(1)  \begin{aligned}F_{n} & = 2^{−n} ∑_{k=0}^{n} \binom{n}{k}(G_{k} + G_{n−k}) + \frac{1}{2}ϕ(a, b) · χ(n) \\ G_{n} & = 2^{−n} ∑_{k=0}^{n} \binom{n}{k}(\frac{F_{k}}{2} + \frac{F_{n−k}}{2} ) + \frac{1}{2}ϕ(b, \frac{a}{2}) · χ(n). \end{aligned}

Here, χ(n) is 1 for n even and 0 for n odd; the \frac{1}{2} is the probability of an “odd–odd split” and ϕ(a, b) is the expected distance between two randomly chosen points one in each of the two halves (with sides \frac{a}{2} and b) of the rectangle with sides a and b. By symmetry, we can write

(2)  \begin{aligned}F_{n} & = 2^{1−n} ∑_{k=0}^{n} \binom{n}{k}G_{k} + \frac{1}{2}ϕ(a, b) · χ(n) \\ G_{n} & = 2^{−n} ∑_{k=0}^{n} \binom{n}{k}F_{k} + \frac{1}{2}ϕ(b, \frac{a}{2}) · χ(n). \end{aligned}

Now, ϕ(a, b) = b · g(\frac{a}{2b}), where g(a) is the expected distance of two points in each of two adjacent rectangles of sides a and 1 each. This quantity is given by^{1}

(3)

g(a) = \frac{1}{a^{2}} ∫_{0}^{1} ∫_{0}^{a} ∫_{0}^{1} ∫_{0}^{a} \sqrt{(x_{1} + x_{2})^{2} + (y_{1} − y_{2})^{2}} , dx_{1} dy_{1} dx_{2} dy_{2} =

\frac{1}{a^{2}} ( − \frac{1}{30} + a^{5} + \frac{\sqrt{1 + a^{2}}}{15} − \frac{a^{2} \sqrt{1 + a^{2}}}{5} + \frac{a^{4} \sqrt{1 + a^{2}}}{15}

− \frac{\sqrt{1 + 4 a^{2}}}{30} + \frac{2 a^{2} \sqrt{1 + 4 a^{2}}}{5} − \frac{8 a^{4} \sqrt{1 + 4 a^{2}}}{15} − \frac{a^{4} log(a)}{12}

+ \frac{5 a^{4} log(2 a)}{4} + \frac{a^{4} log(−1 + \sqrt{1 + a^{2}})}{8} − \frac{a^{4} log(1 + \sqrt{1 + a^{2}})}{24}

− \frac{a log(a + \sqrt{1 + a^{2}})}{6} − \frac{31 a^{4} log(−1 + \sqrt{1 + 4 a^{2}})}{24}

+ \frac{a^{4} log(1 + \sqrt{1 + 4 a^{2}})}{24} + \frac{a log(2 a + \sqrt{1 + 4 a^{2}})}{6}).

1We used Mathematica to compute the four-fold integral

We choose an alternative approach to solve the system of recursions (2). Reingold and Supowit have used the Mellin transform, as demonstrated in Knuth’s book [4,p.131ff], while we use a slightly faster approach, which is called Rice’s method (see [2],[3]).

Let us use two other abbreviations,

(4)  A = \frac{1}{2}ϕ(a, b) = \frac{b}{2} g(\frac{a}{2b})  and  B = \frac{1}{2}ϕ(b, \frac{a}{2}) = \frac{a}{4} g(\frac{b}{a}).

It is useful to deal with the exponential generating functions

F (z) = ∑_{n≥0} F_{n} \frac{z^{n}}{n!} and G(z) = ∑_{n≥0} G_{n} \frac{z^{n}}{n!} .

Then we get, on summing up the recursions for n ≥ 1,

(5)  \begin{aligned}F (z) & = 2e^{z/2}G(\frac{z}{2}) + A\frac{e^{z} + e^{−z}}{2} − A \\ G(z) & = e^{z/2}F (\frac{z}{2}) + B \frac{e^{z} + e^{−z}}{2} − B \end{aligned}

The next step is to introduce

F̃(z) = e^{−z}F (z) = ∑_{n≥0} F̃_{n} \frac{z^{n}}{n!} and G̃(z) = e^{−z}G(z) = ∑_{n≥0} G̃_{n} \frac{z^{n}}{n!} .

Then we obtain

(6)  \begin{aligned}F̃(z) & = 2G̃(\frac{z}{2}) + A\frac{1 + e^{−2z}}{2} − Ae^{−z} \\ G̃(z) & = F̃(\frac{z}{2}) + B \frac{1 + e^{−2z}}{2} − Be^{−z} \end{aligned}

or

F̃(z) = 2F̃(\frac{z}{4}) + B(1 + e^{−z}) − 2Be^{−z/2} + A\frac{1 + e^{−2z}}{2} − Ae^{−z}.

From this we infer that F̃_{0} = 0 and for n ≥ 1

(7)  F̃_{n} = \frac{(−1)^{n}}{1 − 2^{1−2n}} (A(2^{n−1} − 1) + B(1 − 2^{1−n})).

Since F (z) = e^{z}F̃(z), we get for n ≥ 1

(8)  F_{n} = ∑_{k=1}^{n} \binom{n}{k}(−1)^{k} \frac{1}{1 − 2^{1−2k}} (A(2^{k−1} − 1) + B(1 − 2^{1−k})).

As already indicated, we use Rice’s method [3] for the asymptotic evaluation of this (alternating) sum. Set

(9)  ψ(z) = \frac{1}{1 − 2^{1−2z}} (A(2^{z−1} − 1) + B(1 − 2^{1−z})).

The basis of the approach is the integral representation

(10)  ∑_{k=1}^{n} \binom{n}{k}(−1)^{k}ψ(k) = −\frac{1}{2πi} ∫_{C} \frac{Γ(n + 1)Γ(−s)}{Γ(n + 1 − s)} ψ(s) ds,

where C is a positively oriented closed curve that lies in the domain of analyticity of ψ(s), encircles 1, 2, …, n and no other poles of the integrand.

The asymptotic behavior is obtained by enlarging the contour of integration and taking the additional negative residues into account; they are the terms in the asymptotic expansion. In our instance, there are poles at z = \frac{1}{2} + \frac{1}{2}χ_{k}, for k ∈ ℤ, with χ_{k} = \frac{2kπi}{log 2}, since for these values we have 1 − 2^{1−2z} = 0. All these poles are actually simple.

The computation of the negative residues at z = \frac{1}{2} + \frac{1}{2}χ_{k} is standard^{2} and yields

(11)  \frac{Γ(n + 1)Γ(−\frac{1}{2} − \frac{1}{2}χ_{k})}{Γ(n + \frac{1}{2} − \frac{1}{2}χ_{k})} \frac{1}{2 log 2}(A(\frac{(−1)^{k}}{\sqrt{2}} − 1) + B(1 − \sqrt{2}(−1)^{k})).

Note that by [1]

\frac{Γ(n + 1)}{Γ(n + \frac{1}{2} − \frac{1}{2}χ_{k})} = n^{\frac{1}{2} + \frac{1}{2} χ_{k}} (1 + O(\frac{k^{2}}{n} )) = \sqrt{n} · e^{2kπi·log_{4} n} (1 + O(\frac{k^{2}}{n} ))

and |Γ(σ + it)| = \sqrt{2π}|t|^{σ− \frac{1}{2}} e^{− \frac{π}{2} |t|}. It is customary to consider the pole at z = \frac{1}{2} separately. The contribution is

(12)  \sqrt{n} · \frac{Γ(−\frac{1}{2})}{2 log 2}(A + \sqrt{2}B)\frac{1 − \sqrt{2}}{\sqrt{2}}.

Note also that Γ(−\frac{1}{2}) = −2\sqrt{π}, so that we get

\sqrt{n} · \frac{\sqrt{π}(\sqrt{2} − 1)}{\sqrt{2} log 2}(A + \sqrt{2}B).

We collect the residues (11) and (12) and get

Theorem 1. The average cost of a random matching obtained from n points in a rectangle of sides a and b is given by

\sqrt{n} · (\frac{\sqrt{π}(\sqrt{2} − 1)}{\sqrt{2} log 2}(A + \sqrt{2}B) + δ(log_{4} n) + O(\frac{1}{n})),

where δ(t) is a continuous periodic C^{∞}-function of period 1 whose Fourier-expansion is given by

δ(t) = ∑_{k∈ℤ\\{0\}} Γ (−\frac{1}{2} − \frac{1}{2}χ_{k}) \frac{1}{2 log 2}(A(\frac{(−1)^{k}}{\sqrt{2}} − 1) + B(1 − \sqrt{2}(−1)^{k}))e^{2kπit};

2We used Maple to compute the residues

the quantities A and B are given by (4).

For the “classical” case we have B = A/\sqrt{2}, and the factor in front of \sqrt{n} is

\frac{\sqrt{π}(\sqrt{2} − 1)\sqrt{2}A}{log 2},

in accordance with the result in [5], where the quantity D corresponds to 2A.

Now we want to determine which choice of a and b gives the smallest constant. Apart from factors which do not depend on the design parameters, we have to consider A + \sqrt{2}B. As already discussed, we must rescale this quantity by \sqrt{ab}. We thus want to minimize

\frac{A + \sqrt{2}B}{\sqrt{ab}} = \frac{1}{2}\sqrt{\frac{b}{a}} g (\frac{a}{2b}) + \frac{\sqrt{2}}{4}\sqrt{\frac{a}{b}} g (\frac{b}{a}) .

As it is to be expected, this quantity depends only on the ratio x = \frac{a}{b}, and we have to minimize

H(x) := \sqrt{\frac{1}{x}} g (\frac{x}{2}) + \sqrt{\frac{x}{2}} g (\frac{1}{x}) .

It is quite plausible that x = \sqrt{2} minimizes H(x), since then both terms are the same. Maple can differentiate H(x) and confirms that H′(\sqrt{2}) = 0. Since the expression for H′(x) is quite involved, we will use a different approach to confirm that the only minimum is obtained for x = \sqrt{2}.

First, we prove that g(x) is convex by proving that g′(x) is monotonically increasing. For this purpose we compute

g′(x) = \frac{2}{x^{3}} ∫_{0}^{1} ∫_{0}^{1} ∫_{0}^{x} ∫_{0}^{x}(\sqrt{(x_{1} + x)^{2} + (y_{1} − y_{2})^{2}}−
\sqrt{(x_{1} + x_{2})^{2} + (y_{1} − y_{2})^{2}}) dx_{1} dx_{2} dy_{1} dy_{2}.

Let now a < b. Then we get by substituting x_{1} = \frac{a}{b}ξ_{1} and x_{2} = \frac{a}{b}ξ_{2}

g′(a) = \frac{2}{b^{3}} ∫_{0}^{1} ∫_{0}^{1} ∫_{0}^{b} ∫_{0}^{b}(\sqrt{(ξ_{1} + b)^{2} + \frac{b^{2}}{a^{2}} (y_{1} − y_{2})^{2}}−

\sqrt{(ξ_{1} + ξ_{2})^{2} + \frac{b^{2}}{a^{2}} (y_{1} − y_{2})^{2}}) dξ_{1} dξ_{2} dy_{1} , dy_{2}.

Then we can write g′(b) − g′(a) as a four-fold integral with positive integrand and we have proved monotonicity of g′.

The next step in proving H(x) ≥ H(\sqrt{2}) is the simple inequality

g (\frac{x}{2}) g (\frac{1}{x}) ≥ g (\frac{1}{\sqrt{2}})^{2} ,

which is an immediate consequence of

g(x) ≥ g (\frac{1}{\sqrt{2}}) + g′ (\frac{1}{\sqrt{2}}) (x − \frac{1}{\sqrt{2}})

and the fact that \sqrt{2} g(\frac{1}{\sqrt{2}}) − g′(\frac{1}{\sqrt{2}}) > 0. Equality holds, if and only if x = \sqrt{2}. Now we write

H(x) − H(\sqrt{2}) = \sqrt{\frac{1}{x}} g (\frac{x}{2}) + \sqrt{\frac{x}{2}} g (\frac{1}{x}) − \frac{2}{\sqrt[4]{2}} g (\frac{1}{\sqrt{2}}) ≥

\sqrt{\frac{1}{x}} g(\frac{x}{2}) + \sqrt{\frac{x}{2}} g (\frac{1}{x}) − \frac{2}{\sqrt[4]{2}}\sqrt{g (\frac{x}{2}) g (\frac{1}{x})} =

(\sqrt[4]{\frac{1}{x}}\sqrt{g (\frac{x}{2})} − \sqrt[4]{\frac{x}{2}}\sqrt{g (\frac{1}{x})})^{2} ≥ 0,

which proves the desired result.

We note here that we could not find a simpler proof for this fact, which could be read off from the following plot of H(x), which was produced with Mathematica.

Acknowledgment. The second author wants to thank W. Szpankowski for interesting discussions.

References

1. M. Abramovitz and I.A. Stegun, Handbook of Mathematical Functions, Dover Publications, Inc., New York, 1965.

2. P. Flajolet and R. Sedgewick, Digital Search Trees Revisited, SIAM J. Comput. 15 (1986), 748–767.

3. P. Flajolet and R. Sedgewick, Mellin Transforms and Asymptotics: Finite Differences and Rice’s Integrals, Theoret. Comput. Sci. 144 (1995), 101–124.

4. D. E. Knuth, The Art of Computer Programming, Vol. 3., Addison-Wesley, 1973.

5. E. M. Reingold and K. J. Supowit, Probabilistic Analysis of Divide–and–Conquer Heuristics for Minimum Weighted Euclidean Matching, Networks 13 (1983), 49-66.

6. K. J. Supowit and E. M. Reingold, Divide and conquer heuristics for minimum weighted Euclidean matching, SIAM J. Comput. 12 (1983), 118–143.

7. K. J. Supowit, E. M. Reingold and D. A. Plaisted, The traveling salesman problem and minimum matching in the unit square, SIAM J. Comput. 12 (1983), 144–156.

P. Grabner Institut für Mathematik A Technische Universität Graz Steyrergasse 30 8010 Graz, Austria H. Prodinger Institut für Algebra und Diskrete Mathematik Technische Universität Wien Wiedner Hauptstraße 8–10 1040 Wien, Austria

E-mail address: grabner@weyl.math.tu–graz.ac.at, proding@rsmb.tuwien.ac.at