Metaheuristics for the Order Batching Problem in Manual Order Picking Systems

In manual order picking systems, order pickers walk or drive through a distribution warehouse in order to collect items which are requested by (internal or external) customers. In order to perform these operations efficiently, it is usually required that customer orders are combined into (more substantial) picking orders of limited size. The Order Batching Problem considered in this paper deals with the question of how a given set of customer orders should be combined such that the total length of all tours is minimized which are necessary to collect all items. The authors introduce two metaheuristic approaches for the solution of this problem: the first one is based on Iterated Local Search; the second on Ant Colony Optimization. In a series of extensive numerical experiments, the newly developed approaches are benchmarked against classic solution methods. It is demonstrated that the proposed methods are not only superior to existing methods but provide solutions which may allow distribution warehouses to be operated significantly more efficiently.


Introduction
Order picking is a warehouse function dealing with the retrieval of articles (items) from their storage location in order to satisfy a given demand specified by (internal or external) customer orders (Petersen and Schmenner 1999: 481). Order picking arises because incoming articles are received and stored in (large-volume) unit loads while customers order small volumes (less-than-unit loads) of different articles. As a warehouse function, order picking is critical to each supply chain, since underperformance results in an unsatisfactory customer service (long processing and delivery times, incorrect shipments) and high costs (labor cost, cost of additional and/or emergency shipments). The significance of improvements of warehouse operations becomes evident from a joint study of the European Logistics Association and the management consulting company A.T. Kearney (European Logistics Association and A. T. Kearney 2004) which revealed that the total warehousing costs of a company may amount to 25ô of its total logistics costs. Of all warehouse operations, order picking is considered to include the most cost-intensive ones. According to Frazelle (2002) up to 50ô of the total warehousing operating costs can be attributed to order picking. Drury (1988, also see Tompkins, White, Bozer, Frazelle, and Tanchoco 2003) and Coyle, Bardi, and Langley (1996) give estimations of up to 60ô and 65ô, respectively, for these costs. The large proportion of order picking (operations) costs originates from the fact that, like many other repetitive material-handling activities, order picking is still a function, which involves the employment of human operators on a large scale. Even though there have been different attempts to automate the picking process, manual order picking systems are still prevalent in practice. Such manual order picking systems can be differentiated into two categories: Picker-to-parts systems, in which order pickers drive or walk through the warehouse and collect the requested articles; and parts-topicker systems, in which automated storage and retrieval systems deliver the articles to stationary order pickers (Wäscher 2004: 324). With respect to systems of the first kind, three planning problems can be distinguished on the operative level (Caron, Marchet, and Perego 1998: 1), namely the assignment of articles to storage locations (storage location), the grouping of several customer orders into picking orders (order batching), and the routing of order pickers through the warehouse (picker routing). This paper focuses on order batching, which has been proven to be pivotal for the efficiency of warehouse operations (de Koster, Roodbergen, and van Voorden 1999: 232). Improved order batching reduces the (total) lengths of the tours which the order pickers have to cover in order to collect the requested articles. Likewise, the travel time and, consequently, the total picking time (i.e. the total time order pickers spend in the warehouse collecting the articles) will be reduced. According to our experience from practice, this can result in significant savings of labor cost, since it does not only allow for reducing the working time of the pickers, but also for reducing expensive overtime or even for downsizing the (picker) workforce. Since the picking time is also part of the time interval from the moment a customer order is placed until the requested articles are received by the customer, a reduction of the travel time can also immediately be translated into a reduction of the customer order's delivery lead time. In other words, improved order batching also contributes to the advancement of the customer service offered by an order picking warehouse. Both aspects, cost reduction on the one hand, and improvement of customer service on the other, have a positive impact on the competitiveness of the total supply chain. In the past, for the Order Batching Problem predominantly traditional (constructive) heuristics have been suggested as solution methods. The aim of this article is to demonstrate that modern heuristics can lead to substantially improved solutions in order batching. We focus on two metaheuristics, which have demonstrated to provide excellent results for other combinatorial optimization problems. The first one is a variant of Iterated Local Search; the second one is a population-based approach, namely a variant of ant colony optimiza-tion --the Rank-Based Ant System. The remainder of the paper is organized as follows: In Section 2 the Order Batching Problem will be defined and a mathematical formulation will be given. Section 3 contains a literature review of solution approaches for the Order Batching Problem. In the following two sections our implementation of the metaheuristics will be presented, Iterated Local Search (ILS) in Section 4, and Rank-Based Ant System (RBAS) in Section 5. Extensive numerical experiments have been carried out to evaluate the performance of the metaheuristics. The design of these experiments (including the description of warehouse parameters, algorithm parameters, and problem classes) will be presented in Section 6. A comprehensive analysis of the test results will be given in Section 7. The paper will conclude with a summary and an outlook on further research opportunities in Section 8.

Problem Description
When performing his/her tasks, an order picker is guided by a so-called pick list. This list specifies the sequence according to which the requested articles should be collected, as well as the quantities in which they are to be picked. A pick list may contain the articles of a single customer order (pick-byorder) or of a combination of customer orders (pick-by-batch). In practice, the sequence in which the articles are to be picked and the corresponding route of the order picker (which starts at the depot, proceeds to the respective storage locations, and returns to the depot) is usually determined by means of a so-called routing strategy, e.g., by the S-Shape Heuristic or by the Largest-Gap Heuristic. Despite the fact that an optimal polynomial time algorithm for the picker routing problem exists (Ratliff and Rosenthal 1983), it is hardly ever used in practice. Order pickers seem not to accept the optimal routes provided by the algorithm, because of their not-always-straightforward or sometimes even confusing routing schemes (de Koster, Roodbergen, and van Voorden 1999: 230). The heuristic routing schemes on the other hand, are fast to memorize and quite easy to follow. This helps to reduce the risk of missing an article to be picked --an aspect that might be more important than a small reduction of the tour length. The examples of the routing schemes of the S-Shape Heuristic ( Figure 1) and the Largest-Gap Heuristic Figure  2) demonstrate their self-evident character. The black rectangles symbolize the locations of articles to be picked on the respective routes. The S-Shape Heuristic (or traversal strategy) gives a solution in which the order picker only enters an aisle if at least one requested article is located in that aisle and traverses it entirely. Afterwards the order picker proceeds to the next aisle containing a requested article. An exception could be the rightmost aisle containing an article to be picked: if the order picker is positioned in the front cross-aisle, he would pick the articles of that right-most aisle and return to the front aisle, i.e. the cross-aisle which contains the depot. Solutions provided by the Largest-Gap Heuristic are characterized in the following way: The order picker traverses the first and last aisle containing an article to be picked entirely. All the other aisles are entered from the front and back in such a way that the non-traversed distance between two adjacent locations of articles to be picked in the aisle is maximal. Figure 1 shows --for S-Shape Routing --that benefits may arise from collecting the articles of two (or more) customer orders on a single tour instead of two (or more) tours, in particular if the articles are identical or closely located to each other. The upper and middle pictures of Figure 1 depict the tours of two separate customer orders, while the picture at the bottom shows the tour resulting from batching the two customer orders. The length of the resulting tour is obviously smaller than the total length of the two separate tours. The same effect is demonstrated for Largest-Gap Routing in Figure 2. Order pickers are typically using a picking device (a cart or a roll pallet) for the collection of the requested articles. Depending on the size of the picking device, customer orders can be combined until the capacity of the device is exhausted. On the other hand, splitting of customer orders is usually prohibited since it would result in an additional, not-acceptable sorting effort. Order batching can be carried out as proximity batching, where orders are combined considering their locations in the warehouse, or as time window batching, where orders are combined according to their arrival time (Choe and Sharp 1991). Here we consider a situation where we assume that all orders are known beforehand. Therefore, we will concentrate on proximity batching. Due to the high proportion of time-consuming manual operations, order picking is considered to be the most labor-cost-intensive function in a warehouse (Drury 1988). Consequently, the minimization of picking times is of great importance for the efficient control of the picking process. The total order picking time (time spent by order pickers to collect the articles of all customer orders) consists of the setup times for the routes; the travel times that are needed to travel to, from, and between the locations of articles to be picked; the search times for the identification of the articles; and the times needed for picking the articles  (Tompkins, White, Bozer, Frazelle, and Tanchoco 2003). Among these components, the travel time is of outstanding importance, since it consumes the major proportion of the total order picking time, while the other components can be looked upon to be constant (search times and pick times) or neglectable (setup times). Furthermore, given the pickers' travel velocity to be constant, the minimization of the total travel time is equivalent to the minimization of the total length of all picker tours (Jarvis and McDowell 1991: 94). In other words, the Order Batching Problem can be defined as follows: How can a given set of customer orders, with given storage locations, given routing strategy and given capacity of the picking devices, be grouped (batched) into picking orders such that the total length of all necessary picking tours is minimized? (Wäscher 2004: 337)

Model Formulation
According to Gademann and van de Velde (2005), a mathematical formulation of the problem can be given as follows: Let J ã è1; : : : ; né be the set of customer orders, C be the capacity of the picking device, and c j the capacity required for order j (j ∈ J). Furthermore, let each batch of customer orders be described by a vector a i ã äa i1 ; :::; a in å with binary entries a ij stating whether an order j is included in a batch i (a ij ã 1) or not (a ij ã 0). The set I of all feasible batches is characterized by the fact that the capacity of the picking device is not violated, i.e. the following property holds: If we define binary decision variables x i (i ∈ I), which describe if a batch i is chosen (x i ã 1) or not (x i ã 0), and let d i further represent the length of the picking tour --subject to a given routing strategy --on which all orders of a batch i are collected, then the following optimization model can be formulated: The sets of constraints (3) and (4) ensure that a set of batches is chosen in a way that each customer order is included in exactly one of the chosen batches. It is important to keep in mind that the number of possible batches and, therefore, the number of binary variables grows exponentially with the number of orders. The Order Batching Problem as described above is known to be N P-hard (in the strong sense) if the number of orders per batch is greater than two (Gademann and van de Velde 2005).

Review of Solution Approaches
For the Order Batching Problem only a few approaches have been developed which try to solve BuR --Business Research Official Open Access Journal of VHB Verband der Hochschullehrer für Betriebswirtschaft e.V. Volume 3 | Issue 1 | May 10 | 82--105 the problem to or get close to optimality. For the (general) model described in the previous section, Gademann and van de Velde (2005) present a branch-and-price algorithm with column generation that was able to provide optimal solutions for small instances (up to 32 customer orders). Bozer and Kile (2008) present a mixed-integer programming approach, that generates near optimal solutions for small order sizes (up to 25 customer orders). Their problem formulation is limited to S-Shape Routing since they consider complete traversed aisles in their objective function. Furthermore, Chen and Wu (2005) describe an order batching approach based on data mining and integer programming. In this approach, first similarities in customer orders are determined by means of an association rule, then a 0-1 integer programming approach is applied in order to cluster the orders into batches. For larger problems --as they usually occur in practice --the use of heuristics is still suggested. Heuristic solution approaches proposed for the Order Batching Problem are of the constructive type in the first place. They can be distinguished into three groups: priority rule-based algorithms, seed algorithms, and savings algorithms (Wäscher 2004). Priority rule-based algorithms consist of a twostep procedure: In the first step, priorities are assigned to customer orders, which provide the sequence for the second step in which the orders are allocated to batches. Several rules have been described in literature according to which the priorities can be determined. The probably best known and most straightforward way is the application of the First-Come-First-Served (FCFS) rule. Other rules include the application of the two-dimensional and the four-dimensional spacefilling curves by Gibson and Sharp (1992) or the six-dimensional space-filling curve by Pan and Liu (1995). The allocation of the customer orders to batches can either be done by means of the Next-Fit Rule (customer orders are added to a batch until the capacity limit is reached; then a new batch is opened), by the First-Fit Rule (all batches are numbered in the order in which they have been opened; the next customer order is allocated to the first batch which provides sufficient capacity), or the Best-Fit Rule (the next customer order is assigned to the batch with the least remaining capacity) (Wäscher 2004). Elsayed (1981) and Elsayed and Stern (1983), generate batches sequentially. Their procedure can be divided into two phases: seed selection phase and order congruency phase. During the seed selection phase, an initial order --the so-called "seed" --is chosen. A large variety of rules is available for the selection of the seed, e.g., choose the first order, choose the largest order, or choose the order with the longest picking tour, and others. Furthermore, the seed can be determined in a single mode (where only the first order in the batch defines the seed) or in cumulative mode (all orders included in the batch define the seed). Afterwards, unassigned customer orders will be added to the seed according to an order-congruency rule, which measures the "distance" from an order not yet allocated to the seed of the batch. Examples of criteria according to which this "distance" can be determined are the number of additional aisles to be visited, the difference between the centers of gravity of the order and the seed, or the sum of the travel distance between every item of the order to the closest item in the seed. An overview of the various seed selection and order congruency rules is given by de Koster, van der Poort, and Wolters (1999), Ho and Tseng (2006) and Ho, Su, and Shi (2008). Savings algorithms are based on the Clarke-and-Wright Algorithm for the Vehicle Routing Problem (Clarke and Wright 1964) and have been adapted in several ways for the Order Batching Problem. In the initial version for the Order Batching Problem, denoted by CîW(i), for each combination of customer orders i and j the savings sav ij are computed which can be obtained by collecting the articles of the two customer orders on one (large) tour instead of collecting them in two separate tours. Starting with the pair of orders with the highest savings, the pairs are considered for being assigned to a batch in a non-ascending order. This may lead to three situations: in case that none of the two orders have been assigned, a new batch will be opened for them; if one of the orders has already been assigned, the other one is added to the batch if the remaining capacity is sufficient; in the case that not enough capacity is available or that both orders have already been assigned, the next pair of orders will be considered. All orders which are left unallocated at the end of the process will be assigned to an individual batch each. The algorithm can be improved by calculating the savings anew each time orders have been allocated to a batch. This variant of the algorithm is denoted by CîW(ii). In detail, the single orders which have been combined into a batch are deleted from further considerations and the new batch will be interpreted as a new "large" order when the savings are determined again. A savings algorithm that generates batches sequentially is the EQUAL algorithm (Elsayed and Unal 1989). In this method, the first pair of orders which has been allocated to a batch is considered as the initial seed. Then the order which --in combination with the seed --leads to the highest savings and does not exceed the capacity of the picking device, is added to the batch. All orders in the batch form the new seed. A new batch is opened if none of the remaining orders fits into the batch. In the Small-and-Large Algorithm (Elsayed and Unal 1989) two subsets are defined, namely a set of large orders and a set of small ones. In a first step, the large orders are assigned to batches according to the EQUAL algorithm described above. The small orders, in a non-ascending order of their size, are then assigned to the batch where they generate the highest savings without violating the remaining capacity. Remaining orders are again assigned to new batches. Comprehensive numerical experiments (de Koster, van der Poort, and Wolters 1999) have shown that either seed algorithms or savings algorithms provide the best solutions with respect to the total length of all necessary tours, dependent on the warehouse layout and on customer order characteristics. Apart from the constructive heuristics described above, Hsu, Chen, and Chen (2005) have suggested a genetic algorithm for the Order Batching Problem. Their approach includes an aisle-metric for the determination of the tour lengths and, therefore, is limited to S-Shape Routing, only. Tsai, Liou, and Huang (2008) describe an approach for the integrated solution of both the batching and the routing problem by means of a genetic algorithm. Due to these specific conditions under which they can be applied, both approaches will not be considered here, any further. Furthermore, Gademann and van de Velde (2005) describe a simple form of an iterated descent algorithm as a part of their branch-and-price procedure.

Iterated Local Search
Iterated Local Search has been successfully applied to a variety of optimization problems, for instance to the Traveling Salesman Problem (Katayama and Narihisa 1999), Scheduling Problems (Brucker, Hurink, and Werner 1996), or the Quadratic Assignment Problem (Stützle 2006). The main principle can be described as follows (Lorenço, Martin, and Stützle 2003): Let S be the set of feasible solutions of an optimization problem. A solution s äs ∈ Så is called a neighbor of solution s äs ∈ Så if it can be obtained by applying a single local transformation to s. The heuristic consists of two alternating phases: a local search phase and a perturbation phase. In the local search phase one starts from an initial solution s 0 ∈ S and finds a (probably randomized) sequence of feasible solutions s 0 ; s 1 ; : : : ; s m ãŝ. Each element s j äj ∈ è1; : : : ; méå is a neighbor of its predecessor, possessing a smaller (minimization problems) or higher (maximization problems) objective function value than the previous one. Provided the problem is bounded,ŝ is a local optimum. Iterated Local Search aims at exploring the vicinity of this local optimum (used as an incumbent solution) more closely in order to identify a solution with an improved objective function value. Therefore, in the perturbation phase, the incumbent solution is partially modified (perturbed) and a further local search phase is applied to this modified solution. The new solution stemming from the local search phase has to pass an acceptance criterion in order to become the new incumbent solution, otherwise the previous solution remains the incumbent solution for a further perturbation. This criterion may allow for the acceptance of deteriorated solutions. These two phases are repeated until a termination condition is met. The challenge consists in choosing an appropriate perturbation scheme: if the perturbation is too marginal, one will often obtain identical solutions after different local search phases; if the perturbation is too extensive, one might unintentionally leave a good subspace of the solution space. For the Order Batching Problem this general principle has been modified in the following way: An initial solution is generated by means of the FCFS rule. Since all customer orders have the same priority, applying the FCFS rule in this context is equivalent to generating a random sequence of the

BuR --Business Research Official Open Access Journal of VHB Verband der Hochschullehrer für Betriebswirtschaft e.V. Volume 3 | Issue 1 | May 10 | 82--105
customer orders, according to which the orders are assigned to batches afterwards. Two different solutions, i.e. two different sets of batches, are called neighbors if one solution can be achieved from the other by interchanging two orders from different batches (SWAP) or by assigning one order to a different batch (SHIFT). For the execution of these local transformations we consider a list of the customer orders which have been assigned to a batch. The sequence in which the customer orders appear on this list is set up in the chronological order in which the customer orders were assigned to the batch, e.g., the first order assigned to the batch is at position one, the second order assigned to the batch is at position two, etc. A SWAP is applied in a way that the exchanged orders are taking exactly the same positions in the new batches as their counterparts occupied before. In the case of a SHIFT the chosen customer order is assigned to the end of the new batch and is erased from the old batch in a way that all successive customer orders in that batch rise in rank. In order to determine the objective function value of a neighbor, only the tour lengths of the two modified batches have to be recalculated (deltaevaluation). In the local search phase (Function 1), a systematic first improvement strategy is used: We try to reduce the objective function value by a SWAP. If an improvement can be obtained, we proceed from the new solution and search for another improvement by a SWAP. In the case that no SWAP can be identified which improves the objective function value, we try to achieve an improvement by SHIFTs. If no further improvement can be identified, we search again for a SWAP that leads to an improved solution. This is repeated until no improvement by SWAPs and SHIFTs can be identified. The local search phase as described above is a form of a variable neighborhood descent with only two neighborhoods. The perturbation phase (Function 2) has been designed in the following way: We select two different batches k and l randomly and move the first q orders from batch k to batch l, and the first q orders of l to k (q is a random number, which is at most half of the number of orders in k and l, respectively). If this rearrangement is not possible, i.e. capacity constraints would be violated, remaining orders will be assigned to a new batch. The determination of the number of exchanges is a critical step. If the

Function 1 local_search(s)
repeat repeat try to exchange two orders from different batches of s in order to reduce the total length of all tours (SWAP); until no further improvement is possible; repeat try to assign an order from one batch of s to another batch in order to reduce the total length of all tours (SHIFT); until no further improvement is possible; until no further improvement is possible; return s;

Function 2 perturbation(s; )
for i ã 1 to do choose two batches k and l from s randomly; choose q with q > 0 and smaller than half of the number of orders in k and l; remove the first q orders from k and l; insert the removed orders from k into l as long as the capacity constraint is not violated; insert the removed orders from l into k as long as the capacity constraint is not violated; assign the remaining orders to a new batch; end for return s; number is too small, we would probably fall back into the same local optimum. On the other hand, if this number is too high, we would not be able to intensify a good solution subspace. We restrict the number of rearrangements to n ± à 1, where n is the number of batches in the best known solution and (rearrangement parameter) a constant in ae0; 1ç. Regarding the acceptance criterion, a new solution is accepted as an incumbent solution if its objective function value is lower than the best currently known one. Furthermore, we also allow for a few deteriorating steps, i.e. if a sequence of perturbation phases and local search phases applied to a particular incumbent solution does not lead to a new global best solution within a certain time interval t. The solution obtained in the last local search phase will be chosen as the new incumbent solution, if the difference between the length of

BuR --Business Research Official Open Access Journal of VHB Verband der Hochschullehrer für Betriebswirtschaft e.V. Volume 3 | Issue 1 | May 10 | 82--105
all picking tours of the last solution däså and the length of all picking tours of the best known solution däs å is smaller than a threshold ± däs å, where the threshold parameter is in ae0; 1ç. Otherwise the incumbent solution will be perturbated again. Algorithm 1 summarizes the heuristic.

Algorithm 1 Iterated Local Search (ILS)
Parameters: rearrangement parameter , threshold parameter , time interval t for an incumbent solution; let s initial be a solution provided by the FCFS rule; s ñã local_search(s initial ); interval t and däså á däs å < ± däs å then s incumbent ñã s; end if until termination condition is met; s is the solution of ILS;

Rank-Based Ant System
The general idea of ant colony optimization (Ant Systems) is inspired by nature, where it can be observed that ants easily find the shortest path from a nest to a food source by marking their trails with pheromones. Shorter paths will soon get marked with a higher amount of pheromones, such that following ants will choose those marked trails with a higher probability. In the subsequent course this will lead to a self-reinforcing process, ending in a situation where nearly all ants follow the shortest path. Artificial Ant Systems for combinatorial optimization problems were first introduced by Colorni, Dorigo, and Maniezzo (1991) ;Dorigo, Maniezzo, and Colorni (1991) and Dorigo (1992) and have been thoroughly discussed for a multitude of problems, e.g., the Traveling Salesman Problem, Scheduling and Routing Problems (Dorigo and Stützle 2004). Apart from the basic Ant System, several extensions have been developed, e.g., Max-Min Ant Systems (Stützle and Hoos 2000), Ant Systems with elitist strategies or Rank-Based Ant Systems (Bullnheimer, Hartl, and Strauss 1999). In artificial Ant Systems, the pheromone effect taken from nature is combined with another aspect not common for real ants. This is done in order to integrate a greedy behavior into artificial Ant Systems, in addition to the adaptive behavior observed for real ants. For the adaption to the Order Batching Problem, a savings-based Ant System has been chosen which combines greedy (savings) and adaptive (pheromone) aspects when deciding about the combination of two batches. Savings-based Ant Systems in general work very well for those types of problems where the Clarke-and-Wright Algorithm is a good heuristic (e.g., for the Capacitated Vehicle Routing Problem or, as in our case, for the Order Batching Problem). The Ranked-Based Ant System (RBAS) was implemented since it has proven to provide good results in combination with a savings-based approach (Doerner, Gronalt, Hartl, Reimann, Strauss, and Stummer 2002;Reimann, Doerner, and Hartl 2004). The general principle of RBAS has been adopted for the Order Batching Problem as follows: In the beginning each order forms a single batch. In the subsequent steps, the batches are combined as long as the capacity constraint for the picking device is not violated. For each possible combination äk; lå of two batches, we compute the corresponding savings sav kl (consult Section 3) and the pheromone intensity˜ kl . The pheromone intensity of a batch combination is determined as the sum of pheromones of all pairs of customer orders (called order combinations), where one order is in batch k and one order is in batch l , divided by the number (n k ± n l ) of possible order combinations between the batches k and l: pheromone intensity of order combination äi; jå, n k : number of orders assigned to batch k, n l : number of orders assigned to batch l.

BuR --Business Research Official Open Access Journal of VHB Verband der Hochschullehrer für Betriebswirtschaft e.V. Volume 3 | Issue 1 | May 10 | 82--105
The probability p kl for the combination of two batches k and l is defined as pheromone intensity of batch combination äk; lå, : parameter for controlling the influence of the pheromone intensity˜ kl ( OE 0), sav kl : savings obtained by combining the batches k and l, : parameter for controlling the influence of the savings sav kl ( OE 0), Ω: set of all feasible batch combinations.
If no further feasible combinations of batches can be identified, an attempt is made to improve the obtained solution by applying the aforementioned basic local search function (compare Function 1). The process is repeated for each ant that is used. Each ant represents a solution. Afterwards, the pheromones for all order combinations are updated. First, in analogy to the processes in nature, a fraction of the pheromone evaporates. Likewise good solutions will receive an additional amount of pheromones. More precisely, the pheromone increase is calculated according to the ranking position r within the set of the m best solutions in the actual iteration. In addition, the best solution found so far in all iterations will be rewarded with an additional increase of the pheromone intensity. According to Bullnheimer, Hartl, and Strauss (1999) the pheromone intensity is updated as follows: is in a batch of the r-th best solution är ã 1; :::; må; äm à 1å ± ; if order combination (i; j) is in a batch of the best solution so far är ã 0å 0; otherwise; where : uprating parameter ( OE 0).
The updated pheromone intensities provide the basis for the subsequent iteration. The detailed description of the RBAS used here is given in the pseudo-code of Algorithm 2.
Algorithm 2 Rank-Based Ant System (RBAS) Parameters: number of iterations, number of ants, initialization of the ij , number of best solutions m, evaporation parameter , uprating parameter , control parameter , ; for it ã 0 to number of iterations do for a ã 0 to number of ants do create start solution s consisting of a single customer order per batch; repeat for all feasible pairs äk; lå of s do compute sav kl ; determine˜ kl ; end for determine p kl ; choose combination äk;lå according to (6); update s by combiningk andl; until no further combination is possible; s ñãlocal_search(s); update the best m solutions äs 1 ; : : : ; s m å of this iteration; update the best solution s 0 found so far; end for for all order combinations äi; jå do ij ñã ä1 á å ± ij ; end for for r ã 0 to m do for all batches k in s r do for all order combinations äi; jå in k do ij ñã ij à äm à 1 á rå ± ; end for end for end for end for s 0 is the heuristic solution; 6 Design of the Experiments

Warehouse Parameters
In order to evaluate the performance of the proposed heuristics in our numerical experiments we assume a single-block warehouse with two cross

BuR --Business Research Official Open Access Journal of VHB Verband der Hochschullehrer für Betriebswirtschaft e.V. Volume 3 | Issue 1 | May 10 | 82--105
aisles, one in the front and one in the back of the picking area. All aisles are vertically orientated and the depot is located in front of the leftmost aisle. This type of layout was depicted in Figure  1 and Figure 2 already and is identical to one which is frequently used in the literature (compare de Koster, van der Poort, and Wolters 1999, Petersen and Schmenner 1999). The picking area consists of 900 storage locations (cells) and we assume that a different article has been assigned to each storage location. The storage locations are partitioned into 10 picking aisles with 90 storage locations each, i.e. 45 cells on both sides of each aisle. The aisles are numbered from 1 to 10, where aisle no. 1 is the left-most aisle and aisle no. 10 the right-most aisle. Within an aisle two-sided picking is assumed, i.e. being positioned in the center of an aisle, the order picker can pick items from cells on the right as well as from the cells on the left without additional movements. The length of each cell amounts to one length unit (LU). When picking an article, the order picker is positioned in the middle of the cell. Whenever the order picker leaves an aisle, he/she has to move one LU in vertical direction from the first storage location, or from the last storage location respectively, in order to reach the cross aisle. For a transition into the next aisle the order picker has to move 5 LU in horizontal direction, i.e. the center-to-center distance between two aisles amounts to 5 LU. The depot is 1.5 LU away from the first storage location of the left-most aisle, i.e. the distance between cross aisle and depot amounts to 0.5 LU. We assume a class-based storage assignment of articles to storage locations, i.e. the articles are grouped into three classes A, B and C by their demand frequency, where A con- class-based (ABC-storage) tains articles with high, B with medium and C with low demand frequency. Articles of class A are only stored in aisle no. 1, articles of class B in aisles no. 2, no. 3 and no. 4, and articles of class C only in the remaining six aisles. Furthermore, we assume that 52ô of the demand belongs to articles in class A, 36ô to articles in B and 12ô to articles in C. Within a class, the location of each article is determined randomly. Table 1 summarizes all warehouse parameters fixed for the numerical experiments.

Problem Classes
In our numerical experiments we consider five different sizes of customer orders (20,30,40,50,60), where the number of articles per order is uniformly distributed in è5; : : : ; 25é. The capacity of the picking device (maximal number of articles that can be assigned to a batch) is set to 30, 45, 60 and 75 articles; in other words, a batch consists of 2, 3, 4 and 5 customer orders on average. The corresponding problem instances can be considered as typical for real-world applications which we are aware of from our own practical experience (fresh-food or deep-freeze warehouses); they are also in line with other real-world examples from the literature (compare de Koster, Roodbergen, and van Voorden 1999). In combination with the two routing strategies (S-Shape, Largest-Gap) we obtain 40 problem classes, and for each of these classes 40 instances are generated. In total, 1,600 instances have been considered, each with one trial per instance. An overview of the problem classes considered in the experiments is given in Table 2.  (2005) proposed an Iterative Descent approach (ID), which, in accordance with our approach, starts from an initial solution generated by means of the FCFS rule. Major differences can be identified with respect to the local search and the perturbation phase. The authors have implemented a first-improvement strategy in which only SWAPs are applied; i.e. only one type of neighborhood is considered in their approach. A new solution is taken as a new incumbent solution if the objective function value is better than that of the previous one. The perturbation phase consists of three random exchanges of three customer orders from three different batches. Like in the approach of Gademann and van de Velde (2005), in our experiments the number of calls of the local search and the perturbation phase has been limited to 25, too. Furthermore, Gademann and van de Velde (2005) use a different definition of the capacity of the picking device, which restricts the number of customer orders --and not the number of articles like in our approach --that can be served simultaneously by the picking device. This might lead to an infeasible solution, when a SWAP has to be carried out with two customer orders which consist of a different number of articles. In order to avoid such infeasibilities, their approach has been modified in a way that a new batch is opened. To this batch we assign a randomly selected order from the batch which exceeds the capacity constraint.
In order to evaluate the solution quality of the proposed approaches not only in relative terms (i.e. in relation to the solution quality of other existing methods), we have also tried to benchmark their solutions against the corresponding optimal objective function values. For this purpose, it has been attempted to generate the respective Integer Programming Problem (IPP) --as presented in Section 2 --explicitly for each test problem instance as far as it was possible. An optimal solution to it was generated by means of a commercial LP/IP solver. Not unexpectedly, this approach was successful only for a limited number of problem classes, but gives some insights in the absolute solution quality of the proposed methods, nevertheless. Tables 3 and 4 depict the considered problem classes and the average number of feasible batches (number of columns of the mathematical model) which have been generated. This shows the exponential increase of the problem size by an increasing number of customer orders and an increasing capacity of the picking device. Further, the number of in-  For abbreviations see Table 3.

BuR --Business Research Official Open Access Journal of VHB Verband der Hochschullehrer für Betriebswirtschaft e.V. Volume 3 | Issue 1 | May 10 | 82--105
stances is shown where the LP/IP solver was able to generate a solution and prove its optimality. In the remaining instances the LP/IP solver was only able to generate a solution but was unable to prove its optimality, since the computation violated memory restrictions of the used PC. For these instances the average optimality gap as calculated by the LP/IP solver is shown as well. The optimality gap is the difference between the lower bound (objective function value of the best remaining node) and the best integer objective function value found so far. These results demonstrate that the application of exact approaches is limited to Order Batching Problems of small sizes, whereas for larger sizes the application of heuristics is necessary.

Algorithm Parameters
The parameter settings have been determined in a series of pre-tests: For the RBAS we perform 100 iterations, and the number of ants in each iteration has been set to half the number of customer orders. The pheromone intensity of each order combination has been initialized by 2, and the control parameters and have been set to 5. The evaporation parameter has been fixed to 0.1 and the uprating parameter to 0.001. During the process of updating the pheromone intensity, we consider the three best solutions in each iteration, i.e. m ã 3.
In order to provide equivalent conditions for both suggested methods, we allow for an identical computing time, which is determined by RBAS. In other words, each problem instance is solved by RBAS and the resulting computing time is used as the termination condition for ILS. In ILS the rearrangement parameter has been set to 0:3. We further allow 10 deteriorations of maximal ã 0:05 each. In combination with the computing time, this results in t ã t RBAS =10, where t RBAS is the computing time of RBAS. Table 5 summarizes the parameter choice. 0.05 time interval t after which a deterioration is permitted: t RBAS = 10

Implementation and Hardware
The computations for all 1,600 instances have been carried out on a desktop PC with a Pentium processor with 2.21 GHz and 2 GB RAM. The algorithms have been encoded in C++ using the DEV Compiler Version 4.9.9.2. The solver used for the IPP was CPLEX 10.1.

S-Shape Routing Solution Quality
The results from the numerical experiments for S-Shape Routing are depicted in Table 6.

BuR --Business Research Official Open Access Journal of VHB Verband der Hochschullehrer für Betriebswirtschaft e.V. Volume 3 | Issue 1 | May 10 | 82--105
For all problem instances, application of CîW(ii) improved the total length of all picking tours on average by 16.9ô (in comparison to the length of the tours resulting from the application of FCFS only). Differentiated with respect to the problem classes, the (average) improvement ranged from 8.97ô (no. of customer orders: 20; capacity of the picking device: 75) to 21.33ô (n = 60; cap = 45). In general, it can be observed that only small improvements were obtained by CîW(ii) for problem classes characterized by a small number of customer orders. This can be explained by the fact that --for instances from these problem classes --there is just a small probability that the total tour length will be improved since the algorithm opens a new batch for each pair of unassigned orders. We note that our results are in line with those from de Koster, van der Poort, and Wolters (1999), who have run numerical experiments for a similar warehouse and obtained a reduction of the total travel time of the order pickers of 18ô (constant travel velocity assumed) for problem instances with 30 customer orders and a capacity of the picking device of 24 articles. Therefore, we conclude that our implementation of CîW(ii) is credible and that the results obtained by CîW(ii) provide a good benchmark. The addition of a local search phase to CîW(ii) provided even better solutions. CîW(ii) + LS raised the average improvement to 18.05ô, ranging from 11.59ô (n = 20; cap = 75) to 23.21ô (n = 60; cap = 45) for the different problem classes. However, these additional improvements should not be overestimated. The reduction of the total length of all tours achieved by CîW(ii)+LS in comparison to CîW(ii) amounted to only about one percentage point on average. This demonstrates that the solutions provided by CîW(ii) are already close to local minima. Hence, classic local search cannot succeed in improving these solutions significantly. On average, the length of all picking tours provided by the (modified) Iterative Descent method (ID) was by 10.42ô shorter than the length of the tours generated by means of FCFS. The results for the different problem classes varied between an increase of only 6.29ô (n = 20; cap = 75) and up to 14.47ô (n = 50; cap = 45). The results were inferior even to those of the basic savings method CîW(ii), which offers tours that are by 6.48 percentage points shorter on average. This result can be attributed to the fact that ID --as proposed by , the additional improvements amounted to 3.00ô for ILS on average, and 2.94ô for RBAS, respectively. This demonstrates that the solutions obtained by means of the presented metaheuristics were significantly better than those obtained with classic local search.

BuR --Business Research Official Open Access Journal of VHB Verband der Hochschullehrer für Betriebswirtschaft e.V. Volume 3 | Issue 1 | May 10 | 82--105
This result can be attributed to the fact that the metaheuristics were able to reduce the number of batches considerably, e.g., up to 1.4 batches for ILS and up to 1.5 for RBAS (n = 60; cap = 30). Table 7 relates the solution quality of ILS and RBAS to the corresponding optimal objective function values. For all problem classes, ILS and RBAS provided solutions close to the optimum. The average deviation from the optimal objective function value (related to the problem instances in each class which could be solved to a proven optimum) was smaller than 1ô (except for the problem classes (n = 40; cap = 60) and (n = 50; cap = 45), where a deviation of up to 2ô was observed). Average deviations from the optimal objective function value tended to be smaller for problems with a smaller capacity of the picking device than for those with a larger capacity (with the exception of the problem class with n = 20 and cap = 60). RBAS outperformed ILS for small problem instances (n = 20) and for larger instances with cap = 30. In brackets, Table 7 also depicts the average deviation of the total tour lengths obtained by ILS/RBAS from the total tour lengths provided by the LP/IP solver only for those instances in which the optimality of the solutions could not be proven. For some of the instances, where the LP/IP solver was not able to prove the optimality of the obtained solutions, ILS or RBAS even generated solutions with smaller objective function values. These cases are indicated by negative entries (e.g., n = 20; cap = 75). Table 8 depicts the average computing times per problem instance for ILS and RBAS. The times for FCFS, CîW(ii), CîW(ii)+LS and ID have not been listed since these methods consumed less than five seconds per instance. The average computing time for RBAS (which --as has been explained before -sets the limit for the computing time allocated to ILS) varied between 7 seconds (n = 20; cap = 30) and 1088 seconds (n = 60; cap = 75). Basically, two problem parameters determined the computing times, namely the size of the problem instances (i.e. the number of customer orders which has to be dealt with) and the capacity of the picking device. Increasing values of these parameters resulted in increasing computing times. While ILS and RBAS can be considered as equivalent with respect to solution quality, major differences become evident with respect to the computing time which passed (on average per problem instance) until the best solution was found. On average, across all problem instances, the best solution was found by RBAS after 61.3ô of the total computing time. For small problem instances (n = 20) it took 5.5ô (cap = 30) and 31.2ô (cap = 45), and for large instances (n = 60) 73.4ô (cap = 30) and 95.6ô (cap = 75) of the total computing time. ILS, on the other hand, found the best solution much earlier; on average, ILS needed 17.4ô of the total computing time, and even in the worst case only 30.9ô. RBAS is particularly time consuming because the savings have to be recalculated for each ant and, likewise, a new tour length has to be determined for each feasible pair of batches. As for ILS, the only time-consuming part consists of the local search phase, which is also included in RBAS. Given the identical computing times, ILS generates and evaluates more solutions than RBAS does. The numerical experiments have shown that --depending on the problem class --up to 50ô additional solutions can be considered. The development of the solution quality of the two heuristics as a function of the computing time is depicted in Figure 3 for the three problem classes (n = 30; cap = 60), (n = 50; cap = 30) and (n = 60; cap = 75). Each of the functions represents the development of the solution quality for one problem class. The x-axis depicts the computing time of the algorithms, with 100ô representing the termination time of the algorithms. On the y-axis it is represented how far the best solution which has been found until a certain point in time deviates from the best solution obtained after 100ô (final solution) of the computing time. For this representation the deviation from the objective function value of the final solution is measured in 20 intervals of identical size; each point in the graph represents the average deviation computed over all 40 instances in a problem class. The development for ILS is shown in the upper part, and for RBAS in the lower part of Figure 3. The graphs demonstrate that ILS finds the best solution much faster than RBAS does. The functions for ILS are always very similar to each other, independent from the problem class under consideration. In contrast to ILS, the functions for RBAS are different for every problem class. For the problem class with 50 customer orders and a picking device capacity of 30 articles (n = 50; cap = 30) the objective function values are very close to the best ones (deviation only 1ô) within the first 5ô of the allocated computing time. For the problem class (n = 30; cap = 60) already more than 30ô of the computing time is needed to reduce the deviation to 1ô. The problem class (n = 60; cap = 75) requires 30ô of the time to reach a deviation of 2ô and even 90ô to obtain a deviation of 1ô; in order to identify the best solution, the entire computing time needs to be exploited. This development re-

BuR --Business Research Official Open Access Journal of VHB Verband der Hochschullehrer für Betriebswirtschaft e.V. Volume 3 | Issue 1 | May 10 | 82--105
veals that in case of RBAS the size of the picking device has a major influence on the development of the solution quality over time. Table 9 provides information on the average computing times per problem instance for different phases of the optimization process when the LP/IP solver was applied. t gen represents the average time needed to set up the corresponding IPP, and t opt is the time spent from the moment when the optimization run has been started until the moment when the optimal solution has been found. t prov includes t opt plus the time necessary to prove the optimality of the solution. It has to be noted that only those instances are included in the table, for which the LP/IP solver was able to generate an optimal solution and verify its optimality. By comparing the computing times to those of Table 8, another strength of both ILS and RBAS becomes evident. The metaheuristics do not only provide high-quality solutions, they also do so in reasonable computing times, in particular for large problem instances. For problem classes with a small capacity of the picking device (cap = 30) the LP/IP solver is still faster than the metaheuristics for which it finds an optimal solution and proves its optimality in less than two seconds on average. For instances with (n = 20; cap = 45) and (n = 50; cap = 45) the LP/IP solver finds the optimal solution before ILS and RBAS terminate, but takes much additional time for proving its optimality. For all other problem classes the computing times of the metaheuristics are much shorter than the respective times needed by the LP/IP solver to even find the best solution. In general it can be observed that the computing times of the LP/IP solver increase exponentially while those of the metaheuristics exhibit an almost linear behavior. It has to be noted again, however, that the number of optimal solutions which could be determined for each problem class varies quite significantly. As a consequence, Table 9 only gives a biased picture of the actual average computing times and, thus, should be read carefully. For instance, considering the increase of t prov for a move from smaller capacities of the picking device to larger ones (i.e. moving from cap = 30 to cap = 45, cap = 60 and cap = 75) for small problem sizes (n = 20, 30), one would also expect t prov to increase for a transition from problem class (n = 40; cap = 45) to problem class (n = 40; cap = 60). The opposite is true, though, which must be attributed entirely to the fact that the LP/IP solver provided significantly fewer optimal solutions in the case of cap = 60 (16 instances) than in the case of cap = 45 (39). Finally, another interesting aspect is that for some problem classes (e.g. (n = 40; cap = 60)), the generation of the IPP only required almost as much computing time as did the application of the metaheuristics.

Largest-Gap Routing Solution Quality
The results from the numerical experiments for Largest-Gap Routing are presented in Table 10. Application of CîW(ii) provided solutions improving the objective function value of the solutions of FCFS by 16.17ô on average (for all problem instances). The average improvements for the different problem classes ranged from 11.08ô (no. of customer orders: 20; capacity of the picking device: 75) to 19.43ô (n = 60; cap = 60).   As has been shown for S-Shape Routing, it can also be observed here that a large problem size induced an increase of the improvements. Improvements tended to be larger for medium-size capacities of the picking device (cap = 45, 60) than for small (cap = 30) and large capacities (cap = 75). CîW(ii)+LS led to an average improvement of 17.48ô, with a minimum improvement of 13.72ô (n = 20; cap = 30) and a maximum increase of 20.86ô (n = 60; cap = 60). Similar to the S-Shape case, the improvements obtained with the additional classic local search phase are not very convincing: the average improvement amounted to only 1.3 percentage points compared to CîW(ii). Application of the (modified) ID resulted in an average improvement of only 12.17ô with respect to the length of all picking tours. The improvements ranged from 6.76ô (n = 20; cap = 30) to 16.13ô (n = 60; cap = 60). Again, these results are inferior compared to the ones of the basic savings approach. For Largest-Gap Routing the best solutions were provided by means of the newly proposed heuristics ILS and RBAS. Their application improved  Table 11, which depicts the average deviation of the obtained objective function values from the corresponding optimal values. Similar to the case of S-Shape Routing the average deviation is small (ILS: less than 0.85ô; RBAS: less than 0.74ô).   Figure 4, ILS (top) finds the best solutions much earlier than RBAS (bottom) does. Also, the functions of the different problem classes are similar for ILS and different for RBAS. For RBAS good solutions are found much faster for problem classes with smaller order sizes than for problem classes with larger ones. The comparison of the computing times needed for the metaheuristics and the times required for solving the IPP (Table 13) reveals the same result as already observed for the S-Shape case. For the small problem classes with cap = 30 the computing times for solving the IPP and for proving the optimality of the obtained solutions are neglectable, since they are smaller than two seconds. For the other problem classes, again an exponential increase of the computing times can be observed.

Comparison of Routing Strategies
Subject to the parameters of the warehouse and the composition of the customer order data, one routing strategy can be more favorable than the  other. Therefore, we continue our analysis with a comparison of the results obtained for the two routing strategies. Hall (1993) reports that it can be observed that for a large number of articles to be picked on one route S-Shape Routing results in shorter tour lengths, while Largest-Gap Routing should be preferred for a small number of articles within one route. The same results can be observed in our experiments. Independent of the batching method and the problem size (total number of orders), Largest-Gap Routing --on average --outperforms S-Shape Routing for small capacities of the picking device, whereas for high capacities the opposite is true. This can be explained by the fact that an increase in the capacity of the picking device will also increase the picking intensity per aisle (i.e. average number of locations of articles to be picked per aisle). As a consequence, the maximal distance between two adjacent locations of articles to be picked in the same aisle will be reduced. This has no impact on the behavior of the S-Shape Heuristic, but application of the Largest-Gap Heuristic will result in solutions with longer travel distances within the aisles. In a situation, where the order picker has to collect an article from each storage location in an aisle, he/she would have to enter the aisle twice, once from the upper cross-aisle and once from the lower cross-aisle. In each case, the order picker would proceed to the center of the aisle and return from there to the cross-aisle, again. This would provide a tour which is twice as long as the length of the aisle. When the heuristics in our experiment are applied, Largest-Gap Routing outperforms S-Shape Routing for small capacities of the picking device (cap = 30, 45). For a large capacity (cap = 75) the opposite is true. In the case that the picking device has a capacity of 60 articles, no specific precedence ordering of the routing strategies can be given with respect to the solution quality. Instead, which method has to be preferred depends on the size of the problem (number of customer orders).
As mentioned before this result is valid for all batching heuristics. Therefore, it can be said that the decision about the routing strategy should only be based on the characteristics of the warehouse and the customer orders but is independent from the batching heuristics. Nevertheless, our experiments have shown that a better batching heuristic reduces the difference between the lengths of the picker tours resulting from the application of S-Shape on the one hand, and Largest-Gap on the other. With respect to computing times, Largest-Gap Routing always requires more time than S-Shape Routing, due to the fact that one has to search for the largest gap between all picks within each aisle.

Conclusions and Outlook
In this article we considered the Order Batching Problem, a problem which is pivotal for the efficient management and control of manual pickerto-parts order picking systems in distribution warehouses. This problem is of major importance, since --due to the high amount of manual labor --order picking is the most cost-intensive function in a warehouse. We introduced two algorithms, based on Iterated Local Search and on Ant Colony Optimization, for the solution of this problem. By means of extensive numerical experiments it could be demonstrated that --in terms of solution quality --both ap- proaches are superior to existing methods. The results show further that the constructive heuristics as well as the basic local search approaches are not able to generate competitive solutions, whereas ILS and RBAS provide solutions which could improve the efficiency of distribution warehouse operations significantly. When compared to each other, both approaches do not differ very much with respect to the quality of the solutions they provide. However, ILS provides these solutions much faster, namely in only one third of the time required by RBAS. Also the development of the solution quality over time emphasizes how fast ILS provides high-quality solutions compared to RBAS. Additionally, ILS requires less parameter configuration than RBAS. Therefore, Iterated Local Search appears to be more suitable for being implemented in software systems for the solution of practical problems. In detail, the numerical experiments show, that -on average --the application of the presented metaheuristics can reduce the length of the necessary picking tours by more than 20ô in comparison to a FCFS solution. Compared to the best benchmark (CîW(ii) + Local Search) an improvement of up to three percentage points can be observed. Despite the fact, that the magnitude of the obtained improvements varies with respect to the problem size (number of orders), the characteristics of the used equipment (capacity of the picking device) and the chosen routing strategy, the results produced by RBAS and ILS are always superior to all considered benchmarks. Moreover, the experiments have shown that the achieved lengths of the picking tours are very close to the optimal ones (as far as it was possible to generate them). In addition to the superior solution quality it is important to point out that the metaheuristics are able to generate high-quality solutions in reasonable time, even for those problem classes which have proven to be difficult for the LP/IP solver (e.g., large number of orders and large capacity of the picking device).
Our research results are of practical importance with respect to two aspects: Firstly, the extraordinary reduction of the total picker tour lengths and the corresponding reduction of the total picking times establish the opportunity to reduce the regular working hours and/or overtime of the picking staff, and even downsize the workforce. This will turn out to be a powerful measure for improving the low profit margins still prevalent for warehouse operations and secure the existence of picking service providers in the long run. Secondly, reduced tour lengths result in shorter times for picking the customer orders. This gives rise to shorter delivery lead times and, more generally, to improved customer services which could be used by a warehouse operator as a central building block in a strategy for gaining a competitive advantage over its competitors. Future research might follow two roads: From a methodological point of view, research could concentrate on speeding up the suggested methods (e.g., by introducing alternative neighborhoods), or on studying alternative, probably faster metaheuristics. From a practical point of view it would also be desirable to analyze dynamic order batching problems where customer orders come in continuously and batching is carried out under a rolling planning horizon. Situations of this kind evoke additional questions, such as when and how often picking plans should be updated in order to deal with newly arrived customer orders. Moreover --as in many industries customer orders have to be fulfilled at a certain time --it will be important to differentiate the various customer orders with respect to due dates, which necessitates the incorporation of time windows in the analysis. Nevertheless, algorithms will still have to be implemented which solve problems of the kind discussed in this paper as subproblems. We conclude that our solution approaches will also be valuable in dealing with such dynamic situations.