Minimization of acquisition and operational costs in horizontal material handling system design.

J.W. HERRMANN [1]

G. IOANNOU [2]

I. MINIS [3]

J.M. PROTH [4]

This paper considers the problem of minimizing the fixed cost of acquiring material handling transporters and the operational cost of material transfer in a manufacturing system. This decision

problem, which arises during manufacturing facility design, is modeled using an integer programming formulation. Since the problem is NP-complete, two efficient heuristics are developed to solve it. Computational complexity, worst-case performance analysis, and extensive computational tests are provided for both heuristics. The results indicate that the proposed methods are well suited for large-scale manufacturing applications.

1. Introduction

Material handling contributes significantly to overall manufacturing costs. Fixed costs are associated with the investment in material handling equipment during system construction, while variable costs arise from the material transfer between the resources of the manufacturing facility during system operation. These costs are conflicting in nature, since acquiring more equipment may reduce the material handling effort at the expense of increased investment. The goal of this study is to support the material handling system design process by determining the minimum number of transporters required to transfer parts/ batches in the manufacturing facility with minimal material handling effort. Most of the related research discussed below has addressed Automated Guided Vehicle (AGV) systems. However, it is equally applicable to all types of horizontal, carrier-based transportation systems such as rail carts, industrial trucks, forklifts etc. [1,2].

Egbelu [3] has developed a set of formulas to calculate the minimum number of carriers in a manufacturing system, based on the loaded traveling times as well as on some empirical estimates of the unloaded traveling times. This work can be employed in an initial economic justification of AGV systems. To address the same problem, Tanchoco et al. [4] employed a queuing theory-based computer model (CAN-Q), and Wysk et al. [5] used spread-sheet analysis. Their results compared favorably to a simulation-based method (AGVSim [6]). All these approaches provide initial estimates for the number of carriers, which may be further refined by simulation. As such, they do not consider aspects of the problem that may be important during system operation, such as the distribution of moves (loaded and unloaded) among the vehicles. For the problem of determining the AGV fleet size, Sinriech and Tanchoco [7] developed a multi-criteria optimization model that considers the trade-off between investment costs and system throughput. To support the design process they proposed the use of decision tables relating the investment cost, the number of AGVs and their utilization, as well as the trade-off ratio between the corresponding conflicting costs. In this approach, the solution chosen by the designer may or may not be near the optimal one.

The design of efficient horizontal unit-load material handling systems was studied by Maxwell and Muckstadt [8]. They considered the case in which the production rate of each manufacturing resource is constant and known. To determine the minimum number of AGVs, they solved a transportation problem that distributes the unloaded vehicle moves among pairs of resources in a way that minimizes the unloaded vehicle traveling time. Subsequently, they determined sequences of moves that originate and terminate at the same resource (routes) and assigned them to AGVs. Their routing algorithm created a delivery schedule with a near-constant inter-arrival time of material at each resource. The authors, however, did not provide analytical tools to determine vehicle routes, which are critical in assigning moves to vehicles. In large applications, evaluation of these routes may be a complex, and often intractable, problem. Furthermore, the authors did not consider the periodicity of vehicle operations that may lead to addit ional unloaded moves in order to drive the system to its initial state.

In a related research area, the merit of topologically simple flow path designs for better AGV control has been examined by several authors including, Bartholdi and Platzmann [9], Bozer and Srinivasan [10] and Sinriech and Tanchoco [11]. Bozer and Srinivasan [12] have presented a partitioning algorithm for the design of single vehicle loops, in an effort to distribute the workload evenly among the AGVs in the material handling system. Although these designs offer simplicity and allow analytical performance evaluation, no rigorous arguments have been presented to support their advantages over conventional networks. Finally, the review papers by Co and Tanchoco [13] and Ioannou and Minis [14] summarize relevant literature in the AGV system vehicle management area, and in the integration of the material handling system design and facility layout design problems, respectively.

Assuming a fixed shop layout with predetermined material flow paths, the problem of minimizing fixed acquisition and variable transportation costs addresses the following control-related issues: (i) assignment of available transportation equipment to service requests by jobs waiting in the queues of output stations; and (ii) transporter optimal routing from a resource output station to the destination input station. However, since this problem is relevant to the design stage of a manufacturing facility, during which no real-time information is available, we propose a static integer programming formulation. The resulting optimization problem is complex (NP-complete), and two efficient heuristics are developed to solve it. The first is a greedy algorithm similar to the nearest neighbor approach of Clarke and Wright [15] for scheduling vehicles from a central depot to a number of delivery points. The second heuristic is a composite algorithm that solves the assignment relaxation of the integer program, determin es closed sequences of resources with a common origin (routes), and assigns these routes to transporters by solving a two-stage bin-packing problem [16,17]. Both heuristics run in polynomial time and their worst case performance is bounded by ratios of problem parameters. Extensive computational tests against lower bounds also show that they provide satisfactory solutions for applications of practical size.

It is important to note that the optimization problem we consider in this paper is closely related to vehicle routing [18], and it resembles in particular the multiple vehicle routing problem with fixed costs related to vehicle activation, constraints on the amount of work of each vehicle, and lack of a central depot [19]. Laporte [20] and Desrosiers et al. [21] have reviewed formulations and solution approaches to the Vehicle Routing Problem (VRP). As both papers mention, in the traditional version of VRP, the transportation system has one depot, multiple locations, and multiple vehicles. Each vehicle follows a route that starts and ends at the depot. The objective is to find vehicle routes that visit every given location and minimize the total traveling distance. The Distance-constrained VRP (DVRP) adds a constraint to the distance of each route [19]. Our optimization problem is almost identical to the DVRP, but there is no depot: each vehicle can start its route at any location, though travel must end whe re it starts. Because the problem has no central depot, it is difficult to apply directly the formulations and solution approaches reviewed in Laporte [20], which exploit the fact that each route must visit the

central hub. Although an oxymoron, our problem's decentralization complicates both the formulation and the solution technique.

The remainder of this paper is organized as follows. Section 2 introduces our assumptions as well as relevant definitions and notation. Section 3 presents the integer programming formulation of the problem and the assignment model used to compute lower bounds. Section 4 describes the two heuristic solution approaches, and Section 5 presents results on the computational complexity and the worst case performance of these algorithms. Section 6 includes the numerical experiments, and Section 7 summarizes the conclusions of this work.

2. Assumptions, definitions and notation

The development of the mathematical model is based on the following assumptions:

(1) The placement of the manufacturing resources on the shop floor is given, together with the location of the resource pick-up (output) and drop-off (input) stations. (For a review of effective shop layout techniques see Heragu [22] and also Meller and Gau [23].)

(2) The material flow paths between resources are fixed. (For existing methods in flow path design see Herrmann et al. [24].)

(3) The inter-resource material flow rates in terms of loads (trips) per time period are constant and known. They are calculated from the production routings (sequences of operations) of the products to be manufactured, the independent demand of final assemblies over the design horizon and the bill-of-materials. In a periodic production system, the design horizon is equal to the period. In a non-periodic system, the design horizon is divided into equal time periods and the demand is equally distributed over these periods.

(4) Whenever a transporter visits a resource output station, there always exists material to be transferred to subsequent resources. This assumption is necessary since no real-time information is available at the system design stage.

(5) Horizontal material handling transporters are considered (e.g., AGVs, manual or automated rail carts, industrial trucks, and forklifts) with unit load capacity; no sharing of moves between different material flow types (i.e., different batches) is allowed. Note that although several types of material handling equipment with very different characteristics are employed moving material within cells or departments (such as conveyors, cranes, robots, carts, etc.), most inter-department or cell handling moves are performed using some type of unit load horizontal transporters [2].

Consider a set M of manufacturing resources and let [f.sub.rs] be the material flow between resources r, s [epsilon] M, i.e., the number of trips that have to be performed by active transporters between these resources. There exist two types of transporter operations between a pair of manufacturing resources:

(i) A loaded move i is the transporter operation from the output station of the manufacturing resource o(i) to the input station of the destination resource, d(i). The set of loaded moves is denoted by L, and its cardinality (\L\) will be referred to as n throughout the text. Each element of set L is associated with a unit entry of the material flow matrix.

(ii) An unloaded move is the transporter operation from the input station of a manufacturing resource to the output station of another (or the same) resource, during which no load is carried. Due to the unit load assumption, an unloaded move is always packed between two loaded moves. Assuming that there exists a path between each input-output resource pair, it is easy to see that after the completion of a loaded move i [epsilon] L, a transporter can perform an unloaded move to the origin of any other loaded move j [epsilon] L; thus, the set of possible move sequences is L x L, and its cardinality is [n.sup.2] since L x L includes all pairs of loaded moves, with the appropriate unloaded move packed between the loaded ones. An unloaded move can be of zero length (time) only if the input station d(i) and output station o(j) corresponding to the destination of loaded move i and the origin of the subsequent loaded move j, respectively, coincide. Figure 1 illustrates the move types.

A cost [c.sub.ij] is associated with each pair of loaded moves i,j [epsilon] L. Assuming that the costs [[tau].sub.i], [[[tau].sup.p].sub.i], [[[tau].sup.d].sub.i], and [[[tau].sup.j].sub.i] reflect the time needed to perform loaded move i (travel time only), to pick-up the load from o(i), to deliver the load to d(i), and to travel from d(i) to o(j), respectively, then [c.sub.ij] is defined as:

[c.sub.ij] = { [[tau].sub.i] + [[[tau].sup.p].sub.i] + [[[tau].sup.d].sub.i] + [[[tau].sup.j].sub.i] if i [not equal to] j, [infinity] if i = j. (1)

The set of material handling transporters available for material transfer is denoted by V. For each transporter k [epsilon] V the capital investment is denoted by [w.sub.k]. This cost is appropriately scaled to reflect the relative weights of the variable and fixed components of the objective function, and is related to the net present value of each transporter, based on all lifetime costs [25].

A route u is a sequence of moves performed by a transporter that originates and terminates at the same resource output station. This definition is adopted from Maxwell and Muckstadt [8], and will be employed by the second heuristic presented in Section 4. A set of routes that have the same origin is a route set, denoted by [gamma]. The time needed to perform all moves of a route is denoted by [C.sub.u] = [[sigma].sub.ij[epsilon][y.sub.u]] [c.sub.ij], where [y.sub.u] is the set of pairs of loaded moves performed in sequence of route u. In the remainder of the paper we will refer to a route either by its index u or by the corresponding set [y.sub.u]. Appendix A presents an example adopted from the literature [8] to illustrate the definitions of moves, routes, and route sets. Finally, T is the period within which all loaded moves must be performed; T is scaled appropriately to reflect the time costs [c.sub.ij] in (1). Note that in periodic systems T is equal to the design horizon, while in non-periodic systems the design horizon is divided into periods of length T, as discussed above.

3. Mathematical model

To formulate the problem of minimizing the fixed acquisition and the variable operational costs of the material handling system, we use the following additional notation: for all elements of the transporter set V we define a binary variable that indicates which transporters perform at least one loaded move in L, and thus should be acquired; i.e.,

[y.sub.k] = { 1 if transporter k [epsilon] V is employed for some L, 0 otherwise.

Let [[x.sup.k].sub.ij] be a binary variable associated with each pair of loaded moves:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

The design problem can now be expressed as follows:

Problem P

minimize Z = [[sigma].sub.k[epsilon]V][w.sub.k][y.sub.k] + [[sigma].sub.i[epsilon]L][[sigma].sub.j[epsilon]L][[sigma].sub.k[epsi lon]V] [c.sub.ij][[x.sup.k].sub.ij], (2)

subject to:

[[sigma].sub.j[epsilon]L][[sigma].sub.k[epsilon]V][[x.sup.k].sub.ij] = 1 [For all]i [epsilon] L, (3)

[[sigma].sub.i[epsilon]L][[sigma].sub.k[epsilon]V][[x.sup.k].sub.ij]= 1 [For all]j [epsilon] L, (4)

[[sigma].sub.j[epsilon]L][[x.sup.k].sub.ij] - [[sigma].sub.j[epsilon]L][[x.sup.k].sub.ji] = 0 [For all]i [epsilon] L, k [epsilon] V, (5)

[[sigma].sub.i[epsilon]F][[sigma].sub.J[epsilon]F][[x.sup.k].sub.ij] [less than or equal to] [[sigma].sub.i[epsilon]F][[sigma].sub.j[epsilon]L][[x.sup.k].sub.ij] - 1 [For all]F [Subset] L:

2 [less than or equal to] \F\ [less than] [[sigma].sub.i[epsilon]L][[sigma].sub.j[epsilon]L][[x.sup.k].sub.ij], [For all]k [epsilon] V, (6)

[[sigma].sub.i[epsilon]L][[sigma].sub.j[epsilon]L][c.sub.ij][[x.sup.k ].sub.ij] [less than or equal to] T [For all]k [epsilon] V, (7)

[[x.sup.k].sub.ij] [less than or equal to] [y.sub.k] [For all]i [epsilon] L, [For all]j [epsilon] L, [For all]k [epsilon] V, (8)

[[x.sup.k].sub.ij],[y.sub.k] [epsilon] {0, 1} [For all]i [epsilon] L, [For all] j [epsilon] L, [For all] k [epsilon] V. (9)

The objective function Z in (2) accounts for the capital investment to acquire transportation equipment and for the operational cost of material handling. It is important to note that because of assumption (4) of Section 2, the time constants [c.sub.ij] in expression (1) ignore possible waiting times due to congestion or on-line scheduling decisions. As a result, the second term of the objective function (2) represents a lower bound on the total operational costs.

Constraints (3) and (4) ensure that each loaded move in L is performed by exactly one transporter k [epsilon] V. Constraint (5) imposes continuity on the path consisting of loaded moves performed by each transporter, by ensuring material flow conservation. The exponential set of constraints (6) enforces subtour elimination, guaranteeing the existence of a single tour for each transporter k [epsilon] V. A subtour is a sequence of moves of the form ([i.sub.1], [j.sub.2], ..., [i.sub.1]). Constraint (7) limits the time that each transporter operates to the design horizon T. Constraint set (8) prohibits material transfer by non-activated vehicles ([y.sub.k] = 0). Finally, constraint (9) forces the variables [[x.sup.k].sub.ij] and [y.sub.k] to assume binary values.

In order to guarantee a feasible solution to problem P, the following property should be satisfied by each pair of loaded moves: [c.sub.ij] [less than] T/2,[For all]i [epsilon] L,j [epsilon] L. The case in which there exist some loaded moves with [c.sub.ij] [greater than] T/2 may lead to an empty feasible solution space, if these moves have to be performed by the same transporter. Based on this property, and the fact that at least two loaded moves are required to complete a tour, it is easily seen that:

[[sigma].sub.k[epsilon]V][y.sub.k] [less than or equal to] [n/2] [right arrow] \V\ [less than or equal to] [n/2]. (10)

3.1. Assignment lower bounds

The problem that results by considering only the operational cost in the objective function of P, by removing the capacity constraints (7), and by disregarding the transporter indices, is the well known assignment problem (or minimum weight matching problem [26]) which is presented below.

Problem A

minimize [Z.sub.a] = [[sigma].sub.i[epsilon]L][[sigma].sub.j[epsilon]L][c.sub.ij][x.sub.ij ], (11)

subject to:

[[sigma].sub.j[epsilon]L][x.sub.ij] = 1 [For all] i [epsilon] L, (12)

[[sigma].sub.i[epsilon]L][x.sub.ij] = 1 [For all] j [epsilon] L, (13)

[x.sub.ij] [epsilon] {0, 1} [For all] i [epsilon] L, [For all] j [epsilon] L. (14)

The optimal solution [[Z.sup.*].sub.a] of A provides, for every loaded move i [epsilon] L, the loaded move [phi](i) that should follow i, such that the total variable cost is minimized. [[Z.sup.*].sub.a] and the associated [phi](i)'s are derived in polynomial time by the Hungarian algorithm, an application of the primal-dual method [26]. It is clear that [[Z.sup.*].sub.a] bounds from below the variable component [[Z.sup.v].sub.opt] of the optimal solution of P, i.e., [[Z.sup.*].sub.a] [less than or equal to] [[Z.sup.v].sub.opt] = [[sigma].sub.i[epsilon]L] [[sigma].sub.j[epsilon]L] [[sigma].sub.k[epsilon]V] [c.sub.ij][[x.sup.[k.sup.*]].sub.ij], where ([[x.sup.[k.sub.*]].sub.ij], [[y.sup.*].sub.k]) is the optimal variable vector of P.

The above lower bound on the variable cost can be used to obtain a lower bound on the number of transporters given by:

[[R.sup.*].sub.a] = [[[Z.sup.*].sub.a]/T]. (15)

Note that if we assume identical transporters in terms of acquisition costs, i.e., [w.sub.k] = w, [For all] k [epsilon] V, then we obtain a lower bound on the optimal fixed cost in (2): [[R.sup.*].sub.a] x w [less than or equal to] [[Z.sup.f].sub.opt] = [[sigma].sub.k[epsilon]V] [w.sub.k] [[y.sup.*].sub.k] = w[[sigma].sub.k[epsilon]V][[y.sup.*].sub.k]. The lower bounds derived above will be employed for the evaluation of the heuristics during the numerical experiments.

4. Solution algorithms

Problem P is closely related to the well-known vehicle routing problem, and is NP-complete (for a proof of the membership of P to class NP, see Appendix B). Thus, optimal solutions cannot be computed for medium to large-sized problems, and consequently, we develop heuristic approaches that provide near-optimal solutions.

4.1. A greedy heuristic

Let us consider the set of loaded moves. An intuitive approach for minimizing the cost Z is to start matching moves that have the minimum cost coefficients [c.sub.ij], and attempt to allot as many moves to each available transporter as possible. Note that by matching two loaded moves i,j [epsilon] L, an unloaded move from d(i) to o(j) is fixed. Thus, we can proceed by allotting to each transporter moves with minimal cost in a greedy fashion, until capacity constraints are violated. This results in selecting minimum cost moves at the beginning of the procedure. However, as the algorithm proceeds non-favorable selections may be made, as is typical with nearest-neighbor type approaches. The minimization of the number of transporters is implicitly introduced by forcing each transporter to be loaded to near-capacity before another transporter is activated. This greedy algorithm is presented below. Note that [L.sub.g] is an auxiliary set used in the presentation of the algorithm, k is the transporter index, and [D. sub.k] is the associated variable cost.

Algorithm: GREEDY

Step 1. Set k = l, [L.sub.g] = L

Step 2. If [L.sub.g] = [Empty set], go to 8

Step 3. Choose a loaded move p [epsilon] [L.sub.g] at random

Step 4. Set q = p, [D.sub.k] = 0

[L.sub.g] = [L.sub.g] \ {p}

Step 5. Select j [epsilon] [L.sub.g] : [c.sub.pj] = [min.sub.j[epsilon][L.sub.g]] {[c.sub.pj]}

Step 6. If [D.sub.k] + [c.sub.pj] + [c.sub.jq] [less than or equal to] T then:

[D.sub.k] = [D.sub.k] + [c.sub.pj]

[L.sub.g] = [L.sub.g] \ {j}

If [L.sub.g] = [Empty set] then

[D.sub.k] = [D.sub.k] + [c.sub.jq]

Go to 8

P = j

Go to 5

Step 7. If [D.sub.k] + [c.sub.pj] [c.sub.jq] [greater than] T, then

[D.sub.k] = [D.sub.k] + [c.sub.pq]

k = k + l

Go to 2

Step 8. Output number of transporters, assignment of moves, and variable cost

The above algorithm is straightforward. Capacity constraints play an important role in its progress, since they impose the threshold for transporter activation and, thus, for fixed cost allocation. GREEDY attempts to assign pairs of loaded moves of minimum cost to transporters; at the same time it forces each vehicle to perform closed continuous loops. Note that regardless of the value of the fixed costs, the above heuristic will provide the same solution for a given set of variable costs. The routing portion of GREEDY extends the Clarke and Wright [15] algorithm for scheduling vehicles from a central depot to a number of delivery points, to routes that have no central depot.

4.2. An assignment/bin-packing (ABP) composite heuristic

An intuitive first step towards minimization of the variable portion of Z in (2) is to solve the assignment problem A. Starting from any loaded move i [epsilon] L and following the sequence [phi](i), [phi]([phi](i))),. . . , we will return to i, since [Inverted E]j1 [epsilon] L : [phi](i) = j1 and [Inverted E]j2 [epsilon] L : i = [phi](j2). The resulting closed sequences of loaded moves form subtours. If one transporter is allotted to each subtour, then the variable cost in (2) is minimal. However, this ad-hoc allotment of moves to transporters may not be feasible, since capacities may be violated, and is, in general, not economical, since an unnecessarily large number of transporters may be activated. Consequently, an algorithmic approach driven by the preservation of the minimal variable cost is required to translate the solution of A into a near-optimal solution of P. This is the basis of the heuristic ABP presented below.

ABP starts with the subtours of the optimal solution of A and allots moves to transporters in order to minimize the objective of P (see Appendix A for an example of a solution to the assignment problem A). This is accomplished by identifying routes and routes sets and performing a two-stage bin-packing [16,17]. Figure 2 (a and b) illustrates the definition of subtours, routes, and route sets. A subtour derived from the solution of the assignment problem is shown in Fig. 2a. The nodes i = 1,... ,4 represent loaded moves, while the arcs (1,2),...,(4,l) represent the unloaded moves connecting the loaded ones. If some moves of the subtour have a common origin, then the subtour is decomposed into routes which intersect at this origin. This transformation of moves to routes is unique, since it is based on the assignment matchings. In Fig. 2b, two routes are shown, connected at o(l) = o(3). The nodes of a route set represent resource input/output stations, while the arcs represent loaded/unloaded moves. These route s form one route set corresponding to station o(3) = o(l).

In the algorithm given below, [L.sub.b] is an auxiliary set initially equal to the set of loaded moves; pack(i), [For all]i [epsilon] L is the transporter that move i is assigned to; K is the total number of transporters; [D.sub.k] is the total scaled time it would take transporter k to perform the moves allotted to it; [t.sub.u] is a variable that corresponds to each activated transporter; M is the set of resources and [f.sub.rs] is the material flow between resources r, s [epsilon] M, i.e., the number of loaded moves that have to be performed by active transporters between these resources. Furthermore, FFD refers to the first-fit-decreasing algorithm for solving bin-packing problems; the FFD heuristic was employed due to its tight worst-case bounds (11/9 times the optimal [16]). Finally, the load of a transporter is the sum of the cost coefficients associated with the pairs of loaded moves performed by this transporter.

Algorithm: ABP

Step 1. Solve A to obtain [phi](i), [For all] i [epsilon] L

Step 2. Set [L.sub.b] = L

Step 3. Set u = 1, v = 1, [[gamma].sub.u] = [Empty set], [c.sub.u] = 0, [[gamma].sub.v] = [Empty set]

Step 4. Identify r [epsilon] M: [[sigma].sub.s[epsilon]M] [f.sub.rs] = [max.sub.m[epsilon]M]{[[sigma].sub.s[epsilon]M][f.sub.ms]}

Step 5. Select move i [epsilon] [L.sub.b] : o(i) = r

Step 6. do

[L.sub.b] = [L.sub.b] \ {i}

[f.sub.o](i)d(i) = [f.sub.o](i)d(i) - 1

[[gamma].sub.u] = [[gamma].sub.u] [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] {(i,[phi](i))}, [C.sub.u] = [C.sub.u] + [C.sub.i],[phi](i)

i = [phi](i)

until o(i) = r

Step 7. [[gamma].sub.v] = [[gamma].sub.v] [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] {[[gamma].sub.u]}

Step 8. If [[sigma].sub.s[epsilon]M] [f.sub.rs] [greater than] 0

u = u + 1, [[gamma].sub.u] = [Empty set], [C.sub.u] = 0

Go to step 5 (new route)

Step 9. If [L.sub.b] [not equal to] [Empty set]

V = v + 1, [[gamma].sub.v] = [Empty set]

u = u + 1, [[gamma].sub.u] = [Empty set], [C.sub.u] = 0

go to 4 (new route set)

Step 10. For [beta] = 1, ..., v, apply FED to [C.sub.u]'s of [[gamma].sub.u] [epsilon] [gamma] [beta] Get temporary transporters K, loads [D.sub.1], ..., [D.sub.K], and pack(i), [For all] i [epsilon] L

Step 11. Renumber transporters in decreasing order of [D.sub.k]'s, as l = {[t.sub.1], ...[t.sub.k]}

Step 12. Get first element of l

Step 13. Set p = q = 1 and l = l \ {[t.sub.p]}

Step 14. do

Let [t.sub.q] = next element of l, after current [t.sub.q]

If [D.sub.p] + [D.sub.q] [less than or equal to] T then:

[[delta].sub.min] = [min.sub.i[epsilon][t.sub.p]j[epsilon][t.sub.q]][[c.sub.i],[phi](j)] + [C.sub.j], [phi] (i)

- [c.sub.i,[phi](i)] - [c.sub.j[phi](j)]]

If [D.sub.p] + [D.sub.q] + [[delta].sub.min] [less than or equal to] T and [[delta].sub.min] [less than or equal to] [w.sub.k] then:

i) Set l = l \ {[t.sub.q]}, [D.sub.p] = [D.sub.p] + [D.sub.q] + [[delta].sub.min]

ii) [For all] i [epsilon] [t.sub.q], pack(i) = [t.sub.p]

If [t.sub.q] last element of l:

Set [t.sub.p] = first element of l

Set l = l \ {[t.sub.p]}

Set [t.sub.q] = first element of l

until l = [phi]

Step 15. Output transporter number and final assignment of moves to transporters

The solution of A [A in script form] generates the loaded move matchings [phi](i) and the optimal subtours of loaded moves in Step 1 of algorithm ABP. After the initializations in Steps 2 and 3, routes and route sets are identified in Steps 4 through 9. Specifically, Step 4 determines the resource r with the maximum outgoing flow in the current flow matrix. Step 5 selects a loaded move i that originates at resource r. Step 6 follows the sequence of moves i, [phi](i), [phi]([phi](i)), ..., until a new move originating at r is encountered. This sequence forms route [[gamma].sub.u]. Each time the loop of Step 6 is performed, the inter-resource flow matrix and the route cost [c.sub.u], are appropriately updated.

Step 7 adds the route [[gamma].sub.u] that has just been formed to the current route set [[gamma].sub.v], which contains the routes that originate from the output station of resource r. If there remain moves that originate at resource r, Step 8 initiates a new route from this resource. Otherwise, Step 9 initializes a new route set, provided that there exist moves yet unassigned to routes ([L.sub.b] [not equal to] [phi]). Steps 4 through 9 are repeated until all moves are assigned to routes and routes are grouped into route sets.

Routes are allotted to transporters in two stages. Step 10 solves a bin-packing problem for every route set [[gamma].sub.v], in which the routes [[gamma].sub.u] correspond to objects with weights equal to the route costs [c.sub.u], and the bin capacity equals the time horizon T. The First-Fit-Decreasing (FFD) [16] algorithm is employed to determine the minimal number of bins (transporters) required to complete the moves in [[gamma].sub.v]. After this procedure is applied to all route sets [[gamma].sub.v], a temporary assignment of routes (and moves) to K transporters is determined. Note that the optimality of the variable cost is maintained until the end of Step 10. Step 11 renumbers the temporary transporters in list i such that if [t.sub.[k.sub.1]], [t.sub.[k.sub.2]] [epsilon] l : [k.sub.1] [less than] [d.sub.[k.sub.2], then [D.sub.[k.sub.1]] [D.sub.[k.sub.2].

The final stage of the algorithm addresses the fixed component of the objective function. Given the temporary assignment of routes to transporters, Steps 12 through 14 reduce the number of transporters by combining the routes assigned to more than one transporter. The goal is to reassign routes from underutilized transporters to more utilized ones, thus merging several disjoint routes not necessarily from the same route set. To achieve that, a procedure similar to the FFD heuristic is used. Step 13 removes the top element from the list l, i.e., the transporter with maximum load [t.sup.2], and Step 14 examines whether this transporter can be merged with another one, without: (i) violating capacity constraints; or (ii) increasing the variable cost more than the transporter fixed cost. Starting from the transporter that follows [t.sub.p], in l the algorithm iteratively tests each member of l against the two conditions. When a transporter [t.sub.q] that satisfies these conditions is found, its routes are merged with those of [t.sub.p]. In this case, the number of transporters is decreased by one, and [t.sub.q] is removed from l. When all the remaining elements of the list are examined for possible merger with [t.sub.p], the algorithm returns to the top of l, removes its first element, and repeats the procedure of Step 14 until l is exhausted or no more mergers are possible.

To calculate the minimum increase of the variable cost when attempting to merge the routes of two transporters, Step 14 calculates the cost augmentations for all possible connections [[delta].sub.ij] between moves i and j that are allotted to different transporters. The connection that results in the minimum cost increase is implemented. Figure 3 illustrates the process of cost augmentation for moves 2 and 3 that belong to routes {(1, 2), (2, 1)} and {(3, 4), (4, 3)}. In this case there exist two route sets, each consisting of one route. When attempting to merge the two routes, unloaded moves (2,1) and (4,3) are replaced by unloaded moves (2,3) and (4,1). This results in a variable cost increase of [[delta].sub.23] = [C.sub.2,3] + [C.sub.4,1] - [C.sub.2,1] - [C.sub.4,3]. If this is the minimum cost augmentation (i.e., [[delta].sub.23] [less than or equal to] [[delta].sub.24], [[delta].sub.23] [less than or equal to] [[delta].sub.41], [[delta].sub.23] [less than or equal to] [[delta].sub.13]) then the route c onnection is implemented between d(2), o(3) and d(4), o(1). After the route merger, one transporter will perform all moves 1, 2, 3, and 4. Note that in Step 14 optimality with respect to the variable cost may be sacrificed. However, by selecting the minimum augmentation cost at each iteration, the algorithm tries to keep the variable cost as close to the lower bound [[Z.sup.*].sub.a] as possible.

It easy to see that ABP attempts to reduce the fixed cost in Steps 10 and 14, when routes from the same or different route sets are connected. However, it is clear that the main concern of algorithm ABP is the variable cost, the optimality of which is preserved until Step 14. Also, note that in order to apply the FFD in Step 10, we have assumed that [c.sub.u] [less than or equal to] T, ""u. This assumption is well justified in manufacturing applications, since the moving and loading/unloading times are negligible compared to the design horizon T, resulting in small route costs.

5. Evaluation of heuristics

In this section some important properties of the two heuristics are established. The proofs of all theorems are included in Appendix C.

5.1. Computational complexity

The greedy heuristic is particularly fast; its computational complexity is bounded by a low order polynomial. The computational complexity of the composite assignment/ bin-packing heuristic is dominated by the time to compute the solution to the assignment problem.

Theorem 1. The computational complexity of GREEDY is O([n.sup.2]).

Theorem 2. The computational complexity of ABP is O([n.sup.3]).

5.2. Worst case analysis

In this analysis we examine the variable and fixed costs separately. Let [c.sub.max] = [max.sub.i[epsilon]L,j[epsilon]L]{[c.sub.ij]}, and [c.sub.min] = [min.sub.i[epsilon]L,j[epsilon]L]{[c.sub.ij]}, i.e., the maximum and minimum elements of the cost matrix [c.sub.ij]. Also, let [[Z.sup.v].sub.g] and [[Z.sup.v].sub.b] be the variable costs of the solutions derived by applying GREEDY and ABP, respectively, and [R.sub.g] and [R.sub.b] the corresponding numbers of the activated transporters.

The following two theorems establish the worst case performance of the heuristics, with respect to the variable cost.

Theorem 3. For any instance of problem P, algorithm GREEDY provides a solution that satisfies the property

[[Z.sup.v].sub.g]/[[Z.sup.v].sub.opt][less than][c.sub.max]/[c.sub.min]

Thus, if the ratio of the maximum to the minimum entry of matrix [c.sub.ij] is small, this heuristic will perform well even in the worst case, with respect to the variable portion of the cost. Nevertheless, the ratio [c.sub.max]/[c.sub.min] may be arbitrarily bad, thus offering no concrete guarantee on the performance of GREEDY in general problem settings.

Theorem 4. For any instance of problem P, algorithm ABP provides a solution that satisfies the property

[[Z.sup.v].sub.b]/[[Z.sup.*].sub.a][less than or equal to][c.sub.max]/[c.sub.min]

Since [[Z.sup.v].sub.opt][greater or equal to] [[Z.sup.*].sub.a], follows directly from Theorem 4 that:

Corollary 1. For any instance of problem P, algorithm ABP provides a solution that satisfies the property

[[Z.sup.v].sub.b]/[[Z.sup.v].sub.opt][[less than or equal to] [c.sub.max]/[c.sub.min]

Let us now consider the fixed cost, which is proportional to the number of activated transporters if [w.sub.k] = w, "" k [epsilon]V.

Theorem 5. For any instance of problem P, algorithms GREEDY and ABP provide solutions that satisfy the property

[R.sub.h]/[R.sub.opt] [less than or equal to] 1/2 X T/[C.sub.min], where [R.sub.h] represents the number of transporters in the appropriate heuristic (h = g for GREEDY, h = b for ABP) and [R.sub.opt] the number of transporters in the optimal solution.

Thus, although both heuristics are guaranteed to provide solutions for which the fixed part of the cost in (2) is bounded with respect to the optimum, the quality of these solutions may be arbitrarily bad, if the ratio of the time period T over the minimum entry of matrix [[c.sub.ij]] becomes large. However, the results of the next section indicate that on the average the numbers of transporters provided by the heuristics are very close to the optimal.

6. Numerical results

Algorithms GREEDY and ABP were implemented using the C programming language on a Sun Sparc workstation. In order to evaluate the algorithms' effectiveness, we developed an experimental procedure to generate random instances of problem P [P in script form], and we solved each instance using the heuristics. This procedure is described in Section 6.1. Because GREEDY is stochastic, we applied it 20 times to each problem instance and kept the best solution found; on the other hand, since ABP is deterministic, we applied it only once to each problem generated. In addition, we computed the lower bound [[Z.sup.*].sub.a] for each problem instance solved. Section 6.2 compiles the results of the computational experiments.

6.1. Problem generation

The problem generation procedure randomly created 1000 instances, grouped into 10 sets of 100 instances each. Each instance had resources randomly distributed across a square. All instances in the same set used the same size square. Instances in different sets used different squares. In general, for instances in set e, the side of the square measured be, e = 1,..., 10, with [b.sub.1] = 150 units; for e = 2,..., 10, [b.sub.e] = [b.sub.e-1]/[sigma], where [sigma] = 1.5. Figure 4 illustrates the square generation procedure. The transporter velocity [upsilon] = 1 unit and the time horizon T = 500 units were pre-specified and are constant in each example solved. For each instance in a set, the problem generation procedure provided the manufacturing resources, the flow between the resources, and then calculated the cost coefficients for the problem. To generate the resources, the procedure randomly picked a value for \M\, the number of resources, which was normally distributed with a mean of 30. For each of the \M\ resources, the algorithm picked two locations within the specified square, and the resource's input station and output station were placed at these locations. The location probability density function was uniform across the square. To generate the flow, the procedure selected values, from a normal distribution, for the material flow intensity [f.sub.rs] for each ordered resource pair. This is the number of required loaded trips to move material from resource r to resource s. Then [f.sub.rs] loaded moves are added to the set L. Each of these loaded moves (i) started at o(i) = r and ended at d(i) = s. The mean value of the material flow density function for problem set 1 was [f.sub.1] 10. For the remaining problem sets (e = 2,..., l0), [f.sub.e] = [sigma][f.sun.e-1].

To determine the cost coefficients, the algorithm calculated, for each pair of resources, the rectilinear distance between resource r at ([x.sub.r], [y.sub.r]) and resource s at ([x.sub.s], [y.sub.s]): [[micro].sub.rs] = [[micro].sub.sr] = [x.sub.r] - [x.sub.s] + [y.sub.r] - [y.sub.s] Then, for any two loaded moves i [epsilon] L and j [epsilon] L, [c.sub.ij] = (1/[upsilon])([[micro].sub.o(i)d(i)] + [[micro].sub.d(i)o(j)]). The procedure assumed that loading and unloading times are negligible. Note that, by construction, the different example sets have almost invariant average loaded traveling time; i.e., the average value of [[sigma].sub.r[epsilon]M] [[sigma].sub.s[epsilon]M] [f.sub.rs][[micro].sub.rs] for each example set e remains almost constant.

The procedure then calculated a vehicle cost w to balance the fixed and variable cost bounds. The total time required to complete all loaded moves is (l/[upsilon]) [[sigma].sub.i[epsilon]L][[micro].sub.o(i)d(i)] = (1/[upsilon]) [[sigma].sub.r[epsilon]M] [[sigma].sub.s[epsilon]M] [f.sub.rs][[micro].sub.rs], and this is a lower bound on the total time (the variable cost). The minimum number of transporters required to perform the loaded moves of L is [(1/[upsilon]) [[sigma].sub.i[epsilon]L][[micro].sub.o(i)d(i)](1/T)] (since time horizon T limits the amount of work each vehicle can perform). Then, to equate the fixed and variable cost lower bounds, the value of w is given by w = (1/[upsilon]) [[sigma].sub.i[epsilon]L][[micro].sub.o(i)d(i)][(1/[upsilon]) [[sigma].sub.i[epsilon]L][[micro].sub.o(i)d(i)][1/T]].sup.-1]

6.2. Algorithm performance

Table 1 presents the results for the ten example sets. The table includes the average number of moves for each problem instance in the example set, the average value of the variable cost obtained from each heuristic over the 100 instances in the set ([[Z.sup.v].sub.g] for GREEDY and [[Z.sup.v].sub.b] for ABP), the average number of transporters obtained by the heuristics ([R.sub.g] for GREEDY and [R.sub.b] for ABP), and the average percentage of time the transporters are idle within the time period T:

l[c.sub.g,b] = [[sigma].sub.k[epsilon]V] T X yk - [Z.sub.g,b]/[[sigma].sub.k[epsilon]V T X yk X 100.

Excessive idle time indicates that the heuristic solution uses too many transporters. The subscripts g and b refer to the solutions derived by GREEDY and ABP, respectively.

Recall that the resource location range (the square) shrinks as the example set index increases, so the resources are, on average, closer. In addition, the number of moves increases, but they are shorter moves. Consider the trends that Table 1 captures: as the problem set index increases, the average variable cost remains nearly constant; thus, the unloaded traveling time remains nearly constant as well. Furthermore, the average number of transporters decreases, as the number of moves increases, and the inter-resource distances decrease. This is expected, since packing smaller moves into transporters is easier and more efficient for both heuristics; in addition, the return moves that close the loop for each transporter become smaller. Comparing the results for the two heuristics, it is clear that ABP outperforms GREEDY with respect to the variable cost. The opposite occurs in terms of activated transporters. From the last two columns of Table 1, it is evident that as the average length of transporter moves d ecreases, the idle time also decreases. This is due to better assignment of moves to transporters. Finally, the transporters obtained by GREEDY are better utilized than those derived by ABP.

Tables 2 and 3 summarize the heuristics' performance with respect to the lower bounds. Table 2 describes the GREEDY heuristic. The second column lists the average relative deviation between the heuristic variable cost and the lower bound [[Z.sup.*].sub.a]. The third column lists the average relative deviation between the number of transporters and the lower bound [R.sup.*].sub.a]. The fourth column lists the average relative deviation between the total cost [Z.sub.g] and the lower bound Z' = [[Z.sup.*].sub.a] + [[R.sup.*].sub.a]w. (Note that Z' is a lower bound for Z since [[Z.sup.*].sub.a] and [[R.sup.*].sub.a] are lower bounds for the components of Z). Table 3 describes the same performance measures for the ABP heuristic.

Tables 2 and 3 indicate a possible trend: as the problem size (the number of moves) increases while the total loaded travel time remains almost invariant, the heuristics seem to perform better, and the solution values approach the lower bounds. However, due to the construction of the example problems (larger number of moves in smaller areas, thus smaller costs [c.sub.ij] this trend may not capture the true performance of the heuristics in general cases (e.g., when a larger number of moves is associated with longer move lengths). From Tables 2 and 3, it is also evident that ABP finds solutions that are closer to the lower bound for the variable cost, but the GREEDY solutions are closer to the lower bound for the number of transporters.

The solutions obtained by the two heuristics are expected to be more applicable to practical cases than the solutions obtained by the algorithm of Maxwell and Muckstadt [8]. The latter algorithm derives a lower bound for P [P in script form] identical to Z' = [[Z.sup.*].sub.a] + [[R.sup.*].sub.a]w, without considering the assignment of moves to transporters. The second part of the Maxwell and Muckstadt method uses a computationally expensive, simulation-based technique to verify that there exists a feasible assignment (schedule) of routes to transporters. Conversely, both GREEDY and ABP provide in a single step complete solutions to P [P in script form], including the number of necessary transporters as well as the sequence of moves performed by each active transporter.

As a final note, one may need to consider congestion issues when routing transporters over the material flow network. Congestion would affect the solution of P [P in script form] by increasing the transporter travel times (and costs [c.sub.ij] Since no real-time information is available at the design stage to determine when blocking occurs, a static approach to eliminate this potential problem is necessary. Herrmann et al. [24] present such an approach by considering the traffic volume (the number of loaded moves) when designing the material flow network. Since the material flow network segments have capacity constraints, the design approach creates, when necessary, alternative origin-destination paths for some resource pairs r, s [epsilon] M and assigns some loaded moves to these paths, which avoids congestion.

7. Concluding remarks

In this paper we have studied the problem of designing a material handling system that employs the minimum number of transporters to transfer material within a manufacturing facility with minimal handling effort. Fixed acquisition and variable operational costs were explicitly considered. An integer program was formulated to capture the trade-off between these two costs. To solve the resulting optimization problem, we developed two heuristic solution approaches: the first allots in a greedy fashion moves to transporters, while the second starts from the optimal solution of the assignment problem and, after grouping moves to routes, allots them to transporters through a two-stage bin-packing procedure. The heuristics were analyzed in terms of computational complexity and worst-case performance, and extensive computational tests were carried out to evaluate their average performance.

The computational results indicate that both heuristics are efficient and adequate to support the material handling system design process. GREEDY solutions diverge from the lower bound of the number of required transporters (provided by the assignment problem A [A in script form]) less than 2% in large-size problem instances; thus, if the fixed cost is the primary consideration, GREEDY seems to be more appropriate. On the other hand, ABP solutions diverge from the variable cost lower bound less than 5% for large-size problem instances; consequently, this heuristic is more appropriate when material handling cost is the main consideration.

The mathematical formulation accurately models the design-level control problem if the design horizon T is a relatively small period (e.g., a shift). In such cases the same pattern of material flow is repeated in each period, and the unloaded travel of transporters during shop operation is expected to be relatively close to that resulting from the solution of P [P in script form]. However, if the manufacturing system is not expected to demonstrate consistent periodicity, the model may not capture the effects of on-line control. In this case, the assumption of available inventory at each resource may not be satisfied. This may, in turn, lead to an underestimation of the unloaded travel and the number of transporters required. To overcome this drawback, a time window may be introduced to reflect the time within which each loaded move is to be performed [21]. To incorporate time windows in the formulation, minor modifications in the graph of moves are necessary to reflect feasible sequences of moves.

In addition to offering a good initial estimate of the number of transporters required for a manufacturing system, the algorithms presented here may be integrated with facility design methods. Given a machine layout, the material handling flow paths can be optimally designed [24]. The solution of the flow path problem will provide the actual inter-resource distances and, thus, costs [c.sub.ij]. Subsequently, the number of transporters and the unloaded moves can be evaluated by solving P [P in script form]. This will provide a realistic cost for the shop layout under consideration. The procedure can be repeated by generating new layouts, until a "good" global solution, in terms of total material handling investment and operational costs, is found [27]. Simulated annealing or genetic algorithms could be employed in such an integrated facility design methodology.

Acknowledgements

This work was supported in part by the National Science Foundation under Grant number NSF EEC 94-02384. The authors would like to thank Michael J. Vukovich for programming the Hungarian and FFD algorithms on the CIM Lab Sun Sparc Workstations, as well as Dr. Ronald G. Askin, focused area Editor of the IIE Transactions, and the anonymous referees for their useful comments and suggestions that helped improve the content and the presentation of this work.

Biographies

Jeffrey W. Herrmann is an Assistant Professor at the University of Maryland, where he holds a joint appointment with the Department of Mechanical Engineering and the Institute for Systems Research. He earned his B.S. in Applied Mathematics from Georgia Institute of Technology, and as a National Science Foundation Graduate Research Fellow from 1990 to 1993, he received his Ph.D. in Industrial and Systems Engineering from the University of Florida. His dissertation investigated production scheduling problems motivated by semiconductor manufacturing. He held a post-doctoral research position in the Institute for Systems Research from 1993 to 1995. He has worked on applied research projects supported by NSF, NIST, ONR, the Semiconductor Research Corporation, the US Army Tank-Automotive Command, Harris Semiconductor, Northrop Grumman Electronic Sensors and Systems Division, Black & Decker, and other manufacturers in the state of Maryland. His publications cover topics in production scheduling, manufacturing facili ty design, and design evaluation and partner selection for agile manufacturing. His current research interests include the design and control of manufacturing systems and the integration of product design and manufacturing system design. He is a member of INFORMS and ASME.

George Ioannou is an Assistant Professor at the Department of Industrial and Systems Engineering of Virginia Polytechnic Institute and State University. He is also the Director of the Manufacturing Systems Integration Laboratory at Virginia Tech, and holds a faculty position at the Department of Operations Research and Marketing of Athens University of Economics and Business. He received his undergraduate diploma in Mechanical Engineering from the National Technical University of Athens in 1990, and his M.Sc./D.I.C. in Industrial Robotics and Manufacturing Automation from Imperial College, London, UK in 1991. Between 1992 and 1995 he was a Graduate Research Assistant at the Institute for Systems Research of the University of Maryland at College Park, where he received his Ph.D. in Mechanical Engineering. His research concentrates on the quantitative and analytical study of production systems, and merges operations research tools with modern information technology to address open problems faced by today's comp lex enterprises and supply chain networks. His work is sponsored by Virginia's Center for Innovative Technology, NSF, the US Navy, and several manufacturing companies. His publications have appeared in various archival journals, and cover topics ranging from facility and material handling system design and operation, to environmentally conscious manufacturing. He is a member of IIE and INFORMS.

Joannis Minis is a manager in Planning SA., a Greek Industrial Management consulting firm. His areas of expertise include production planning and scheduling, industrial facility design, integrated product and process development and business planning. Dr. Minis holds a Ph.D. degree from the University of Maryland in Mechanical Engineering, and has held the positions of Assistant and Associate Professor of Mechanical Engineering at the same University (1988-1997). He has consulted with major US and Greek firms, has conducted research in Manufacturing Systems for over 15 years and has authored several book chapters, over 30 journal articles and over 40 articles in conference proceedings. Dr. Minis has been awarded the Best Paper Award in Database Management in the 1992 ASME CIE conference, and he is the recipient of the 1993 Earl E. Walker Outstanding Young Manufacturing Engineer Award from the Society of Manufacturing Engineers (SME).

Jean-Marie Proth received the Ph.D. degree in Computer Sciences from the University of Nancy (France) and the Doctorat d'Etat in Management Science from the University of Paris-Dauphine (France). He was a Professor in several French Universities from 1969 to 1979. Between 1979 and 1982, he served as a Professor at the European Institute for Advanced Studies in Management in Brussels (Belgium), and then joined INRIA (Institut National de Recherche en Informatique et en Automatique) in France. Since 1986, he has been a Research Director at INRIA and head of the SAGEP team located in Metz (France). He is also Associate Member of the Institute for Systems Research of the University of Maryland (USA). His research interests include: preliminary design of manufacturing systems, design and use of hierarchical production management systems, scheduling, logistics, and modeling/evaluation of flexible manufacturing systems. He is also a consultant in the fields of manufacturing system design and management. He is an aut hor and co-author of several books, Associate Editor of the Journal of Applied Stochastic Models and Data Analysis, a member of the Editorial Board of the Journal of Intelligent Manufacturing and he was Associate Editor of the journal IEEE Transactions on Robotics and Automation until December 1998.

Contributed by the Manufacturing Systems Modeling Department.

(1.) Department of Mechanical Engineering and Institute for Systems Research, University of Maryland, College Park, MD 20742, USA

(2.) Department of Industrial and Systems Engineering, Virginia Tech, Blacksburg, VA 24061-0118, USA E-mail: ioannou@vt.edu

(3.) Planning S.A., Leoforos Pentelis 77-79, Chalandri-Athens, 152-33, Greece

(4.) INRIA-Lorraine, 4 Rue Marconi, Metz 2000, 57070 Metz, France

References

(1.) Jacobs-Blecha, C. and Goetschalckx, M. (1992) Analysis of horizontal transportation systems with multiload carriers, in Proceedings of the 1992 NSF Design and Manufacturing Systems Conference, Atlanta, GA, pp. 871-874.

(2.) Tompkins, J. A., White, JA., Bozer, Y. A., Frazelle, E. H., Tanchoco, J.M.A. and Trevino, J. (1996) Facilities Planning, 2nd edn., John Wiley, New York.

(3.) Egbelu, P.J. (1987) The use of non-simulation approaches in estimating vehicle requirements in an automated guided vehicle based transport system. Material Flow, 4, 17-32.

(4.) Tanchoco, J.M.A., Egbelu, P.J. and Taghaboni, F. (1987) Determination of the total number of vehicles in an AGV-based material transport system. Material Flow, 4, 33-51.

(5.) Wysk, R.A., Egbelu, P.J., Zhou, C. and Ghosh, B.K. (1987) Use of spread sheet analysis for evaluating AGV systems. Material Flow, 4, 53-64.

(6.) Egbelu, P.J. and Tanchoco, J.M.A. (1982) AGVSim user's manual. Technical Report 8204, Department of Industrial Engineering and Operations Research, Virginia Tech, Blacksburg, VA.

(7.) Sinriech, D. and Tanchoco, J.M.A. (1992) An economic model for determining AGV fleet size. International Journal of Production Research, 30, 1255-1268.

(8.) Maxwell, W.L. and Muckstadt, J.A. (1982) Design of automated guided vehicle systems. IIE Transactions, 14, 114-124.

(9.) Bartholdi, III J.J. and Platzmann, L.K. (1988) Decentralized control of a fixed route automated guided vehicle system. lIE Transactions, 21, 76-81.

(10.) Boxer, Y.A. and Srinivasan, M.M. (1991) Tandem configurations for automated guided vehicle systems and the analysis of single vehicle loops. lIE Transactions, 23, 23-27.

(11.) Sinriech, D. and Tanchoco, J.M.A. (1992) Impact of unloaded vehicle flow on performance of single loop AGV systems. International Journal of Production Research, 30, 2237-2252.

(12.) Bozer, Y.A. and Srinivasan, M.M. (1992) Tandem AGV systems: a partitioning algorithm and performance comparison with conventional AGV systems. European Journal of Operational Research, 63, 173-191.

(13.) Ca, C.G. and Tanchoco, J.M.A. (1991) A review of research on AGVS vehicle management. Engineering Costs and Production Economics, 21, 35-42.

(14.) Ioannou, G. and Minis, I. (1998) A review of current research in manufacturing shop design integration. Journal of Intelligent Manufacturing, 9, 57-72.

(15.) Clarke, G. and Wright, J.W. (1964) Scheduling of vehicles from a central depot to a number of delivery points. Operations Research, 12, 568-581.

(16.) Johnson, D., Demers, A., Ullman, J., Garey, M. and Graham, R. (1974) Worst-case performance bounds for simple one-dimensional packing algorithms. SIAM Journal of Computing, 3, 297-325.

(17.) Lewis, R.T. and Parker, R.G. (1982) On a generalized bin-packing problem. Naval Research Logistics Quarterly, 29, 119-145.

(18.) Golden, B.L. and Assad, A.A. (1988) Vehicle Routing: Methods and Studies, Elsevier, Amsterdam.

(19.) Li, C.-L., Simchi-Levi, D. and Desrochers, M. (1992) On the distance constrained vehicle routing problem. Operations Research, 40, 790-799.

(20.) Laporte, G. (1992) The vehicle routing problem: an overview of exact and approximate algorithms. European Journal of Operational Research, 59, 345-358.

(21.) Desrosiers, J., Dumas, Y., Solomon, M.M. and Soumis, F. (1994) Time constrained routing and scheduling, in Handbook in Operations Research and Management Science, Ball, ME., Magnanti, T.L., Monma, C. and Nembauser, G.L. (eds), Elsevier Science, Amsterdam, The Netherlands, pp. 553-578.

(22.) Heragu, S.S. (1992) Recent models and techniques for solving the layout problem. European Journal of Operational Research, 57, 1-13.

(23.) Meller, R.D. and Gau, K.-Y. (1996) The facility layout problem: recent and emerging trends and perspectives. Journal of Manufacturing Systems, 15, 351-366.

(24.) Herrmann, J.W., Ioannou, G., Minis, I., Nagi, R. and Proth, J.M. (1995) Design of material flow networks in manufacturing facilities. Journal of Manufacturing Systems, 14, 277-289.

(25.) Johnson, M.E. and Brandeau, M.L. (1993) An analytic model for the design of a multi-vehicle automated guided vehicle system. Management Science, 39, 1477-1489.

(26.) Papadimitriou, C.H. and Steiglitz, K. (1982) Combinatorial Optimization: Algorithms and Complexity, Prentice-Hall, Englewood Cliffs, NJ, pp. 247-255.

(27.) Ioannou, G. (1996) Integrating shop layout and material handling system design decisions, in Progress in Material Handling Research, Graves, R.J., McGinnis, L.F., Medeiros, D.J., Ward, R.E. and Wilhelm, M.R. (eds), Braun-Brumfield, Inc., Ann Arbor, MI, pp. 227-244.

(28.) Lenstra, J. and Rinnooy Kan, A.H.G.(1981) Complexity of vehicle routing and scheduling problems. Networks, 11, 221-228.

(29.) Carmen, T.H., Leiserson, C.E. and Rivest, R.L. (1990) Introduction to Algorithms, MacGraw-Hill, New York.

                 Results of heuristics for 10 example sets
                (average values over 100 examples per set)
Example set Number of moves       Variable cost
                            [[Z.sup.v].sub.g]
 1              107.45               4035.38
 2              141.68               4012.27
 3              205.31               4001.55
 4              258.42               3967.48
 5              315.57               3943.94
 6              412.28               3917.32
 7              547.32               3985.47
 8              693.67               4005.12
 9              901.55               4020.24
10             1157.42               4043.12
Example set                   Transporters
            [[Z.sup.v].sub.b]  [R.sub.g]   [R.sub.b]
 1                   3837.45     9.28       9.85
 2                   3803.58     9.17       9.36
 3                   3785.13     8.92       9.18
 4                   3746.76     8.34       9.03
 5                   3728.31     8.25       8.78
 6                   3705.67     8.16       8.69
 7                   3721.83     8.04       8.55
 8                   3736.55     7.95       8.60
 9                   3785.18     8.02       8.49
10                   3772.26     7.98       8.32
Example set percent of unutilized
                  capacity
                 [lc.sub.g]       [lc.sub.b]
 1                  21.15           24.32
 2                  18.24           21.37
 3                  14.32           17.56
 4                  11.56           14.28
 5                   9.12           12.03
 6                   7.23           10.64
 7                   6.58            9.47
 8                   6.04            8.71
 9                   5.79            8.25
10                   5.30            7.98
                     Performance of GREEDY: Solutions
                            versus lower bounds
Example set ([[Z.sup.v].sub.           ([R.sub.g] - [[R.sup.*].
             g] - [[Z.sup.*].sub.a]/  sub.a]/[[R.sup.*].sub.a])
            [[Z.sup.*].sub.a]) x 100%          x 100%
 1                    14.32                     36.42
 2                    14.03                     31.72
 3                    13.95                     25.74
 4                    13.48                     20.86
 5                    11.97                     17.38
 6                    11.08                     14.63
 7                    10.76                      9.54
 8                    10.12                      6.97
 9                     9.75                      3.76
10                     9.26                      1.98
Example set ([Z.sub.g] - Z'/Z')
                  x 100%
 1                 22.67
 2                 20.76
 3                 18.49
 4                 16.33
 5                 14.07
 6                 12.46
 7                 10.28
 8                  8.87
 9                  7.36
10                  6.33
             Performance of ABP: solutionsversus lower bounds
Example set ([[2.sup.v]          ([R.sub.b] - [[R.sup.*] ([Z.sub.b] - Z'/Z')
            .sub.b] - [[2.sup.*] .sub.a]/[[R.sup.*]      X 100%
            .sub.a]/[[Z.sup.*]   .sub.a]) X 100%
            .sub.a]) X 100%
 1                  9.85                  41.27                  21.73
 2                  9.13                  35.37                  19.11
 3                  8.74                  30.18                  17.03
 4                  8.12                  24.32                  14.38
 5                  7.31                  20.71                  12.51
 6                  6.25                  19.68                  11.36
 7                  6.14                  14.20                   9.30
 8                  5.67                  12.15                   8.22
 9                  5.28                   8.56                   6.59
10                  4.56                   3.74                   4.23

Appendix A

This Appendix provides a simple example from the literature [8] to demonstrate some key attributes of problem "" and facilitate the discussion of the models and solution algorithms. Consider a manufacturing shop comprising five resource groups, arranged as shown in Fig. A1. The locations of the input and output stations of the resource groups are shown in Fig. A1. Any segment of the thick lines surrounding the resource groups may serve as part of the material flow network. The material flow matrix, provided in Table A1, includes 32 loaded moves per time period.

Table A2 provides, for each origin and destination pair, the index of each loaded move that is part of the associated material flow. In Table A2, [O.sub.r] and [I.sub.r] denote the input and output station of resource group r. Note that one loaded move is equivalent to one unit of a nonzero entry of the material flow matrix.

A tour comprises a closed sequence of loaded moves. The optimal solution of the assignment problem "", for example, provides such tours. The tour that results by solving the assignment problem for our example, using the costs provided in Maxwell and Muckstadt [8], includes the following sequence of loaded moves: {l, 27, 2, 28, 3, 29, 4, 5, 6, 7, 8, 9, 10, 11, 12, 30, 13, 31, 14, 32, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, l}. This tour can be transformed to the following sequence of resource groups: {5, 1, 1, 5, 5, 1, 1, 5, 5, 1, 1, 5, 5, 2, 5, 2, 5, 2, 5, 2, 5, 2, 5, 2, 5, 2, 5, 2, 5, 3, 3, 5, 5, 3, 3, 5, 5, 3, 3, 5, 5, 4, 5, 4, 5, 4, 5, 4, 5, 4, 5, 4, 5, 4, 5, 4, 2, 5, 2, 5, 2, 5, 2, 5}, using the index relationships of Table A2. Note that the tour starts at the output station of resource group 5, and includes both loaded and unloaded moves. For example, the second and third elements of this set imply an unloaded move between the input and the output station of resource group 1. If the above tour r epresents the sequence of loaded moves performed by transporters, the unloaded moves populate the matrix of Table A3; in this table, the elements of the first row and column represent resource groups.

We defined a route [[gamma].sub.u] to be a sequence of moves that start and end at the same input or output station. A collection of routes with a common starting point is a route set [[gamma].sub.v] According to those definitions and the route determination procedure employed by the ABP heuristic, the above sequence of resource groups leads to the following two route sets: [[gamma].sub.1] = {(5,l,1,5), (5,2,5), (5,3,3,5), (5,4,5)} and [[gamma].sub.2] = {(2,5,2)}, and the set of route sets is simply [gamma] = {[[gamma].sub.1], [[gamma].sub.2]}. Note that routes within each route set are enclosed in between parentheses.

Appendix B

In this Appendix we provide a formal proof of the membership of problem P [P in script form] to the class NP.

Theorem 6. Problem P [P in script form] is NP-complete.

Proof. We know that the Distance-constraint Vehicle Routing Problem (DVRP) is NP-complete [28]. To show that P [P in script form] is also in NP, we need to provide a polynomial transformation from an instance of DVRP to an instance of P [P in script form] and show that DVRP [less than or equal to]P P [P in script form]. This can be accomplished as follows. Consider a DVRP instance defined on the complete graph G' = (L', E'), where L' = L P [P in script form] {[xi]}, [xi] is a node corresponding to the depot, and

E' = {(i,j) : i [epsilon] L' [lambda] j [epsilon] L'}. Define the cost coefficients of this DVRP as follows:

[c.sub.ij]' = {c if i [epsilon] L [lambda] j [epsilon] L, 0 if i = [xi] [xi] j = [xi],

i.e., the moves to and from the depot are of zero length for all nodes of G'. Suppose also that the path (i,[xi],j) of DVRP on G' = (L',E') is substituted by the path (i,j), [For all]i, j [epsilon] L when we consider problem P on G = (L,E), where E = E' \ {(i,j) [epsilon] E' : i = [xi] [vee] j = [xi]}, and vice-versa. Furthermore, assume that for DVRP the distance bound is T - c, while the distance bound for P is T. Then, we need to show that DVRP has a solution with k vehicles and [Z.sub.[upsilon]] routing cost on G', if and only if P has a solution with k vehicles and a [Z.sub.[upsilon]] + kc routing cost on G. This is straightforward by the construction of the problem instances. If DVRP has a solution with k vehicles and [Z.sub.[upsilon]] routing cost, then substitute the path (i,[xi],j) on G' with the path (i,j) on G; the equivalent P instance has a solution with k transporters and [Z.sub.[upsilon]] + kc routing cost on G, and the capacity constraints are satisfied since the cost augmentation is exactly c for e ach transporter. Conversely, if P has a solution with k transporters and [Z.sub.[upsilon]] + kc routing cost on G, substitute one path (i,j) for each transporter by the corresponding path (i,[xi],j) on G'; the equivalent DVRP instance on G' has a solution with k vehicles and routing cost [Z.sub.[upsilon]], and the capacity constraints are satisfied. Thus, DVRP [less than or equal to]P P, and consequently, P is NP-complete.

Appendix C

This Appendix provides the proofs of the theorems of Section 5.

Proof of Theorem 1. At each iteration of GREEDY exactly one loaded move is assigned to a transporter. However, when the capacity of a transporter is exceeded the last move of the current iteration must be reexamined in a subsequent iteration, since it cannot be assigned to the current transporter. This will occur at most once for each transporter, and if I is the total number of iterations, we obtain: I [less than or equal to] n + [[sigma].sub.k[epsilon]V] [y.sub.k]. From inequality (10), [[sigma].sub.k[epsilon]V] [y.sub.k] [lesser than] [n/2]. Thus, I [less than or equal to] n + [n/2] [less than or equal to] 3n/2+ 1.

Identifying the minimum row element of the matrix [[c.sub.ij]] in Step 5 of GREEDY, requires at most n comparisons at each iteration; actually, we need n - 1 at the first iteration, n - 2 at the second one, etc. As a result, the maximum total number of operations, [OP.sub.g], of this heuristic is bounded as follows:

O[P.sub.g] [lesser than] n x I [less than or equal to] n x (3n/2 + 1) = [3n.sup.2] + 2n/2. (A1)

Consequently, the computational complexity of GREEDY is O([n.sup.2]).

Proof of Theorem 2. The primal-dual algorithm for the assignment problem in Step 1 of ABP requires 0([n.sup.3]) operations, since the n loaded moves in L must be matched [26]. The identification of routes and route sets in Steps 4 through 9 requires inspection of each loaded move in the sequences given by the optimal solution of A. Thus, exactly n operations are required before the criterion of Step 9 is no longer satisfied.

Since the maximum number of routes is n/2, Step 10 requires at most: (i) (n/2) x log(n/2) operations to sort these routes in decreasing order of length [29]; and (ii) n/2 operations to pack them into transporters of size T [16]. Finally, the enumeration of the best connecting points in Step 14 requires at most O([n.sup.2]) operations, as shown below. Consider K transporters to be merged, and let [n.sub.1],...,[n.sub.K] be the number of moves assigned to each of them. Then, in the worst case, Step 14 would examine each possible pair. This would require [[[sigma].sup.K].sub.i=1] [[[sigma].sup.K].sub.j=1,j[not equal to]i] [n.sub.i][n.sub.j] operations. Since [[[sigma].sup.K].sub.i=1] [n.sub.i] = n, the maximum number of operations at this Step is [n.sup.2].

Combining the results of the arguments, the total number of operations, [OP.sub.abp], is bounded by:

[OP.sub.abp] [less than or equal to] [n.sup.3] + n + n/2 x log(n/2) + n/2 + [n.sup.2]. (A2)

Consequently, the computational complexity of ABP is O([n.sup.3]).

Proof of Theorem 3. The GREEDY algorithm selects the minimum element of the cost matrix in the first iteration. However, in the remaining n - 1 iterations the algorithm may select elements different than those of the optimal solution. Since the maximum deviation of the cost elements is [c.sub.max] - [c.sub.min], it follows that:

[[Z.sup.v].sub.g] - [[Z.sup.v].sub.opt] [less than or equal to] (n - 1) x ([c.sub.max] - [c.sub.min]) [less than] n x ([c.sub.max] - [c.sub.min]). (A3)

The variable cost of the optimal solution to problem P is bounded from below by [nc.sub.min]. Thus, from inequality (A3), we conclude that:

[[Z.sup.v].sub.g] - [[Z.sup.v].sub.opt]/[[Z.sup.v].sub.opt] [less than] n([c.sub.max] - [c.sub.min])/[nc.sub.min] [Left/Right arrow] [[Z.sup.v].sub.g]/[[Z.sup.v].sub.opt] [less than] [c.sub.max]/[c.sub.min]

Lemma 1 shows that the bound provided by the first of inequalities (A3) is tight, i.e., there exist problem instances of P, for which (A3) holds at equality.

Lemma 1. There exists an n x n matrix [[c.sub.ij]] of cost coefficients, a time horizon T, and a sequence Q = {[p.sub.1],[p.sub.2],...} of loaded moves selected in Step 2 of GREEDY for which:

[[Z.sup.v].sub.g] - [[Z.sup.v].sub.opt] = (n - 1) x ([c.sub.max] - [c.sub.min])

Proof. Consider an instance of problem P with [c.sub.max] = 2[c.sub.min] and T = 5[c.sub.min]. Let [[c.sub.ij]] be the 5 x 5 matrix with entries [c.sub.max] and [c.sub.min] only, as shown below:

               [infinity]  [c.sub.min] [c.sub.min] [c.sub.min] [c.sub.min]
               [c.sub.min] [infinity]  [c.sub.max] [c.sub.max] [c.sub.max]
[[C.sub.ij]] = [c.sub.max] [c.sub.min] [infinity]  [c.sub.min] [c.sub.min]
               [c.sub.min] [c.sub.min] [c.sub.min] [infinity]  [c.sub.max]
               [c.sub.min] [c.sub.min] [c.sub.min] [c.sub.max] [infinity]

We claim that for this problem instance (A3) holds at equality.

Let Q = {1, 4} be a sequence that could be followed by GREEDY. This sequence could provide the following matchings (assigned to the two transporters, [l.sub.1] and [l.sub.2]), depending on the tie-breaking rules:(1) [l.sub.1] = {(1,2), (2, 3), (3, l)} with subtour length equal to [C.sub.min] + 2[c.sub.max] = T;(2) [l.sub.2] = {(4, 5), (5, 4)} with subtour length equal to 2[C.sub.max] [less than] T. The resulting variable cost is [[Z.sup.v].sub.g] [C.sub.min] + 4[C.sub.max]. The optimal solution to this instance of P gives the following assignment matchings for a single transporter:[l.sub.1] = {(1,5),(5,3),(3,4),(4,2),(2,1)}, and the variable cost is [[Z.sup.v].sub.opt] = 5[c.sub.min]. Thus, [[Z.sup.v].sub.g] - [[Z.sup.v].sub.opt] = 4([c.sub.max] - [c.sub.min]).

Proof of Theorem 4. Consider the subtours derived by solving the assignment problem corresponding to p. If each subtour comprises two moves of length [C.sub.min], and since the maximum number of subtours is n/2, the total variable cost is [[Z.sup.*].sub.a] = n x [C.sub.min]. When ABP attempts to connect route sets, the total cost is augmented at most by 2 x ([C.sub.max] - [C.sub.min]). This connection removes two complete moves from the solution and introduces two new complete moves. In the worst case, [n/2] - 1 connections will be performed, and, consequently

[[Z.sup.v].sub.b] - [[Z.sup.*].sub.a] [less than or equal to]([n/2] -1) x (2[c.sub.max] - 2[c.sub.min])

[less than or equal to] n X ([c.sub.max] - [C.sub.min]). (A4)

From (A4) and the assignment cost, we conclude that:

[[Z.sup.v].sub.b] - [[Z.sup.*].sub.a]/[[Z.sup.*].sub.a] [less than or equal to] ([c.sub.max] - [c.sub.min])/[c.sub.min] [Left/Right arrow] [[Z.sup.v].sub.b]/[[Z.sup.*].sub.a] [less than or equal to] [c.sub.max]/[c.sub.min].

Proof of Theorem 5. The optimal number of transporters is bounded from below by the ratio of the variable cost over the design horizon T, i.e., [R.sub.opt] [greater than or equal to] [[[Z.sup.v].sub.opt]/T][greater than or equal to] [[Z.sup.u].sub.opt]/T. From inequality (10) we know that [R.sub.h] [less than or equal to] n/2. Also [[Z.sup.v].sub.opt] [greater than or equal to]n[c.sub.min], and consequently [R.sub.h] [less than or equal to] [[Z.sup.v].sub.opt]/2[c.sub.min]. Thus,

[R.sub.h]/[R.sub.opt] [less than or equal to] 1/2 X T/[C.sub.min].

                       Material flow per time period
From \ To 1 2 3 4 5
1         0 0 0 0 3
2         0 0 0 0 4
3         0 0 0 0 3
4         0 0 0 0 0
5         3 8 3 8 0
                          Move-flow relationships
Origin    Destination Move index
[O.sub.5]  [I.sub.1]     1-3
[O.sub.5]  [I.sub.2]     4-11
[O.sub.5]  [I.sub.3]    12-14
[O.sub.5]  [I.sub.4]    15-22
[O.sub.2]  [I.sub.5]    23-26
[O.sub.1]  [I.sub.5]    27-29
[O.sub.3]  [I.sub.5]    30-32
                 Matrix of unloaded moves per time period
From \ To 1 2 3 4 5
1         3 0 0 0 0
2         0 0 0 0 8
3         0 0 3 0 0
4         0 1 0 0 7
5         0 3 0 0 7