Mixed Integer Linear Programming Models to Solve a Real-Life Vehicle Routing Problem with Pickup and Delivery

: This paper presents multiple readings to solve a vehicle routing problem with pickup and delivery (VRPPD) based on a real-life case study. Compared to theoretical problems, real-life ones are more difﬁcult to address due to their richness and complexity. To handle multiple points of view in modeling our problem, we developed three different Mixed Integer Linear Programming (MILP) models, where each model covers particular constraints. The suggested models are designed for a mega poultry company in Tunisia, called CHAHIA. Our mission was to develop a prototype for CHAHIA that helps decision-makers ﬁnd the best path for simultaneously delivering the company’s products and collecting the empty boxes. Based on data provided by CHAHIA, we conducted computational experiments, which have shown interesting and promising results.


Introduction
Transportation studies can be considered as an intersection between several domains, methods, and techniques to propose, solve, and develop solutions for real-life problems. Artificial intelligence [1][2][3][4] and optimization [5][6][7] stand behind the majority of developed solutions. Among transportation and mobility problems, we recall the Vehicle Routing Problem (VRP). VRP is a well-known problem studied in Operational Research and Computational and Decision sciences. It consists of routing products from depots to customers by drawing adequate distribution plans. The routing plans, also called circuits or tours, need to respect multiple constraints, such as vehicle loading capacities and customer demands. The main objective of VRP is to provide significant savings in terms of transportation costs while ensuring the lowest delivery delays. A variant of VRP is called VRP with Pickup and Delivery (VRPPD). VRPPD joins the processes of collecting and delivering either simultaneously or separately. To deal with a real-life VRPPD, several requirements need to be considered, such as time window restriction, travel cost, vehicle capacity, heterogeneous vehicle fleet, vehicle travel time, vehicle speed, multi-dimensional capacity, route restrictions, and uncertain decisions to choose between picking up first and delivering second, or the opposite, or even performing both simultaneously. In addition, there might be different objective functions that could be minimized, such as those related to travel cost, satisfying customer requirements, loss of funds, etc. This study is motivated by an industrial context, which is a problem related to a mega poultry company in Tunisia called CHAHIA. CHAHIA supplies a large number of

Related Literature
The VRPPD could be classified into three different categories [8,9], which are as follows: • Delivery first and pickup second: vehicles perform a pickup operation after the delivery process. This category is called Vehicle Routing Problems with Backhauls (VRPB). • Mixed pickup and delivery: vehicles deliver or pickup in any sequence along their routes. This category is called Mixed Vehicle Routing Problem with Backhauls (MVRPB). • Simultaneous Pickup and delivery: vehicles simultaneously perform the delivery and pickup. This category is called Vehicle Routing Problems with Simultaneous Pickup and Delivery (VRPSPD).
The majority of published VRPPD models assume that each customer is visited only once [10]. Some exceptions have considered multi-customer visits, such as in [11] where the authors developed a large-scale model, allowing two visits per customer for a many-tomany pickup and delivery routing problem. In [12,13], the authors considered a VRPSPD involving the delivery of full bottles and the collection of empty ones. The authors assumed that such a problem is classified among reverse logistics operations.
Several real-life applications encountered in the beverage industry are described in [14]. For example, in [15], authors considered a truck scheduling problem for container transportation in a local area with multiple depots and terminals. They proposed an approach based on an integer programming heuristic that determines pickup and delivery sequences for daily drayage operations while ensuring a minimum transportation cost. In the same context, the authors in [16] demonstrated that drayage operations could be considered as a multi-stop VRPB. More details regarding VRPB and VRPPD can be found in two comprehensive overviews provided in [17,18], respectively. The problem of container drayage has been considered by other research work. For example, in [19], authors have developed models capable of minimizing both present and future operating costs. The authors in [20] modeled the container drayage problem as a pickup and delivery problem and proposed Lagrangian relaxation to solve the problem. In [21], the authors proposed three approaches for the container movement problem with time windows at origins and destinations. These approaches are based on an asymmetric multiple traveling salesman problem.
Some researchers have considered a single vehicle while modeling the VRPPD. For example, the single-vehicle VRPPD with deterministic demands and predetermined customer visit sequence has been considered in [22]. Authors in [22] have developed the routing of a single-vehicle that delivers multiple products under stochastic demand. Other researchers, such as [23,24], have considered a VRP with mixed pickup and delivery, called (VRPM), where pickups are made before deliveries.
Recently, Zhang et al. [25] introduced a new VRP variant which encourages the reuse of collected items. They developed a segment-based evaluation procedure to reduce the computation time. They provided a mathematical formulation and a metaheuristic algorithm, and analyzed several features of the problem. In [26], a hybrid heuristic algorithm for the 3L-PDP problem is extended by using two key improvements: the first is the usage of a tabu strategy for enlarging the local search space, the second is the employment of complex block generation and depth-first heuristics for incrementally finding one proper box at a time in the packing phase. The experimental results show that the improved hybrid heuristic algorithm outperforms its origin regarding total travel distance on benchmark state-of-the-art instances. More variants of VRPPD including MRVRP, VRPP, SDVRP, SVRP, VRPSF, and VRPTZ can be found detailed in [27,28].
Finally, we recall the main contributions of this paper. We have been motivated by the large number of studies that tackled the VRPPD. It is worth noting that the majority of possible constraints have already been suggested in the last two decades. However, our work provides prototypes to solve the VRPPD based on CHAHIA's needs. These needs include three different variants of VRPPD. The first variant is a VRP with a Simultaneous Pickup and Delivery (VRPSPD) without leaving any empty boxes at any customer. The second variant is a flexible VRPPD that ensures separate and/or simultaneous pickup and delivery. The third one is a simultaneous pickup and delivery VRP allowing the flexibility to drop and leave some boxes, which will be consequently penalized. More specificities regarding the proposed models can be found in Section 3.

Problem Ingredients and Approach of the Proposed MILP Models
In order to state the problem more formally and clearly, we present our ideas which are based on a real-life context. We have taken into account a set of assumptions, which are as follows: • Time window: Includes the service time at each customer, loading time for each vehicle, travel time between each couple of customers, time necessary to park the vehicle, and the workday duration. • Combined capacity; volume and weight: The majority of VRP formulations consider the physical capacity of vehicles, which is generally expressed in terms of weight. In our problem, goods are delivered using boxes. Therefore, beside the capacity in terms of weight, we have to make sure that each vehicle has the capacity in terms of number of boxes to be fitted into the vehicle. Thus, we defined a new volume capacity that includes both loaded weight and box numbers. We denote that we have considered the volume occupation in each vehicle to be discrete. Motivated by this assumption, which complicates the use of vehicle space, we assume that the combined capacity has not been addressed in the literature. •

Customer balance (stock) in terms of boxes:
We suppose, initially, that each customer has boxes ready for pickup. These boxes are considered by our models before performing any delivery. • Optimization choice: Three different assumptions are considered. The first one guarantees that the pickup and delivery are performed simultaneously, without leaving any empty boxes for any customer. The second assumption is more flexible. It provides two choices of pickup and delivery operation, separately and/or simultaneously. The third one ensures that the pickup and delivery operation are performed simultaneously. However, it provides a flexibility to drop and leave some boxes, which will be consequently penalized. • Heterogeneity of the fleet Multiple vehicle types exist, having different sizes.

Mathematical Formulations
Three different formulations of the real-life VRPPD are detailed in this section.

Simultaneous Pickup and Delivery Vehicle Routing Problem: MILP 1
This section models the Simultaneous Pickup and Delivery Vehicle Routing Problem variant (SPDVRP). In SPDVRP, the operations of pickup and delivery are performed simultaneously at each customer node. Each customer is visited once, their order (in terms of goods) is all delivered, and all boxes are picked up. The variables and parameters of the proposed model are as follows: -Indices and parameters: • n: number of customers (i = 1, . . . , n).  i : volume demand of customer i in terms of number of boxes. This parameter is deduced from the following relation: d α i : average weight of a full loaded box dedicated for a customer i. • Sc ini i : initial balance of boxes that exists at the store of a customer i.
formulation is as follows: Knowing that: According to the MILP1 formulation, the objective function (1a) seeks to minimize the total transportation cost. Constraints (1b) and (1c) ensure that the first vertex (node) associated with each customer is visited once for a simultaneous pickup and delivery. Constraint (1d) is related to the time service. The latter includes the travel time between two customers and the time necessary for loading goods into the vehicle. It prevent exceeding the predefined time capacity. Constraint (1e) ensures that the weight of the vehicle load before leaving a customer j does not exceed the vehicle capacity if there is a path between customers i and j. Constraint (1f) ensures that the weight of goods required by all customers does not exceed the vehicle capacity in terms of weight. Constraint (1g) guarantees that the number of boxes loaded into the vehicle before leaving a customer j does not exceed the vehicle volume whenever it travels from j to another customer i. Constraint (1h) ensures that the number of boxes required by all customers does not exceed the vehicle volume capacity. Constraints (1i) and (1j) express the relation between two successive customers in terms of weight. We considered here two types of weight; one is associated with the delivery boxes, while the other one is associated with the initial boxes existing at each customer store. Constraint (1k) defines the vehicle load before leaving the depot. Constraints (1l) and (1m) present the weight of goods loaded into the vehicle before traveling from the depot to a customer j. Constraints (1n) and (1o) represent the relation between two successive customers in terms of volume (number of boxes). Constraint (1p) defines the boxes loaded into the vehicle before leaving the depot. Constraints (1q) and (1r) show the relation between the depot and a customer j in terms of volume. Finally, (1s) and (1t) represent the sub-tour elimination constraints.

Flexible VRPPD: MILP 2
MILP2 is more flexible than MILP1. It considers a different representation of customer services. We associate two vertices i and i + n with each customer i, where i and i + n are the same, but have different types of service operations.
In reality, customers can be classified into two types according to the way of required services. Indeed, MILP2 allows two types of services for customer i. The first type assumes that the pickup and delivery are performed simultaneously. In such cases, vertex i is visited and vertex i + n is left. The second type ensures that a customer i is visited twice, where the delivery is performed at vertex i, while pickup is done at vertex i + n. Figure 2 illustrates the representation of customers based on a virtual replica.  Figure 3 represents two feasible tours for the same vehicle trip. The solution includes nodes visited once with simultaneous pickup and delivery (black nodes), and nodes visited twice due to two separated operations (encircled black nodes). The originality in the proposed configuration is that certain feasible subtours must be allowed in appearance. However, thanks to the duplication of nodes, it became easy to model the problem by conserving classical subtour elimination principles, while using a subtle manner of formulation that takes into account the non-conservative progression of vehicle load. For MILP2 formulation, we define additional parameters and variables, which are as follows. Additional parameters for MILP2: (∀j = 1, ..., 2n; ∀i = 1, ..., n; i = j) (2j) (∀j = 1, ..., 2n; ∀i = 1, ..., n; i = j) (2k) x k ij + y j−n )M 12 (∀j = n + 1, ..., 2n; ∀i = 1, ..., 2n; i = j) (2l) j−n ) * β + (1 − x k 0j + y j−n )M 17 (∀j = n + 1, ..., 2n; ∀k = 1, ..., m) (2q) x k ij − 1)M19 (∀j = 1, ..., n; ∀i = 1, ..., 2n; i = j) (2s) x k ij )M20 (∀j = 1, ..., n; ∀i = 1, ..., 2n; i = j) (2t) x k ij + y j−n )M21 (∀j = 1, ..., n; ∀i = 1, ..., 2n; i = j) (2u) (∀j = n + 1, ..., 2n; ∀i = 1, ..., 2n; i = j; k = 1, ..., m (2v) )M27 (∀j = 1, ..., n; ∀j = n + 1, ..., 2n; ∀k = 1, ..., m) (2aa) x k ij )M 28 (∀j = 1, ..., n; ∀i, j = 1, ..., 2n; i = j) (2ab) x k ij − 1)M 29 (∀j = 1, ..., n; ∀i, j = 1, ..., 2n; i = j) (2ac) Knowing that: • TC j = TPc u * d shares the same objective function as MILP1. However, constraints (2b) and (2c) mean that the first vertex associated with each customer must be visited once, either for a single operation of delivery or for a simultaneous pickup and delivery. Constraint (2d) ensures that the second vertex associated with the customer is visited only if a combined pickup and delivery do not occur at the first vertex. Constraint (2e) defines the time capacity. It contains three parts which are; the time capacity at the depot, the time capacity when the pickup and delivery are performed simultaneously, and the time capacity when the pickup are performed separately. Both Constraints (2l) and (2m) define W j in terms of Wi whenever j is visited immediately after i, knowing that i and j are two virtual customers. Constraints (2q) and (2r) express the relation between a customer j and the depot in terms of weight. Both constraints define W j in terms of W(0, k) whenever j is visited immediately after the depot, knowing that j represents a virtual customer. Constraints (2z) and (2aa) express the link between a customer j and the depot in terms of volume (number of boxes). Constraints (2ab) and (2ac) define Tj in terms of Ti whenever j is visited immediately after i. Both constraints represent the time left between i and j. Constraints (2ad) and (2ae) express the time relation between a customer j and the depot. Constraints (2af) and (2ag) represent the time relation at a customer j. Both constraints reveal whether the pickup and delivery are performed simultaneously or not. Finally, Constraints (2ah) and (2ai) represent the subtour elimination constraints. We note that yj is considered in MILP1 to make sure of whether the pickup and delivery are performed simultaneously or not.

Flexible VRPSPD: MILP 3
MILP3 represents a simultaneous pickup and delivery operation such as the one in MILP1. However, the novelty in MILP3 is its flexibility, which allows the dropping/releasing of some boxes. Thus, a partial collection of boxes is acceptable due to either a limited vehicle capacity and/or the high cost of collecting all boxes. Despite MILP3 being more realistic in terms of collecting flexibility, several trips become necessary to collect all boxes, which could be high in terms of cost. To deal with this issue, we added a penalty cost to the objective function that penalizes the process of leaving boxes. Figure 4 represents an overview of MILP3. We define the following additional variables for MILP3.

•
Sc end j : boxes left at customer j. • c ca : penalizing cost associated to the left boxes.
MILP3 formulation is as follows: (∀i, j = 1, ..., n; i = j) (3j) (∀i, j = 1, ..., n; i = j) (3k) (∀i, j = 1, ..., n; i = j) (3p) Knowing that: It is worth noting that MILP3 and MILP1 are quite similar. As mentioned earlier in this section, MILP3 allows leaving some boxes due to either vehicle capacity or trip cost. To satisfy this flexibility, we added constraint (3s), which ensures that boxes related to customer demand plus initial balance must not exceed total required boxes.

Computational Study
We adopted a configuration quite similar to the one considered in [29], where authors have implemented their models in the IBM ILOG CPLEX Optimization Studio (Version: 12.6). All the experiments were conducted on a computer with an Intel(R) Core (TM) i7-7700 CPU@3.6 GHz and 8 GB memory under the Windows 10 Pro system. As we highlighted earlier in the introduction section, our work is dedicated to a poultry company in Tunisia, called CHAHIA. Due to a confidentiality issue, CHAHIA has fed our research work with limited data. Based on the collected data, we created a set of instances. To make sure that data features are diversified, we created 44 instances. Several factors have been used to calibrate the quality of solutions, enhance the solving process, and diversify the benchmark of instances. These factors are as follows: • The types of demands, i.e., the quantities of required goods (low/medium/high). Results show that MILP2 achieves the lowest objective function values compared to both MILP1 and MILP3. This achievement is mainly due to its flexibility and high number of specifications that it covers.
In case of low demand or low number of boxes to be picked up, both MILP1 and MILP3 provide similar costs. However, if the demand or box numbers increase, MILP3 outperforms MILP1 and provides the lowest travel cost. This fact shows that sometimes it would be more profitable to leave some boxes instead of collecting them. We recall that the left boxes could be collected in future trips, otherwise, a penalization is considered.
In case we reduce the cost of left boxes by relying on the possibility to recover them in the future with better conditions, MILP3 will provide the lowest cost for the majority of instances. Moreover, if we change the capacity of vehicles from small to big, such as in instances "inst5.2.004" and "inst5.2.005", which are the same in fact, MILP1 and MILP3 diversify and provide different value costs, however, MILP2 remains unchangeable. Furthermore, if we increase both, the demand to be between 151 Kg and 500 Kg per customer, and the number of customers to exceed 18, MILP3 becomes the most efficient model. With regard to the computational time performance indicator, MILP3 converges quickly. Contrariwise, MILP2 consumes more time to provide optimal solutions. In case the number of customers exceeds six, MILP2 remains incapable of retrieving an optimum solution during a reasonable period of time (i.e., less than 3600 min). With regard to the number of served customers, MILP1 and MILP3 provide a better performance compared to MILP2.
It is noticeable that when the number of vehicle stops for pickup or delivery increase, it will affect the quality of poultry products from CHAHIA. This affection is related to the unstable air conditioning of the vehicle refrigeration cabin. Therefore, it is recommended to reduce the number of load breaks.
The different readings of the same problem led to an interesting diversification in terms of travel cost, computational time, and number of served customers.

Discussion of Limitations
Our contribution develops a 'proof of concept' useful to establish better awareness and foresight in decision-making processes. It provides flexibility in reacting against miscellaneous realistic situations in transportation and goods delivery activities. However, the major limitations of this work could be summarized as follows: • Due to a confidentiality issue, CHAHIA has fed our research work with limited data. Therefore, bigger data is required for a more relevant benchmark. • CHAHIA has refused to reveal all serve stores, customers, and quantities, due to confidentiality. • CHAHIA was looking for a prototype that could be enhanced in the future by its engineers. • Due to lack of equipment, we were not able to perform any experiments with more than 23 customers for MILP1 and MILP3, and 6 customers for MILP2. • Exact methods considered in this work are incapable of solving large instances. Therefore, integrating heuristics is highly recommended.

Conclusions and Future Work
In this paper, the Vehicle Routing Problem (VRP) with Pickup and Delivery (VRPPD) was addressed based on real-life perceptions. Three different VRPPD models have been developed based on Mixed Integer Linear Programming approaches. The first model is designed as a VRP with simultaneous pickup and delivery. The second one is modeled as a VRP with either simultaneous and/or separated pickup and delivery. The last one is modeled as a VRP with simultaneous and flexible delivery and pickup allowing both partial pickups and penalized losses. In addition, this study proposes a proof of concept based on a rich benchmark. The data adopted in this work were provided by a poultry company in Tunisia, called CHAHIA. The developed models have been compared against each other using multiple state-of-the-art key performance indicators. Based on the presented richness in term of modeling, interesting perspectives could be considered in the future, such as matching operational issues with tactical planning elements in order to make the models more dynamic and robust facing various demand scenarios. Finally, various heuristics could be proposed in order to solve large instances and enable quicker convergence towards the optimum.