The consistent vehicle routing problem under uncertain environment

Abstract

In the increasingly competitive logistics industry, much more emphasis has been put on customer satisfaction by firms to distinguish themselves from their competitors. Motivated by the practices, the consistent vehicle routing problem (ConVRP) incorporates service consistency into the vehicle routing problem to improve customer satisfaction. However, the majority of the existing research considers the ConVRP in a deterministic environment while the uncertainties are not fully studied. Therefore, this paper solves the consistent vehicle routing problem under uncertain environment (UnConVRP) taking into account uncertain customer demands, travel times, and service times. Over a multi-day planning horizon, customers may have multi-day or single-day service requirements. The same driver is assigned to each customer almost at the same time each day when customers require service. The objective is to design the routes for vehicles over the planning horizon under uncertain environment while maintaining service consistency. Two uncertain programming models are established based on different decision criteria and the crisp equivalents are proposed using uncertainty theory. An efficient template-based solution framework is designed to solve the models where the artificial bee colony algorithm is embedded. Initially, the template route is constructed to guarantee the service consistency of frequent customers. Then the final daily routes can be derived from the template route. Finally, numerical experiments are performed to show the effectiveness of the proposed algorithm.

Keywords

Vehicle routing consistent service uncertain programming artificial bee colony

1 Introduction

With the intensified competition and the increasing demand of customers, customer satisfaction has received substantial interest in recent years, as this has a direct impact on customer loyalty. In many industries, customer satisfaction has been seen as a key differentiator of business to gain a competitive advantage. Consistency with respect to the service starting time and the service provider can enhance brand loyalty and customer satisfaction by allowing the employees to develop relations and create bonds with the customers. Consistency in customer service also boosts high-quality service and, hence, high customer satisfaction. According to Wong [46], customer schedule consistency is of great importance in small packages delivery. By assigning the same driver to each customer, drivers can gain valuable familiarity with customers, which will improve the service and contribute to additional business profits. The practice of UPS also demonstrates the value of service consistency [45]. As mentioned in Woodward et al. [47], home care clients whose care provision with high consistency of personnel and timing reported the greatest satisfaction with service delivery. As a consequence, in recent years, more and more companies begin to focus on customer relationship management. By providing consistent service, the relationship between customers and service personnel can be enhanced, which in turn can positively increase sales. Furthermore, the efficiency of service can be improved by the familiarity of customers and drivers. Therefore, the company can enjoy the benefits of both increased sales and efficiency from the provision of consistent service.

Vehicle routing problem (VRP) with consistency considerations has received substantial interest in recent years because of the practical importance of providing consistent service. However, providing a high level of service consistency will cause an increase in the routing cost. In practice, if a driver is assigned to specific several customers, the total demand of the customers may surpass the vehicle capacity. To deal with this problem, one way is to decrease the number of customers assigned to a driver, which in turn will increase the manpower cost. But if customers can be allowed to serve by different drivers over the planning horizon, the benefits of service consistency will be violated. Therefore the tradeoff between the service consistency and the routing cost propels companies to explore the most appropriate routes for drivers to provide consistent service for customers at reasonable costs, in other words, to deal with the consistent vehicle routing problem (ConVRP).

To the best of our knowledge, a series of papers have studied the ConVRP. However, few of them involve uncertainty, which, indeed, challenges routes planning and execution. Actually, uncertainty arises from several aspects in the practical operations of small package delivery. For example, the demands are difficult to calculate before the requests arrive since there is no historical data for new customers. Besides, on a specific day, there may be some unpredictable conditions happening to drivers or vehicles or roads, which will lead to the uncertainty of the travel time between two locations. In addition, due to the unknown customer demands or indeterminacy related to customers, the precise service times of customers are not available until the service is completed. Therefore, in order to improve the robustness of the ConVRP, uncertainty needs to be taken into account. In this paper, we put forward the uncertain consistent vehicle routing problem (UnConVRP) with uncertain customer demands, travel times, and service times.

As we mentioned above, due to the lack of historical data or insufficient information, uncertainty will appear and is different from randomness, which can be tackled employing probability theory. While faced with this kind of uncertainty, probability theory may not be applicable in the absence of sufficient observed data. To deal with this type of uncertainty, uncertainty theory established by Liu [24] is a useful tool to model indeterminacy. This paper deals with uncertain variables utilizing uncertainty theory and establishes two types of uncertain programming models for the UnConVRP. According to the operational law of independent uncertain variables, the proposed models can be respectively transformed into the corresponding crisp equivalent models. Due to the NP-hardness of UnConVRP, we design a template-based artificial bee colony (TABC) algorithm to address this problem. A series of numerical experiments are conducted to show the performance of the proposed algorithm.

The remainder of this paper is organized as follows. Section 2 provides a review of the related literature about the consistent vehicle routing problem and uncertainty theory. Section 3 formally defines the problem and presents the formulations of mathematical models. The algorithm designed for the problem is portrayed in Section 4. In Section 5, some numerical experiments are conducted to show the performance of the proposed algorithm. Finally, the conclusion is drawn in Section 6.

2 Related literature

Recently, the importance of service consistency has gained extensive attention [18 , 46]. Examples of real-life applications where such consistent policies have been applied include, among others, courier services [14, 43], vendor-managed inventory distribution [5 , 33], home care and nursing services for the elderly [1, 11].

The first attempt in the literature towards taking consistent service into vehicle routing problem is presented by Groër et al. [14], introducing a new VRP variant named the ConVRP. In this paper, service consistency is imposed by the same driver (driver consistency) at almost the same time (arrival time consistency). A formulation of ConVRP as a mixed-integer program is provided, and a ConVRP record-to-record travel (ConRTR) heuristic is employed to obtain near-optimal solutions quickly. Based on the template route, the solution framework can guarantee driver consistency. Two additional papers adopt the concept of template route to get a better solution for the ConVRP [19, 44], the template-based tabu search algorithm and the template-based adaptive large neighborhood search, respectively. Xu and Cai [48] propose a variable neighborhood search algorithm with a new shaking method outperforming template-based approaches. The first exact solution for the ConVRP is offered by Goeke et al. [13], the proposed clustered column generation is demonstrated to be able to obtain the exact solution of medium-sized instances with 30 customers over a five-day horizon.

In the basic ConVRP, customers are served by the same driver whenever they need service. This definition of driver consistency has been adopted and modeled by many researchers [2 , 39]. In terms of generalized driver consistency, some scholars extend the service consistency to a broader level, in which each customer can be visited by a limited number of different drivers [1 , 27], in which a large neighborhood search heuristic is applied by the first two papers. Luo et al. [27] propose a three-stage heuristic called the decomposition, repair, and distance reduction approach. Recently, customers are divided into two groups according to their requirements of driver consistency: hard driver consistency (the same driver) and relaxed driver consistency (limited number of different drivers) [10]. Most of the existing studies define arrival time consistency as the maximum arrival time difference over a planning horizon [20 , 40–42]. Feillet et al. [11] introduce the concept of time-classes to refine the arrival time difference, minimizing the total number of different time classes where each customer is located, and the maximum number of time classes is bound for each customer.

Service consistency has been modelled with different methods in prior literature. Papers mentioned above model the ConVRP as a single-objective problem, in which consistency is considered as constraints. As the ConVRP involves several conflicting objectives, a multi-objective optimization approach facilitates the trade-off analysis. In some related literature, service consistency is also modelled as weighted objective [20, 33] or independent objective [21, 23]. In multi-objective ConVRP, improving driver consistency and arrival time consistency, and minimizing routing cost are modeled as independent objectives of the problem. In the context of multi-objective optimization, two solutions can’t be compared directly, thus Pareto dominance is employed to compare the quality of different solutions. Kovacs et al. [21] propose two exact approaches based on the ε-constraint framework and one heuristic. Exact approaches can be valid when dealing with small problem instances. And for large-scale problems, the multi-directional large neighborhood search combining the multi-directional local search and the large neighborhood search can provide good solutions. The trade-off analysis of total travel time and service consistency has also been performed on large problem instances. Subsequently, an improved multi-directional local search algorithm is designed by Lian et al. [23] to solve the multi-objective model.

There are also some extensions involving the service consistency in vehicle routing problems. Coelho et al. [5] propose delivery quantity consistency in inventory routing problems to improve the quality of service, and the consistent inventory routing problem is solved by adaptive large neighborhood search heuristic. A unified inventory routing model with a homogeneous and heterogeneous fleet, transshipment options, and consistency features is provided, and a branch-and-cut algorithm is designed to solve the several above classes of multi-vehicle inventory routing problems in Coelho et al. [6]. Furthermore, the multi-product multi-vehicle inventory-routing problem is introduced with consistency features in Coelho et al. [7]. Other variants have also been incorporated in the ConVRP context, including ConVRP with profits [40], the collaborative ConVRP with workload balance [28], ConVRP with shift length [22] and ConVRP with simultaneous pickup and delivery [9, 50]. To solve these complicated ConVRPs, an adaptive tabu search [40], an iterated local search [28], a variable neighborhood search [22], and three heuristics including the record-to-record travel algorithm-based solution method, the local search with variable neighbourhood search-based solution method, and the tabu search-based method [50] are developed, respectively.

Uncertainty theory, proposed in Liu [24] and refined in Liu [26], provides a new non-probabilistic method to cope with the empirical data in optimization models. Subsequently, it has been widely studied by scholars and applied in social network problem [29 –31], transportation problem [3, 12], pricing decision [4 , 17], product configuration [36 –38] etc. Liu [26] firstly proposes an uncertain vehicle routing problem and models it by uncertain programming in which the travel times are assumed to be uncertain variables with known uncertainty distributions. Based on the above proposed model, Zhang et al. [49] introduce two models: maximum the measure that each customer is visited within its time window and minimize the expected value of the arrival time. These models then are transformed into their deterministic form and solved by the genetic algorithm. Furthermore, Ning et al. [32] assume that the travel times between two customers are uncertain variables and propose an uncertain multilevel programming model for a vehicle routing problem, of which the leader’s objective is to minimize the total waiting times of the customers, and the follower’s objective is to minimize the waiting times of the vehicles for the beginning moments of the customers’ time windows.

By reviewing the related literature, we can conclude that the literature in recent years has flourished about the ConVRP but little focus on uncertainty, which definitely appears in real life operations. As mentioned in Kovacs et al. [18], uncertain travel times have not yet been incorporated into the ConVRP. Thus, in this paper, we investigate the ConVRP where uncertain customer demands, travel times, and service times are taken into account.

3 Problem description

In this paper, we consider a logistics company serving customers in a certain period. In the planning horizon, the company arranges the drivers to visit the customers who need service on some days of the service period. In terms of cultivating customer relationships, the customers who require service frequently must be visited by the same driver at the approximately same time over the given horizon. Therefore, the firm must plan the most appropriate routes which can meet the requirement of service consistency and achieve the shortest total routing time. In our mathematical model, we will show how to find the optimal routes, with service consistency modelled as constraints.

A route must be developed for each vehicle so that all customers are serviced and the total travelling time of the fleet is minimized. The UnConVRP is defined on a complete directed graph G = (N, A), where N ={ 0, 1, 2, . . . , |N|} is the nodes set, node 0 denotes the origin and the destination depot and the other nodes correspond to the customers, and A = {(i, j) |i, j ∈ N, i ≠ j} is the set of arcs. Customers are visited by a fleet of vehicles in the set of K ={ 1, 2, . . . , |K|}, each vehicle k ∈ K has an upper boundary of capacity Q and the maximum operating time T on each day. D ={ 1, 2, . . . , |D|} represents the set of days in the planning period. Each customer i ∈ N ∖ {0} is associated with a demand ${\tilde{q}}_{id}$ and a service time ${\tilde{s}}_{id}$ , which are uncertain in our model and characterized as uncertain variables. We also assume the travel times ${\tilde{t}}_{ijd}$ of each arc are uncertain variables. A binary variable w_id denotes whether customer i ∈ N ∖ {0} require service (w_id = 1) or not (w_id = 0) on day d. L is the arrival time variation (between the latest and the earliest arrival time) of each customer over |D| days.

The objective is to determine a set of routes over |D| days such that the total traveling time is minimized, satisfying the following constraints:

each vehicle starts from the depot at time 0 and ends at the depot on all days

vehicles are homogeneous and capacitated, subject to a certain capacity boundary Q

a maximum operating time is imposed on all vehicles, at most T time units

customers must be visited exactly once on each day they need service

frequent customers must be visited by the same driver over the given planning period

the maximum difference between the earliest and latest vehicle arrival times over the planning horizon does not exceed L time units

In this paper, a four-index formulation is introduced. Two groups of decision variables associated with the customer visiting sequence, the customers’ assignment to vehicle routes are utilized. The binary decision variable x_ijkd equals 1 if customer j is serviced after customer i by vehicle k on day d and equals 0 otherwise. Binary variable y_ikd equals 1 if customer i is visited by vehicle k on day d and equals 0 otherwise. Continuous variable a_id is the arrival time at customer i ∈ N ∖ {0} on day d.

The notations of input parameters and decision variables are provided in the following Table 1.

Table 1
Notations

Indices and sets

i, j index of node (customer or depot).

N set of nodes and the depot {0}.

d, c index of day.

D set of service days.

k index of vehicle.

K set of vehicles.

Parameters

Q capacity of each vehicle.

T maximum travel time of each vehicle.

L maximum arrival time variation between the latest and the earliest arrival time).

w _id equals 1 if customer i asks for service on day d, otherwise 0.

${\tilde{t}}_{ijd}$ travel time between any two nodes i and j on day d.

${\tilde{q}}_{id}$ demand at customer i on day d.

${\tilde{s}}_{id}$ the service times at customer i on day d.

Decision variables

x _ijkd binary variable, equals 1 if customer j is serviced after customer i by vehicle k on day d and equals 0 otherwise.

y _ikd binary variable, equals 1 if customer i is visited by vehicle k on day d and equals 0 otherwise.

a _id arrival time at customer i on day d, equals 0 if customer i asks no service on day d.

Indices and sets
i, j	index of node (customer or depot).
N	set of nodes and the depot {0}.
d, c	index of day.
D	set of service days.
k	index of vehicle.
K	set of vehicles.
Parameters
Q	capacity of each vehicle.
T	maximum travel time of each vehicle.
L	maximum arrival time variation between the latest and the earliest arrival time).
w _id	equals 1 if customer i asks for service on day d, otherwise 0.
${\tilde{t}}_{ijd}$	travel time between any two nodes i and j on day d.
${\tilde{q}}_{id}$	demand at customer i on day d.
${\tilde{s}}_{id}$	the service times at customer i on day d.
Decision variables
x _ijkd	binary variable, equals 1 if customer j is serviced after customer i by vehicle k on day d and equals 0 otherwise.
y _ikd	binary variable, equals 1 if customer i is visited by vehicle k on day d and equals 0 otherwise.
a _id	arrival time at customer i on day d, equals 0 if customer i asks no service on day d.

We first deal with the arrival times a_id. Obviously, a_id is the arrival time function of vehicles at customer i on day d and it is determined by the decision variable x_ijkd and uncertain variables ${\tilde{s}}_{id}$ and ${\tilde{t}}_{ijd}$ . Because the arrival time at each customer heavily depends on the routing plan. Here we assume a delivery always starts at the time 0 from the depot every day, a_0d = 0. If customer i is served before customer j by vehicle k, then we have $a_{jd} = a_{id} + {\tilde{s}}_{id} + {\tilde{t}}_{ijd} .$ This recursive process may produce the arrival time at each customer each day. Dependent on the operational decision x_ijkd, arrival times can be obtained from the following equation $a_{jd} = \sum_{i \in N} x_{ijkd} (a_{id} + {\tilde{s}}_{id} + {\tilde{t}}_{ijd}) .$ Assume that ${\tilde{t}}_{ijd}$ and ${\tilde{s}}_{id}$ are independent uncertain variables with regular uncertainty distributions Φ_{t
_ijd} and Φ_{s
_id}, respectively. According to Theorem 1, the arrival time a_jd is an uncertain variable whose inverse uncertainty distribution is $Ψ_{a_{jd}}^{- 1} (α) = \sum_{i \in N} x_{ijkd} (Ψ_{a_{id}}^{- 1} (α) + Φ_{s_{id}}^{- 1} (α) + Φ_{t_{ijd}}^{- 1} (α)) .$ Where α is the predefined confidence level.

Since an uncertain objective function cannot be directly minimized and the uncertain constraints do not define a crisp feasible set, the model involving the uncertain variables, the travel time ${\tilde{t}}_{ijd}$ , the service time ${\tilde{s}}_{id}$ , and the customer demand ${\tilde{q}}_{id}$ , is not well defined. Next, we adopt two decision criteria and propose the corresponding uncertain programming models to solve the problem proposed in this paper.

3.1 Expected value model

The first model is the expected value model (EVM), which is to optimize the expected value of uncertain variables in an indeterminate environment. We also propose some crisp equivalents when assuming that the uncertain variables have particular distributions.

The objective of the model is to minimize the expected total travelling time, given constraints mentioned above, we formulate the EVM for the UnConVRP as follows: $min E [\sum_{d \in D} \sum_{k \in K} \sum_{i \in N} \sum_{j \in N} {\tilde{t}}_{ijd} x_{ijkd}]$ (1) $\begin{matrix} s . t . \\ y_{0 kd} = 1, \forall k \in K, d \in D, \end{matrix}$ (2) $\sum_{k \in K} y_{ikd} = w_{id}, \forall i \in N ∖ {0}, d \in D,$ (3) $\begin{matrix} \sum_{j \in N} x_{ijkd} = \sum_{j \in N} x_{jikd} = y_{ikd}, \\ \forall i \in N, k \in K, d \in D, \end{matrix}$ (4) $\begin{matrix} M {\sum_{i \in N ∖ {0}} {\tilde{q}}_{id} y_{ikd} \leq Q} \geq β_{kd}, \\ \forall k \in K, d \in D, \end{matrix}$ (5) $\begin{matrix} M {a_{id} + {\tilde{s}}_{id} + {\tilde{t}}_{i 0 d} \leq T} \geq γ_{id}, \\ \forall i \in N ∖ {0}, d \in D, \end{matrix}$ (6) $\begin{matrix} M {(a_{id} - a_{ic}) w_{id} w_{ic} \leq L} \geq δ_{i}, \\ \forall i \in N ∖ {0}, d, c \in D, d \neq c, \end{matrix}$ (7) $\begin{matrix} | y_{ikd} - y_{ikc} | \leq 2 - w_{id} - w_{ic}, \\ \forall i \in N ∖ {0}, k \in K, d, c \in D, \end{matrix}$ (8) $\begin{matrix} \sum_{i, j \in S} x_{ijkd} \leq | S | - 1, \\ \forall S \subseteq N, 2 \leq | S | \leq n, k \in K, d \in D, \end{matrix}$ (9) $\begin{matrix} x_{ijkd} \in {0, 1}, \forall i, j \in N, k \in K, d \in D, \end{matrix}$ (10) $y_{ikd} \in {0, 1}, \forall i \in N, k \in K, d \in D,$ (11) $a_{id} \geq 0, \forall i \in N, d \in D,$ (12) where β_kd, γ_id, and δ_i represent the predetermined confidence levels for all i ∈ N, k ∈ K, d ∈ D.

The objective function (1) minimizes total travel time of all vehicles over the planning horizon. Equation (2) ensures that all vehicles departs from the depot. Equation (3) guarantees that customers only can be visited exactly once when they require service. Equation (4) is flow conservation constraint at each node. Equation (5) makes sure that the total loads of each vehicle should not be larger than its capacity with the predetermined confidence level β_kd. Equation (6) defines the maximum operating time of each vehicle everyday with the predetermined confidence level γ_kd. Arrival time consistency constraint is ensured in Equation (7) with the predetermined confidence level δ_i. Driver consistency is enforced in Equation (8). Equation (9) prevents the subtours. The domain of decision variables is defined in Equation (10), Equation (11) and Equation (12).

Assume that ${\tilde{t}}_{ijd}$ , ${\tilde{a}}_{id}$ , and ${\tilde{q}}_{id}$ are independent uncertain variables with regular uncertainty distributions Φ_{t
_ijd}, Φ_{s
_id}, and Φ_{q
_id}, respectively. The uncertain programming model can be converted to the corresponding deterministic equivalents.

According to Theorem 3, $\begin{matrix} E [\sum_{d \in D} \sum_{k \in K} \sum_{i \in N} \sum_{j \in N} {\tilde{t}}_{ijd} x_{ijkd}] \\ = \sum_{d \in D} \sum_{k \in K} \sum_{i \in N} \sum_{j \in N} x_{ijkd} E [{\tilde{t}}_{ijd}] . \end{matrix}$ It follows from Theorem 2 that the expected value is $E [{\tilde{t}}_{ijd}] = \int_{0}^{1} Φ_{t_{ijd}}^{- 1} (α) d α .$ Then we have $\begin{matrix} E [\sum_{d \in D} \sum_{k \in K} \sum_{i \in N} \sum_{j \in N} {\tilde{t}}_{ijd} x_{ijkd}] \\ = \sum_{d \in D} \sum_{k \in K} \sum_{i \in N} \sum_{j \in N} x_{ijkd} \int_{0}^{1} Φ_{t_{ijd}}^{- 1} (α) d α . \end{matrix}$ According to Theorem 4, we can transform the constraint $M {\sum_{i \in N ∖ {0}} {\tilde{q}}_{id} y_{ikd} \leq Q} \geq β_{kd}$ into a deterministic form $\sum_{i \in N ∖ {0}} y_{ikd} Φ_{q_{id}}^{- 1} (β_{kd}) \leq Q .$ In a similar manner, the constraint $M {a_{id} + {\tilde{s}}_{id} + {\tilde{t}}_{i 0 d} \leq T} \geq γ_{id}$ can be converted into $Ψ_{a_{id}}^{- 1} (γ_{id}) + Φ_{sid}^{- 1} (γ_{id}) + Φ_{ti 0 d}^{- 1} (γ_{id}) \leq T .$ And the constraint $M {(a_{id} - a_{ic}) w_{id} w_{ic} \leq L} \geq δ_{i}$ is equivalent to $(Ψ_{a_{id}}^{- 1} (δ_{i}) - Ψ_{a_{ic}}^{- 1} (1 - δ_{i})) w_{id} w_{ic} \leq L .$

Then the deterministic transformations of EVM are as follows: $min \sum_{d \in D} \sum_{k \in K} \sum_{i \in N} \sum_{j \in N} x_{ijkd} \int_{0}^{1} Φ_{t_{ijd}}^{- 1} (α) d α$ (13) $\begin{matrix} s . t . \\ y_{0 kd} = 1, \forall k \in K, d \in D, \end{matrix}$ (14) $\sum_{k \in K} y_{ikd} = w_{id}, \forall i \in N ∖ {0}, d \in D,$ (15) $\begin{matrix} \sum_{j \in N} x_{ijkd} = \sum_{j \in N} x_{jikd} = y_{ikd}, \\ \forall i \in N, k \in K, d \in D, \end{matrix}$ (16) $\begin{matrix} \sum_{i \in N ∖ {0}} y_{ikd} Φ_{qid}^{- 1} (β_{kd}) \leq Q, \\ \forall k \in K, d \in D, \end{matrix}$ (17) $\begin{matrix} Ψ_{a_{id}}^{- 1} (γ_{id}) + Φ_{sid}^{- 1} (γ_{id}) + Φ_{ti 0 d}^{- 1} (γ_{id}) \leq T, \\ \forall i \in N ∖ {0}, d \in D, \end{matrix}$ (18) $\begin{matrix} (Ψ_{a_{id}}^{- 1} (δ_{i}) - Ψ_{a_{ic}}^{- 1} (1 - δ_{i})) w_{id} w_{ic} \leq L, \\ \forall i \in N ∖ {0}, d, c \in D, d \neq c, \end{matrix}$ (19) $\begin{matrix} | y_{ikd} - y_{ikc} | \leq 2 - w_{id} - w_{ic}, \\ \forall i \in N ∖ {0}, k \in K, d, c \in D, \end{matrix}$ (20) $\begin{matrix} \sum_{i, j \in S} x_{ijkd} \leq | S | - 1, \\ \forall S \subseteq N, 2 \leq | S | \leq n, k \in K, d \in D, \end{matrix}$ (21) $\begin{matrix} x_{ijkd} \in {0, 1}, \forall i, j \in N, k \in K, d \in D, \end{matrix}$ (22) $y_{ikd} \in {0, 1}, \forall i \in N, k \in K, d \in D,$ (23) $a_{id} \geq 0, \forall i \in N, d \in D,$ (24) An illustrative example of the crisp equivalents for the model is given in the following if uncertain variables are of particular uncertainty distributions. Assume that ${\tilde{t}}_{ijd}$ , ${\tilde{s}}_{id}$ , and ${\tilde{q}}_{id}$ are independent normal uncertain variables $N (e_{t_{ijd}}, σ_{t_{ijd}})$ , $N (e_{s_{id}}, σ_{s_{id}})$ , $N (e_{q_{id}}, σ_{q_{id}})$ , respectively, then the objective function is equivalent to $min \sum_{d \in D} \sum_{k \in K} \sum_{i \in N} \sum_{j \in N} x_{ijkd} \int_{0}^{1} (e_{t_{ijd}} + \frac{σ_{t_{ijd}} \sqrt{3}}{π} ln \frac{α}{1 - α}) d α,$ and the chance constraints are equivalent to $\sum_{i \in N ∖ {0}} y_{ikd} (e_{q_{id}} + \frac{σ_{q_{id}} \sqrt{3}}{π} ln \frac{β_{kd}}{1 - β_{kd}}) \leq Q,$ $\begin{matrix} Ψ_{a_{id}}^{- 1} (γ_{id}) + (e_{s_{id}} + \frac{σ_{s_{id}} \sqrt{3}}{π} ln \frac{γ_{id}}{1 - γ_{id}}) \\ + (e_{t_{i 0 d}} + \frac{σ_{t_{i 0 d}} \sqrt{3}}{π} ln \frac{γ_{id}}{1 - γ_{id}}) \leq T, \end{matrix}$

$(Ψ_{a_{id}}^{- 1} (δ_{i}) - Ψ_{a_{ic}}^{- 1} (1 - δ_{i})) w_{id} w_{ic} \leq L .$

And $Ψ_{a_{id}}^{- 1} (δ_{i})$ can be deduced by

$\begin{matrix} Ψ_{a_{id}}^{- 1} (δ_{i}) = & \sum_{j \in N} x_{jikd} Ψ_{a_{jd}}^{- 1} (δ_{j}) + e_{s_{id}} + e_{t_{ijd}} \\ + \frac{(σ_{s_{id}} + σ_{t_{jid}}) \sqrt{3}}{π} ln \frac{δ_{i}}{1 - δ_{i}} . \end{matrix}$

3.2 Chance-constrained programming

The second model is the chance-constrained programming (CCP) model. Under uncertain environment, the decision maker hopes to minimize the pessimistic value with a predefined confidence level α. So the objective function in CCP can be established as follows: $min \bar{f} .$

We have the chance constraint: $M {\sum_{d \in D} \sum_{k \in K} \sum_{i \in N} \sum_{j \in N} {\tilde{t}}_{ijd} x_{ijkd} \leq \bar{f}} \geq α .$ Other constraints are the same as those shown in EVM. α is the predetermined confidence level, and $\bar{f}$ is the α-pessimistic value of the objective value.

We formulate the following chance-constrained programming model for the UnConVRP: $min \bar{f}$ (25) $\begin{matrix} s . t . \\ M {\sum_{d \in D} \sum_{k \in K} \sum_{i \in N} \sum_{j \in N} x_{ijkd} {\tilde{t}}_{ijd} \leq \bar{f}} \geq α, \end{matrix}$ (26) $y_{0 kd} = 1, \forall k \in K, d \in D,$ (27) $\sum_{k \in K} y_{ikd} = w_{id}, \forall i \in N ∖ {0}, d \in D,$ (28) $\begin{matrix} \sum_{j \in N} x_{ijkd} = \sum_{j \in N} x_{jikd} = y_{ikd}, \\ \forall i \in N, k \in K, d \in D, \end{matrix}$ (29) $\begin{matrix} M {\sum_{i \in N ∖ {0}} {\tilde{q}}_{id} y_{ikd} \leq Q} \geq β_{kd}, \\ \forall k \in K, d \in D, \end{matrix}$ (30) $\begin{matrix} M {a_{id} + {\tilde{s}}_{id} + {\tilde{t}}_{i 0 d} \leq T} \geq γ_{id}, \\ \forall i \in N ∖ {0}, d \in D, \end{matrix}$ (31) $\begin{matrix} M {(a_{id} - a_{ic}) w_{id} w_{ic} \leq L} \geq δ_{i}, \\ \forall i \in N ∖ {0}, d, c \in D, d \neq c, \end{matrix}$ (32) $\begin{matrix} | y_{ikd} - y_{ikc} | \leq 2 - w_{id} - w_{ic}, \\ \forall i \in N ∖ {0}, k \in K, d, c \in D, \end{matrix}$ (33) $\begin{matrix} \sum_{i, j \in S} x_{ijkd} \leq | S | - 1, \\ \forall S \subseteq N, 2 \leq | S | \leq n, k \in K, d \in D, \end{matrix}$ (34) $\begin{matrix} x_{ijkd} \in {0, 1}, \forall i, j \in N, k \in K, d \in D, \end{matrix}$ (35) $y_{ikd} \in {0, 1}, \forall i \in N, d \in D,$ (36) $a_{id} \geq 0, \forall i, j \in N, k \in K, d \in D,$ (37) where α, β_kd, γ_id, and δ_i are the predetermined confidence levels ∀i ∈ N, k ∈ K, and d ∈ D.

Based on uncertainty theory, we can also convert the constraint $M {\sum_{d \in D} \sum_{k \in K} \sum_{i \in N} \sum_{j \in N} x_{ijkd} {\tilde{t}}_{ijd} \leq \bar{f}} \geq α$ into a deterministic constraint: $\sum_{d \in D} \sum_{k \in K} \sum_{i \in N} \sum_{j \in N} x_{ijkd} Φ_{t_{ijd}}^{- 1} (α) d α \leq \bar{f} .$

In the same way, assume that ${\tilde{t}}_{ijd}$ are independent normal uncertain variables $N (e_{t_{ijd}}, σ_{t_{ijd}})$ , then the constraint is equivalent to $\sum_{d \in D} \sum_{k \in K} \sum_{i \in N} \sum_{j \in N} x_{ijkd} (e_{t_{ijd}} + \frac{σ_{t_{ijd}} \sqrt{3}}{π} \ln \frac{α}{1 - α}) \leq \bar{f} .$

4 Solution framework

As we have mentioned, traditional VPR is an NP-hard problem, and some algorithms designed to solve it are very complicated. As the ConVRP is more constrained in terms of service consistency, adding the driver consistency and time consistency constraints to the traditional VRP, we need to design an effective algorithm to solve it.

In this section, we design a template-based solution method to work out our problem. The concept of template route introduced by Groër et al. [14] features by a precedence principle to guarantee the service consistency: if customer a and b ask service on the same day, and customer b is serviced after customer a by the same vehicle, then customer b will always be serviced after customer a on all days they both require service. Following this precedence principle, the customers who require service on multiple days will always be assigned to the same vehicle. While for customers who have a single-day service requirement, driver consistency can always be achieved as they are only visited once. Therefore, the template route contains only customers who ask for service more than one day over the whole planning days. We refer them to frequent customers and define customers who require service only once as non-frequent customers. Using the template route, daily routes adhering to service consistency for both frequent and non-frequent customers can be derived.

On the basis of what we described above, our algorithm adopts a two-stage algorithmic scheme to solve the UnConVRP. The first stage involves the generation of the template route for those frequent customers. The artificial bee colony (ABC) algorithm is utilized to improve the template route. Although the neighborhood operators operate directly on the template route, we evaluate the feasibility and quality of the neighborhood solution according to the final daily routes. In the second stage, daily routes are derived from the generated template route through deleting the customers who do not ask for service on a specific day and inserting customers who require service on some day. Efficiency of this algorithmic structure greatly depends on the structure of the template route in the first stage, so ABC is utilized to improve the template routes. A Greedy method is employed in the second stage to construct daily routes, each non-frequent customer is inserted at his best position, which causes the least increase in total travel time.

The framework of the template-based artificial bee colony algorithm is provided below.

Algorithm 1 Template-based ABC framework

Input: t_ijd, s_id, q_id, w_id, Q, T, L and the parameters of ABC algorithm

Output: the set of daily routes R^best and its objective value f (R^best)

1: Initialization phase

2: Generate initial template routes;

3: Template routes improvement phase using ABC algorithm

4: return The best template route $R_{tem}^{best}$

5: while d≤ Day do

6: Construct an initial daily route R_d

7: Deletion and insertion operation

8: Daily routes improvement by a greedy method

9: end while

10: return the best daily routes R^best

4.1 Initialization phase

As we mentioned above, the template route only contains frequent customers so the first step is to divide customers into frequent customers and non-frequent customers and define the relevant input parameters associated with the frequent customers in the generation of template route.

An example is illustrated in the following Table 2, there are 10 customers over a 3-days planning horizon. The first three rows denote whether the service requirement of a customer occurs someday (value 1) or not (value 0). The last but one row gives the sum of the number of service days over the planning horizon and the last row indicates the customers considered in the template construction. Customers who ask for service more than once are contained in the template (value 1), the others are not considered (value 0). Then we get the set of customers (1, 2, 3, 4, 6, 8, 10) in the template.

Table 2
Customers division

Customer 1 2 3 4 5 6 7 8 9 10

Day 1 1 0 1 1 0 1 1 1 0 1

Day 2 0 1 1 1 1 0 0 1 1 1

Day 3 1 1 0 1 0 1 0 1 0 1

Number of service days 2 2 2 3 1 2 1 3 1 3

Template 1 1 1 1 0 1 0 1 0 1

Customer	1	2	3	4	5	6	7	8	9	10
Day 1	1	0	1	1	0	1	1	1	0	1
Day 2	0	1	1	1	1	0	0	1	1	1
Day 3	1	1	0	1	0	1	0	1	0	1
Number of service days	2	2	2	3	1	2	1	3	1	3
Template	1	1	1	1	0	1	0	1	0	1

Next, for each template customer, we compute a template demand ${\tilde{q}}_{i}^{tem}$ and a template service time ${\tilde{s}}_{i}^{tem}$ . For the travel times, we compute a template value ${\tilde{t}}_{ij}^{tem}$ . These values are used throughout the template construction. They are set to the mean values of frequent customers across all days. For example, we deal with the template demand ${\tilde{q}}_{i}^{tem}$ , let ${\tilde{q}}_{i}^{tem} = \frac{\sum_{d \in D} {\tilde{q}}_{id}}{| D |} .$ Applying Theorem 1, then $Φ_{q_{i}^{tem}}^{- 1} (α) = \frac{\sum_{d \in D} Φ_{q_{id}}^{- 1} (α)}{| D |} .$ Similarly, we deal with ${\tilde{s}}_{i}^{tem}$ and ${\tilde{t}}_{ij}^{tem}$ , then we have $Φ_{s_{i}^{tem}}^{- 1} (α) = \frac{\sum_{d \in D} Φ_{s_{id}}^{- 1} (α)}{| D |},$ $Φ_{t_{ij}^{tem}}^{- 1} (α) = \frac{\sum_{d \in D} Φ_{t_{ijd}}^{- 1} (α)}{| D |},$ where α is the predefined confidence level. Next, we can put the processed variables in the model and deal with it.

In the template constraints, the maximum capacity and travel time for template routes also need to declare. We estimate the value of template maximum travel time T_tem and maximum capacity Q_tem based on the original values of Q and T. The way we define the values follows the approach used in Groër [14]. An expansion factor is defined by the formula $e = n_{f} / \bar{s}$ , where n_f is the number of frequent customers and $\bar{s}$ is the mean value of all the services on each day. Then we get $T_{tem} = T / \sqrt{e}$ and Q_tem = Q × e.

4.2 Construction of initial solution

In this section, we illustrate the procedure of constructing an initial template route for frequent customers under the constraints of vehicle capacity Q_tem and total operating time T_tem. Input data is designed following the method mentioned in 4.1. An insertion heuristic is designed to obtain feasible template routes with high-quality solutions. The feasibility of template route is checked through two dimensions: vehicle capacity Q_tem, and the tour length T_tem constraints in the template; vehicle capacity Q, tour length T and arrival time consistency L constraints in the daily routes obtained from the template route, since driver consistency has been guaranteed by the template concept. The insertion heuristic is used to insert customers one by one into a route at their best position. At first, a customer is selected as the starting point of one route according to the roulette wheel principle where the probability is the reciprocal of travel time between its position and the depot 0. The customer with a shorter travel time to the depot has a greater probability to be chosen as the first customer in a route. In the next step, for unvisited customers, the feasibility of each position is checked, which must satisfy the constraints of vehicle capacity Q_tem and total operating time T_tem. Then the customer who leads to the least insertion cost is selected and inserted in his best position. If there is no feasible position for unvisited customers, a new route is initialized. New routes will be added until all the customers are served. The initial template route will be checked feasibility again in the context of daily routes, which refers to vehicle capacity Q, tour length T, and arrival time consistency L constraints in the daily routes.

4.3 ABC implementation

ABC is one of the recently introduced swarm-based algorithms, which simulates the intelligent foraging behavior of a honeybee swarm [16]. It consists of three phases: employed bees’ phase, onlooker bees’ phase, and scout bees’ phase. Firstly, each employed bee finds a food source and searches for a new one in the neighborhood. Then each onlooker bee chooses to follow one employed bee according to the roulette wheel principle. And they also explore a new food source in the neighborhood. In the two stages, if the new food source is better than the old one, the old food source will be replaced by the new one. However, if a food source is not improved through the above stages, it will be abandoned and its corresponding employed bee will transform into the scout bee to find a new food source. When a new food source is found, the scout bee changes into the employed bee, then one iteration is completed. The whole process will repeat until it meets a stopping condition.

In the ABC algorithm, the position of a food source represents a possible solution to the optimization problem and the nectar amount of a food source corresponds to the quality of the associated solution. At the first step, NP initial solutions are constructed adopting the insertion heuristic mentioned in 1. Then each initial solution is selected randomly to generate a new solution by performing some neighborhood search methods. In each neighbor search, the solution with a better fitness value will be chosen. In onlooker bees’ phase, a solution is selected according to the wheel roulette principle, the probability of being selected is calculated by the formula $p_{i} = \frac{{fit}_{i}}{\sum_{i = 1}^{NP} {fit}_{i}},$ solutions with higher quality will be selected and updated again in this phase. To avoid getting a local optimal solution, if there is no improvement for limit times for a solution, the solution will be replaced by a regenerated initial solution. The process will iterate for iter times.

In our algorithm, we apply traditional local moves, namely reverse, swap, and insert. Each time we need to obtain a new solution, one neighborhood operator is selected randomly with equal probability. The operators are illustrated as follows:

(a) Random reverse: Randomly select two points and the subsequence of consecutive customers between the two points is selected. Then the order of the sequence is reversed. See Fig. 1a.

Fig. 1

Three neighbourhood operators.

(b) Random swap: This operator randomly selects two points and swaps them. See Fig. 1b.

The algorithm is provided in the following:

Algorithm 2 ABC implemention

Input: t_ijd, s_id, q_id, w_id, Q, T, L, the number of bee colony NP, the maximum limit for the repeated solution limit, the iteration number iter

Output: the best template route $R_{tem}^{best}$

generate NP initial solutions;

set t = 1 and l₁ = l₁ = . . . = l_NP = 0.

while t ≤ iter do

// employed bees¡¯ phase;

for m in NP do

For each initial solution R_tem, apply a neighbourhood operator to generate a new solution $R_{tem}^{'}$ ;

if $R_{tem}^{'}$ is better than R_tem then

set $R_{tem} = R_{tem}^{'}$ ;

else

l_m=l_m+1;

end if

end for

calculate the selection possibility of each solution;

// onlooker bees¡¯ phase;

for m in NP do

select a solution R_tem according to roulette wheel principle and apply a neighbourhood operator to generate a new solution $R_{tem}^{'}$ ;

if $R_{tem}^{'}$ is better than R_tem then

set $R_{tem} = R_{tem}^{'}$ ;

else

l_m=l_m+1;

end if

end for

// scout bees¡¯ phase;

for m in NP do

if l_m = limit then

randomly generate a new solution $R_{tem}^{'}$ and replace R_tem;

end if

end for

t = t + 1;

end while

return the best template route $R_{tem}^{best}$

4.4 Generation of daily routes

Based on the generated template route in the first stage, we can construct daily routes quickly. The procedure mainly consists of deletion and insertion operation and then routes improvement. Customers who do not ask for service on a specific day can be removed directly from the template route. For customers requiring service while not being contained in the template route on one day, we choose the insertion position by checking feasibility: vehicle capacity, route duration and arrival time consistency. Then for a feasible daily route, we decide the best insertion position using the greedy operator. The procedure of the generation of daily routes is provided below:

Algorithm 3 Generation of daily routes

Input: the template route $R_{tem}^{best}$

Output: the set of daily routes R = (R₁, R₂, R₃, . . . , R_D)

while d ≤ D do

set the initial daily route $R_{d} = R_{tem}^{best}$

end while

while d ≤ D do

delete the customers who do not require service on day d and get R_d

end while

while d ≤ D do

obtain the set of customers to be inserted into the daily route on day d, C_d

for i in C_d

insert the customer i into the R_d, check the feasibility and select the position that has the least travel time increment by greedy method do

end for

end while

return the best daily routes set R

5 Numerical experiments

In this section, we conduct several computational experiments to show the effectiveness of the proposed algorithm. The number of employed bees and the number of onlookers are set to be equal to the number of food sources, which is 25. This is based on Karaboga and Basturk [16], which set the number of employed bees to be equal to the number of onlookers to reduce the number of parameters, and find that the colony size of 50 can provide an acceptable convergence speed for search. We set the maximum cycle number to 1000 and the limit number to 75 in the experiments.

A benchmark data set has been provided in Groër et al. [14]. However, none of the problem instances involve uncertainty parameters. Thus, we design new test problem instances for our proposed models. In our experiments, these uncertain parameters are assumed to be normal uncertain variables. Thus the new problem instance is generated by adding standard deviances of the uncertain variables and setting the primal values as the expected values. The new problem instance contains 50 customers with a five-day planning horizon. Vehicle capacity and the travel time limit are 160 and 200 units. The maximum arrival time differential is not concluded in the solution process of the ConRTR method. It is calculated after obtaining daily routes and reported in the literature, thus in our algorithm, we test the value of the parameter and use the best value as the maximum time variation for our experiment. The uncertain parameters are assumed to be with normal uncertainty distributions. For simplicity, we assign the standard variance for each parameter according to its expected value. The distributions of the parameters are provided in Table 3.

Table 3
Distributions of uncertain variables

Parameters t _ijd s _id q _id

Uncertainty $N (t_{ijd}, t_{ijd} / 5)$ $N (s_{id}, s_{id} / 5)$ $N (q_{id}, q_{id} / 5)$

distribution

Parameters	t _ijd	s _id	q _id
Uncertainty	$N (t_{ijd}, t_{ijd} / 5)$	$N (s_{id}, s_{id} / 5)$	$N (q_{id}, q_{id} / 5)$
distribution

In the first group of experiments, we evaluate the impact of the three problem parameters, vehicle capacity Q, travel time T, and the maximum arrival time difference L. At first, we set L=70, and evaluate the impact vehicle capacity Q, and travel time T. We set vehicle capacity varying between 160 and 220 and maximum travel time of the vehicle between 160 and 220. To solve EVM, we assume the confidence levels are set in the following way: β_kd=γ_id=δ_i=0.9 for i ∈ N, k ∈ K, and d ∈ D. Results of total travel time with different combination of the two parameters are presented in Table 4. Obviously, these two parameters have significant impacts on the total travel time. We fix the vehicle capacity Q and increase the maximum travel time of the vehicle T, then the solution obtained is getting much worse. However, further increase of T can not improve the total travel time for some cases, due to the constraint of the fixed Q. Similarly, we can draw the same conclusion by fixing Q and varying T.

Table 4

Sensitivity analysis of vehicle capacity and working time

Travel time T	Vehicle capacity Q
	160	180	200	220
160	4628.10	4612.86	4612.86	4612.86
180	3903.34	3901.29	3901.29	3901.29
200	3570.66	3541.21	3541.21	3541.21
220	3398.56	3358.20	3337.36	3299.31

In the benchmark data set, the original maximum arrival time difference is reported as 63.47. In this experiment, we set β_kd=γ_id =δ_i=0.9, then vary the time difference limit between 40 and 80. In the CCP model, α is set as 0.9. Results of varying this parameter are presented in Table 5.

Table 5

Sensitivity analysis of maximum arrival time difference

Arrival time difference	40	50	60	70	80
Total travel time of EVM	–	3732.66	3580.88	3570.66	3586.98
Total travel time of CCP	–	4640.99	4480.28	4477.16	4434.64

We can see that decreasing the maximum arrival time difference will tighten the constraint thus slightly increase the obtained total travel time. While when we decrease this value further, there will be fewer solutions in the solution space that satisfy this constraint. When we decrease the time difference limit to 40, it is difficult to obtain feasible solutions in the experiment. On the other hand, increasing the maximum arrival time difference will relax the constraint but the influence on total travel time is not obvious.

In the second group of experiments, we focus on the influence of the confidence levels on the results of the models. As shown in Section 3, the results of solving EVM and CCP are dependent on the values of confidence levels β_kd, γ_id, and δ_i for i ∈ N, j ∈ N, k ∈ K, and d ∈ D, and the CCP model is also relevant to the parameter α, it is meaningful to conduct a sensitivity analysis to evaluate the impact of these parameters. Thus, we first focus on solving EVM to study the effect of the confidence levels β_kd, γ_id, and δ_i. In the experiments, we fix γ_id=0.9 and δ_i=0.9 then change the value of confidence level β_kd between 0.6 and 0.9 to show the impact of β_kd. To show the influence of different values of γ_id, we fix β_kd =0.6 and δ_i=0.8. The same operation is imposed on the parameter δ_i, where β_kd =0.6 and γ_id=0.7. In the following Fig. 2a, we show the results and it can be seen that the objective value gets worse with the confidence level increasing. Because increasing confidence level will tighten the associated constraint, the experiment conforms to theoretical inference. Then we solve the CCP model and study how the confidence levels α, β_kd, γ_id, and δ_i impact the result of the model. The parameters β_kd, γ_id, and δ_i are set the same as those in EVM when validating the effects of each one and α is set to 0.9. The results are shown in the following Fig. 2b and support theoretical observation similarly.

Fig. 2

Comparison of solutions for different confidence levels.

In Fig. 3, we change the value of confidence level α to show the performance of the CCP models. We study the impact of parameter α in CCP model together with the three problem parameters. Under different combinations of parameters, we validate their impacts. It can be found that the value of $\bar{f}$ increases with α. The results support the theoretical inference that a larger confidence level will lead to a more strict constraint, in turn a worse solution. Sensitivity analysis of vehicle capacity, travel time limit, and the maximum arrival time difference is also conducted in this part. As shown in Fig. 3, these parameters impact significantly the results of CCP model. We also see that solutions are getting better with the increase of vehicle capacity, travel time limit, and the maximum arrival time difference, which coincides with our intuition.

Fig. 3

Sensitivity analysis of parameters.

In the last group of experiments, we evaluate the effect of uncertainty on the total travel time in the two models. By varying the uncertain standard variance of parameters and fixing the expected values, the experiments aim to examine how the uncertainty of the uncertain parameters ${\tilde{t}}_{ijd}$ , ${\tilde{q}}_{id}$ , and ${\tilde{s}}_{id}$ influence the results. The solutions of EVM and CCP are shown in Table 6. The results in Table 6 reveal that the uncertainty of the three parameters has different degrees of influence on the problem. We can see that ${\tilde{q}}_{id}$ has a slight impact on the results, followed by ${\tilde{s}}_{id}$ . However, the total travel time greatly increases as the standard variance of ${\tilde{t}}_{ijd}$ increases. Thus for decision-makers, it’s wise to invest in the acquisition of travel times-related data to reduce uncertainty.

Table 6

Comparisons of different variances

	Distribution	Total travel time of EVM	Total travel time of CCP
${\tilde{t}}_{ijd}$	$N (t_{ijd}, t_{ijd} / 20)$	3295.62	3487.90
	$N (t_{ijd}, t_{ijd} / 10)$	3359.93	3771.78
	$N (t_{ijd}, t_{ijd} / 5)$	3570.66	4477.16
${\tilde{q}}_{id}$	$N (q_{id}, q_{id} / 20)$	3541.21	4465.19
	$N (q_{id}, q_{id} / 10)$	3541.21	4465.19
	$N (q_{id}, q_{id} / 5)$	3570.66	4477.16
${\tilde{s}}_{id}$	$N (s_{id}, s_{id} / 20)$	3420.02	4285.29
	$N (s_{id}, s_{id} / 10)$	3343.97	4322.42
	$N (s_{id}, s_{id} / 5)$	3570.66	4477.16

6 Conclusions and future research

In this paper, we present the consistent vehicle routing problem under uncertain environment, where the travel times between customers, the service times, and the demands of customers are assumed to be uncertain variables. Specifically, uncertain travel times, service times, and customer demands are defined and tackled by uncertainty theory. Following different decision criteria, the expected value model and the chance-constrained programming model are established. Taking advantage of uncertainty theory, the proposed models could be respectively transformed into the corresponding deterministic equivalents.

The proposed models are solved by the artificial bee colony algorithm based on the template concept, in which the template route is first generated to guarantee the service consistency and the daily routes are created from the template route. The ABC algorithm is employed to optimize the template route and the daily routes. Finally, numerical experiments are given to demonstrate the effectiveness of the proposed solution approach.

The study in this paper points out several potential directions of future research. First, we can extend driver consistency to coordinate with more realistic situations. The presented problem in this paper is defined for the same driver. In the future, we can extend the problem by allowing more than one driver for each customer and design a new measurement method for driver consistency. Furthermore, for customers with different consumption habits, we can assign different numbers of the different drivers for them. On the other hand, in our research we study the basic ConVRP under uncertain environment. As many variants of ConVRP develop, in the future we can extend our work to more complex consistent vehicle routing problems, for instance, the values of different customers can be taken into consideration.

Footnotes

Acknowledgments

This work was supported by National Natural Science Foundation of China (No. 71471038), Program for Huiyuan Distinguished Young Scholars, UIBE (No. 17JQ09), “the Fundamental Research Funds for the Central Universities” in UIBE (No. CXTD10-05).

Appendix: Uncertainty theory

In this section, some basic concepts and theorems in uncertainty theory are introduced.

Let Γ be a nonempty set, and ℒ be a σ-algebra over Γ. Each element Λ ∈ ℒ is called an event. A number M {Λ} indicates the possibility that Λ will occur. Uncertain measure M is introduced as a set function satisfying the following axioms [24]:

Axiom 1. (Normality Axiom) M {Γ} =1 for the universal set Γ.

Axiom 2. (Duality Axiom) M {Λ} + M {Λ^c} =1 for any event Λ.

Axiom 3. (Subadditivity Axiom) For every countable sequence of events {Λ_i}, we have $M {⋃_{i = 1}^{\infty} Λ_{i}} \leq \sum_{i = 1}^{\infty} M {Λ_{i}} .$

The triplet (Γ, ℒ, ℳ) is called an uncertainty space. In addition, the product uncertain measure [25] was defined as following.

Axiom 4. (Product Axiom) Let (Γ_k, ℒ_k, _Mk) be uncertainty spaces for k = 1, 2, ⋯ , n. The product uncertain measure M is an uncertain measure satisfying $M {\prod_{k = 1}^{\infty} Λ_{k}} = ⋀_{k = 1}^{\infty} M_{k} {Λ_{k}} .$

Definition 1. (Liu [24]) An uncertain variables is a function ξ from an uncertainty space (Γ, ℒ, ℳ) to the set of real numbers such that ξ ∈ B is an event for any Borel set B.

Definition 2. (Liu [24]) The uncertainty distribution Φ of an uncertain variable ξ is defined by $Φ (x) = M {ξ \leq x}, \forall x \in R$ for any real number x.

Definition 3. (Liu [26]) An uncertainty distribution Φ (x) is said to be regular if it is a continuous and strictly increasing function with respect to x at which 0 < Φ (x) <1, and $lim_{x \to - \infty} Φ (x) = 0, lim_{x \to + \infty} Φ (x) = 1 .$

Definition 4. (Liu [24]) Let ξ be an uncertain variable. Then the expected value of ξ is defined by $E [ξ] = \int_{0}^{+ \infty} M {ξ \geq r} dr - \int_{- \infty}^{0} M {ξ \leq r} dr$ provided that at least one of the two integrals is finite.

Theorem 1. (Liu [26]) |IT Let ξ₁, ξ₂, ⋯ , ξ_n are independent uncertain variables with regular uncertainty distributions Φ₁, Φ₂, ⋯ , Φ_n, respectively. If f (ξ₁, ξ₂, ⋯ , ξ_n) is strictly increasing with respect to ξ₁, ξ₂, ⋯ , ξ_m and strictly decreasing with respect to ξ_m+1, ξ_m+2, ⋯ , ξ_n, then $ξ = f (ξ_{1}, ξ_{2}, \dots, ξ_{n})$ has an inverse uncertainty distribution

$\begin{matrix} Ψ^{- 1} (α) = f (Φ_{1}^{- 1} (α), \dots, Φ_{m}^{- 1} (α), \\ Φ_{m + 1}^{- 1} (1 - α), \dots, Φ_{n}^{- 1} (1 - α)) . \end{matrix}$

Theorem 2. (Liu [26]) Let ξ be an uncertain variable with regular uncertainty distribution Φ. Then $E [ξ] = \int_{0}^{1} Φ^{- 1} (α) d α .$

Theorem 3. (Liu [26]) Assume ξ₁, ξ₂, ⋯ , ξ_n are independent uncertain variables with regular uncertainty distributions Φ₁, Φ₂, ⋯ , Φ_n, respectively. If f (ξ₁, ξ₂, ⋯ , ξ_n) is strictly increasing with respect to ξ₁, ξ₂, ⋯ , ξ_m and strictly decreasing with respect to ξ_m+1, ξ_m+2, ⋯ , ξ_n, then $ξ = f (ξ_{1}, ξ_{2}, \dots, ξ_{n})$ has an excepted value

$\begin{matrix} E [ξ] = \int_{0}^{1} f (Φ_{1}^{- 1} (α), \dots, Φ_{m}^{- 1} (α), \\ Φ_{m + 1}^{- 1} (1 - α), \dots, Φ_{n}^{- 1} (1 - α)) d α . \end{matrix}$

Theorem 4. (Liu [26]) Assume the constraint function g (x, ξ₁, ξ₂, ⋯ , ξ_n) is strictly increasing with respect to ξ₁, ξ₂, ⋯ , ξ_m and strictly decreasing with respect to ξ_m+1, ξ_m+2, ⋯ , ξ_n. If ξ₁, ξ₂, ⋯ , ξ_n are independent uncertain variables with uncertainty distributions Φ₁, Φ₂, ⋯ , Φ_n respectively, then the chance constraint $M {g (x, ξ_{1}, ξ_{2}, \dots, ξ_{n}) \leq 0} \geq α$ holds if and only if

$\begin{matrix} g (x, Φ_{1}^{- 1} (α), \dots, Φ_{m}^{- 1} (α), Φ_{m + 1}^{- 1} (1 - α), \\ \dots, Φ_{n}^{- 1} (1 - α)) \leq 0 . \end{matrix}$

References

Braekers

and Kovacs

A.A.

, A multi-period dial-a-ride problem with driver consistency, Transportation Research Part B: Methodological 94 (2016), 355–377.

Campelo

, Neves-Moreira

, Amorim

and Almada-Lobo

, Consistent vehicle routing problem with service level agreements: A case study in the pharmaceutical distribution sector, European Journal of Operational Research 273(1) (2019), 131–145.

Chen

, Peng

and Zhang

, Uncertain goal programming models for bicriteria solid transportation problem, Applied Soft Computing 51 (2017a), 49–59.

Chen

, Peng

, Liu

and Zhao

, Pricing and effort decisions for a supply chain with uncertain information, International Journal of Production Research 55(1) (2017b), 264–284.

Coelho

L.C.

, Consistency in multi-vehicle inventory-routing, Transportation Research Part C: Emerging Technologies 24 (2012), 270–287.

Coelho

L.C.

and Laporte

, A branch-and-cut algorithm for the multi-product multi-vehicle inventory-routing problem, International Journal of Production Research 51(23-24) (2013a), 7156–7169.

Coelho

L.C.

and Laporte

, The exact solution of several classes of inventory-routing problems, Computers and Operations Research 40(2) (2013b), 558–565.

Day

J.M.

, Wright

P.D.

, Schoenherr

, Venkataramanan

and Gaudette

, Improving routing and scheduling decisions at a distributor of industrial gasses, Omega 37(1) (2009), 227–237.

Emadikhiav

, Bergman

and Day

, Consistent routing and scheduling with simultaneous pickups and deliveries, Production and Operations Management 29(8) (2020), 1937–1955.

10.

Esmaili

M.H.

and Mousavi

S.M.

, An Integrated Perishable Inventory Routing Problem with Consistent Driver Services and Deliver Fresh Products using Possibility and Necessity Measures, International Journal of Industrial Engineering and Production Research 31(2) (2020), 231–242.

11.

Feillet

, Garaix

, Lehu’ed’e

, P’eton

and Quadri

, A new consistent vehicle routing problem for the transportation of people with disabilities, Networks 63(3) (2014), 211–224.

12.

Gao

, Yang

and Li

, Uncertain models on railway transportation planning problem, Applied Mathematical Modelling 40(7-8) (2016), 4921–4934.

13.

Goeke

, Roberti

and Schneider

, Exact and heuristic solution of the consistent vehicle-routing problem, Transportation Science 53(4) (2019), 1023–1042.

14.

Groër

, Golden

and Wasil

, The consistent vehicle routing problem, Manufacturing and Service Operations Management 11(4) (2009), 630–643.

15.

Huang

and Ke

, Pricing decision problem for substitutable products based on uncertainty theory, Journal of Intelligent Manufacturing 28(3) (2017), 503–514.

16.

Karaboga

and Basturk

, On the performance of artificial bee colony (ABC) algorithm, Applied Soft Computing Journal 8(1) (2008), 687–697.

17.

, Li

and Huang

, Uncertain pricing decision problem in closed-loop supply chain with risk-averse retailer, Journal of Uncertainty Analysis and Applications 3(1), 14, 2015.

18.

Kovacs

A.A.

, Golden

B.L.

, Hartl

R.F.

and Parragh

S.N.

, Vehicle routing problems in which consistency considerations are important: A survey, Networks 64(3) (2014a), 192–213.

19.

Kovacs

A.A.

, Parragh

S.N.

and Hartl

R.F.

, A template-based adaptive large neighborhood search for the consistent vehicle routing problem, Networks 63(1) (2014b), 60–81.

20.

Kovacs

A.A.

, Golden

B.L.

, Hartl

R.F.

and Parragh

S.N.

, The generalized consistent vehicle routing problem, Transportation Science 49(4) (2014c), 796–816.

21.

Kovacs

A.A.

, Parragh

S.N.

and Hartl

R.F.

, The multi-objective generalized consistent vehicle routing problem, European Journal of Operational Research 247(2) (2015), 441–458.

22.

Kulachenko

I.N.

, Kononova

P.A.

and Kochetov

Y.A.

, The variable neighborhood search for a consistent vehicle routing problem under the shift length constraints, IFAC-PapersOnLine 52(13) (2019), 2314–2319.

23.

Lian

, Milburn

A.B.

and Rardin

R.L.

, An improved multi-directional local search algorithm for the multi-objective consistent vehicle routing problem, IIE Transactions 48(10) (2016), 975–992.

24.

Liu

, Uncertainty theory, Springer-Verlag, Berlin, 2007.

25.

Liu

, Some research problems in uncertainty theory, Journal of Uncertain Systems 3(1) (2009), 3–10.

26.

Liu

, Uncertainty theory: a branch of mathematics for modeling human uncertainty, Springer-Verlag, Berlin, 2010.

27.

Luo

, Qin

and Che

C.H.

, On service consistency in multi-period vehicle routing, European Journal of Operational Research 243(3) (2015), 731–744.

28.

Mancini

, Gansterer

and Hartl

R.F.

, The collaborative consistent vehicle routing problem with workload balance, European Journal of Operational Research, 2021.

29.

, Shi

and Wei

, Optimizing influence diffusion in a social network with fuzzy costs for targeting nodes, Journal of Ambient Intelligence and Humanized Computing 8(5) (2017a), 819–826.

30.

, Ning

and Ke

, Modeling and minimizing information distortion in information diffusion through a social network, Soft Computing 21(18) (2017b), 5281–5293.

31.

, Sequential seeding to optimize influence diffusion in a social network, Applied Soft Computing 56 (2017c), 730–737.

32.

Ning

and Su

, A multilevel approach for modelling vehicle routing problem with uncertain travelling time, Journal of Intelligent Manufacturing 28(3) (2017), 683–688.

33.

Ortega

E.J.A.

, Schilde

and Doerner

K.F.

, Matheuristic search techniques for the consistent inventory routing problem with time windows and split deliveries, Operations Research Perspectives, 100152, 2020.

34.

Paquette

, Bellavance

, Cordeau

J.F.

and Laporte

, Measuring quality of service in dial-a-ride operations: The case of a Canadian city, Transportation 39(3) (2012), 539–564.

35.

Rodr’iguez-Mart’in

, Salazar-Gonz’alez

J.J.

and Yaman

, The periodic vehicle routing problem with driver consistency, European Journal of Operational Research 273(2) (2019), 575–584.

36.

Song

and Ni

, Optimal platform design with modularity strategy under fuzzy environment, Soft Computing 23(3) (2019), 1059–1070.

37.

Song

, Ni

and Ralescu

D.A.

, Product configuration using redundancy and standardisation in an uncertain environment, International Journal of Production Research (2020), 1–20.

38.

Song

, Ni

and Ralescu

D.A.

, The impact of lead-time uncertainty in product configuration, International Journal of Production Research 59 (2020), 1–23.

39.

Song

, Ulmer

M.W.

, Thomas

B.W.

and Wallace

S.W.

, Building trust in home services-stochastic team-orienteering with consistency constraints, Transportation Science 54(3) (2020), 823–838.

40.

Stavropoulou

, Repoussis

P.P.

and Tarantilis

C.D.

, The vehicle routing problem with profits and consistency constraints, European Journal of Operational Research 274(1) (2019), 340–356.

41.

Subramanyam

and Gounaris

C.E.

, A branch-and-cut framework for the consistent traveling salesman problem, European Journal of Operational Research 248(2) (2016), 384–395.

42.

Subramanyam

and Gounaris

C.E.

, A decomposition algorithm for the consistent traveling salesman problem with vehicle idling, Transportation Science 52(2) (2017), 386–401.

43.

Sungur

, Ren

, Ord’oñez

, Dessouky

and Zhong

, A model and algorithm for the courier delivery problem with uncertainty, Transportation Science 44(2) (2010), 193–205.

44.

Tarantilis

C.D.

, Stavropoulou

and Repoussis

P.P.

, A template-based tabu search algorithm for the consistent vehicle routing problem, Expert Systems with Applications 39(4) (2012), 4233–4239.

45.

UPS, Driving success: Why the UPS model for managing 103,500 drivers is a competitive advantage. http://pressroom.ups.com/mediakits/popups/factsheet/0,1889,1201,00.html, 2007.

46.

Wong

R.T.

, Vehicle routing for small package delivery and pickup services, The Vehicle Routing Problem: Latest Advances and New Challenges 43, 475–485, 2008.

47.

Woodward

C.A.

, Abelson

and Tedford

, What is important to continuity in home care? Perspectives of key stakeholders, Social Science And Medicine 58(1) (2004), 177–192.

48.

and Cai

, Variable neighborhood search for consistent vehicle routing problem, Expert Systems with Applications 113 (2018), 66–76.

49.

Zhang

and Zhang

, Vehicle Routing Problem Based on Uncertainty Theory, in: Proceedings of Proceedings of the First International Conference on Uncertainty Theory, Urumchi, China, 2010.

50.

Zhen

, Lv

, Wang

, Ma

and Xu

, Consistent vehicle routing problem with simultaneous distribution and collection, Journal of the Operational Research Society 71(5) (2020), 813–830.

The consistent vehicle routing problem under uncertain environment

Abstract

Keywords

1 Introduction

2 Related literature

3 Problem description

4.1 Initialization phase

Table 2 Customers division Customer 1 2 3 4 5 6 7 8 9 10 Day 1 1 0 1 1 0 1 1 1 0 1 Day 2 0 1 1 1 1 0 0 1 1 1 Day 3 1 1 0 1 0 1 0 1 0 1 Number of service days 2 2 2 3 1 2 1 3 1 3 Template 1 1 1 1 0 1 0 1 0 1

4.3 ABC implementation

5 Numerical experiments

Table 3 Distributions of uncertain variables Parameters t ijd s id q id Uncertainty N ( t ijd , t ijd / 5 ) N ( s id , s id / 5 ) N ( q id , q id / 5 ) distribution

Footnotes

Acknowledgments

Appendix: Uncertainty theory

References

Table 2
Customers division

Customer 1 2 3 4 5 6 7 8 9 10

Day 1 1 0 1 1 0 1 1 1 0 1

Day 2 0 1 1 1 1 0 0 1 1 1

Day 3 1 1 0 1 0 1 0 1 0 1

Number of service days 2 2 2 3 1 2 1 3 1 3

Template 1 1 1 1 0 1 0 1 0 1

Table 3
Distributions of uncertain variables

Parameters t _ijd s _id q _id

Uncertainty $N (t_{ijd}, t_{ijd} / 5)$ $N (s_{id}, s_{id} / 5)$ $N (q_{id}, q_{id} / 5)$

distribution