Computers & Industrial Engineering 178 (2023) 109143

Contents lists available at ScienceDirect

Computers & Industrial Engineering

journal homepage: www.elsevier.com/locate/caie

A stabilized branch-and-price-and-cut algorithm for the waste

transportation problem with split transportation
Li Zhang a, Zhongshan Liu a, Wenxuan Shan a, Bin Yu a, b, *
a
Sch Transportat Sci & Engn, Beihang University, Beijing 100191, PR China
b
Beijing Key Lab Cooperat Vehicle Infrastruct Syst, Beijing 100191, PR China

A R T I C L E I N F O A B S T R A C T

Keywords: In the context of waste classification, different categories of waste are collected and stored separately in each
Waste transportation collection location, which makes the split transportation implement readily in the waste transportation problem.
Split transportation The split transportation refers to that trailers are allowed to visit each location several times to transport different
Branch-and-price-and-cut algorithm
categories of waste from collection locations to transfer stations. Since split transportation relaxes the basic
Stabilized column generation
assumption that each collection locationshould be visited exactly once, it would lead to flexible combinations of
routes and help waste transportation departments or companies save on operational costs. In this paper, we
present a novel waste transportation problem that extends the traditional waste transportation problem by
introducing the concept of split transportation. We develop a mixed integer programming formulation and a
branch-and-price-and-cut algorithm on an extended network to solve the investigated problem. The extended
network simplifies the pricing subproblem, but renders high degeneration in column generation. Then, the
stabilized technique is devised to speed up the convergence speed of the proposed algorithm. The efficiency of
the proposed algorithm is tested on a real-world instance and a set of benchmark instances. Then, we discuss the
impact of split transporation, demand intervals, and trailer capacity to give managerial insights for waste
transportation departments or companies.

1. Introduction municipal waste. From 2000 to 2021, the waste generation in China
increases from 118 million tons to 248 million tons (Statista, 2023b). To
The population growth, urbanization, and economic development dispose waste properly and protect environment, Chinese government
have significantly increased the worldwide waste generation recently. launched mandatory policy of waste classification in some pilot cities in
By 2050, the worldwide waste generation is expected to increase by 2017, such as Beijing and Shanghai, and then adopted it nationwide. The
about 70% compared to that in 2016 and reach 3.4 billion tons (Statista, waste classification requires residents to sort waste based on predefined
2023a). The increasing waste generation may not only lead to heavy standards and discard the sorted waste into specific containers. The
environment issues but also bring considerable economic investment in waste classification can bring great convenience in waste collection and
waste management. In 2020, it may cost 383 billion dollars (Mordor transportation.
Intelligence, 2022) in waste management worldwide. The waste collection and transportation are core procedures in urban
The waste management is typically segmented by waste categories, waste management system, which is to collect waste from waste
such as municipal waste, hazardous waste, recyclable waste. Due to the collection locations in communities or at roadsides, and then transport
significant amount of generated waste, it is essential to apply reasonable the collected waste to dump sites. The introduction of waste classifica-
and distinct treatment techniques for various waste categories to tion can bring independent collection and transportation procedure of
improve the efficiency of waste management. As a result, the collection different categories of waste, which brings a novel waste transportation
and transportation of different categories of waste should be indepen- strategy, the split transportation. In the traditional waste collection and
dent, which introduces the context of waste classification. transportation applications, when trailers visit waste collection loca-
China is regarded as one of the countries associated with the largest tions, the tailers are required to carry all waste. This follows from the
waste generation in the world, and accounts for more than 15% fact that all waste is mixed and stored in the same waste containers

* Corresponding author at: Sch Transportat Sci & Engn, Beihang University, Beijing 100191, PR China.
E-mail address: ybzhyb@163.com (B. Yu).

https://doi.org/10.1016/j.cie.2023.109143
Received 18 December 2021; Received in revised form 12 February 2023; Accepted 2 March 2023
Available online 8 March 2023
0360-8352/© 2023 Elsevier Ltd. All rights reserved.
L. Zhang et al. Computers & Industrial Engineering 178 (2023) 109143

before applying waste classification. Collecting a part of mixed waste is waste transportation problem.
impractical. The split transportation strategy relaxes this restriction and To the best of our knowledge, it is the first time that the aspect of split
only let trailers carry specific categories of waste in one visit. The un- transportation is introduced in the waste transportation problem. The
collected waste can be collected by the other trailers. The introduction of waste transportation problem with split transportation can be regarded
split transportation strategy can make trailer routes flexible and save as a variant of split delivery VRP. The split delivery VRP was proposed
operational cost. by Dror and Trudeau (1990), and a local search heuristic algorithm was
The procedures of waste collection and transportation are mainly designed. The split delivery VRP allows multiple times of service at each
studied as the waste transportation problem in the current literature. customer site. For the SDVRP, Dror et al. (1994) proposed a branch-and-
Researchers and practitioners have shown great interest in this area. Due cut algorithm. Then, Desaulniers (2010) designed a branch-and-price-
to some special characteristics of waste transportation issues, a variety and-cut algorithm to address the split delivery VRP with time win-
of vehicle routing problems (VRPs) have been proposed to model the dows, in which the linear relaxation of a bound knapsack problem was
waste transportation problem. Waste trailers can stop at transfer stations included in the column generation to find a convex combination of
(the transfer station hereinafter is called depot) to empty loads during extreme delivery patterns. The algorithm proposed by Desaulniers
their routes. In the current literature, the routing optimization of waste (2010) was applied by Luo et al. (2017), where the branch-and-price-
transportation are mainly studied as a VRP with inter depot. The VRP and-cut algorithm to solve the split delivery VRP with time windows
with inter depot is similar to the multi-trip VRP. The difference is that and load-related travel cost. Salani and Vacca (2011) proposed a
multi-trip VRP requires that trailers can only visit their original depots, discrete split delivery VRP with time windows, in which the demand
while the VRP with inter depot does not. Kim et al. (2006) studied the could be split into feasible combinations. They designed a branch-and-
waste transportation problem with inter depot, where the lunch breaks price algorithm to solve the discrete split delivery VRP. The
of drivers were considered. They used an insertion algorithm inspired by commodity-constrained split delivery VRP was introduced by Archetti
Solomon (1987) to solve the proposed model. A similar problem was et al. (2015 and 2016), where each customer site can be visited more
studied by Buhrkal et al. (2012), and they applied an adaptive large than once, but each category of commodity should be delivered within
neighborhood search algorithm. Benjamin and Beasley (2010) intro- one visit. To solve the commodity-constrained split delivery VRP,
duced the multiple inter depots, unlimited vehicle fleets, and rest pe- Archetti et al. (2015) developed a branch-and-price-and-cut algorithm.
riods of drivers, and a variable neighborhood search-based algorithm Based on the research of Archetti et al. (2015), Gschwind et al. (2019)
was designed to address the problem. Markov et al. (2016) investigated proposed a modified branch-and-price-and-cut algorithm, which
a waste transportation problem with inter depot and heterogeneous involved a stabilization strategy of dual inequalities, to address the
vehicles. In their research, the trailers were not obligated to return to the commodity-constrained split delivery VRP. The split delivery VRP has
source depots at the end of planning horizon, which inducing a flexible been studied widely (Yu et al.,2019; Zhang et al., 2019). A detailed re-
assignment of destination depots. To solve their realistic instances, a view of split delivery VRP and its heuristic algorithms can be found in
multiple neighborhood search algorithm is designed. The difference Archetti and Speranza (2012). The current literature of the commodity-
between source and sink depot was also studied by Ramos et al. (2013). constrained split delivery VRP, such as Gschwind et al. (2019), did not
Some other researchers have focused on several other features of the take the inter depot into account. When the inter depot is involved in the
waste transportation problem, such as periodic routing, multiple com- commodity-constrained split delivery VRP, their proposed branch-and-
modities, and rollon-rolloff routing. The selected papers for the waste price-and-cut algorithm would not be applied directly.
transportation problem are summarized in Table 1. Beliën et al. (2014) In this paper, we study the waste transportation problem with split
and Al-Salem et al. (2009) provide detailed literature reviews of the transportation, where multi-commodity and inter depot are also

Table 1
Selected references about the waste transportation problem.
Reference Time Inter Periodic Multiple Rollon- Rest Exact Heuristic Split
Windows Depot Routing Commodities Rolloff period Method Method Delivery

Archetti and Speranza ✓ ✓ ✓ ✓ ✓

(2004)
Baldacci et al. (2006) ✓ ✓ ✓
Benjamin and Beasley ✓ ✓ ✓ ✓
(2012)
Buhrkal et al. (2012) ✓ ✓ ✓ ✓
Cárdenas-Barrón et al. ✓ ✓
(2019)
De Bruecker et al. (2018) ✓ ✓ ✓
Expósito-Márquez et al. ✓ ✓ ✓
(2019)
Faccio et al. (2011) ✓ ✓
Hauge et al. (2014) ✓ ✓ ✓
Kim et al. (2006) ✓ ✓ ✓ ✓
Li et al. (2018) ✓ ✓ ✓ ✓
Markov et al. (2016) ✓ ✓ ✓ ✓
Markov et al. (2020) ✓ ✓ ✓
Miranda et al. (2015) ✓ ✓ ✓
Muter et al. (2014) ✓ ✓
Nuortio et al. (2006) ✓ ✓ ✓ ✓ ✓
Ramos et al. (2013) ✓
Ramos et al. (2014) ✓ ✓
Sahoo et al. (2005) ✓ ✓ ✓
Teixeira et al. (2004) ✓ ✓ ✓
Tirkolaee et al. (2018) ✓ ✓
Wy et al. (2013) ✓ ✓ ✓ ✓ ✓
This paper ✓ ✓ ✓ ✓ ✓

2
L. Zhang et al. Computers & Industrial Engineering 178 (2023) 109143

involved. That is, trailers collect specific categories of waste in each our study, which results in the flexible combinations of routes. Customer
waste collection location within time windows. During the planning sites and depots both have predefined time windows [ai , bi ], and trailers
horizon, trailers can visit depots to empty the carried waste. Trailers are must serve these customer sites and depots within their time windows.
allowed to visit each collection location several times to collect all The travel time tij and travel cost cij along arc (i, j) ∈ A are given in
waste, which is named as split transportation. To bridge the gap between advance and keep static. The objective is to minimize the total trans-
the current literature and practical applications, we give a mixed integer portation cost. Since the investigated problem relaxes the assumption
programming formulation to model the investigated problem. A stabi- that each customer site is only visited once, the optimal solution of the
lized branch-and-price-and-cut algorithm is developed to solve the proposed problem gives a lower bound of the traditional waste trans-
model, where the column generation procedure based on an extended portation problem without split transportation.
network is applied at each branch node to obtain the linear relaxation An explanatory example of the proposed problem is given in Fig. 1. In
solution. The main contributions of this paper can be summarized as this example, we assume that there are two depots in the planning area.
follows: One trailer, which can carry four units of waste, is available in depot A.
There are two categories of waste that should be collected from three
(1) We introduce the waste transportation problem with split trans- customer sites. These categories of waste are named waste a and waste b
portation, multi-commodity, and inter depot to seek reasonable in this example and represented by squares and circles, respectively. The
combinations of routes and save on operational cost for the waste amount of waste is given in parenthesis. The trailer, which leaves depot
transportation departments or companies. The investigated A with empty load, should serve these three customer sites, and finally
problem is modeled as a variant of the VRPs, and a four-index return to depot A or depot B. When the load of trailer reaches its ca-
compact formulation is built. pacity, the trailer is allowed to dump waste in depots. Time windows are
(2) In this paper, we design a stabilized branch-and-price-and-cut ignored in this example. A feasible route is given in Fig. 1, which in-
algorithm to solve the investigated problem. An extended cludes two trips (trip one is represented by bold lines, and trip two is
network is proposed to simplify the label setting algorithm that is represented by dashed lines). The trailer leaves depot A and serves
used to address the pricing subproblem, and a stabilization customers 1 and 2 before it reaches depot B. Note that at customer site 2,
technique is applied to decrease the degeneration of the column one unit of waste a and two units of waste b should be collected.
generation. However, since only two units of capacity are available, the trailer only
(3) We test our algorithm on a set of modified benchmark instances collects all waste b in trip one. Then, the trailer dumps the waste
and a real-world instance. The computational results reveal that collected from customer sites 1 and 2 in depot B. Note that in this
the proposed algorithm can find acceptable results even in large example, the source depot and sink depot of the trailer are the same, but
instances. Some managerial analyses are given to study the it is not mandatory. In this example, the split transportation occurs at
impact of split transportation, demand intervals, and trailer ca- customer site 2. If the split transportation is prohibited in this instance,
pacity in the waste transportation problem. at least two trailers are required to serve all of customer sites.

This paper is organized as follows. Section 2 gives the description 2.2. Mathematical formulation
and model formulation of the proposed waste transportation problem. In
Section 3, Dantzig-Wolfe (1960) decomposition is applied in an In this subsection, a compact formulation is given as a four-indexed
extended network to obtain the master problem and pricing problem. integer programming formulation. Necessary notations used in the
The label setting algorithm is developed in Section 4 to solve the pricing proposed formulation are defined in Table 2.
subproblem. Cutting planes, branching strategies, and stabilized column With above notations, the proposed model can be formulated as
generation are also given in Section 4. Section 5 reports the computa- follows:
tional results of the modified benchmark instances and a real-world ∑∑ ∑
instance. Finally, Section 6 concludes this paper. min cij xklij (1)
k∈K l∈L (i,j)∈A

2. Problem description and formulation Subject to

2.1. Problem description

The proposed waste transportation problem can be described as

follows: A waste transportation company entails to collect in waste
collection locations of living units (these locations hereinafter are called
customer sites). The overall problem is defined on a complete network
G(V, A), where V = {0, 1, ..., nc , ..., nc + nd } indicates the vertex set and
A = {(i, j)|i, j ∈ V } denotes the arc set. In vertex set V, vertex 0 stipulates
a dummy source depot, which is used to design the proposed algorithm.
Let Vc = {1, 2, ..., nc } and Vd = {nc + 1, ..., nc + nd } depict the vertex sets
of customer sites and depots, respectively. A fleet of homogenous trailers
K has a capacity Q. Let M denotes the set of waste categories. Parameter
Dmi indicates that whether category m ∈ M of waste is required to be
collected at customer site i ∈ Vc , and the amount of waste m ∈ M at
customer site i ∈ Vc is dm i . Trailers can visit each customer site several
times to collect all waste. Suppose that the amount of each category of
waste at each customer site does not exceed the capacity of tailers This
assumption is reasonable and can help us simplify our model and al-
gorithm. During the planning horizon, trailers can dump carried waste
to one of the depots. Analogous to Markov et al. (2016), trailers are not
obligated to return to their source depots after serving all customers in
Fig. 1. A small example for the proposed waste transportation problem.

3
L. Zhang et al. Computers & Industrial Engineering 178 (2023) 109143

Table 2
Notations of sets, parameters, and variables. qkli = 0 ∀i ∈ Vd , k ∈ K, l ∈ L (14)

Sets
xklij ∈ {0, 1} ∀(i, j) ∈ A, k ∈ K, l ∈ L (15)
V Vertex set of customer sites and depots, V = Vc ∪ Vd , indexed by
letters i and j.
Vc Vertex set of customer sites, Vc = {1, 2, ..., nc }.
yklim ∈ {0, 1} ∀i ∈ Vc , m ∈ M, k ∈ K, l ∈ L (16)
Vd Vertex set of depots, Vd = {nc + 1, ..., nc + nd }.
Objective function (1) aims to minimize the total cost. Constraints
A Arc set, which includes all arcs that connect vertices in vertex set V,
A = {(i, j)|i, j ∈ V }. (2) indicate that each customer should be served at least once. Since the
Γ− (i),Γ+ (i) Preceding and succeeding vertex sets of vertex i ∈ V. split delivery is allowed, each customer can be visited by different trips
K Set of available trailers, indexed by k. more than once. Constraints (3) are flow balance constraints that ensure
L Set of trips for each trailer, indexed by l. the solution is composed of feasible trips. Constraints (4) and (5) denote
M Set of waste categories, indexed by m.
that trailers should start and end their trips in one of the depots. Con-
Parameters
cij Route cost between two vertices i ∈ V and j ∈ V. straints (6) ensure that trailers can collect waste from customer sites
tij Travel time between two vertices i ∈ V and j ∈ V. only when visiting such customer sites. Constraints (7) indicate that
dmi Amount of waste m ∈ M at customer site i ∈ Vc . each category of waste should be collected within one visit. The time
Dm i Binary parameter which indicates whether waste m ∈ M is required to consistency and load consistency are assured in constraints (8) and (10).
be collected at customer site i ∈ Vc . When waste m ∈ M is required to
The service start time and vehicle load are bounded in constraints (9)
be collected from customer site i ∈ Vc , Dm m
i = 1, otherwise Di = 0.
Q Capacity of each trailer. and (11). Constraints (12) guarantee that each trip should be started in
T Overall planning horizon. the depot that is the ending depot of its preceding trip. Constraints (13)
σ Unloading time at depot site. indicate that the unloading time should be reserved after each trip.
si Service time at customer site i ∈ Vc When trailers visit any depot, all load should be empty (constraints
ai , bi Earliest and latest service time at customer or depot site i ∈ V.
Variables
(14)). Domains of decision variables are given in constraints (15) and
xkl Binary variable which indicates whether arc (i, j) ∈ A(V) is traversed (16).
ij
by trip l ∈ L of trailer k ∈ K. When arc (i, j) ∈ A(V) is traversed by trip
l ∈ L of trailer k ∈ K, xkl kl
ij = 1, otherwise xij = 0. 3. Column generation
ykl
im
Binary variable which denotes whether waste m ∈ M of customer site
i ∈ Vc is transported by trip l ∈ L of trailer k ∈ K. When waste m ∈ M of
The formulation given in subsection 2.2 can be solved by commercial
customer site i ∈ Vc is transported by trip l ∈ L of trailer k ∈ K, ykl
im = 1,
otherwise ykl
optimizers directly, such as CPLEX. However, solving this formulation
im = 0.
τkl Nonnegative variable which denotes the service starting time at by commercial optimizers is not practical in the medium or large scale
i
customer site i ∈ Vc of trip l ∈ L of trailer k ∈ K instances. The branch-and-price-and-cut algorithm has shown great
qkl
i
Nonnegative variable which denotes current load of trailer k ∈ K when ability in achieving great results in the literature of SDVRP (Desaulniers,
it visits customer site i ∈ Vc in the trip l ∈ L.
2010; Luo et al., 2017; Salani and Vacca, 2011; Archetti et al., 2015;
Gschwind et al., 2019). In this paper, we devise a branch-and-price-and-
∑∑ ∑ cut algorithm where the column generation procedure is applied at each
xklij ⩾1 ∀i ∈ Vc (2)
k∈K l∈L j∈Γ+ (i) branch node to access the linear relaxation of the restricted master
problem. To simplify the column generation, we build an extended
∑ ∑
xhikl = xklij ∀i ∈ Vc , k ∈ K, l ∈ L (3) customer-demand network so that the pricing subproblem can be
h∈Γ− (i) j∈Γ+ (i) viewed as an elementary shortest path problem with resource con-
straints, which can be solved by the label setting algorithm.
∑∑
xklij = 1 ∀k ∈ K, l ∈ L (4)
i∈Vd j∈Γ+ (i)
3.1. Extended customer-demand network
∑ ∑
xhikl = 1 ∀k ∈ K, l ∈ L (5)
An extended network defined in a complete directed graph G (V , A )
′ ′ ′
i∈Vd h∈Γ− (i)
facilitates the implementation of the branch-and-price-and-cut algo-
∑
yklim ⩽ xklij ∀i ∈ Vd , k ∈ K, l ∈ L, m ∈ M (6) rithm. The extended vertex set V = Vc ∪ Vd includes all extended
′ ′

j∈Γ+ (i)
vertices and depot vertices, where the vertex set Vc = {(i, m)|i ∈Vc , m ∈
′

∑∑ Φ } denotes the set of extended customer-demand vertices. Note that if

yklim = Dmi ∀i ∈ Vc , m ∈ M (7) Dmi = 0, the extended vertex (i, m) is not involved in vertex set Vc . Arc set
′

k∈K l∈L
A = {((i, m), (j, l) )|(i, m), (j, l) ∈ V } depicts arcs between each pair of
′ ′

( )
τkli + si + tij − T 1 − xklij ⩽τklj ∀(i, j) ∈ A, k ∈ K, l ∈ L (8) vertices in the extended network. The travel cost ̃cim,ju and travel time
̃t im,ju along arc ((i, m), (j, l) ) are defined in equations (17) and (18):
{
ai ⩽τkli ⩽bi ∀i ∈ V, k ∈ K, l ∈ L (9) cij i ∕
=j
cim,ju =
̃ (17)
( ) 0 i=j
∑
qkli + dim yklim − Q 1 − xklij ⩽qklj ∀(i, j) ∈ A, k ∈ K, l ∈ L (10) {
m∈M tij + si i∕
=j
̃tim,ju = (18)
0 i=j
0⩽qkl
i ⩽Q ∀i ∈ V, k ∈ K, l ∈ L (11)
An illustration of the applied extension procedure is given in Fig. 2.
∑ ∑ The example given in Fig. 2 has the same parameters as the example
xhikl = xk,l+1
ij ∀i ∈ Vd , k ∈ K, l ∈ L\{|L|} (12)
h∈Γ− (i) j∈Γ+ (i) shown in Fig. 1. In Fig. 2, all customer sites are extended. Since there are
two categories of waste in this instance, each customer site is extended
τk,l+1 > τkli + σ ∀i ∈ Vd , k ∈ K, l ∈ L\{|L|} (13) to two dummy vertices. In the extended network, each extended vertex
i
of customer is visited exactly once. When the extended vertex (i, m) is
visited by a trailer, all of waste m should be collected. The travel time

4
L. Zhang et al. Computers & Industrial Engineering 178 (2023) 109143

linear relaxation solution of the restricted master problem. If the feasible

routes with negative reduced cost are found by the pricing subproblem,
these routes are added in subset Ω . If no such route can be found, the
′

optimal solution of linear master problem is obtained.

3.3. Pricing subproblem

Let αmi denotes the non-negative dual variable associated with con-
straints (20), and let β denotes the dual variable of constraint (21). The
(∑ )
m im
reduced cost cr associated with route r is cr = cr − (i,m)∈V ′ αi δr + β .

The pricing subproblem can be formulated as follows:

∑ ∑ ∑ ∑ ∑
minimize l
̃cim,ju xim,ju − αmi xlim,ju − β (23)
l∈L (i,m)∈V ′ (j,u)∈V ′ ′
(i,m)∈V (j,l)∈V
′

Subject to
∑
xlim,ju ⩽1∀(i, m) ∈ V , l ∈ L (24)
′

Fig. 2. An illustrational example of the network extension. (j,u)∈V

′

∑ ∑
and cost between different dummy vertices of the same customer site are xlim,ju = xlju,im ∀(i, m) ∈ V , l ∈ L
′
(25)
both zero. The Dantzig-Wolfe decomposition can be applied directly in (j,u)∈V
′
(j,u)∈V
′

the extended network since each vertex should be visited exactly once.
∑ ∑
The label setting algorithm can be applied to find a feasible route for xli0,ju = 1∀l ∈ L (26)
the pricing subproblem because each dummy vertex is visited exactly i∈Vd (j,u)∈V ′

once. However, the network extension enlarges the scale of network. ∑ ∑

More vertices in the network induce high degeneration in the master xlju,i0 = 1∀l ∈ L (27)
problem and decrease the convergence speed of label setting algorithm. i∈Vd (j,u)∈V ′

In this paper, two accelerating strategies are implemented to speed up ( )

the label setting algorithm, and a stabilized column generation is τlju ⩾τlim + ̃tim,ju − T 1 − xlim,ju ∀(i, m), (j, m) ∈ V , l ∈ L
′
(28)
designed to reduce the degeneration.
ai ⩽τlim ⩽bi ∀(i, m) ∈ V , l ∈ L (29)
′

3.2. Master problem

( )
qlju ⩾qlim + dim − Q 1 − xlim,ju ∀(i, m), (j, m) ∈ V \Vd , l ∈ L (30)
′ ′

To give the master problem, several additional notations are defined

as follows: Let route set Ω, which is indexed by letter r, denote all
qli0 = 0∀i ∈ Vd (31)
feasible routes. The cost of route r ∈ Ω is represented by cr . The
parameter δim
r indicates whether route r ∈ Ω visits extended vertex (i, m). 0 < qlim < Q∀(i, m) ∈ V (32)
′

When route r visits extended vertex (i, m), δim im

r = 1, otherwise δr = 0.
The binary variable λr = 1 when route r ∈ Ω is chosen in the solution, ∑ ∑
xlju,i0 = xl+1
i0,ju ∀i ∈ Vd , l ∈ L\{|L|} (33)
otherwise λr = 0. The master problem is modeled as follows: (j,u)∈V
′
(j,u)∈V
′

∑
minimize cr λr (19)
r∈Ω τl+1
i0 ⩾τi0 + σ ∀i ∈ Vd , l ∈ L\{|L|}
l
(34)

Subject to
xlim,ju ∈ {0, 1}∀(i, m), (j, u) ∈ V , l ∈ L (35)
′

∑
δim
r λr = 1 ∀i ∈ Vc , m ∈ M (20)
r∈Ω
Where binary variable xlim,ju indicates whether the arc between
extended vertices (i, m) and (j, u) is traversed by trip l ∈ L. Nonnegative
∑
λr ⩽|K| (21) variables τlim and qlim denote the service start time and vehicle load when
r∈Ω the trailer visits extended vertex (i, m), respectively.
Objective function (23) aims to find the shortest route with the
λr ∈ {0, 1} ∀r ∈ Ω (22)
minimum reduced cost. Constraints (24) guarantee that each extended
Objective function (19) aims to minimize the total cost. Set partition vertex can be visited at most once. Flow balance constraints (25) - (27)
constraints (20) stipulate that each extended vertex should be visited ensure that the optimal solution is a continuous route. The arrival time
exactly once. The number of required trailers is bounded by constraint of each extended vertex is calculated with constraints (28), and con-
(21). Constraints (22) denote that λr is a binary variable. straints (29) indicate that time windows should not be violated. Con-
The column generation procedure is applied to find the optimal so- straints (30) - (32) ensure that capacity of each trailer should be
lution of the linear relaxed master problem, which is obtained by respected. Constraints (33) and (34) state that the source depot of trip
relaxing constraints (22) to 0 < λr < 1 for any feasible route. However, l +1 ∈ L is the sink depot of trip l ∈ L, and the unloading time should be
route set Ω is difficult to be obtained since enumerating all feasible preserved. The domain of variable xlim,ju is defined in constraints (35).
routes is impractical even in medium scale instances. Instead of using The pricing subproblem can be regarded as a variant of elementary
the route set Ω, a subset of routes Ω ⊆ Ω is used to build a restricted shortest path problem with resource constraints, where resource con-
′

master problem. The optimal solution of the linear relaxation of the straints include capacity constraints and time window constraints
restricted master problem gives an upper bound of linear relaxed master (Hernandez et al., 2016; Feillet et al., 2004). We have recourse to the
problem. A pricing subproblem is needed to test the optimality of the label setting algorithm (Irnich and Desaulniers 2005; Righini and Salani

5
L. Zhang et al. Computers & Industrial Engineering 178 (2023) 109143

2008) to solve the pricing subproblem. However, due to the existence of by the branching strategies. The overall algorithm terminates when
inter depot, the load extension function is not monotonic. The tradi- there is no unexplored branch node, and the optimal solution is the
tional dominance strategy of the label setting algorithm should be current upper bound.
modified to ensure that the label extension procedure is feasible. The
label setting algorithm is introduced in detail in the next section. 4.1. Label setting algorithm

4. Stabilized branch-and-price-and-cut algorithm In this subsection, a label setting algorithm based on the extended
network is developed, and the dominance strategy is redesigned because
In this section, we introduce the stabilized branch-and-price-and-cut of the involvement of inter depot. In this section, two dummy vertices
algorithm to find optimal combinations of routes. The framework of the 0 and N are introduced to denote the source and sink vertices, respec-
proposed stabilized branch-and-price-and-cut algorithm is given in tively. In the label setting algorithm, let label Lim denotes the multiple
Fig. 3, where the stabilized column generation is embedded in the dimensional label that corresponds to the partial route from the source
branch-and-bound algorithm. At each branch node, the branch-and- vertex to extended vertex (i, m). There are several components in the
bound procedure is applied to find the optimal integer solution. At label Lim , which can be denoted as Lim = {(i, m), N(i, m), V(i, m),
each node of the branch-and-bound, the linear relaxation relaxed solu- q im , ̂c im }. Additional notations used in label Lim are defined as fol-
̂τ im , ̂
tion is obtained by the stabilized column generation. Specifically, a lows:
stabilized master problem and a pricing subproblem is solved interac-
tively to generate linear relaxation solution. The stabilized master • (i, m) denotes the current ending vertex of this partial route;
problem is solved by CPLEX, and pricing subproblem is solved by a
• N(i, m) ⊆ V is the set of all visited vertices of this partial route;
′

designed label setting algorithm. When there is no route with negative

• V(i, m) ⊆ V is the set of all vertices that can be extended by this
′

reduced cost and the introduced parameters of stabilized column gen-

partial route;
eration procedure are zero, the stabilized column generation terminates.
• ̂τ im is the earliest service start time in vertex (i, m);
Then, we test whether any cutting planes can be separated. When any
cutting planes is separated, the obtained cutting planes are added in the • ̂q im is the current load of the trailer when trailer visits vertex (i, m);
stabilized master problem, and the stabilized column generation pro- • ̂c im is the current reduced cost of this partial route.
cedure is performed again. When no cutting plane can be separated, the
optimal linear relaxation solution associated with current branch node is In the source vertex, the multi-dimension label is initialized to L0 =
obtained. If the obtained solution is an integer solution, the upper bound {(0, 0), N(0, 0), V(0, 0), 0, 0, 0 }. Feasible labels are extended to sink
will be updated. Otherwise, the children branch nodes will be generated vertex N with resource extension functions (36) – (40). Note that, the
initial label should be firstly extended to depot i ∈ Vd . The optimal

Fig. 3. Framework of proposed stabilized pricing-and-price-and-cut algorithm.

6
L. Zhang et al. Computers & Industrial Engineering 178 (2023) 109143

solution can be achieved by identifying the label with the minimum reach the crucial resource bound, are combined to obtain complete
reduced cost in the sink vertex. If label Lim = {(i, m), N(i, m), V(i, m), τim , routes. According to Righini and Salani (2006), the bounded bidirec-
q im , ̂c im } can be extended to vertex (j, u), components of new label can
̂ tional search strategy can be applied in the large scale instances to speed
be extended based on the following resource extension functions: up the convergence speed of the label setting algorithm. We use the
{ b b }
N(j, u) = N(i, m) ∪ {(j, u) } (36) notation Lbim = (i, m), Nb (i, m), V b (i, m), ̂τ bim , ̂
q im , ̂c im to depict the
backward label of the partial route between extended vertex (i, m) and
( { } { })
V(j, u) = V(i, m)\ {(j, u) } ∪ (k, h)|τju + tjk > bk ∪ (k, h)|̂
q im + dkh ⩾Q dummy sink vertex N. In the backward label notation, sets Nb (i, m),
b b
(37) V b (i, m) and parameters ̂c im , ̂
q im share similar definitions with that in the
{ } forward label. ̂τ bim
denotes the latest service start time in vertex (i, m).
̂τ ju = max ̂τ im + ̃tim,ju , aj (38) The resource extension function and dominance rule of the backward
searching are similar to that of the forward searching. In this paper, the
{ time resource is selected as the most crucial resource to build the
q im + dim , i ∕
∈ Vd
searching bound, namely, ̂τ im ⩽T/2 and ̂τ bim ⩾T/2.
̂
q ju =
̂ (39)
0, i ∈ Vd K-cycle elimination: K-cycle elimination of the shortest path
problem with resource constraints with k⩾3 is proposed by Irnich and
c ju = ̂c im + ̃cim,ju − αmi
̂ (40) Villeneuve (2006) to forbid cycles with short lengths. Namely, only
Within the extension procedure, the number of labels increases cycles with at least k +1 length are allowed to be generated. According
exponentially, which would affect the convergence speed of the label to numerical results given by Irnich and Villeneuve (2006), k-cycle
setting algorithm. The dominance rules can improve the convergence elimination with k⩾3 cases would improve the lower bound and accel-
speed by eliminating useless labels. The useless label denotes the erate the convergence speed of the label setting algorithm. The case k =
dominated label, whose extended labels can also be extended by other 2, which forbids sequences with the form i −j −i, is applied in this paper,
existing labels with lower reduced cost. The dominance rule presented since this case can give a strong lower bound with a practical complexity
below is a feasible dominance rule without eliminating any solution that increase in the label setting algorithm.
could be the possible optimal solution.
Dominance rule: Given two labels Lim = {(i, m), N(i, m), V(i, m), ̂τ im , 4.2. Branching strategies
{ }
q im , ̂c im } and Lju = (j, u), N(j, u), V(j, u), ̂τ ju , ̂
̂ q ju , ̂c ju associated with
The column generation is solved to obtain the linear relaxation so-
two partial routes from source vertex 0 to vertices (i, m) and (j, u), lution of the restricted mater problem at each branch node. Then,
respectively. If the following conditions hold, label Lim can dominate Lju branching strategies are applied if the linear relaxation solution at this
safely. node is lower than the current best upper bound. Note that the master
(i, m) = (j, u) (41) problem is solved by CPLEX directly to find an integer solution in this
paper. The integer solution obtained by CPLEX can give a good upper
V(j, u) ⊆ V(i, m) (42) bound in the early age of the branch-and-bound procedure and help us
cut branches. When a better integer solution is achieved, the upper
̂τ im ⩽̂τ ju (43) bound is updated. In this research, new nodes are generated based on
four branching strategies, which are branching on the total number of
q im ⩽̂
̂ q ju (44) required trailers, branching on the numbers of each customer be visited,
branching on the flow on each arc, and branching on the flows on
c im ⩽̂
̂ c ju (45) consecutive arcs. Branching nodes are searched based on the best-first
policy.
Proof of dominance rule: When conditions (41) and (42) hold, all
feasible extensions of L2 must be feasible extensions of L1 (without
consideration of time window and capacity constraints). Assuming that 4.3. Cutting planes
vertex (k, h) is a feasible extension of L1 and L2 , conditions (43) and (45)
1 2 1 Two types of cuttings, which are k-path inequality (Desaulniers,
still hold because ̂τ kh = ̂τ im + ̃t im,kh ⩽̂τ kh = ̂τ ju + ̃t ju,kh , ̂c kh =
2 2010, Kohl et al., 1999) and subset row cuts (Jepsen et al., 2008, Pecin
q im +̃cim,kh −α
̂ m
c kh
i ⩽̂ = ̂c ju +̃cju,kh −α and (i, m) = (j, u). With regard to
u
j et al., 2016), are implemented. K-path cut is a robust cut and does not
1
condition (44), when k ∈ Vd the condition (44) still holds since ̂
q kh = change the label setting algorithm. However, the sub-row cut is a non-
2
q kh = 0. When (k, h)∈ V\Vd , expressions ̂
1
q kh = ̂
q im +dm
2
q kh = ̂
q ju +duj robust cut, and the label setting algorithm should be modified to
i ⩽̂
̂
and (i, m) = (j, n) make condition (44) hold. Therefore, any feasible include the dual variables of sub-row inequalities.
extension of L2 is also feasible to L1 , and the cost of every feasible label K-path inequality: The basic formulation of the k-path inequality is
x(S)⩾k(S), where S is a subset of extended vertices, that is S ⊆ V . In the
′
extended by L1 is lower than that of the label extended by L2 . □.
The convergence speed of the label setting algorithm is slow, when integer solution, subset S should be served by at least k(S) trailers.
there are many extended vertices. Therefore, in this paper, several Namely, there are at least k(S) paths that should enter subset S in the
accelerating strategies, which are bounded bidirectional search and k- integer feasible solution. The notation x(S) denotes the total entering
cycle elimination, are applied to speed up the convergence speed of the flow of subset S. According to the previous notations, the formulation of
label setting algorithm. the k-path inequality is presented as following:
Bidirectional search: The bidirectional search (Righini and Salani, ∑ ∑
δim (46)
′
r λr ⩾k(S) ∀S ⊆ V
2006, 2008) is an accelerated strategy for the label setting algorithm, r∈Ω (i,j)∈Γ− (S)
which extends labels in both forward and backward directions. The
typical procedure can be summarized as follows: Firstly, we set the Let πS denotes the dual variable of the k-path inequalities (46). The
∑ ∑ ∑
searching bound of the bidirectional search with the most crucial reduced cost of route r should ∑be∑
modified to minimize
∑ ∑ l∈L (i,m)∈V
′

(j, u) ∈ V ̃cim,ju xlim,ju − αmi xlim,ju − πS xlim,ju − β

′

resource; The forward and backward partial routes are generated based
S⊆V i∈Γ (S)
−
on the resource extension function and eliminated based on the domi-
′ ′ ′
(i,m)∈V (j,l)∈V

nance rules; Finally, the forward and backward partial routes, which . The most difficult work for applying k-path inequality is to identify the

7
L. Zhang et al. Computers & Industrial Engineering 178 (2023) 109143

violation of inequalities (46). Several separation algorithms were pro- Table 3

posed in the literature. For instance, Desaulniers (2010) gave a heuristic Additional notations for the stabilized column generation.
recursive algorithm. In this paper, we use the k-path inequality with the Parameters
case k = 2, which is called 2-path cuts (Kohl et al., 1999). To identify the
Left and right bounds of the trust region of extended vertex (i, m) in
im ,ηim
η−,l + ,l
violation, we apply the greedy heuristic algorithm proposed by Kohl major iteration l.
et al. (1999) to recursively find subset S with x(S)⩽2. Once a subset S is Left and right bounds of the inner penalty region of extended vertex (i,
im ,μim
μ−,l + ,l

obtained, the label setting algorithm proposed in Subsection 4.1 is m) in major iteration l.
applied in S to find a feasible route to cover all of the vertices. If the label Penalty parameters in the left (or right) inner penalty region of
im ,εim
ε−,l + ,l

extended vertex (i, m) in major iteration l.

setting algorithm cannot give a feasible route, subset S cannot be served Penalty parameters in the left (or right) outer penalty region of
im ,ζim
ζ−,l + ,l
by one trailer. Then, a violation of 2-path inequality is identified. extended vertex (i, m) in major iteration l.
Limited-memory subset row cut: Jepsen et al. (2008) proposed the Variables
subset row cut to strengthen the lower bound at each branch node. The y−
im ,yim Surplus and slack variables of extended vertex (i, m) in major iteration l
+

form of subset row cut can be written as with the bound of inner region penalty parameters ε−,l im ,εim .
+ ,l

im ,zim
z− + Surplus and slack variables of extended vertex (i, m) in major iteration l
∑ 1 ∑ with the bound of outer region penalty parameters ζ−,l
im ,ζim .
|C| + ,l
δim ⌋λr ⩽⌊ ⌋ (47)
′
⌊ ∀C ⊆ V
r∈Ω
p (i,m)∈C r p

Where C is a subset of the extended vertex set. There are various ∑

λr ⩽|K| (50)
subset row cuts with different |C| and p. In this paper, we choose the 3- r∈Ω
subset row cut with |C| = 3 and p = 2. With the application of subset
row cut, the reduced cost of route r should be modified. The modified (51)
′
0⩽z−
im ⩽ζ im
−,l
∀(i, m) ∈ V
reduced cost of the pricing subproblem would change the label setting
algorithm proposed in Subsection 4.1. A new resource should be added im ⩽εim
0⩽y− −,l
∀(i, m) ∈ V
′
(52)
to represent every k times visit the subset C. Since this modification
makes the label setting algorithm difficult to be solved, we implement a 0⩽z+
′
(53)
im ⩽ζ im ∀(i, m) ∈ V
+,l

weaker variant of the subset row cut, called the limited-memory subset
row cut (Pecin et al., 2016). The limited-memory subset row cut can be (54)
im ⩽εim
′

∑ 0⩽y+ +,l
∀(i, m) ∈ V
defined as r∈Ω α(C, M, k, r)λr ⩽⌊|C|
k
⌋, where subset M is the memory set
and satisfies the relaxationship C ⊆ M ⊆ V . The separation algorithm of
′
λr ⩾0 ∀r ∈ Ω (55)
the limited-memory subset row cut and the calculation algorithm of Note that in this research only set partitioning constraints (20) are
memory set M and α(C, M, k, r) are given in Pecin et al. (2016). stabilized. Penalties of surplus and slack variables are added into
objective function (48). Constraints (49) give stabilized set partitioning
4.4. Stabilized column generation constraints. Note that this stabilization strategy would enlarge the
domain of the master problem, and feasible solutions of the master
As noted above, the extended network leads a simplified label setting problem are obtained when surplus and slack variables are all zero.
algorithm, since each extended vertex can only be visited once. How- Constraints (50) are as same as constraints (21). The upper bound of
ever, this technique enlarges the network and leads to higher degener- surplus and slack variables in major iteration l are defined in constraints
ation in the master problem. When we use the standard column (51) – (54). In this subsection, we use αlim to define the dual variable of
generation algorithm, the convergence rate is unsatisfactory. Therefore,
constraints (49). Let u−im , vim , uim and vim be non-negative dual variables of
− + +
the stabilization strategy proposed by Amor and Desrosiers (2006) is
constraints (51) - (54). The stabilized dual master problem can be
used to modify the standard column generation for acquiring the faster
written as follows:
convergence speed. Detailed theory of stabilization strategy can be
∑ ( )
found in Amor and Desrosiers (2006). The basic idea of applied stabi- maximize αlim − ζ−,l
im uim − εim vim − ζ im uim − εim vim − β
− −,l − +,l + +,l +
(56)
lization strategy can be described as follows: When dual variables are far (i,m)∈V
′

away from the stability center, a penalty based on the penalty function ∑
should be added into the objective function. In the following para- αlim δim
r + β⩽cr ∀r ∈ Ω (57)
graphs, we will describe the stability center, penalty function, stabilized (i,m)∈V
′

master and dual problem, and update strategy that are used in this
research. η−,l
im − uim ⩽αim ⩽ηim + uim
l
(58)
′
− +,l +
∀(i, m) ∈ V
l
The stability center ̂π is initialized to the zero vector and updated to
μ−,l
im − vim ⩽αim ⩽μim + vim
l
(59)
′
current dual variables if no feasible route with negative reduced cost can
− +,l +
∀(i, m) ∈ V
be found. A five-piecewise linear penalty function gl with a non-
(60)
′
u−
im , uim , vim , vim ⩾0
+ − +
∀(i, m) ∈ V
penalized trust region around stability center ̂ π l is applied. In this no-
tation, index l denotes the major iteration, in which the stability center
β⩾0 (61)
and parameters in the linear penalty function are updated. Additional
notations for building the penalty function and the stabilized master From constraints (58) and (59), we can find that the dual variable αlim
problem (also called the stabilized primal problem) are presented in is restricted within the trust region and inner penalty region with a soft
Table 3, and penalty function gl is illustrated in Fig. 4. bound. If dual variables do not fall in the trust region, a penalty is added.
With above additional notations, the stabilized master problem In this paper, the stability center is initialized to the zero vector. The
(stabilized primal problem) can be formulated as follows: half-width of the trust region is initialized to Δη = 5, and the half-width
∑ ∑ ( ) of the inner region is initialized to Δμ = 50. Penalty parameters within
minimize cr λr + − η−,l
im yim − μim zim + ηim yim + ηim zim
− −,l − +,l + +,l +
(48)
r∈Ω ′
the inner region are ε−,0 = ε+,0 = 0.3. Note that the relationship be-
(i,m)∈V
tween inner and outer penalty parameters ε +ζ⩾1 should be satisfied.
∑ We set initial values of outer penalty parameters to ζ−,0 = ζ+,0 = 1.
δim (49)
′
r λr − yim − zim + yim + zim = 1 ∀(i, m) ∈ V
− − + −

r∈Ω Once the optimal solution of model (48) – (55) is obtained, the stability

8
L. Zhang et al. Computers & Industrial Engineering 178 (2023) 109143

,l
,l

gl gl
l
ˆ

,l ,l ,l ,l

,l ,l
,l ,l
gl gl

Fig. 4. Five-piecewise linear penalty function in major iteration l.

center is updated to new optimal dual variables. Parameters of the penalties ζ− and ζ+ are constant. However, the maximum width of the
penalty function are updated dynamically depending on the relationship trust region, the maximum width of the inner region and the penalty in
between current dual variables and the trust region (or inner region). A the inner region are 20, 500 and 1, respectively. The minimum width of
detailed description of the updating strategy is described as follows: the trust region, the minimum width of the inner region and the penalty
[ ]
If πlim ∈ η−,l in the inner region are 1, 10 and 0.001, respectively.
im , ηim , we decrease the length of the trust region by half,
+,l

and double the penalty parameter of the inner penalty region. If πlim ⩽η−,l
im , 5. Computational results
we double the left trust region length between the left bound and the
stability center and decrease the penalty in the left inner penalty region In this section, the proposed stabilized branch-and-price-and-cut
ε− by half. If πlim ⩾η+,l
im , we double the right trust region length between algorithm is tested on a set of modified benchmark instances and a
the right bound and the stability center and decrease the penalty in the real-world instance collected in Beijing. All of the computations are
right inner penalty region ε+ by half. When πlim ⩽μ−,l l
im (or π im ⩾μim ), the
+,l
performed on a Windows PC with an Intel Core i7 8700 CPU, 8 GB RAM.
width of the left (right) outer penalty region is enlarged similarly, but The stabilized branch-and-price-and-cut algorithm is coded in C++ with
the Visual Studio 2015 compiler, and the master problem is solved with

Fig. 5. Customer and depot sites of the practical instance.

9
L. Zhang et al. Computers & Industrial Engineering 178 (2023) 109143

IBM CPLEX 12.8. customer sites 3, 13, and 18. In route 1, waste categories 3 and 4 of
customer site 3 are collected, and the remaining two categories of waste
are collected by route 4. The services of customer sites 13 and 18 are
5.1. A real-world case study split to routes 2 and 3, respectively.
We also test the case that the split transportation is not allowed. The
In this subsection, the proposed stabilized branch-and-price-and-cut minimum travel time is 1264 min, and five trailers are required. The
algorithm is applied to design the routes for a waste transportation case transportation strategy with split transportation can save 4.03% on
study in Beijing, China. A fleet of trailers collects four categories of waste travel time. The detailed routes associated with the case without split
(i.e., recyclables, toxic substances, compostable organic matter, and transportation is presented in Table 5 and Fig. 7. It can be found that
soiled waste) from 31 customer sites. Three depots (i.e., transfer sta- without split transportation, trailers tend to return to depots more
tions) are available. The geographical distributions of given customer frequently than that in the case with split transportation.
sites and depots are illustrated in Fig. 5, where triangles and circles In conclusion, the case with 31 customer sites, three depots, and four
represent depots and customer sites, respectively. The diameter of each waste categories can be solved to the optimum by the proposed algo-
circle denotes the total waste demand at this customer site. rithm within several seconds. In the context of waste classification, split
Within the planning horizon, each trailers starts route from its transportation leads to flexible routes and helps operators save on
parking depot and return an arbitrary depot. That is, trailers are not operational cost.
obligated to return to their source depots. For the convenience of
calculation, we build dummy source and sink vertices to represent the
source and sink depot of trailers. The travel time and travel cost between 5.2. Algorithm performance
two dummy vertices and three depots are set as zero. The travel time and
travel cost between dummy vertices and customer sites are set as the We test the performance of proposed algorithm on several modified
minimum values of travel time and travel cost between three depots and benchmark instances based on those proposed by Solomon (1987). As
{ }
customer sites, namely di,n+nd +1 = minj∈Vd dij . In most real-world Archetti et al. (2015) did, only benchmark instances C101 and R101 are
considered, where “C” and “R” denote the instances with clustered and
municipal waste managerial applications, some waste, such as com-
random customer sites, respectively. The instances with different scales
postable organic matter and hazardous materials (Ghaderi and Burdett,
are generated, including small scale instances with 20 customers, me-
2019), is expected to be collected and transported early because of their
dium scale instances with 40 and 60 customers and large scale instances
specific environmental impacts. In this instance, the time windows of
with 80 and 100 customers. The instances with two and three categories
customer sites are restricted within several hours so that the accumu-
of waste are generated. Each customer sites have the demand of different
lation of waste would not affect the surrounding residents. The total
categories of waste with probability p (p = 1 or p = 0.6). The demand at
travel time is minimized. The computational results of applying stabi-
each customer site is randomly generated from intervals [1,100] or
lized branch-and-price-and-cut algorithm to solve this problem is 5.78 s,
and the details of obtained routes are given in Table 4. [40,60] (Archetti et al., 2015). Two types of trailers with capacity Q =
( ∑ )
The total travel time of all trailers is 1215 min. In Table 4, the column α maxi∈n m∈Φ dm i are used, where parameter α can take values 1 or 0.6..
“Routes” shows the detailed routes, where the letter di (i = 1, 2, 3) in- For each instance, the locations of depots are generated randomly within
{ [ ]}
dicates the visiting depots. The numbers in parentheses denotes the the region (x, y)|x ∈[xmin , xmax ], y ∈ ymin , ymax , where [xmin , xmax ]
[ ]
categories of collected waste. The column “Split transportation” gives ( ymin , ymax ) denotes the minimum and maximum x-coordinate (y-co-
the customer sites that are visited more than once. The column “Inter- ordinate) of all customer sites in each instance.
mediate depot index” denotes depots visited in the routes.
According to the routes presented in Table 4, the total demands can 5.2.1. Performance of stabilized column generation
be satisfied by five trailers. Four trailers start their routes in depot 2, and We firstly test the performance of the introduced stabilization col-
only two trailers return to this depot. Details of these routes are reported umn generation technique on the instance C101 with 60 customers and
in Fig. 6. From Fig. 6, we can find that trailers mainly serve customer three demand categories. Each customer site requires all different
sites that are close to their starting depot, except for trailers 4 and 5. commodities (p = 1). The demand is generated within the demand in-
Trailers 4 and 5 start their routes in depot 2, and finally return to depots terval [1,100]. The capacity coefficient α of the trailers is 1.5.
1 and 3, respectively. The choice of parking depot is affected mainly by We only apply the stabilized column generation technique at the root
the last customer sites of routes (such as customer sites 3 and 21 in node, since the column generation procedure converges quickly at other
routes 4 and 5, respectively). Three split transportation occurs at nodes. Detailed computational results are reported in Fig. 8, and pa-
rameters of the stabilization technique are set as that presented in sub
Table 4 *** Section 4.4. Since the surplus and slack variables are involved in
The details of obtained routes with split transportation. constraints (49), solutions of stabilized column generation are infea-
Index Routes Split Intermediate sible. Thus, the objective function values of the stabilized column gen-
transportation depot index eration are lower than those of the standard column generation within
1 d1 −25(1,2,3,4)-6(1–2-3–4)-d1 -3 3 1 the computational process. In the standard column generation without
(2,4)-7(1,2,3,4)-27(1–2-3–4)-26 the stabilization technique, the optimal solution is accessed in the 721st
(1–2-3–4)-d1 iteration. However, the stabilized column generation can access the
2 d2 −13(3,4)-14(1,2,3,4)-d2 -11 13 2, 2 optimal solution within 300 iterations. Therefore, we can conclude that
(1,2,3,4)-12(1,2,3,4)-d2 -13(1,2)-
29(1,2,3,4)-28(1,2,3,4)-d2
the stabilization technique can reduce the degeneration of the proposed
3 d2 −15(1,2,3,4)-18(1,2)-17 18 3, 2 algorithm.
(1,2,3,4)–23(1,2,3,4)-d3 -19
(1,2,3,4)-20(1,2,3,4)-18(3,4)- 5.2.2. Performance on small scale instances
d2 -10(1,2,3,4)-30(1,2,3,4)-d2
Then, we perform the proposed algorithm to solve the small scale
4 d2 −31(1,2,3,4)-2(1,2,3,4)-1 3 1
(1,2,3,4)-d1 -5(1,2,3,4)-4(1,2,3,4)- instances with 20 customers, and the computational results are given in
3(1,3)-d1 Table 6. The first column shows the names of instances. The column “m”
5 d2 −9(1,2,3,4)-8(1,2,3,4)-16 – 2 denotes the number of categories. The column “α” depicts the parameter
(1,2,3,4)-d2 -24(1,2,3,4)–22 of trailer capacity, and the column “p” gives the probability that cus-
(1,2,3,4)-21(1,2,3,4)-d3
tomers require each category. Values of demand are generated within

10
L. Zhang et al. Computers & Industrial Engineering 178 (2023) 109143

Fig. 6. Illustration of routes with split transportation.

with 20 customers, eight medium scale instances with 40 or 60 cus-

Table 5
tomers and eight large scale instances with 80 or 100 customers are
The details of obtained routes without split transportation.
studied in this part. The numbers of categories are two or three.
Index Routes Source Sink Intermediate depot The overall computational results are presented in Table 7. The
depot depot index
columns “z* ” “z* ”, and “Gap(%)” share the same meaning with those in
1 d1 −25–6-d1 -3–5-7-d1 -4- d1 d1 1, 1 Table 6. The computational time is reported in column “t(s)”. Except for
d1
2 d2 −11–12-d2 -29–28- d2 d2 2, 2
C101 and R101 instances with 100 customers and three demand cate-
d2 -10–30-d2 gories, all other instances can be solved within one hour. For C101 in-
3 d2 −15–13-d2 -18–17-23- d2 d3 2, 3 stances, when numbers of customers and categories increase, gaps
d3 -19–20-d3 shown in the column “Gap(%)” increase accordingly. This is expected
4 d2 −14-d2 -31–2-1- d2 d1 2, 1
because the increasing of customers and categories leads to a larger
d1 -27–26 d1
5 d2 −9–8-16-d2 -24–22- d2 d3 2 extended network. In general, the proposed algorithm can find relatively
21-d3 good solutions within the acceptable computational time, even in large
scale instances.

the parameters given in the column “Δ”. The obtained integer solutions 5.3. Analysis and managerial implications
are presented in the column “z* ”, and the linear relaxation solutions are
given in the column “z* ”. The column “Gap(%)” shows percentage gaps 5.3.1. The impact of split transportation
between integer solutions and linear relaxation solutions, which are To discuss the managerial meaning of split transportation, we
calculated with Gap = (z* − z* )/z* × 100%. Finally, the column “t(s)” compare the computational results of the modified benchmark instances
reports the computational time in seconds. with and without split transportation. When split transportation is not
The optimal solutions of twelve instances are obtained by the pro- allowed, the demand values of each customer site are the summation of
posed algorithm within the time limitation, and the other fourteen in- demand values over all categories. The number of categories, capacity
stances can be solved within 1% gaps. Gaps larger than 1% arise in six coefficient, and demand interval are the same as those in Table 7. In
instances with three commodities and probability 1.0. Due to the Table 8, columns “z*NSP ” and “z*NSP ” give integer solutions and linear
network extension, more extended customer vertices in these instances relaxation solutions without split transportation. The results shown in
make the pricing subproblem more cumbersome to solve. In general, all ( )/
the column “Gap1 (%)” are calculated with Gap1 = z*NSP − z*NSP z*NSP ×
instances can be solved to be optimal or near-optimal within several
100%. The column “z*SP ” denotes integer solutions of instances with split
dozens of seconds, which means that our algorithm has good conver-
gence efficiency in small scale instances. transportation, which are as same as results shown in the column “z* ” in
Table 7. The results shown in the column “Gap2 (%)” are calculated with
( )/
5.2.3. Performance on medium and large scale instances Gap2 = z*NSP − z*SP z*NSP × 100%. Columns “|KNSP |” and “|KSP |” report
We test the proposed algorithm on large scale instances in this sub- the number of required trailers with and without split transportation,
section.. The demand interval is set as [1,100], and the capacity coef- respectively.
ficient α is set as 1.5. There are twenty different instances with different From Table 8, we can find that allowing split transportation in all
numbers of categories and customers, in which four small scale instances instances can achieve better solutions than those in the instances
without split transportation, except for the C101 instance with 100

11
L. Zhang et al. Computers & Industrial Engineering 178 (2023) 109143

Fig. 7. Illustration of routes without split transportation.

Fig. 8. Objective function values of the stabilized and standard column generations.

customers and two categories and the R101 instance with 20 customers R101 instances, introducing split transportation can only save 0.89%
and two categories. The optimal solutions of these instances with split operational cost. Thus, it is more efficient to introduce split trans-
transportation are equal to those without split transportation. portation in the instances with cluster customer sites than in those with
Table 8 also shows that the number of required trailers with split random customer sites. Moreover, in the most real-world applications,
transportation may be smaller than that of the case without split the customer sites are typically clustered, such as the waste collection
transportation. For example, 21 trailers are required in the case with sites in residential areas. Therefore, the waste transportation companies
split transportation in the C101 instance with 80 customers and three will benefit from applying split transportation.
categories, while more trailers are required without split transportation.
The last column “Gap2 (%)” shows the improvement of introducing 5.3.2. The impact of demand intervals
split transportation. It can be found that for C101 instances, introducing In this subsection, the impact of demand intervals is discussed.
split transportation can lead to average 3.28% operational cost save. For Specifically, the instances with different demand intervals [1,100],

12
L. Zhang et al. Computers & Industrial Engineering 178 (2023) 109143

Table 6
Computational results of small scale instances.
Instance z* z* Gap(%) t(s)

m α p Δ

C101 2 1.1 0.6 [0,100] 353.48 353.48 0.00 4.13

[40,60] 552.48 552.48 0.00 5.47
1 [0,100] 564.13 560.79 0.59 12.78
[40,60] 754.58 754.58 0.00 13.87
1.5 0.6 [0,100] 257.35 255.83 0.59 17.78
[40,60] 446.56 446.56 0.00 9.78
1 [0,100] 335.69 333.83 0.55 17.52
[40,60] 684.35 682.66 0.25 21.33
3 1.1 0.6 [0,100] 582.03 578.91 0.54 83.71
[40,60] 494.13 494.13 0.00 17.01
1 [0,100] 550.78 540.33 1.90 206.96
[40,60] 908.30 903.72 0.50 150.38
1.5 0.6 [0,100] 267.16 267.16 0.00 281.61
[40,60] 409.04 406.02 0.74 80.68
1 [0,100] 494.00 479.39 2.96 315.36
[40,60] 621.21 611.68 1.53 620.77
R101 2 1.1 0.6 [0,100] 457.21 457.21 0.00 4.10
[40,60] 555.53 555.53 0.00 7.73
1 [0,100] 663.10 663.10 0.00 4.71
[40,60] 953.29 949.18 0.43 3.59
1.5 0.6 [0,100] 536.40 536.40 0.00 2.34
[40,60] 514.22 510.75 0.67 1.84
1 [0,100] 609.96 609.96 0.00 15.99
[40,60] 753.49 751.28 0.29 14.37
3 1.1 0.6 [0,100] 511.03 506.03 0.98 15.79
[40,60] 680.88 677.47 0.50 13.05
1 [0,100] 866.26 849.71 1.91 235.64
[40,60] 991.70 983.77 0.80 69.52
1.5 0.6 [0,100] 548.76 546.86 0.35 16.28
[40,60] 537.13 537.13 0.00 15.61
1 [0,100] 866.95 849.71 1.99 233.75
[40,60] 870.80 855.36 1.77 336.96

Moreover, the average gaps over instances with [1,100], [20,80],

Table 7
and [40,60] are 2.84%, 4.93%, and 5.24%, respectively. It can be found
Computational results of medium and large scale instances.
that the gaps increase when demand interval changes from [1,100] to
Instance z* z* Gap(%) t(s) [40, 60]. Namely, introducing split transportation in the instance with
n m
centralized demand values (i.e., the demand values of different customer
C101 20 2 300.39 300.39 0.00 4.24 sites are similar) can help operators implement more flexible routes and
3 388.44 381.62 1.76 96.21
save more operational costs.
40 2 584.32 576.94 1.26 49.01
3 844.04 815.22 3.41 241.24
60 2 926.79 914.14 1.37 221.07 5.3.3. The impact of trailer capacity
3 1486.66 1419.79 4.50 2629.56 This subsection shows the impact of trailer capacity on the instances
80 2 1397.37 1347.98 3.53 547.75 with different trailer capacity, ranging from 150 units to 350 units. The
3 2255.02 2143.46 4.95 1467.60
numbers of customers and demand categories are 20 and 3, respectively.
100 2 2082.12 1957.56 5.98 1947.56
3 2968.38 2817.46 5.08 3945.21 The demand of each category is generated from the ranges [1,100] or
R101 20 2 444.90 444.90 0.00 4.78 [40, 60], and the probability of requiring each category is 1.0. Table 10
3 514.34 513.19 0.22 16.78 presents the computational results. It can be found from Table 10 that
40 2 830.23 830.23 0.00 19.06
with the increase of trailer capacity, the operational costs of the cases
3 972.41 971.78 0.06 129.18
60 2 1222.71 1222.71 0.00 78.98 with and without split transportation both decrease. This is expected
3 1459.27 1443.12 1.11 434.67 since the tailers with larger capacity can perform more flexible routes
80 2 1658.09 1647.54 0.64 193.38 and yield smaller operational cost. Moreover, with the increase of trailer
3 1934.80 1870.27 3.34 1934.35 capacity, the gaps between the operational costs of cases with and
100 2 1998.47 1984.89 0.68 432.46
without split transportation decrease correspondingly. It can be
3 2184.57 2122.57 2.84 4834.14
concluded that the impact of introducing split transportation is more
significant if the tailers with smaller capacity are applied.
[20,80] and [40,60] are solved. The number of customers, number of
demand categories, parameter of trailer capacity, and probability of 5.3.4. Summaries
requiring each category are 20, 2, 1.1, and 1.0, respectively. The According to the analysis on the impact of split transportation, de-
computational results are given in Table 9. mand intervals, and trailer capacity, we can obtain the following
The average gaps over C101 instances with two categories, C101 implications:
instances with three categories, R101 instances with two categories, and Implication 1: Introducing split transportation can save operational
R101 instances with three categories are 3.32%, 6.68%, 2.71%, and cost and may help operator decrease the fleet of trailers.
4.64%, respectively. It can be concluded that, for the instances with Implication 2: The impact of introducing split transportation is more
more categories, split transportation can help operators implement more significant in the instance with cluster customer sites than in the
flexible routes and save more operational cost. instance with random customer sites.

13
L. Zhang et al. Computers & Industrial Engineering 178 (2023) 109143

Table 8
Computational results of instances with and without split transportation.
Instance z*NSP z*NSP |KNSP | Gap1 (%) z*SP |KSP | Gap2 (%)

n m

C101 20 2 310.01 310.01 3 0.00 300.39 3 3.10

3 417.75 417.75 4 0.00 388.44 4 7.02
40 2 585.26 585.26 5 0.00 584.32 5 0.16
3 904.69 904.69 7 0.00 844.04 7 6.70
60 2 936.40 927.72 11 0.93 948.40 11 1.03
3 1560.06 1559.58 14 0.03 1486.66 14 4.70
80 2 1413.43 1380.43 17 2.33 1397.37 17 1.14
3 2379.35 2367.35 22 0.50 2255.02 21 5.23
100 2 2082.12 2069.18 24 0.62 2082.12 24 0.00
3 3084.18 3070.54 30 0.44 2968.38 31 3.75
R101 20 2 444.90 444.90 5 0.00 444.9 5 0.00
3 527.73 527.73 7 0.00 514.34 7 2.54
40 2 836.13 836.13 12 0.00 830.231 12 0.71
3 975.45 975.45 14 0.00 975.446 14 0.31
60 2 1227.93 1225.86 15 0.17 1225.86 15 0.43
3 1496.68 1487.73 23 0.60 1459.27 23 2.50
80 2 1658.12 1635.66 22 1.35 1658.09 23 0.01
3 1937.70 1931.23 26 0.33 1934.8 26 0.15
100 2 2000.92 1998.72 27 0.11 1998.47 27 0.12
3 2231.48 2223.94 31 0.34 2184.57 31 2.10

similar, designing the routes with split transportation is beneficial for

Table 9
operations.
The computational results of instances with different demand intervals.
Implication 5: It is better to introduce split transportation in the case
Instance z*SP z*NSP Gap(%) where the tailers with smaller capacity are applied.
m Δ

C101 2 [1,100] 415.22 428.70 3.14 6. Conclusions

[20,80] 639.29 658.19 2.87
[40,60] 612.92 638.08 3.94
3 [1,100] 538.08 555.11 3.07 In this paper, we study a waste transportation problem with split
[20,80] 633.43 690.29 8.24 transportation and inter depot. This problem is motivated by the waste
[40,60] 574.63 629.63 8.74 transportation application with the background of waste classification.
R101 2 [1,100] 866.70 885.21 2.09 One key feature distinguishing the investigated problem from the
[20,80] 824.87 858.21 3.88
traditional waste transportation problems is that the split transportation
[40,60] 904.37 924.36 2.16
3 [1,100] 723.31 746.12 3.06 is allowed, which can save on operational cost for waste transportation
[20,80] 834.29 875.81 4.74 companies. We propose a stabilized branch-and-price-and-cut algorithm
[40,60] 898.40 957.06 6.13 based on an extended network. With the network extension, the inves-
tigated problem can be simplified to a transportation problem with the
single visit. A stabilization technique is applied to improve the conver-
Table 10 gence speed of the proposed algorithm. The proposed algorithm is tested
The computational results of instances with different trailer capacities. on the modified benchmark instances. The results of modified bench-
Instance z*SP z*NSP Gap(%) mark instances show that the proposed algorithm can solve the inves-
Δ Q tigated problem efficiently. The managerial importance of split
C101 [1,100] 150 639.47 606.79 5.11
transportation is studied based on modified benchmark instances. Ac-
200 571.28 543.97 4.78 cording to computational results, the variance of demand distribution
250 524.10 507.75 3.12 and the number of commodities would affect the effect of split trans-
300 497.32 483.49 2.78 portation. The case study of a real-world instance is used to show the
350 488.66 476.25 2.54
practical importance of the investigated problem and proposed algo-
[40,60] 150 711.82 646.12 9.23
200 602.11 547.62 9.05 rithm. The optimal solution of this instance can be approached within
250 533.72 489.90 8.21 several seconds. When allowing split transportation, a combination of
300 527.28 486.78 7.68 routes with lower total travel time can be provided to help operators
350 506.95 469.79 7.33
save on operational cost compared with that in the case without split
R101 [1,100] 150 873.80 813.77 6.87
200 769.23 737.54 4.12
transportation.
250 692.10 670.58 3.11
300 672.45 653.22 2.86 CRediT authorship contribution statement
350 631.28 613.79 2.77
[40,60] 150 1123.49 1026.53 8.63
200 934.86 860.26 7.98
Li Zhang: Conceptualization, Methodology, Writing – original draft.
250 873.55 820.96 6.02 Zhongshan Liu: Conceptualization, Software. Wenxuan Shan:
300 832.15 786.88 5.44 Conceptualization, Validation. Bin Yu: Conceptualization, Supervision.
350 824.06 779.73 5.38

Declaration of Competing Interest

Implication 3: Introducing split transportation is more significantly
when there are more demand categories. The authors declare that they have no known competing financial
Implication 4: When demand values of different customer sites are interests or personal relationships that could have appeared to influence
the work reported in this paper.

14
L. Zhang et al. Computers & Industrial Engineering 178 (2023) 109143

Data availability Irnich, S., & Desaulniers, G. (2005). In Shortest path problems with resource constraints (pp.
33–65). New York): Column Generation (Springer.
Irnich, S., & Villeneuve, D. (2006). The shortest-path problem with resource constraints
Data will be made available on request. and k-cycle elimination for k ≥ 3. INFORMS Journal on Computing, 18(3), 391–406.
Jepsen, M., Petersen, B., Spoorendonk, S., & Pisinger, D. (2008). Subset-row inequalities
Acknowledgment applied to the vehicle routing problem with time windows. Operations Research, 56
(2), 497–511.
Kim, B. I., Kim, S., & Sahoo, S. (2006). Waste collection vehicle routing problem with
This work was supported by National Natural Science Foundation of time windows. Computers and Operations Research, 33(12), 3624–3642.
China (72201017) and the Fundamental Research Funds for the Central Kohl, N., Desrosiers, J., Madsen, O. B. G., Solomon, M. M., & Soumis, F. (1999). 2-Path
Cuts for the Vehicle Routing Problem with Time Windows. Transportation Science, 33
Universities. (1), 101–116.
Li, H., Jian, X., Chang, X., & Lu, Y. (2018). The generalized rollon-rolloff vehicle routing
References problem and savings-based algorithm. Transportation Research Part B: Methodological,
113, 1–23.
Luo, Z., Qin, H., Zhu, W., & Lim, A. (2017). Branch and Price and Cut for the Split-
Al-Salem, S. M., Lettieri, P., & Baeyens, J. (2009). Recycling and recovery routes of
Delivery Vehicle Routing Problem with Time Windows and Linear Weight-Related
plastic solid waste (PSW): A review. Waste Management, 29(10), 2625–2643.
Cost. Transportation Science, 51(2), 668–687.
Amor, H. B., & Desrosiers, J. (2006). A proximal trust-region algorithm for column
Markov, I., Bierlaire, M., Cordeau, J.-F., Maknoon, Y., & Varone, S. (2020). Waste
generation stabilization. Computers & Operations Research, 33(4), 910–927.
Collection Inventory Routing with Non-stationary Stochastic Demands. Computers &
Archetti, C., Bianchessi, N., & Speranza, M. G. (2015). A branch-price-and-cut algorithm
Operations Research, 113, Article 104798.
for the commodity constrained split delivery vehicle routing problem. Computers &
Markov, I., Varone, S., & Bierlaire, M. (2016). Integrating a heterogeneous fixed fleet and
Operations Research, 64, 1–10.
a flexible assignment of destination depots in the waste collection VRP with
Archetti, C., Campbell, A. M., & Speranza, M. G. (2016). Multicommodity vs.
intermediate facilities. Transportation Research Part B: Methodological, 84, 256–273.
single–commodity routing. Transportation Science, 50(2), 461–472.
Miranda, P. A., Blazquez, C. A., Vergara, R., & Weitzler, S. (2015). A novel methodology
Archetti, C., & Speranza, M. G. (2004). Vehicle routing in the 1-skip collection problem.
for designing a household waste collection system for insular zones. Transportation
Journal of the Operational Research Society, 55, 717–727.
Research Part E: Logistics and Transportation Review, 77, 227–247.
Archetti, C., & Speranza, M. G. (2012). Vehicle routing problems with split deliveries.
Mordor Intelligence, 2022. Global waste management market - growth, trends, COVID-
International Transactions in Operational Research, 19(1–2), 3–22.
19 impact, and forecasts (2023 - 2028). https://www.mordorintelligence.com/
Baldacci, R., Bodin, L., & Mingozzi, A. (2006). The multiple disposal facilities and
industry-reports/global-waste-management-market.
multiple inventory locations rollon–rolloff vehicle routing problem. Computers &
Muter, I., Cordeau, J.-F., & Laporte, G. (2014). A branch-and-price algorithm for the
Operations Research, 33(9), 2667–2702.
multidepot vehicle routing problem with interdepot routes. Transportation Science,
Beliën, J., De Boeck, L., & Van Ackere, J. (2014). Municipal Solid Waste Collection and
48(3), 425–441.
Management Problems: A Literature Review. Transportation Science, 48(1), 78–102.
Nuortio, T., Kytojoki, J., Niska, H., & Braysy, O. (2006). Improved route planning and
Benjamin, A. M., & Beasley, J. E. (2010). Metaheuristics for the waste collection vehicle
scheduling of waste collection and transport. Expert Systems with Applications, 30(2),
routing problem with time windows, driver rest period and multiple disposal
223–232.
facilities. Computers & Operations Research, 37(12), 2270–2280.
Pecin, D., Pessoa, A., Poggi, M., & Uchoa, E. (2016). Improved branch-cut-and-price for
Buhrkal, K., Larsen, A., & Ropke, S. (2012). The Waste Collection Vehicle Routing
capacitated vehicle routing. Mathematical Programming Computation, 9(1), 61–100.
Problem with Time Windows in a City Logistics Context. Procedia - Social and
Ramos, T. R. P., Gomes, M. I., & Barbosa-Póvoa, A. P. (2013). Planning waste cooking oil
Behavioral Sciences, 39, 241–254.
collection systems. Waste Management, 33(8), 1691–1703.
Cárdenas-Barrón, L. E., González-Velarde, J. L., Treviño-Garza, G., & Garza-Nuñez, D.
Ramos, T. R. P., Gomes, M. I., & Barbosa-Póvoa, A. P. (2014). Economic and
(2019). Heuristic algorithm based on reduce and optimize approach for a selective
environmental concerns in planning recyclable waste collection systems.
and periodic inventory routing problem in a waste vegetable oil collection
Transportation Research Part E: Logistics and Transportation Review, 62(2), 34–54.
environment. International Journal of Production Economics, 211, 44–59.
Righini, G., & Salani, M. (2006). Symmetry helps: Bounded bi-directional dynamic
Dantzig, G. B., & Wolfe, P. (1960). Decomposition principle for linear programs.
programming for the elementary shortest path problem with resource constraints.
Operations Research, 8(1), 101–111.
Discrete Optimization, 3(3), 255–273.
De Bruecker, P., Beliën, J., De Boeck, L., De Jaeger, S., & Demeulemeester, E. (2018).
Righini, G., & Salani, M. (2008). New dynamic programming algorithms for the resource
A model enhancement approach for optimizing the integrated shift scheduling and
constrained elementary shortest path problem. Networks, 51(3), 155–170.
vehicle routing problem in waste collection. European Journal of Operational
Sahoo, S., Kim, S., & Kim, B. I. (2005). Routing optimization for waste management.
Research, 266(1), 278–290.
Interfaces, 35(1), 24–36.
Desaulniers, G. (2010). Branch-and-price-and-cut for the split-delivery vehicle routing
Salani, M., & Vacca, I. (2011). Branch and price for the vehicle routing problem with
problem with time windows. Operations Research, 58(1), 179–192.
discrete split deliveries and time windows. European Journal of Operational Research,
Dror, M., Laporte, G., & Trudeau, P. (1994). Vehicle routing with split deliveries. Discrete
213(3), 470–477.
Applied Mathematics, 50(3), 239–254.
Solomon, M. (1987). Algorithms for the vehicle routing and scheduling problems with
Dror, M., & Trudeau, P. (1990). Split delivery routing. Naval Research Logistics, 37(3),
time window constraints. Operations Research, 35(2), 254–265.
383–402.
Statista, 2023a. Global waste generation - statistics & facts. https://www.statista.com/
Expósito-Márquez, A., Expósito-Izquierdo, C., Brito-Santana, J., & Moreno-Pérez, J. A.
topics/4983/waste-generation-worldwide/#topicOverview.
(2019). Greedy Randomized Adaptive Search Procedure to Design Waste Collection
Statista, 2023b. Solid waste in China - statistics & facts. https://www.statista.com/
Routes in La Palma. Computers & Industrial Engineering, 137, Article 106047.
topics/5655/solid-waste-in-china/.
Faccio, M., Persona, A., & Zanin, G. (2011). Waste collection multi objective model with
Teixeira, J., Antunes, A. P., & de Sousa, J. P. (2004). Recyclable waste collection
real time traceability data. Waste Management, 31(12), 2391–2405.
planning––a case study. European Journal of Operational Research, 158(3), 543–554.
Feillet, D., Dejax, P., Gendreau, M., & Gueguen, C. (2004). An exact algorithm for the
Tirkolaee, E. B., Mahdavi, I., Esfahani, M. S., & M.,. (2018). A robust periodic capacitated
elementary shortest path problem with resource constraints: Application to some
arc routing problem for urban waste collection considering drivers and crew’s
vehicle routing problems. Networks, 44(3), 216–229.
working time. Waste Management, 76, 138–146.
Ghaderi, A., & Burdett, R. L. (2019). An integrated location and routing approach for
Wy, J., Kim, B. I., & Kim, S. (2013). The rollon–rolloff waste collection vehicle routing
transporting hazardous materials in a bi-modal transportation network.
problem with time windows. European Journal of Operational Research, 224(3),
Transportation Research Part E: Logistics and Transportation Review, 127, 49–65.
466–476.
Gschwind, T., Bianchessi, N., & Irnich, S. (2019). Stabilized Branch-Price-and-Cut for the
Yu, M., Jin, X., Zhang, Z., Qin, H., & Lai, Q. (2019). The split-delivery mixed capacitated
Commodity-constrained Split Delivery Vehicle Routing Problem. European Journal of
arc-routing problem: Applications and a forest-based tabu search approach.
Operational Research, 278(1), 91–104.
Transportation Research Part E: Logistics and Transportation Review, 132, 141–162.
Hauge, K., Larsen, J., Lusby, R. M., & Krapper, E. (2014). A hybrid column generation
Zhang, Y., Sun, L., Hu, X., & Zhao, C. (2019). Order consolidation for the last-mile split
approach for an industrial waste collection routing problem. Computers & Industrial
delivery in online retailing. Transportation Research Part E: Logistics and
Engineering, 71, 10–20.
Transportation Review, 122, 309–327.
Hernandez, F., Feillet, D., Giroudeau, R., & Naud, O. (2016). Branch-and-price
algorithms for the solution of the multi-trip vehicle routing problem with time
windows. European Journal of Operational Research, 249(2), 551–559.

15