The vehicle routing problem with heterogeneous locker boxes

To achieve logistic efficiency and customer convenience in last-mile delivery processes, a system with alternative delivery points in the form of locker box stations can be used. In such a system, customers can be served either at their home address within a certain time window, or at a locker box station where parcels can be picked up at any time. Customers can get a compensation payment when being served at a locker box. They can have a request of more than one parcel and the parcels can be of different sizes. At a locker box station, a limited number of slots of different sizes is available; we assume that parcels of one customer can be stored together in a slot. We consider the vehicle routing problem with heterogeneous locker boxes, where the total cost—consisting of routing and compensation costs—has to be minimized while taking into account the packing of parcels into locker boxes. We provide a mathematical formulation of the problem and propose a metaheuristic solution method. Instances and results from the literature for the problem with a single parcel and a single slot size are used to benchmark our metaheuristic solution method. For the problem with different sizes, we compare a unit-size model to a multi-size model, packing being considered in the latter. Finally, we analyze how different configurations of locker box stations work for different demand scenarios.


Introduction
As a consequence of ongoing trends in e-commerce, direct-to-consumer deliveries have become a meaningful part of today's business. Customer convenience and logistic efficiency play an important role in developing competitive strategies. This brings new challenges for the related logistic process, in particular in the last mile of the delivery chain. When delivery companies have to deliver parcels to the customers, one of the problems they face is that of unsuccessful delivery attempts when the customer is not at home to receive the parcel(s) in person. This involves, above all, several unnecessarily driven kilometers for the deliverer. Moreover, it can also pose inconveniences for the customer, for instance, when the parcels are brought to an arbitrary shop, that can be accessed only during its opening hours, or to a neighbor who is not at home at the right time, or when they have just been dropped in front of the door having the inherent risk of being stolen. One option to overcome these problems is to provide the customers-in addition to home delivery-the possibility of picking up their parcels at what are known as locker box stations. Such locker box stations can be found frequently in practice, located, for example, at post offices, shopping malls, railway stations, residential areas, etc. They are accessible at any time, which means that customers can decide when it fits best for them to pick up the order. The slot where the customer's parcels are stored can only be opened with the electronic private code that the customer receives from the delivery company, which minimizes the risk of the parcels being stolen.
In the context of vehicle routing, the possibility of locker box delivery brings new aspects. In contrast to the classical vehicle routing problem, there can be more than one location associated with a customer's request. A customer can be served at the home address during the preferred time window to receive the parcels, which can be also referred to as attended home delivery. Even though this delivery location can be any specific location where the customer can be met in person, we call this option home delivery in the following. In addition, there is the possibility to serve a customer at a locker box station. At a locker box station, there are several slots of different sizes; if a customer receives more than one parcel, they can be packed together into a larger slot or into several separate smaller slots. We assume that the customer indicates in advance which locker box stations are acceptable for delivery, since some may be more favorable than others. For example, there may be stations closer to the home address or on the way home from work; alternatively, there may be stations that are totally out of the route and the customer may not be willing to travel all the way just to pick up the parcels.
When customers are served at a locker box station, they receive a compensation payment which is predefined and independent of the number of parcels or the chosen station. Paying a compensation value should provide an incentive to the customers for using locker box delivery, since it allows for a more efficient route planning on the side of the delivery company.
The compensation payment is taken into account in the objective function of the total cost which consists of routing costs, measured by traveled distance, and the total compensation costs. A delivery plan has to be constructed such that the total cost is minimized. The delivery plan includes the decision regarding which customers are served at home (during their time window) and which are served at a locker box station, in which case it also has to be decided to which station a customer is assigned to and how the different sized parcels are packed into the different sized slots. In addition, the delivery plan includes the routing decision for the used delivery locations, that is, the used locker box stations and the used home addresses.
The contribution of this article is threefold. First, we introduce the vehicle routing problem with heterogeneous locker boxes (VRPHLB) as a relevant problem in the context of sustainable city logistics and provide a mathematical model in order to offer a specific description of the problem. Second, we propose a metaheuristic solution to solve realistically sized instances. Using the metaheuristic method to solve existing instances from the literature for the simpler problem with only one slot and parcel size, we can show that the metaheuristic method is competitive with respect to the solution quality and time. Third, we use the proposed metaheuristic method to solve the newly introduced VRPHLB in order to obtain valuable insights on the structure of the problem and derive managerial implications.
The remainder of the article is organized as follows: In Sect. 2, we provide an overview of the relevant existing literature. In Sect. 3, we define the VRPHLB formally and introduce a corresponding mathematical model. In Sect. 4, we describe the proposed metaheuristic method in detail. In Sect. 5, we conduct an experimental study.

Literature review
A general idea in modern city logistics is to provide the possibility of alternative or flexible delivery locations to respond to the lifestyle of today's customers. This is to some extent already well established in practice and has also become of interest in research in the last years. Savelsbergh and Van Woensel (2016) discuss the challenges and opportunities of city logistics for the present and future. They mention pickup-point systems (which is the category to which locker box stations belong) as one promising concept to efficiently design last-mile delivery systems. One of the early works in the field of alternative delivery locations is provided by Reyes et al. (2017). The authors introduce the vehicle routing problem with roaming delivery locations, where, for each customer, different delivery locations with non-overlapping time windows are given. This model was inspired by the idea that customers get their parcels delivered to the trunk of their car which is parked at different locations throughout the day. The authors show that roaming delivery provides significant cost advantages compared to a traditional home delivery system. A further study which deals with alternative delivery locations is provided by Zhou et al. (2018). They consider a problem where the parcel can be delivered to each customer directly or through a pickup point and formulate the resulting problem as a multi-depot two-echelon vehicle routing problem. However, they do not consider capacity restrictions at the pickup points. Orenstein et al. (2019) investigate a distribution system where customers are served at capacitated locker box stations. Home deliveries, and thus time windows for customers, are not considered in this problem. They consider different sizes for the locker boxes. However, at most one parcel can be assigned to a slot; thus, packing decisions do not have to be taken. In the work of Sitek and Wikarek (2019), customers can be served at home or at an alternative delivery location. Capacity constraints at the alternative delivery locations are considered. However, time windows at the customers' home locations are not part of the problem. Mancini and Gansterer (2020b) present the vehicle routing problem with private and shared delivery locations, which is closely related to our problem. The private locations are the home delivery locations and the shared locations are the locker box stations in our problem. However, they do not consider different sizes for the parcels and slots. The authors build a mathematical formulation of the problem and propose a matheuristic based on large neighborhood search, where a repair step is done using an MIP model. Iterated local search is used to prevent the matheuristic from getting stuck in a local minimum. A recent work in the field of locker boxes has been produced by Schwerdfeger and Boysen (2020). However, they deal with mobile locker box stations whose locations can change over time. They address the problem of optimally locating the mobile locker boxes such that customers are within a certain range of their assigned station at some time during the planning horizon; the locker fleet is minimized.
The vehicle routing problem with heterogeneous locker boxes, which is introduced in this paper, has to deal with a routing problem as well as with a packing problem. A combination of packing and routing can be also found in the works of Iori et al. (2007) and Fuellerer et al. (2009), addressing the two-dimensional loading capacitated vehicle routing problem, where rectangular weighted demand items have to be loaded onto capacitated vehicles with a rectangular loading surface. Iori et al. (2007) develop an approach based on a branch-and-cut approach to solve the problem and integrate lower bounds, heuristics and a branch-and-bound algorithm to detect infeasible loading patterns that were allowed in the first place. Fuellerer et al. (2009) propose a heuristic solution method based on ant colony optimization. The feasibility of loading patterns is checked by using lower bounds, heuristics and a truncated branch-and-bound.

Problem definition and mathematical model
There is one depot represented by node 0, where every tour starts and ends. An unlimited fleet of homogeneous vehicles v ∈ V can be used. We assume that the capacity of the vehicles is not restrictive, following the idea that parcels are small enough to not constrain a tour.
There is a set of m customers C = {1, . . . , m}, where each customer requests parcels of different sizes. The number of parcels of size l ∈ L which a customer i gets is denoted by q il .
There is a set of n locker box stations B = {m + 1, . . . , m + n}. At each locker box station k ∈ B there is a limited number Q kl of slots of size l ∈ L available.
The whole set of nodes is denoted by N = {0} ∪ C ∪ B. Travel distance between two nodes i, j ∈ N is denoted by d i j . In this study we assume that travel time between two nodes is equal to distance.
Each customer i accepts his home address, denoted by {i}, as a delivery location. Delivery there can happen only within a time window [E i , L i ], where E i is the earliest possible start of service and L i the latest. Service at a customer's home address takes time s i . The time window at the depot, denoted by [0, T ], is set such that it represents the maximal duration of a tour. The same holds for the time windows at the locker box stations, since we assume that such stations can be entered by delivery personnel at any point in time within the planning horizon. Visiting a locker box station requires a service time in general longer than that of the customers.
In addition to home delivery, a customer can accept locker box delivery. However, he does not need to accept all available stations but only a subset of them; he may accept none when he does not want to have locker box delivery. The accepted locker box stations of customer i are denoted by B i ⊆ B. Serving a customer at a locker box station causes compensation cost c.
The whole demand of a customer has to be brought to the same locker box station. Parcels of one customer can be packed together into one slot if the slot space allows it. However, not more than one customer can be assigned to one slot. This is for practical security reasons, since only the lawful customer is allowed to have access to the parcels.
The decision variable x v i j is 1 if node j is visited directly after node i using vehicle v; it is 0 otherwise. The binary decision variable y v ik is 1 if customer i is served at locker box station k using vehicle v. The binary decision variable z v i is 1 when a customer is served at home and 0 when a customer is served at a locker box station. The continuous decision variable S i gives the service start time at node i.

Relaxed model with unit-size locker boxes and parcels
In a first modeling approach, we transform the different sizes to unit sizes and add up the parcels for a customer and thus get an aggregated demand q i = l∈L q il for each customer i. Similarly, this is done for the different slot sizes at a locker box station; in other words, we get an aggregated capacity Q k = l∈L Q kl for each locker box station k. With this we obtain a relaxed version of the original problem, since we do not consider how the parcels are packed into the slots but only check if the total demand of a customer can be brought to a locker box station.
The objective function of the problem is to minimize the total cost that consists of total traveled distance plus total compensation costs for serving customers at locker boxes.
Constraints (2) ensure route continuity. Constraints (3) ensure that a vehicle is used at most once. Constraints (4) connect decision variables y to decision variables x and constraints (5) connect decision variables z to x. By (6) it is guaranteed that every customer is served and service can happen either at home or at a locker box station. Constraints (7) ensure that capacity constraints of the locker box stations are respected. Constraints (8) define the time window constraints. Constraints (9) care for the continuity of service time variables. Constraints (10)-(13) define the scope of the decision variables.

Model with multi-size locker box stations and parcels
We now give a formulation of the complete problem, which also contains the decision regarding how the different sized parcels are packed into the different sized slots. Hence, we are required to integrate a packing part into the previous modeling formulation. We refer to this as the model with multi-size locker box stations and parcels. The packing problem in this model belongs on the one hand to the class of bin packing problems with conflicts (Jansen 1999;Gendreau et al. 2004) and on the other hand to the class of variable-sized bin packing problems (Friesen and Langston 1986;Kang and Park 2003;Haouari and Serairi 2009). For our case, a bin is a slot of a locker box station. A conflict appears when two parcels do not belong to the same customer, as such parcels are not allowed to be packed together into a slot. We therefore have to deal with a capacitated bin packing problem with conflicts and heterogeneous bins. Parcels have to be assigned to slots such that the capacity of a locker box station and the space restriction (size) of a slot are respected, no conflicts appear and no split deliveries occur.
For the bin packing part of our problem, we use the following notation. There are in total i∈C,l∈L q il parcels and each parcel g ∈ G is characterized by size a g and the customer i g which it belongs to. There are in total k∈B,l∈L Q kl slots and each slot h ∈ H is characterized by size b h and the locker box station k h which it belongs to.
The decision variable u gh is 1 if parcel g is assigned to slot h and 0 otherwise. The decision variable w h is 1 if slot h is used and 0 otherwise. Constraint (7) can be taken out of the model (1)-(13) and replaced by the following constraints: Constraints (14) ensure that the space restriction of a slot is met. Constraints (15) guarantee that each parcel is assigned to a slot. By constraints (16), we forbid that parcels, that are not of the same customer, are packed together into one slot. Constraints (17) ensure that parcels of one customer are brought to the same station. Constraints

Metaheuristic solution method
When implementing the model of Sect. 3, the limits concerning solvable instance size in reasonable time are reached quickly. Hence, an alternative approach is needed to solve realistically sized problems. We propose a metaheuristic solution where the general idea is as follows: -Solve the problem as a VRPTW where all customers are served at home during their time window. We call this a pure home delivery plan. -Then we do the following for a number of iterations: -Decide which customers are taken out of the pure home delivery plan and served at locker box stations. -Re-route the remaining home delivery customers plus the used locker box stations. -Solve bin packing problem to reduce capacity usage.
-Re-optimize the selection of home delivery and locker box customers.
In the following, the different steps of the method are described in more detail. Subsequently, we provide the corresponding pseudocode of the complete method in Sect. 4.6.

Adaptive large neighborhood search
To solve the routing problem, we use adaptive large neighborhood search (ALNS) which was first introduced by Ropke and Pisinger (2006). Since then, ALNS has been applied and adapted to a broad range of routing problems and has shown to work very well with respect to solution quality and computational time. Kovacs et al. (2012) developed an ALNS to solve the service technician routing and scheduling problem, where a given number of technicians have to complete a given number of service tasks. They implemented a new adaptive mechanism based on pairwise evaluation of the performance of destroy-repair operator pairs. A recent work related to ALNS has been undertaken by Ortega et al. (2020), where a matheuristic solution approach based on ALNS is presented to solve the consistent inventory routing problem with time windows and split deliveries.
An outline of the basic structure of the ALNS based on the work of Ropke and Pisinger (2006) can be seen in Algorithm 1. (Ropke and Pisinger 2006) 1: Generate a starting solution s 0 2: s...incumbent solution, s best ...best solution 3: s ← s 0 , s best ← s 0 4: repeat 5:

Algorithm 1 Adaptive Large Neighborhood Search
Choose a destroy operator d and a repair operator r based on adaptive weights 6: Apply d and r to s yielding s 7: if s is better than s best then 8: set s best = s 9: set s = s 10: else if s meets the acceptance criteria then 11: set s = s 12: end if 13: Update the scores and weights of the operators 14: until some stopping criteria is met 15: return s best To construct an initial feasible solution, we use sequential cheapest insertion. Weights and scores are evaluated pairwise for each pair of destroy/repair operator (Kovacs et al. 2012).
For destroying a solution, the number of customers to remove in each iteration is chosen randomly from the interval [min(0.1m, 30), min(0.4m, 60)], where m is the total number of customers. As destroy operators we use the following: random removal, worst removal with/without perturbation, related removal seed (a seed customer plus a number of nodes that are most related to this seed customer are removed) and related removal chain (a chain of customers that are most related to each other is removed).
When repairing a solution, we have to insert the removed customers back into the routes; the feasibility of such an insertion with respect to time windows has to be checked. This is done as suggested by Campbell and Savelsbergh (2004). As repair operators we use the following: cheapest insertion parallel, cheapest insertion sequential, arbitrary insertion, m-regret insertion (m = 2, 3, 4) with two variations. As in the work of Hemmelmayr et al. (2012), the m best insertions can be within one route or, as in the work of Ropke and Pisinger (2006) and Kovacs et al. (2012), the m best insertions for a customer can be in different routes.
Every time a new global best solution is found, we apply the Or-Opt operator to this solution. Or-opt is the moving of a sequence of up to three consecutive nodes to another position within the same route (intra-route) or to another route (inter-route). Deteriorating solutions are accepted based on simulated annealing (Ropke and Pisinger 2006). The whole algorithm stops when the time limit, limit on the number of total iterations, or limit on the number of iterations without improvement is reached.

Decide which customers are served at locker box stations
The decision regarding which customers are taken out from the pure home delivery plan and served at locker box stations is made by the following operators.

Reduce distance by assigning customers to locker boxes
We shift those customers to locker boxes that provide the biggest savings in distance when they are removed from the routing solution. This is done as long as moving customers to locker box delivery is beneficial with respect to solution quality. One has to consider here that every time we move a customer from home delivery to locker box delivery, on the one hand, we save traveled distance, but on the other hand, we have to add the compensation cost (for locker box delivery) to the objective value. Thus, the compensation value has to be subtracted from the distance savings in order to evaluate total savings. Total savings are multiplied with a random factor from the range [0.8, 1.2] in order to allow for more possible solutions. Moreover, deterioration is allowed when the benefit of moving a customer from home delivery to locker box delivery is evaluated. How far deterioration can go is decided by a threshold value determined randomly and dependent on the total cost of the solution. In our computational experiments, we take a threshold value from the range [0, f * 0.05], where f is the total cost of a solution.
If a customer is decided to be served at a locker box station, we have to decide which station this customer should be assigned to. The "reduce distance" operator is applied with two different ways of opening locker box stations: sequential and parallel. For the sequential approach, a new station is opened only when none of the already open ones have enough free capacity for the currently considered customer. For the parallel approach, we assume that all locker box stations are open and we try the stations in increasing order based on the remaining capacity, that is, the station with the smallest remaining capacity is considered first for taking the customer.

Fill up locker box stations
In this operator, locker box stations are opened sequentially and each new opened locker box station is filled up with appropriate customers until capacity limit is reached. More precisely, home delivery customers are sorted such that the first customer in the list L is the one that provides the biggest improvement in the objective function when this customer is not in the routing solution but served at locker box stations. The last customer in the list is the one that is least beneficial for serving at locker boxes. When evaluating how the objective function changes when a certain customer is shifted from home delivery to locker box delivery, the corresponding change value is multiplied with a random factor, as it has been done in the previous operator "reduce distance".
When a new station has to be opened, we choose the station that is on the one hand accepted by the first customer in list L and on the other hand closest (with respect to distance) to the last customer in L. The latter is included as a decision criterion since we can assume that a customer at the end of list L will rather remain in the routing solution and therefore additional routing costs incurred by an opened station may not be too high in the end.
Home delivery customers are assigned to the currently open station following the order of list L, taking into account that only those customers that accept the currently open station can be assigned to it. When no suitable customer can be feasibly moved to the currently open station with respect to locker box capacity, we check if the overall solution has improved. If it has, we open the next station. If it has not, we stop opening stations; in other words we stop moving customers from home delivery to locker box delivery.

Remove tour
In case there are more than two tours in the pure home delivery solution, we try removing one tour by moving all customers from that tour to locker box delivery or reassigning them to another tour if locker box delivery is not possible. First, we open the locker box station that is accepted by the largest number of customers from the closed tour. From the customers of the closed tour, we first move those that have the lowest number of accepted stations, which means that there is less flexibility in assigning them to stations. When there is a customer that does not accept locker box delivery, we route it in another tour that is not closed. We evaluate all tours and close the one that provides the biggest improvement in the objective value when it is closed.

Checking locker box capacity feasibility
When assigning customers to locker box stations, we have to check feasibility with respect to the available number of slots. This step is not straightforward due to the following characteristics of the problem. First, customers have a demand for more than one parcel and each parcel can have a different size. Moreover, the slots are of different sizes. Furthermore, the whole demand of a customer has to be brought to the same locker box station. Parcels of one customer can be packed together into a slot, but no other customer's parcel(s) can be packed into a slot that is already occupied by another customer.
In order to check if a customer's demand can be feasibly packed into the available slots of a locker box station, we have to evaluate the different possibilities regarding how the parcels can be packed in general, assuming in this step that there is an unlimited number of slots available.
Therefore, a promising (sub-)set of feasible packing solutions for each customer is generated in a pre-processing step by the Iterative First-Fit Decreasing (IFFD) algorithm. IFFD was introduced by Kang and Park (2003) to solve the variable-sized bin packing problem. In IFFD, all the items (parcels) are first assigned to the largest size bins (slots) using First-Fit Decreasing (FFD). In FFD, an item is assigned to the first slot that has enough capacity left. Then another solution is obtained by repacking the items from the last bin of the solution to the next smaller bins using FFD. This is continued until repacking of the items is impossible (either because the next bin size is not big enough to take the items of the last bin or the smallest bin size is reached). The first solution in the list of feasible solutions obtained in this way is the one that uses only bins of the largest size.
After having obtained a set of feasible packing solutions for every customer, when checking if a customer can be feasibly assigned to a locker box station, we can use this set to check if one of those solutions is feasible with respect to the available slot capacities.
If there is sufficient capacity left at a locker box station, several or even all packing solutions may fit. However it could make a difference for later decisions concerning feasible assignments of other customers, which of the packing solutions is used, that is, which and how many slots are used, and which ones are left free for other customers. Thus, how we go through the list of packing solutions when checking feasibility makes a difference. Two possibilities are the following.

Start with the first feasible solution
One strategy is to always start at the beginning of the list, which is the solution that only takes slots of the largest size. If this solution is not feasible concerning the available capacities, we consider the next solution in the list and so on. Starting in this manner, we may have locally minimum capacity requirements. However, if we do this for every customer, we soon reach the limit on bins of the largest size and customers that would need a bin of the largest size later (since a large parcel may be packed only into a large bin) cannot be moved to the locker box station anymore.

Start with a random feasible solution
We consider another strategy, where we start at a random point in the list, in order to not always check the solution with only the largest bins first. With this, a more balanced usage concerning the different sizes of slots should be attained.

Bin packing
When assigning customers to locker box stations in the previous step, an assignment of parcels to locker box slots was also obtained as a result of the capacity check. However, there, the decision regarding which packing solution is used was taken separately for each customer when moving them to locker boxes. After deciding which customers are served at locker boxes, we have the full set of locker box customers and their respective parcels; the capacity requirements may be optimized by solving it as a simultaneous bin packing problem. This means that we have to decide how to pack the parcels of the locker box customers into the locker box slots such that capacity usage is minimized.

Heuristic set covering formulation
We solve the bin packing problem with a heuristic set covering formulation. A (sub-)set of feasible packing solutions for each locker box customer is generated as described in Sect. 4.3 according to the Iterative First-Fit Decreasing (IFFD) algorithm. Since IFFD yields in general only a subset of feasible solutions, the according set covering formulation is heuristic.
We want to mention here that there would be other algorithms, as for example the work of Delorme and Iori (2020), to solve the bin packing part of our problem. However, since the overall problem is already quite complex, we want to keep the computational effort of the packing part tractable.
For the set covering formulation the set of feasible solutions for a customer i is denoted by Ω i . The cost, that is, the number of used slots, of a packing solution p ∈ Ω i is denoted by c i p . The number of slots of size l needed by a packing solution p ∈ Ω i is denoted by q l i p . Decision variable X k i p is 1 if, for customer i, the packing solution p ∈ Ω i is used and assigned to locker box station k. With this, the following model is set up: The objective function (22) is to minimize the number of used slots. Constraints (23) ensure that the parcels of each customer are packed. By (24), the capacity restriction at the locker box stations is taken into account. Constraints (25) define the decision variables to be binary.

Re-optimize the selection of home delivery and locker box customers
In this step, the solution which has been found so far is improved by moving customers from home delivery to locker box delivery and vice versa. For this, we use several operators: - When a locker box station becomes empty in any of the operators (because the assigned customers are shifted to home delivery or to another station), it is closed, which in particular means that it is removed from the routing solution.

The complete metaheuristic method
A pseudocode of the complete metaheuristic is shown in Algorithm 2.
Algorithm 2 Metaheuristic solution method 1: Generate a pure home delivery starting solution s by ALNS 2: for m iterations do 3: Move customers in s from home delivery to locker box delivery; yields s 4: Solve routing problem with nodes from s by ALNS 5: Solve bin packing 6: while limit on number of iterations without improvement not reached do 7: Choose a re-optimize operator randomly and apply it to s 8: end while 9: if f (s ) < f (s best ) then 10: s best ← s 11: end if 12: end for 13: return s best In the following, some additional information for Algorithm 2 is provided: -Line 2: The number of iterations m is chosen to be 500 throughout the computational experiments. -Lines 4-11 are only executed if the set of moved home delivery customers (including the respective assignment of parcels to slots) has never appeared before. -Line 6: The number of iterations without improvement is chosen to be 10 in the computational experiments.

Computational study
We now provide an experimental analysis of the previously introduced models and the metaheuristic solution method. All algorithms and models are implemented using the C++ programming language. CPLEX 12.9 is used for solving the exact models (MIP). Multithreading is deactivated. All programs are run on an Intel Xeon Processor E5-2670 v2 (25M Cache, 2.50 GHz) with a 3 GB RAM limit. The operating system is Linux.

VRP with unit-size locker boxes and parcels
First, we test the proposed metaheuristic method for the relaxed problem with unit-size locker box and parcels where no packing decisions have to be taken. This means in particular that the capacity feasibility check, when assigning a customer to a locker box station, boils down to just checking if the available capacity of that station is big enough to take the whole demand of that customer. The step of reducing capacity usage by solving the bin packing set covering problem can be omitted.

Comparison between exact model and metaheuristic method
We use self-generated instances in order to compare the results of the exact model (MIP) for the problem with unit-size locker boxes and parcels of Sect. 3.1 with the results of the metaheuristic method.
The coordinates of the customer locations are generated randomly in a predefined square area within range [0, 100]. The locations of the locker box stations are chosen based on k-means clustering (Lloyd 1982). This means that, depending on the number of locker box stations n we want to have, n clusters are generated and a locker box station is located in the centroid of the corresponding cluster. However, customers need not be assigned to stations according to the clusters. The k-means is just used for defining the locations of the locker box stations.
Travel distance d i j between two points i and j are calculated based on the Euclidean distance.
Each customer requests up to 5 parcels and all parcels have the same size. We assume that a request of 1 parcel is the most likely case (40%), a request of 2 or 3 parcels is medium likely (20%) and 4 or 5 parcels is not so likely (10%).
The length of a customer's time window is between 2 and 4 h. The time window at the depot has a length of 12 h, which represents the maximal route duration. The time window of a locker box station is the same as the one from the depot. A service time of 9 min is assumed for every customer, following the setting of the VRPTW Solomon instances (Solomon 1987). We assume that service at a locker box station takes about twice as long as service at a customer's home address. Thus, the service time for a locker box station is set to 20 min.
At each locker box station there is a number of homogeneous slots available calculated as m * 1.5, where m is the total number of customers.
A customer always accepts home delivery. The number of locker box stations that a customer accepts for delivery is determined randomly. It may also happen that a customer does not accept any locker box delivery. When a customer accepts locker box delivery, he or she always accepts the locker box station that is closest to the home address. Depending on the number of accepted stations, a customer further accepts stations within a certain radius ρ of the home address. This parameter ρ is chosen to be 15 in this study, following the conclusions of Mancini and Gansterer (2020b) suggesting that it represents a reasonable travel distance (in relation to the predefined area) that customers are willing to accept to visit locker box stations for pickup and to study the problem on hand. Moreover, a customer may accept other stations that are chosen randomly and represent stations either close to his or her workplace or a railway station or on the way home.
The compensation cost c for serving a customer at a locker box station is set to 5 for the following experiments. The value of this parameter is again chosen based on the study of Mancini and Gansterer (2020b) providing a reasonable value based on the predefined area and expected total route lengths.
The used set of instances can be found online under Grabenschweiger et al. (2020). The results are given in Table 1. The name of an instance is given by "X_Y_Z", where X, Y and Z correspond to the number of customers, number of locker box stations and instance number respectively. For solving the problem with the MIP, we set a time limit of 12 h and report the best known solution (BKS) that was found in that time. t opt shows the time in seconds that was needed to solve the MIP. "time limit" indicates that no proven optimal solution could be found within the time limit. In that case, the optimality gap is shown in the column opt gap. For the metaheuristic method, we did 5 runs for every instance. f best , f avg and f worst report the total cost of the best, average and worst solution of those 5 runs respectively. The objective values are compared to the BKS from the MIP and the corresponding gaps are given in the column gap best , gap avg and gap worst respectively. t avg gives the average computation time in seconds over the 5 runs.
When we compare the results with respect to solution quality, we can see that the metaheuristic could find the optimal solution for every instance. In some cases, when the MIP did not provide the proven optimal solution, we are able to find new best known solutions.
Concerning computational time, we can see that for the instances with 10 and 14 customers, the MIP is quite fast. However, with increasing instance size, computational time increases considerably. Solving some of the instances with more than 22 customers already takes a few hours and for some instances no optimal solution could be found within the given time limit of 12 h. In comparison, the metaheuristic method needs on average only a couple of seconds to provide comparably good solutions.

Comparison with Mancini/Gansterer
There is a publicly available instance set by Mancini and Gansterer (2020a) which the authors used for solving the "vehicle routing problem with private and shared delivery locations" (see Literature Review in Sect. 2). The characteristics of these instances are such that we can use them for solving the problem with unit-size locker boxes and parcels with our metaheuristic method. We then compare the obtained results to those of Mancini and Gansterer (2020b), for which they used an exact solution approach (MIP) and a matheuristic based on Iterated Local Search (ILS) for solving the given instances.
The instance set comprises of 10 instances with 25 customers, 10 instances with 50 customers and 10 instances with 75 customers. In each instance, 5 locker box stations are available. Customers are randomly distributed in a [10 × 10] area and the depot is located in the southern part of the customer area. The locker box stations are in every instance at the same locations which are in the South-east, South-west, North-east, North-west and centrally. Travel distances between two points are calculated based on the Euclidean distance and multiplied by a factor 3. The maximal duration of a tour (which corresponds to the time window of the depot) is 12 h. The length of the customers' time windows is 60 min.
At each locker box station there is a number of homogeneous slots available which is always slightly larger than the total number of parcels (30, 55 and 80 slots for instances with 25, 50 and 75 customers respectively). There is always one bigger station that has more capacity than the others. The other 4 stations have equal capacity. For example, for the instances with 25 customers, there is 1 station with 10 slots and 4 stations with 5 slots each. For 50 customers, it is 1 station with 15 slots and 4 stations with 10 slots, and for 75 customers, it is 1 station with 20 slots and 4 stations with 15 slots.
We assume that a customer accepts all locker box stations within a radius ρ of the home address. The radius ρ and the compensation value c are taken as before in Sect. 5.1.1, that is, ρ is 15 and c is 5.
In the routing problem of Mancini and Gansterer (2020b), a cost is charged for every vehicle that is used, which is part of the total costs that are to be minimized. However, we do not consider this in our model when comparing the results; we add posteriori the value for each used vehicle to the total cost of our solution.
The results can be seen in Table 2. An instance is denoted as rX_Y_Z, where X, Y, and Z are the number of customers, the number of locker box stations and the instance number respectively.
The results from Mancini and Gansterer (2020b) are shown in the columns MIP and ILS, where for the ILS average results are reported.
Following that, the results that we obtained with the proposed metaheuristic method are shown. We did 5 runs for every instance. f best , f avg , and f worst report the best, average , and worst solution across the performed runs respectively. t avg shows the average runtime in seconds. gap best , gap avg , and gap worst report the percentage gap from the corresponding f to the best known solution from the MIP. We average the values over all instances from one size category (i.e., the 25, 50, and 75 customer instances) and over all instances together.
It can be seen that the metaheuristic method outperforms the ILS with respect to solution quality and is competitive with respect to runtime. For the small instances, we can find the optimal solution for every instance almost every time. For the medium instances, we can find in the best case the optimal solution or BKS for every instance and have on average and, in the worst case, only a very small gap. For the large instances, we can report for 9 out of 10 instances a new BKS.

Extended Mancini/Gansterer instances
In order to perform computational experiments for the VRPHLB as a complete problem, that is, using the multi-size model where packing is considered, we extend the Mancini/Gansterer instances. This means that we take the spatial and time information of those instances and add the needed information for demand and capacity which have to be given now for different sizes.
The size categories of the parcels are small (S), medium (M), and large (L), where S corresponds to a 1 = 1, M to a 2 = 2, and L to a 3 = 4. Each customer has a demand of up to 5 parcels in total, where it is now decided with a uniform probability distribution if total demand is 1, 2, 3, 4 or 5 parcels. The size of each parcel is determined randomly and each size is equally likely to appear.
The locker box stations have slots of size S, M, and L equivalent to the parcel sizes. Thus, slot sizes are b 1 = 1, b 2 = 2, and b 3 = 4. At each station the number of available slots is equal for each size; we take for each size just the capacity that was given in the unit-size setting. So if there were 10 slots for one size, now there are 10 slots of every size.
Everything else remains as described before in Sect. 5.1.2. The instances that were used in the following sections for the VRPHLB can be found under Grabenschweiger et al. (2020).

The cost of packing
We now want to get an insight on what we call the cost of packing. We compare the two different models for tackling the problem with heterogeneous locker boxes, which we formulated in Sect. 3.1 with unit-size locker boxes and parcels and in Sect. 3.2 with multi-size locker box stations and parcels taking into account packing of the different sized parcels.
When we solve the problem with unit-size locker boxes and parcels, we use the instance that was generated with different sizes and add up parcels and capacities. In more detail, we add up the demands of the different sizes of a customer to get the unified, total demand of that customer; we also add up the capacities of the different sizes of a locker box station to get the total capacity of that locker box station.
Then we also solve the problem according to the multi-size model, where packing decisions are included when feasibility of locker box station capacity is checked or when capacity usage is reduced. For the multi-size model, we just take the instance as given. Table 3 shows the obtained results; "f" gives the total cost, "HD" the number of customers that are served at home, "stations" the number of used locker box stations and "vehicles" the number of used vehicles in the solution. All values are averaged over 5 runs.
The percentage gap reports the difference of the multi-size model to the unit-size model in terms of total cost. We can see that, on average, over all instances, the multisize model is 1.57% more expensive than the unit-size model. This value is not as high as one would expect. A reason for that may lie in the compensation cost which may induce in both models equally well balanced solutions with respect to routing and compensation costs, since it is rather the cost for compensating locker box customers than the restriction on locker box capacity that prevents serving more customers at locker box stations. Thus, with a lower compensation value c the cost of packing may increase since capacity constraints may become binding earlier in the multi-size model where feasible packing is required.
Moreover, we observe that the cost of packing is 2.06% on average higher for the small instances compared to 1.32% and 1.34% for the medium and large instances. This could be due to the fact that, for small instances, it is more significant that the different sized parcels are packed properly and feasibly into the different sized slots, while for larger instances this effect evens out over the larger number of locker box customers.
Another observation is that, for some of the large instances, the average number of customers that are served at home in the multi-size model is lower than that in the unit-size model. One may think that, on the contrary, more customers have to be served at home when packing is taken into account, since more capacity may be used, leaving less space for other customers to be served at locker boxes. This unexpected result can be explained with the already described characteristics of the problem, that is, the compensation cost and the lower importance of proper packing with an increasing number of locker box customers. Table 3 The cost of packing

Configuration of locker box stations
We now analyze the relation between different demand scenarios and different locker box station configurations and their impact on total costs. For this purpose, we generate four different demand scenarios varying in the percentage regarding the number of parcels that are of small, medium, and large size.
-demand333: This is the base case of Sect. 5.1.2, where each size is equally likely to appear. Thus, across all customers there will be in total around 1/3 small parcels, 1/3 medium parcels, and 1/3 large parcels. -demand622: 60% S, 20% M, 20% L -demand262: 20% S, 60% M, 20% L -demand226: 20% S, 20% M, 60% L We further generate four different configuration scenarios of locker box stations. For that we assume that the total size of a station remains unchanged across all configuration scenarios. The total size of a station is the one of the base configuration where there are the same number of slots for each size (as described in Sect. 5.2.1). Each slot adds to the total size correspondent to its slot size. Therefore, total size of a station k can be calculated by 1 * Q k1 + 2 * Q k2 + 4 * Q k3 . The four configuration types vary in the share of total size that corresponds to a certain slot size. A different total number of slots arises in every scenario, because, for example, when most of the space is used for small slots, a higher number of slots can become available in total compared to when most of the space is used for large slots.
The four configuration types are as follows: -config333: 1/3 of the total size corresponds to slots of size S, 1/3 to slots of size M, and 1/3 to slots of size L. -config622: 60% S, 20% M, and 20% L -config262: 20% S, 60% M, and 20% L -config226: 20% S, 20% M, and 60% L Each demand scenario is then compared to each configuration. The proposed metaheuristic is used to solve the problems according to the multi-size model with packing.
The results can be found in Table 4. The values shown are the total costs averaged over all instances of the group of small, medium, large instances. Total average values over all instances are also given for each demand scenario and each configuration.
One can observe that across all demand scenarios, config226, that is, the one where most of the space is used for large slots, always gives the best solutions. The second best configuration is config333 where the space is shared equally between the different slot sizes. The configurations config622 and config262 are worst in every demand scenario. The cost difference between the different configurations becomes bigger when more parcels are of medium or large size.
The observation that config226 works best for demand scenarios with a lot of large parcels is not surprising, since large parcels can only be packed into large slots. However, the fact that it also performs well for scenario demand622, where most of the parcels are small, is somewhat unexpected. One possible explanation can be that in instances with many small parcels also, there will likely be many customers that have many small parcels and these can be stored together in a large box. Furthermore, even when 60% of the space is used for large slots and 20% for small and medium slots, a comparably high number of small and medium slots is still available since they do not need to use so much space.
Moreover, when comparing a configuration across the four demand scenarios, it can be seen that total costs increase when more parcels of larger size are present in the problem. This can be, on the one hand, due to the fact that more locker box stations may be needed to serve a comparable number of customers at locker boxes. On the other hand, capacity at locker box stations may be used up earlier, thus allowing for less customers to be served at locker boxes, which may lead to a worse trade-off between routing and compensation costs.

Conclusion
We presented the vehicle routing problem with heterogeneous locker boxes. A mathematical formulation was given for the problem with unit-size locker boxes and parcels where no packing decisions have to be taken and for the complete problem with multisize locker box stations and parcels where packing is part of the problem. Furthermore, a metaheuristic solution method was proposed. The numerical experiments for small self-generated instances with unit-size locker boxes and parcels show that the metaheuristic provides for every instances the optimal solution or a new BKS (in case the exact model did not find the proven optimal solution within the given time limit of 12 h). We solved instances from the literature for the problem with one request and one size and compared the results obtained with our proposed metaheuristic method to those from the literature. Our metaheuristic method performs better with respect to solution quality and it is competitive with respect to runtime. It can find the optimal solution for every instance or can give a new BKS. For the problem with different locker box and parcel sizes, we extended the instances from the literature and added information about demand for different sized parcels and capacity information for different sized slots at the locker box stations. We then compared the unit-size model (without packing) to the multi-size model (with packing) and found that the cost of packing is 1.57% which is not as high as one would expect. However, we concluded that this behavior depends on the choice of the compensation cost value c; with a lower compensation value the cost of packing may increase, since capacity constraints may become more binding in the multi-size model where feasible packing is required. As the last part of our computational study, we generated different configuration settings for the locker box stations, which vary in terms of how the total space of a station is shared among slots of the different size categories and different demand scenarios which vary in the percentage of the number of parcels that are of small, medium, and large size. We then compared each configuration to each demand scenario and observed that across all demand scenarios, the configuration where most of the space is used for large slots gives on average the best solution. The second best configuration is one where the space is shared equally between the different slot sizes. The configurations where most of the slots are small or medium sized should not be chosen in the given setting. Moreover, the results show that total costs are highest when we face a demand scenario with many large parcels.