Exact Solutions to the Symmetric and Asymmetric Vehicle Routing Problem with Simultaneous Delivery and Pick-Up

In reverse logistics networks, products (e.g., bottles or containers) have to be transported from a depot to customer locations and, after use, from customer locations back to the depot. In order to operate economically beneficial, companies prefer a simultaneous delivery and pick-up service. The resulting Vehicle Routing Problem with Simultaneous Delivery and Pick-up (VRPSDP) is an operational problem, which has to be solved daily by many companies. We present two mixed-integer linear model formulations for the VRPSDP, namely a vehicle-flow and a commodity-flow model. In order to strengthen the models, domain-reducing preprocessing techniques, and effective cutting planes are outlined. Symmetric benchmark instances known from the literature as well as new asymmetric instances derived from real-world problems are solved to optimality using CPLEX 12.1.


Introduction
The vehicle routing problem with simultaneous delivery and pick-up is an extension of the Capacitated Vehicle Routing Problem (CVRP) where products have to be transported from the depot to customer locations and other products have to be taken from customer locations back to the depot. Each customer requires two types of service á a delivery and a pick-up á and both activities have to be carried out simultaneously by the same vehicle. Traditionally, the objective is to find a set of routes that minimizes the transportation costs or the distances traveled, where the customer demands are satisfied and vehicle load capacities are not exceeded at any time.
The VRPSDP is an operational problem which typically occurs in reverse logistics. In addition to the distribution process to customers, used or obsolete products have to be transported in the opposite di-rection. Depending on the type of product and its previous usage, the reverse material flow may enter the supply chain of the original manufacturer or a different supply chain. In closed-loop supply chains, companies reintegrate returned products into their own processes in order to reap the maximum economic benefit from end-of-use products and to utilize all materials already integrated in the value added chain. Additionally, a reduction of waste and a decrease of natural resource extraction may be achieved. In order to handle the arising flow of used products, companies need to set up new logistics infrastructures, e.g., sorting and disassembly centers. Therefore, many authors consider network design and architecture issues. Less attention has currently been paid to tactical and operational aspects concerning routing and handling in reverse logistics networks. For that reason, this paper considers the VRPSDP, specifies an application area, and provides exact solutions to symmetric and asymmetric benchmark instances. The remainder of the paper is organized as follows: In section 2 we investigate an interesting practical example of the VRPSDP that occurs in closed-loop supply chains. Section 3 is devoted to a literature review on simultaneous delivery and pick-up problems. In section 4 two mathematical formulations for the VRPSDP are presented which can be posted to a standard solver (e.g., CPLEX). Based on the formulations, we proceed to describe methods for improving the quality of the models in terms of computation time and solution gap. The results of a comprehensive performance analysis, where symmetric and asymmetric instances are considered, are given in section 5. Finally, conclusions are presented in section 6.

Practical application and closed-loop supply chains
Vehicle routing with simultaneous delivery and pick-up typically arises in closed-loop supply chains. Figure 1 shows the structure of a general closed-loop supply chain, where raw materials enter the system at the manufacturing facility. After the production process, finished goods are delivered to customers or retailers (forward flow). From customer/retailer locations, a flow of used or obsolete products (e.g., spare parts or commodities which are unsold) runs to the place of recovery (reverse flow). There, inspection and sorting activities are carried out in order to differentiate between components that may be part of new production processes as well as recyclable materials and disposal. A closed-loop supply chain for security contain- Raw ers may be found in the field of data destruction. The constant evolution of data processing and the increase in volume of sensitive data has brought about the need for responsible service providers to destroy data in compliance with international standards. Documents (or disks, magnetic tapes, videos, and microfilms) must be deleted once the legal retention period is met. In order to ensure that data goes directly and unseen into the destruction plant, customers appoint security containers that provide protection from unauthorized access. Standard container capacities are 240 liters and 350 liters. Special vehicles are required to deliver empty containers to customers (e.g., companies, public bodies, or private households) and to pickup full containers at customers. The filled containers are depleted in a safe place (i.e., at the depot location). The documents are sorted into various paper types and then shredded as well as baled to reduce the amount of space required. The paper bales are temporarily stored in the depot and then transported to paper mills, where they are recycled into new paper products. The security containers themselves are cleaned after the depletion process and then reused directly at new customer locations. The material flows to and from customers are usually lower than the vehicle capacity, and the number of delivery and pick-up containers at customer locations generally differ from each other. Therefore, a VRPSDP has to be solved in order to generate a feasible set of routes for the vehicles. Since a data destruction service provider usually operates on a regional level (i.e., within a 100 km radius), varying numbers of urban and country customers have to be considered during the planning phase.

Literature review
In the last twenty years a number of articles appeared in the domains of reverse logistics and closed-loop supply chains (see e.g., the survey of Bostel, Dejax, and Lu 2005). However, most papers deal with strategic and tactical problems, since they directly affect the companies' profitability. For a comprehensive review on strategic problems, we refer to Akçal , Çetinkaya, and Üster (2009). Tactical applications can be found in Sheu, Chou, and Hu (2005), Gu and Ji (2008) or Novais (2009 and. Operational problems in closed-loop supply chains deal with vehicle routing and scheduling decisions. If customers require both a delivery and a pick-up service, companies have to decide whether to simultaneously pick-up returns and distribute new products. The trade-off occurs between the exploitation of existing vehicle capacities and the reduction of flexibility due to mixed pick-up and delivery operations (De Koster, De Brito, and Van De Vendel 2002). A survey on different variants of pick-up and delivery problems was presented in Parragh, Doerner, and Hartl (2008). For many companies (particularly for data destruction service providers), it is economically beneficial to simultaneously deliver and pick-up. Therefore, in what follows we will concentrate on articles that are related to the VRPSDP. As an extension of the CVRP, the VRPSDP is N P-hard in the strong sense (Garey and Johnson 1979;Toth and Vigo 2002). This fact and the necessity to solve real-world instances have pushed researchers to develop heuristic algorithms. In contrast, exact algorithms have received little attention in the literature. The VRPSDP is introduced by a case study that deals with a public library delivery and pick-up system in Ohio. To solve the problem, Min (1989) presented a solution approach which is analogous to a``cluster-first route-second'' method, i.e., the customers are initially clustered into groups and afterwards a traveling salesman problem is solved for each group. Similarly, Halse (1990) devised an improved``cluster-first route-second'' method. Dethloff (2001) used the``cheapest-insertion'' concept to solve the VRPSDP. Based on a chosen single customer route, customers are inserted into the current route according to three criteria; extra travel distance, capacity remaining on the route, and distance to the depot and back. Nagy and Salhi (2005) introduced the idea of``weakly feasible'' and``strongly feasible'' routes. If the length of a route does not exceed a maximum distance, the route is called weakly feasible. If, in addition to the weak feasibility, the total load on the routes is below the capacity of the vehicles, the route is called strongly feasible. The algorithm relies on a route construction procedure to provide an initial (weakly feasible) solution and on a variety of routines to improve the routes and to reach strong feasibility. Mitra (2005) considered the VRPSDP, where the delivery and pick-up quantities at customer locations are not limited to the capacity of the vehicles. Therefore, it may be possible that a customer must be visited more than once to fulfill the requirements. As a result, the next visits can be made by the same vehicle or by other ones. The author developed a mixed-integer linear program and a route construction heuristic. Most of the metaheuristics developed for the VRPSDP are based on pure or hybrid versions of the tabu search method. Montané and Galvão (2006) presented a tabu search algorithm where initial solutions are found using constructive heuristic procedures. The neighborhood is built on the moves relocation, interchange, crossover, and 2-opt. Chen and Wu (2006) introduced a hybrid heuristic based on record-to-record travel, tabu lists, and route improvement procedures such as 2-exchange, swap, 2-opt, and or-opt. Aǹ`i nsertion-based'' algorithm is used to generate initial solutions. Privé, Renaud, Boctor, and Laporte (2006) studied a problem that arises from the distribution of soft drinks and collection of recyclable containers of a Quebec-based company. The problem is modeled as a multi-product vehicle routing problem with a heterogeneous vehicle fleet, time windows, and simultaneous delivery and pick-up at customer locations, where the weight of the material collected is always less than the weight of the material delivered. The objective is the minimization of routing costs, lower revenues resulting from the sale of recyclable material. Three construction heuristics, one inspired from the``nearest-neighbor'' concept and two derived from petal algorithms (see, e.g., Ryan, Hjorring, and Glover 1993), are developed for the problem. A subsequent improvement method uses a combination of 3-opt and 2-interchange moves as well as route merges. Bianchessi and Righini (2007) presented and compared several heuristic algorithms for the VRPSDP, in particular greedy insertion, local search, and tabu search procedures. Wassan, Wassan, and Nagy (2008) designed a reactive tabu search metaheuristic where an initial solution is constructed using a modification of the sweep algorithm (Gillett and Miller 1974). The core of the algorithm is a fast feasibility check for neighboring moves (e.g., shift, swap, and reverse). Zachariadis, Tarantilis, and Kiranoudis (2009) proposed a hybrid framework based on two metaheuristics, namely tabu search and guided local search. The algorithm is aimed at achieving a vast exploration of the search space by escaping from local optima

BuR --Business Research Official Open Access Journal of VHB German Academic Association for Business Research (VHB) Volume 6 | Issue 1 | May 2013 | 77--92
and intensifying at promising solution areas. To guarantee the robustness of the method, the number of search parameters is kept to a minimum. Subramanian, Drummond, Bentes, Ochi, and Farias (2010) presented a multi-start heuristic that consecutively generates a solution by means of a greedy procedure. In order to improve this solution, a local search phase as well as some perturbation mechanisms are used. Recently, some population based metaheuristics were developed for the VRPSDP. Ai and Kachitvichyanukul (2009) proposed a particle swarm optimization procedure based on a solution representation which consists of a priority list of customers and its preferred vehicle. The particle is converted into the problem specific solution through a decoding procedure. Gajpal and Abad (2009) presented an ant colony system algorithm that uses a construction phase as well as three local search schemes based upon interchange and insertion operations as well as an exchange of sub-paths. An initial solution is constructed using the nearest neighbor heuristic. Zachariadis, Tarantilis, and Kiranoudis (2010) introduced an adaptive memory programming methodology for the VRPSDP. The proposed adaptive memory collects solution characteristics obtained through the search process. These characteristics are combined to produce new solutions that are subsequently improved by a tabu search method. Most of the proposed metaheuristics are tested on instances with up to 400 customers. Since the authors usually neglected the consideration of lower bounds, no solution gap values are determined. Exact solution procedures presented for the VRPSDP may be separated in: • set-partitioning or set-covering formulations containing a binary variable for every potential vehicle route as well as • two-index or three-index formulations containing binary variables indicating whether an arc in the underlying digraph is selected (by a specified vehicle) or not. Angelelli and Mansini (2002) considered the VRPSDP, where the delivery and pick-up at customer locations must be carried out within given time windows. They implemented a branch-and-price approach based on a set-covering formulation for the master problem.
A relaxation of the elementary shortest path problem with capacity constraints and time windows is used as pricing problem. Dell'Amico, Righini, and Salani (2006) examined the application of branch-and-price techniques, too. Thereby, the pricing problem is solved using dynamic programming. In order to improve the effectiveness of the pricing algorithm, the authors considered two new strategies: bidirectional search and bounded number of steps. Furthermore, different branching strategies are implemented, namely branching on arcs, branching on resources, and branching on cycles.
Both Mitra (2005) as well as Chen and Wu (2006) used two-/three-index formulations in combination with the standard solver CPLEX to validate results obtained by heuristic algorithms. Subramanian, Uchoa, and Ochi (2010) presented an undirected and a directed two-commodity-flow formulation that make use of variables indicating the pick-up, the delivery, as well as the simultaneous pick-up and delivery flows. The models are embedded into a branch-and-cut approach containing the``CVRPSEP package'' (Lysgaard 2003). Furthermore, Subramanian, Uchoa, Pessoa, and Ochi (2011) proposed a branch-and-cut method over a symmetric formulation with only edge variables.
All exact methods introduced in this section are constructed to solve symmetric problem instances. For the asymmetric case, only the one-commodityflow model presented by Dell'Amico, Righini, and Salani (2006) and the directed two-commodityflow model presented by Subramanian, Uchoa, and Ochi (2010) can be considered. Due to the lack of exact solution algorithms designed for asymmetric instances, we propose below exact solution procedures for both the symmetric and the asymmetric VRPSDP based on branch-and-cut.

Model formulations
In this section, all relevant notation and terminology used throughout the paper are introduced (subsection 4.1). Furthermore, we present two model formulations for the VRPSDP, namely a vehicle-flow formulation (subsection 4.2) and a commodity-flow formulation (subsection 4.4), which can be posted to a standard solver. Subsections 4.3, 4.5, and 4.6 are devoted to efficient modeling techniques that can help to ensure better

BuR --Business Research Official Open Access Journal of VHB German Academic Association for Business Research (VHB) Volume 6 | Issue 1 | May 2013 | 77--92
performance results.

General model components
The VRPSDP is defined on a complete digraph G ã äV ; Aå, where V ã è0; : : : ; né ã C ∪ è0é is the node set and A ã è i; j | i; j ∈ V é is the arc set. Each node i ∈ C represents a customer, while node 0 corresponds to the depot. A non-negative weight c ij is associated with each arc i; j ∈ A that represents the transportation costs between nodes i and j. We assume that the cost matrix satisfies the triangle inequality; this condition is always fulfilled in practical situations. A set of K identical vehicles, each with capacity Q, is available at the depot. The capacity of the vehicles cannot be exceeded in any route. Each customer i ∈ C is associated with a demand d i OE 0 that should be delivered and a demand p i OE 0 that should be picked-up, where d i à p i > 0 and d i ; p i ‹ Q. The demands of customers are measured in abstract transport units, which can be calculated from the dimension or weight of the individual product. For the depot, we set d 0 ñã p 0 ñã 0. Both the vehicle-flow and the commodity-flow formulation contain binary variables indicating whether an arc in the underlying digraph G is selected or not, i.e., they make use of decision variables (1) x ij ñã 1; if arc i; j belongs to the solution 0; otherwise with i; j ∈ V :

Vehicle-flow model
A vehicle-flow model describes the vehicle's motion and filling level in the road network. In order to construct a suitable vehicle-flow model, we therefore use the idea of``Resource Extension Functions'' that determine the updates of resources along an arc in the underlying network (Irnich 2008). With binary variables x ij (1) and auxiliary variables: l i OE 0 amount of load after visiting customer i ∈ C ld i OE 0 amount of load that has to be delivered to node i ∈ V and to all other following nodes the vehicle-flow model for the VRPSDP has the form: Objective function (2) represents the transportation costs which are to be minimized. The indegree and outdegree constraints (3) and (4) ensure that each customer is visited exactly once by a vehicle. Inequality (5) states that no more than K routes are created. With constraints (6), we specify the delivery quantity that has to be loaded at the depot. Additionally, constraints (6) force an order for customer visits in the routes, which ensure that no subtours without the depot are generated. Inequalities (7) and (8) indicate the amount of load in the vehicles after the visit of the first customer and the other customers in the routes, respectively. Constraints (6) and (8) are disjunctive constraints that are linearized by using large multipliers (``big-M values''). In order to create valid inequalities, we set M 1 = M 2 ñã Q. Constraints (9) and (10) guarantee that the capacity of the vehicles cannot be exceeded.

Efficient modeling techniques for the vehicle-flow model
The vehicle-flow model contains | V | 2 binary as well as | C | à | V | real auxiliary variables. Furthermore, the model has no special structural characteristics and, therefore, as in generic binary models, only instances with a small, or medium, number of customers can be solved to optimality using a standard solver (e.g., CPLEX). Since

BuR --Business Research Official Open Access Journal of VHB German Academic Association for Business Research (VHB) Volume 6 | Issue 1 | May 2013 | 77--92
the run times are quite large for medium-scale instances, we want to give the solver supplementary knowledge such as preprocessing information and additional constraints. Under preliminary tests, we have realized that the following modeling techniques result in significant improvements in terms of computation time and memory space. In order to facilitate the solution process, it is appropriate to restrict the domains of auxiliary variables. Considering model (2)á(11), the bounds of auxiliary variables l i and ld i ; i ∈ C; can be strengthened by investigating the differences between delivery and pick-up quantities.  Figure 2 shows a vehicle route with two customers i and j. The solid arcs are weighted with the quantities that flow along the corresponding road segment. The quantities given in the boxes indicate the maximum amount of load after visiting customer i and before visiting customer j. The latter quantity is always larger than, or equal to, the quantity that has to be delivered to j and to all other following nodes. Using these information, we obtain the conditions: The vehicle-flow model involves two big-M values, which have to be chosen smallest possible while maintaining feasibility of all potentially optimal solutions. We therefore replace the general``big-M values'' (M 1 and M 2 ) by specified M ij 1 and M ij 2 for the corresponding node pairs èi; jé. In order to ensure that inequalities (6) are redundant, if x ij ã 0, M ij 1 must be larger than, or equal to, variable ld j . Thus, we are able to utilize conditions (13) and replace the value M 1 by Constraints (8) are always fulfilled, if we choose M ij 2 as larger than, or equal to, the difference between l i and d j . Using conditions (12), the``big-M value'' M 2 can then be replaced by In addition to the auxiliary variables used in the vehicle-flow model, information on pick-up loadings can be collected and exploited. To that end, we introduce new auxiliary variables that are used to formulate the following additional constraints: With inequalities (14), we determine the pick-up quantity that has to be offloaded at the depot. Constraints (15) guarantee that the capacity of the vehicles is not exceeded, and constraints (16) are structural equalities that help to correlate the load variables with each other. In general case, the``big-M value'' in disjunctive constraints (14) can be set to M 3 ñã Q. However, to choose``big-M'' as small as possible, constants M ij 3 are introduced that are larger than, or equal to, variable lp i . Since the amount of pick-up load after visiting customer i is always larger than, or equal to, the total quantity of products after visiting customer i, conditions (12) can be used to formulate:

Commodity-flow model
In addition to binary variables x ij (1), the commodity-flow model makes use of two additional sets of continuous variables which represent the amounts of demand that flow along the associated arcs (Dell'Amico, Righini, and Salani 2006). We identify the decision variable y ij OE 0 with the amount of delivery load carried along arc i; j and the decision variable z ij OE 0 with the amount of collected load carried along i; j , i; j ∈ V : Then,

BuR --Business Research Official Open Access Journal of VHB German Academic Association for Business Research (VHB) Volume 6 | Issue 1 | May 2013 | 77--92
the mixed-integer linear program for the VRPSDP can be formulated as follows: As previously reported, objective function (17) expresses the transportation costs. Constraints (18) ensure that each customer is reached precisely once. Constraints (19) characterize the flow on the path to be followed by a vehicle. Inequality (20) guarantees that no more than K routes are created. Constraints (21) and (22) are the balance equations to satisfy the delivery and pick-up demand at customer j, respectively. It is easy to see that inequalities (21) and (22) eliminate subtours, since they sequence customers in the routes beginning and ending at the depot (source and sink of products). Disjunctive constraints (23) guarantee that the capacity of the vehicles is not exceeded in any route; the``big-M value'' can be set to M 4 ñã Q.

4.5
Efficient modeling techniques for the commodity-flow model The commodity-flow formulation requires | V | 2 binary variables as well as 2 | V | 2 continuous variables. In comparison with the vehicle-flow model, the number of continuous variables is relatively high. However, the formulation uses only one set of disjunctive (``big-M'') constraints. Generally, the reduction or elimination of big-M constraints improves the linear programming lower bounds of an integer program. We will investigate this statement in subsection 5.2. In order to restrict the domains of continues variables in the model, Subramanian, Uchoa, and Ochi (2010) replaced constraints (24) by the following conditions: Constraints (26) ensure that the quantity on a used road segment i; j is at least equal to the delivery demand of customer j. Since node i is visited before j, the maximum quantity that flows along the arc is less than, or equal to, the difference between capacity Q and d i . Considering the pickup demands of nodes i and j, constraints (27) can be specified in an analogous manner. Disjunctive constraints (23) are linearized by using a big-M value. The general value M 4 can be replaced by specified M ij 4 if the differences between delivery and pick-up quantities are considered (see also Fig. 2). We then obtain Under preliminary tests, we have determined that the consideration of additional constraints improves the solver solution process. The following inequalities: set lower bounds on the total quantities of products entering and leaving a customer location.

4.6
Valid inequalities for both models In our computational experiments, each VRPSDP problem instance is solved using CPLEX 12.1 and the Concert interface for communications with the solver. CPLEX uses a branch-and-cut approach, and the features of CPLEX allow general cuts to be generated during optimization. We identified clique, implied-bound, and flow-cover cuts as particularly suitable for the vehicle-flow model. For the commodity-flow model, the integration of implied-bound, flow-cover, and flow-path cuts improves the results significantly. Furthermore, we

BuR --Business Research Official Open Access Journal of VHB German Academic Association for Business Research (VHB) Volume 6 | Issue 1 | May 2013 | 77--92
allow a priority branching on special ordered sets of type one. Models and cutting planes are implemented using C++ code compiled with MS Visual Studio 2008. In order to further strengthen both models, we include rounded capacity cuts which are special inequalities for the symmetric CVRP. To insert rounded capacity cuts, we use an adjusted variant of the``CVRPSEP package'' implemented by Lysgaard (2003). If we transfer the linear programming solution x to the separation algorithm, we substitute the sum äx ij àx ji ) for the decision variable x äi;jå (specifying if edge äi; jå is in the solution). The rounded capacity inequalities i∈S j∈V \S impose the vehicle capacity restrictions, where räSå is the minimum number of vehicles needed to serve a customer set S ⊆ C. For the VRPSDP, we obtain räSå ã maxè i∈S d i ; i∈S p i é=Q .
The separation problem of rounded capacity cuts is known to be N P-hard. For the computational tractability, we only add rounded capacity cuts at the first 400 nodes of the branch-and-bound tree. Since all decision and auxiliary variables are bounded, the convex hull H of the set of feasible solutions to the VRPSDP has either the form of a convex polytope or it is empty. This implies that there always exists an optimal solution if H ã ∅ is satisfied.

Computational results
The computational tests have been performed on benchmark instances known from the literature as well as on instances that are derived from realworld data (cf. the application in section 2). In total, we used 206 symmetric and asymmetric test instances to evaluate the performance of our two models. During a series of preliminary tests, we have investigated the effectiveness of the directed twocommodity-flow model presented by Subramanian, Uchoa, and Ochi (2010). A comparison of run times showed that the model formulation is indeed able to cope with small-scale asymmetric instances that are, in general, easy to solve. For more difficult instances, our computing environment was not able to complete the enumeration within the time limit. In the remaining part of this paper, we therefore exclude the two-commodityflow model from the computational study. We start off by describing the composition and generation of the problem instances used for testing the models (subsection 5.1). In subsection 5.2 a comprehensive performance study is presented.

Benchmark instances
In order to investigate the increases in computation time caused by different problem structures, we subdivided the benchmark set in different classes. Particularly, class 1 contains asymmetric real-world instances generated in accordance with our cooperation partner (a data destruction service provider). Classes 2 to 5 cover symmetric instances that were proposed in the literature. Finally, classes 6 and 7 incorporate asymmetric instances constructed on the basis of benchmark sets from the literature, which are usually considered to evaluate heuristic solution methods. Class 1 is the largest test set and it includes instances that are generated using real-world data. Thereby, the customers are positioned in a special geographic area in the North of Germany (region around Hannover with a radius of 100 km). Each instance contains urban customers in and near Hannover and country customers usually placed in clusters. Transportation costs are determined by the Euclidean distance (rounded downward to the second decimal place) multiplied by a route factor for the different route segments (motorway, highway, and street). Delivery and pick-up quantities are measured in abstract transport units, where the different filling levels of security containers are considered, i.e., we set d i ∈ è0; : : : ; 20é and p i ∈ è0; : : : ; 66é for all customers i ∈ C. Furthermore, each instance involves two to six vehicles with capacity Q ∈ è100; 120é. In the second class, a subset of problem instances introduced by Mitra (2005) is considered, where c ij ∈ è10; : : : ; 28é for all i; j ∈ V ; c ii ã ‰; d i ; p i ∈ è1; 5é for all i ∈ C, and Q ã 10 transport units. Class 3 contains the benchmark of Chen and Wu (2006) that is derived from the extended Solomon instances R121, R141, R161, R181, and R1101 (Homberger and Gehring 2002). The benchmark is composed of 25 instances with 15, 17, and 20 customers. Transportation costs are determined by the Euclidean distance rounded to

BuR --Business Research Official Open Access Journal of VHB German Academic Association for Business Research (VHB) Volume 6 | Issue 1 | May 2013 | 77--92
the third decimal place. In order to incorporate simultaneous deliveries and pick-ups, a ratio r i ñã minèäx i =y i å; äy i =x i åé is calculated for each customer i ∈ C; x i ; y i are the x-and y-coordinate values of i. The integer delivery demand of customer i is equal to the rounded solution of d i ñã r i i , where i is the original demand of the underlying problem instance and r i is rounded to the second decimal place. The pick-up demand p i is set to p i ñã i á d i . Classes 4 and 5 are introduced by Dell'Amico, Righini, and Salani (2006). Both benchmarks consist in equal parts of instances with 20 and 40 customers. Transportation costs are given by the Euclidean distance rounded downward to the next integer value. The numbers of vehicles are assumed to be fixed; the numbers come from the solution of a multiple bin packing problem. Class 4 contains instances that are derived from the Solomon (1987) instances C101, R101, and RC101, where the delivery quantity d i equals the original demand i of customer i ∈ C. The pick-up quantity p i is calculated by p i ñã ä1 á åd i if i ∈ C is even, and by p i ñã ä1 à åd i if i ∈ C is odd; ∈ è0:2; 0:8é. Class 5 is obtained from the LIBrary of capacitated Vehicle Routing Problems (VRPLIB) 1 . The original demands in the instances are considered as delivery demands. Furthermore, the pick-up demand of each customer i is computed as p i ã ä0:5 à råd i , where r is taken from an uniform distribution in the interval ae0; 1ç. Classes 6 and 7 contain instances with 50 customers that were introduced by Dethloff (2001) as well as Salhi and Nagy (1999). Previous heuristic results for these instances are based on nonrounded input data. Since the benchmarks are constructed using well-conceived rules, we used them as a basis for generating asymmetric instances. In order to avoid numerical problems during presolve and cut generation, we decided to use a rounding scheme for the transportation costs, the vehicle capacity, as well as the pick-up and delivery quantities. Due to the rounding scheme, the route composition obtained by aforementioned heuristics and our solution procedure varies widely. For example, quantities of 9:6; 5:6, and 4:6 units can be combined in a vehicle with a capacity of 20 units, whereas rounded quantities 9:6 š 9:7; 5:6 š 5:7, and 4:6 š 4:7 cannot be combined. Further, we 1 The VRPLIB is available at http://elib.zib.de/pub/Packages/ mp-testdata/vehicle-rout/vrplib/index.html. introduced a matrix Λ ã ä ij å; i; j ∈ V ; to modify all instances of classes 6 and 7; ij ∈ [1:0; 1:5] for i ã 0; j ∈ C; i; j ∈ C; i < j and ij ñã 1:0 otherwise. Then, an asymmetric cost matrix Γ ã ä ij å; i; j ∈ V ; is obtained by ij ñã c ij ij , where c ij is the Euclidean distance between pairs of nodes rounded to the fourth decimal place. The costs that do not satisfy the triangle inequality are appropriately decreased by using the Floyd-Warshall algorithm. We calculate the delivery and pick-up demands with a similar technique as used in class 3, i.e., d i ñã r i i ; p i ñã i á d i and both values as well as the capacity of the vehicles are rounded to the second decimal place. It should be pointed out that the specification of the exact number of digits after the decimal point is crucial in order to use the results for experimental performance comparison of different solution procedures. If the number of digits is not specified, rounding differences in the objective function values as well as different route compositions or customer sequences can occur and may lead to a misinterpretation of the optimality. The consideration of a rounding scheme is non-critical for real-world VRPSDP instances, because the number of digits can individually be adjusted to the underlying problem.

5.2
Performance study The following section will give an overview of and a discussion on the results obtained for our two model formulations. All tests are performed on an Intel 6-core processor, 3.46 GHz, 24 GB RAM, running Windows 7. We begin our analysis with the first test set which consists of real-world asymmetric instances. In order to show the effects of modeling techniques described in subsections 4.3, 4.5, and 4.6, we solved each instance with 20 and 30 customers six times: • test run m 1 -I: using the vehicle-flow model (m 1 ).
• test run m 1 -II: using the vehicle-flow model (m 1 ) and valid inequalities given in subsection 4.6.
• test run m 1 -III: using the vehicle-flow model (m 1 ), valid inequalities, and efficient modeling techniques given in subsection 4.3.
• test run m 2 -I: using the commodity-flow model (m 2 ). • test run m 2 -II: using the commodity-flow model (m 2 ) and valid inequalities.

BuR --Business Research Official Open Access Journal of VHB German Academic Association for Business Research (VHB) Volume 6 | Issue 1 | May 2013 | 77--92
• test run m 2 -III: using the commodity-flow model (m 2 ), valid inequalities, and efficient modeling techniques given in subsection 4.5. Table 1 shows the average computation times for small-scale instances with up to 30 customers. In column``K=Q'' the number of vehicles used with respective capacity and in column``no'' the total number of instances with the same |C|=K=Q values is given. Since the considered VRPSDP is an N P-hard optimization problem, we cannot expect that a branchand-cut approach will terminate within a reasonable time limit. That is why we allow a maximum computation time of twelve hours after which the best solution found is returned. In column``íopt'' the number of optimal solutions found and in column``gap'' the gap (in percentage) between the best integer objective and the current lower bound of all active nodes is given. Column``t cpu '' displays the average computation time (in seconds) for all optimally solved instances. The results highlight that the average computation times depend on the number of customers and on the number of vehicles. Furthermore, the consideration of valid inequalities is highly important for the performance of the models. Without valid inequalities, only 18 (39) out of 40 instances could be solved to optimality using the vehicleflow (commodity-flow) model. In addition, test runs m i -II and m i -III, i ã 1; 2; are efficient. All instances are solved to optimality and columǹ`t cpu '' shows the average CPU time over the respective ten instances. The resulting average run times are lower than six minutes, but it is apparent that the commodity-flow model outperforms the vehicle-flow model. Considering the computation time differences between related test runs (m 1 -II and m 1 -III as well as m 2 -II and m 2 -III), the additional value of modeling techniques cannot be clearly realized. However, this result will change by analyzing medium-scale instances with more than 30 customers. In order to improve the bounding accuracy and to speed up the solver, particularly for medium-scale instances, an upper bound on the objective function can be supplied. We used the solution obtained after passing 10,000 iterations of a simplified version of the multi-start procedure introduced in Rieck and Zimmermann (2009) to construct an upper bound and a start solution. Then, to investigate the effects of modeling techniques in combination with start solutions, we solved each instance with 40, 50, and 60 customers eight times. Thereby, we used the test runs: m i -I, m i -II, m i -II-start, and m i -III-start, i ã 1; 2, where the suffix``start'' imply that a start solution is inserted. Table 2 summarizes the results for the different test runs. As expected, the solution process without valid inequalities is characterized by many instances that could not be solved to optimality. In particular, the vehicle-flow model is never able to terminate the enumeration. If we concentrate on test runs m i -III-start, i ã 1; 2, the commodity-flow model produces the best results in terms of computation time and solution gap. All instances with 40 customers and most of the larger instances are solved to optimality. Additionally, the positive performance effect of the described modeling techniques (subsections 4.3 and 4.5) can now be realized in its whole extent. In test run m 2 -IIstart, for example, only 46 instances are solved to optimality, whereas in test run m 2 -III-start 51 optimal solutions are generated. With additional constraints, the value of the linear programming relaxation at the root node increases on average by 0.1ô for the vehicle-flow model and by 0.2ô for the commodity-flow model.

BuR --Business Research Official Open Access Journal of VHB German Academic Association for Business Research (VHB) Volume 6 | Issue 1 | May 2013 | 77--92
Since the results of test runs m i -I, m i -II, and m i -II-start, i ã 1; 2, are sobering, we skip pursuing further details for them. We continue our analysis using small-scale symmetric problem instances with up to 22 customers. Table 3 shows the optimal solutions and run times. The instances are denoted by the same descriptors that are used in the literature. In column``costs'' the lowest transportation costs (or distances traveled) and in column``t cpu '' the corresponding CPU time (in seconds) are given. With column``model'' we demonstrate the model formulation that was used to create the best run time; m 1 is taken for the vehicle-flow model (test run m 1 -III-start) and m 2 for the commodity-flow model (test run m 2 -III-start). If the costs are written in bold type, the underlying instance is solved to optimality and the optimality is proven for the first time. The optimal solution for the Min (1989) problem instance was found and verified by Montané and Galvão (2006) within two hours. Thereby, the authors compared the heuristic solution to a lower bound obtained by CPLEX. With our approach, the optimality could be proven in two seconds. Class 2 contains instances that are generated using deterministic rules. All instances are solved to optimality with low computational effort. Instances that are derived from the Solomon benchmark (here classes 3 and 4) stand out with the composition of geographical data. In group``C'' the customers are placed in clusters, in group``R'' the geographical data is randomly generated, and group``RC'' contains a mix of clustered and randomly generated customer locations. For these instances the vehicle-flow model performs reasonably well. We are able to solve five problem instances to optimality for the first time. Class 5 contains vehicle routing instances from real-world projects of the University of Cologne. In contrast to the other classes, each instance considers a large number of vehicles with high capacity. All instances with 20 customers are solved to optimality, but for 3C_20_80_1 the number K of vehicles used has to be increased to 11 to obtain the result indicated in the literature. Additionally, we discover that the optimal objective function value of instance 3C_20_50_2 (written in italic type) is lower than those published in the literature. The deviation of 10.5ô might occur from insufficient rounding during the optimization process in Dell'Amico, Righini, and Salani (2006). The Min (1989) problem is obtained from a real-world data set. It is aimed at studying the effects and capabilities of the solution procedure in a more realistic way. Here, the procedure based on the vehicle-flow model provides the fastest run time. The best results for instances considered in Table 3 are created with the vehicle-flow model m 1 . In spite of the large number of``big-M constraints'', m 1 seems to be the best choice to solve small-scale symmetric instances of the VRPSDP. Tables 4 and 5 show the results concerning run times and solution gap for the remaining 57 instances with 40 customers (classes 4 and 5) and with 50 customers (classes 6 and 7). If the optimal solution is not reached within the time limit, column``costs'' shows the lowest transportation costs of a feasible solution (best integer solution) and  Subramanian, Uchoa, and Ochi (2010) or Subramanian, Uchoa, Pessoa, and Ochi (2011). In order to construct asymmetric instances for the VRPSDP, the instances are modified (subsection 5.1). Table 5 depicts the corresponding run times and solution gaps for the asymmetric benchmark sets under consideration. For all instances,  the vehicle-flow model provides the best solution quality. All instances with three or four vehicles are solved to optimality. For instances with nine or ten vehicles the enumeration could not be completed within the time limit. This means that the number of vehicles is a significant parameter in order to identify the difficulty of the underlying problem. The analytical study of both models shows that the CPU time needed to obtain an optimal solution depends on the number of customers, the number of vehicles, and the structure of the considered instances. The approach based on the vehicle-flow formulation is effective in finding optimal solutions for small-scale symmetric instances, and in finding good upper bounds for asymmetric instances that contain either clustered customers or randomly generated customers (instances of type``C'' or``R''; see classes 6 and 7). In contrast, the approach based on the commodity-flow model turns out to be more robust in generating upper bounds as well as optimal solutions for medium-scale symmetric problem instances. Moreover, model m 2 performs well for real-world asymmetric instances that contain a mix of clustered and randomly generated customer locations (instances of type``RC''; see class 1). Table 6 shows the average numbers of cuts added during the optimization process. Thereby, the cuts of the individual models, i.e., m 1 and m 2 , are separated into different columns. We constitute the rounded capacity cuts as``user cuts'' just like CPLEX does.

BuR --Business Research Official Open Access Journal of VHB German Academic Association for Business Research (VHB) Volume 6 | Issue 1 | May 2013 | 77--92
Considering the average numbers of cuts, it is obvious that more cuts are generated for the vehicleflow model than for the commodity-flow model. For the vehicle-flow model, user cuts and impliedbound cuts are particularly created in the subproblems of the branch-and-bound tree. Rounded capacity cuts are known to be effective for routing problems (see, e.g., Fukasawa, Longo, Lysgaard, De Aragao, Reis, Uchoa et al. 2006) and impliedbound cuts are useful if the binary variables imply bounds on continuous variables (relevant for load constraints (6)á(8), and (14)á(16)). A clique can be introduced if at most one variable in the group of binary variables can be positive in any integer feasible solution (constraints (3) and (4)). For the commodity-flow model, mainly user cuts are generated. Furthermore, implied-bound cuts are only added for symmetric instances (classes 2 to 5), this means that inequalities (23), (26), and (27) are stronger for asymmetric instances. Flow-cover and flow-path cuts treat the continuous variables as nodes and apply binary variables to indicate the flow entering and leaving a node.

Conclusion
In this paper, we have presented the vehicle routing with simultaneous delivery and pick-up that typically occurs in closed-loop supply chains. In order to solve the VRPSDP, a vehicle-flow and a commodity-flow formulation are proposed.
Considering the models, optimal solutions for small-scale and medium-scale problem instances can efficiently be generated by a standard solver (e.g., CPLEX). To facilitate the solution process, we add problem-specific preprocessing techniques and cutting planes. The computational tests have

BuR --Business Research Official Open Access Journal of VHB German Academic Association for Business Research (VHB) Volume 6 | Issue 1 | May 2013 | 77--92
been performed on symmetric and asymmetric benchmark instances with up to 60 customers. Most instances are solved to optimality and for the remaining instances good upper bounds could be generated. Future work will include the explicit consideration of multiple objective functions in the mathematical models. In addition, combining the VRPSDP with the problem of optimizing transport flows to solve the joint problem of determining the optimal connections between pairs of nodes and finding the optimal set of vehicle routes may be an interesting issue.