A novel model for sustainable waste collection arc routing problem: Pareto-based algorithms

Abstract

Municipal solid waste (MSW) management is known as one of the most crucial activities in municipalities that requires large amounts of fixed/variable and investment costs. The operational processes of collection, transportation and disposal include the major part of these costs. On the other hand, greenhouse gas (GHG) emission as environmental aspect and citizenship satisfaction as social aspect are also of particular importance, which are inevitable requirements for MSW management. This study tries to develop a novel mixed-integer linear programming (MILP) model to formulate the sustainable periodic capacitated arc routing problem (PCARP) for MSW management. The objectives are to simultaneously minimize the total cost, total environmental emission, maximize citizenship satisfaction and minimize the workload deviation. To treat the problem efficiently, a hybrid multi-objective optimization algorithm, namely, MOSA-MOIWOA is designed based on multi-objective simulated annealing algorithm (MOSA) and multi-objective invasive weed optimization algorithm (MOIWOA). To increase the algorithm performance, the Taguchi design technique is employed to set the parameters optimally. The validation of the proposed methodology is evaluated using several problem instances in the literature. Finally, the obtained results reveal the high efficiency of the suggested model and algorithm to solve the problem.

1 Introduction

Sustainable waste management (SWM) mainly tries to address three pillars of economy, environment and social in operational (collection, transportation and disposal) activities and minimize the amount of municipal solid waste (MSW) that is generated and disposed. According to reports from the United Nations, about 11.2 billion tonnes of MSW is generated and collected throughout the world such that almost all of which originates from humans (Bolingbroke et al., 2021). Hence, managers should concurrently manage this waste and come up with efficient strategies to tackle this issue sustainably. SWM consists of any process involved in the organization of MSW, from generation to final disposal. At a higher level, the main goals of SWM systems are to minimize the amounts of natural resources consumed and maximize reusing of the materials extracted from nature. Operational activities include the most critical elements of SWM systems in which economic, environmental and social aspects should be scrutinized (Seadon, 2010; Xu & Yeh, 2017).

Recently, a few researchers have been trying to address to sustainable MSW collection through vehicle routing problems (VRPs) (Ghannadpour et al., 2021; Inghels et al., 2016; Lin & Kernighan, 1973; Lou et al., 2020; Paletta & Triki, 2004), location-routing problems (LRPs) (Tirkolaee et al., 2021) and inventory-routing problems (IRPs) (Cárdenas-Barrón et al., 2019) within a given network. On the other hand, the application of arc routing problems (ARPs) in SWM was neglected even though it is a well-known type of routing problems that has been specifically studied to formulate the vehicle routing operations in waste collection, snow removal of streets, street painting and other street-related services (de Armas et al., 2019; Foulds et al., 2015). This type of routing was firstly introduced by Golden and Wang (1981), which is a specific area of VRP. Furthermore, ARP is known as an extension of the Chinese Postman Problem (CPP) which can be found in many studies such as Pearn (1994), Filipiak et al. (2009) and Çodur and Yılmaz (2020).

The most important research works on ARP and in the field of municipal solid waste management are summarized as follows. Sniezek and Bodin (2006) employed two mixed integer programming (MIP) techniques to formulate the capacitated arc routing problem (CARP) for residential sanitation collection. They incorporated the vehicle/site dependencies into the model and then investigated several examples. Bautista et al. (2008) modified the ARP by applying changes to the node routing to solve the waste collection problem in the municipality of a town of Barcelona. These changes were due to road constraints such as no U-turn to return from the edge. They solved the problem using the ACS algorithm, which was based on the nearest-neighborhood and closest-placement methods. Rodrigues and Soeiro Ferreira (2015) developed a mixed CARP for the collection and transportation of municipal waste. They considered a heterogeneous fleet of vehicles and multiple disposal sites in the proposed problem and solved the problem by CPLEX solver using benchmark problems. Tirkolaee et al. (2016) applied an SA algorithm to solve the robust CARP under fuzzy demand for urban waste collection. To improve the performance of the proposed algorithm, they employed the Taguchi design method to tune the algorithm parameters and evaluated the performance of the algorithm compared to the CPLEX solver. A comprehensive review of ARP was conducted by Mourão and Pinto (2017) which considers the most important researches from 2010 to 2017. Tirkoaee et al. (2018a) developed a robust CARP to solve the urban waste collection problem by considering driver and crew shifts. They designed random examples and solved them by an exact method and SA algorithm in order to validate the applicability of the proposed mathematical model. A hybrid genetic algorithm (GA) was developed by Tirkolaee et al. (2018b) to solve the multi-trip capacitated green ARP to collect MSW. To test the applicability of the proposed algorithm, they generated several random examples in different sizes and evaluated the results compared to the CPLEX solver. Tirkoaee et al. (2019a) offered a MILP formulation for the problem of multi-trip CARP for urban waste collection. To solve the problem, they developed an improved max-mini ant system (MMAS) algorithm and evaluated the performance of the algorithm in comparison with the heuristic algorithms and the simplified version of the max–min ant colony optimization algorithm by benchmark samples. Recently, several studies have been done in the area of uncertainty. A robust bi-objective MILP model was introduced by Tirkoaee et al. (2019b) for multi-period multi-trip CARP under fuzzy demand. The model aimed to minimize the total cost and the longest travelled distance of vehicles simultaneously. They implemented multi-objective invasive weed optimization (MOIWO) algorithm and ε-constraint method to treat the problem. Tirkoaee et al. (2020) developed an improved hybrid ant colony optimization (ACO) algorithm to solve the multi-trip CARP under fuzzy demand for urban waste collection. They applied their proposed model on a real case study problem in Iran. Keenan et al. (2021) offered a strategic oscillation simheuristic to solve a time-capacitated ARP under stochastic demand which is applicable to waste collection problem. They tested the efficiency of the algorithm using a set of numerical experiments and against existing algorithms. Recently, a comprehensive review was conducted by De Maio et al. (2021) to survey the most important studies in the literature addressing different applications of ARP under uncertainty, particularly waste collection problem.

To the best of our knowledge and according to the above survey, there is no study in the literature addressing sustainable waste collection ARP. Furthermore, there are only a few number of studies that incorporated the sustainable development concept in their waste collection problems. In this study, providing a sustainable waste collection and transportation system requires the study of conflicts and trade-offs between economic, environmental and social objectives subject to the main operating constraints. To this end, total cost minimization, total pollution minimization, total job opportunity maximization and workload deviation minimization are taken into account to address these objectives, respectively. An efficient methodology including mathematical model and solution techniques is designed treat the problem. Accordingly, a tri-objective novel mixed-integer linear programming (MILP) model is developed to formulate a sustainable periodic capacitated arc routing problem (PCARP). To tackle the complexity of the problem efficiently, a hybrid multi-objective optimization algorithm, namely, MOSA-MOIWOA is designed and then tested against the ε-constraint method using different well-known indexes. Moreover, the Taguchi design method is employed to provide a higher efficiency for the proposed algorithm.

The remaining sections of the paper are structured as follows. Section 2 explains the problem and the main assumptions of the model. Moreover, the proposed mathematical model is given in this section. Section 3 introduces the proposed solution techniques. The results and comparative analyses are given in Sect. 4. Finally, Sect. 5 concludes the study and suggest a useful outlook for future studies.

2 Problem description

This section introduces the sustainable PCARP and its specific characteristics. The operational objective is to find the optimal number of required waste-collection vehicles among a fleet of homogeneous vehicles as well as constructing the optimal planning of routes. The sustainability is addressed by studying its three main aspects (economic, environmental and social) as separate objective functions in the proposed model. The first objective function seeks to minimize the total cost including the traversing costs and usage costs of vehicles. Total pollution emission is taken into account as the second objective function to be minimized. The third objective function is to maximize job opportunities through increasing the number of hired labor. Finally, total workload deviation of drivers and crew is minimized as the fourth objective function. Accordingly, the sustainability aspect of the problem can be studied.

The proposed network can be displayed by a graph G = (V, E), where V = {1, 2… n} represents the set of nodes and E represents the set of edges. Each pair of nodes constitutes an edge of the network and each edge contains two arcs in opposite directions. Here, each edge in the set E may be a waste edge (required edge) or non-waste edge (traversing edges). Thus, ${E}_{R}\subseteq E$ displays the set of required edges so that ${E}_{R}=\left\{\left(i,j\right)|{d}_{ijt}>0\right\}$, where ${d}_{ijt}$ represents the demand of edge $(i,j)$ in tth period. So, each non-waste edge has a demand equals to 0. Moreover, node number 1 stands for the depot and node number n shows the disposal site.

Figure 1 provides a schematic view to the problem.

Fig. 1

Schematic view of the proposed problem

As can be seen, the suggested MSW collection network in Fig. 1 is composed of seven required edges (1–2, 2–3, 3–4, 4–8, 5–6, 6–7, 7–8), one depot node and one disposal site. As a possible solution for a single planning period, all required edge are covered exactly by one vehicle and two trips. The required and traversing edges are specified by trash bins and dashed lines, respectively. The optimal routing for the first and second trips are 1–2 → 2–3 → 3–4 → 4–8 and 1–8 → 8–7 → 7–6 → 6–5 → 5–8, respectively. Finally, the vehicle moves back to the depot for completing its tour.

The main assumptions of the model are given in the following:

• I. Each required edge is covered only by one vehicle.

• II. Separate locations are considered for the depot and disposal site in the network.

• III. Each vehicle begins its first trip from the depot and end it at the disposal site. Then start its possible 2nd, 3rd,…, pth trips from the disposal site and end at there again.

• IV. A fleet of homogeneous vehicles is taken into account.

• V. Vehicles have a maximum available service time.

• VI. Usage cost of vehicles include drivers’ and crew’s wage cost, fuel cost, hiring cost, etc.

• VII. A set of planning periods is regarded to serve the required edges.

• VIII. Environmental pollution [greenhouse gas (GHG) emission] is considered for the transportation of waste in urban districts.

• IX. Social aspect is studied by maximizing the number of hired labor and minimizing workload deviation.

Now, the main components of the suggested mathematical model including sets and indices, parameters and variables are described in the following.

Sets and Indices

$V$ :

Set of nodes; $i,j\in V$

$E$ :

Set of edges; $(i,j)\in E$

${E}_{R}$ :

Set of required edges; ${E}_{R}=\left\{\left(i,j\right)|{d}_{ijt}>0\right\},$

$K$ :

Set of homogeneous vehicles; $k\in K$

P :

Set of vehicle trips; $p\in P$

T :

Set of planning periods; $t\in T$

$S$ :

Each possible subset of edges

$V[S]$ :

Set of nodes constituting S

Parameters

${c}_{ij}$ :

Distance of edge $(i,j)$

$W$ :

Available capacity for each vehicle

${d}_{ijt}$ :

Demand of edge $(i,j)$ in tth period

${T}_{max}$ :

Maximum available time for vehicles

$M$ :

A large number,

${t}_{ij}$ :

Traversing time of edge $(i,j)$

${cv}_{k}$ :

Usage cost of $k$th vehicle

$\theta$ :

Conversion factor of distance to cost

${G}_{ij}$ :

Amount of pollution emission released by traversing edge $(i,j)$

${\sigma }$ :

Number of required workforce (drivers and crew) for each vehicle

Variables

${x}_{ijkt}^{p}$ :

Number of traversing the edge $(i,j)\in E$ by $k$th vehicle in $p$th trip and tth period

${y}_{ijkt}^{p}$ :

1 if edge $(i,j)\in {E}_{R}$ is served by $k$th vehicle in $p$th trip and tth period, otherwise 0

${u}_{kt}$ :

1 if $k$th vehicle is employed in tth period, otherwise 0

${LT}_{k}^{p}$ :

Total loading time of $k$th vehicle in $p$th trip and tth period

${UT}_{k}^{p}$ :

Total unloading time of $k$th vehicle in $p$th trip and tth period

${WT}_{kt}$ :

Total service time of $k$th vehicle in tth period;

$${WT}_{kt}=\sum_{p\in P}{LT}_{kt}^{p}+\sum_{p\in P}{UT}_{kt}^{p}+\sum_{p\in P}\sum_{(i,j)\in {E}}{t}_{ij} {x}_{ijkt}^{p}\quad \forall k\in K,\forall t\in T,$$

Now, the developed MILP formulation of the problem is as follows:

$$\mathrm{minimize}\;{Z}_{1}=\theta \left(\sum_{\left(i,j\right)\in E}\sum_{p\in P}\sum_{t\in T}\sum_{k\in K}{c}_{ij }{x}_{ijkt}^{p}\right)+\sum_{t\in T}\sum_{k\in K}{{cv}_{k}}{u}_{kt}$$

(1)

$$ {\text{minimize}}\;Z_{2} = \mathop \sum \limits_{{\left( {i,j} \right) \in E}} \mathop \sum \limits_{p \in P}^{ } \mathop \sum \limits_{t \in T} \mathop \sum \limits_{k \in K}^{ } G_{ij } x_{ijkt}^{p} $$

(2)

$$ {\text{maximize}}\; Z_{3} = \mathop \sum \limits_{t \in T} \mathop \sum \limits_{k \in K} \sigma u_{kt} $$

(3)

$$\mathrm{minimize}\;{Z}_{4}=\sum_{t\in T}\sum_{k\in K}\frac{{T}_{max}-{WT}_{kt}}{{T}_{max}}$$

(4)

subject to

$$\sum_{i\in V[E]}{x}_{ijkt}^{p}=\sum_{j\in V[E]}{x}_{jikt}^{p}\quad \forall i\in V\left[E\right];\left(i,j\right)\in E, \forall k\in K,\forall p\in P,\forall t\in T,$$

(5)

$$\sum_{p\in P}\sum_{k\in K}({y}_{ijkt}^{p}+{y}_{jikt}^{p}) =1\quad \forall \left(i, j\right) or \left(j,i\right)\in {E}_{R},\forall t\in T,$$

(6)

$$\sum_{\left(i,j\right)\in {E}_{R}}{d}_{ijt} {y}_{ijkt}^{p}\le W\quad \forall k\in K,\forall p\in P,\forall t\in T,$$

(7)

$${y}_{ijkt}^{p}\le {x}_{ijkt}^{p} \forall \left(i,j\right)\in E ,\quad \forall k\in K,\forall p\in P,\forall t\in T,$$

(8)

$$\sum_{p\in P}\sum_{(i,j)\in {E}}{x}_{ijk}^{p}\le M {u}_{kt}\quad \forall k\in K,\forall t\in T,$$

(9)

$${LT}_{kt}^{p}=ul \sum_{\left(i,j\right)\in {E}_{R}}{d}_{ij} {y}_{ijkt}^{p}\quad \forall k\in K,\forall p\in P,\forall t\in T,$$

(10)

$${UT}_{kt}^{p}=uu \sum_{\left(i,j\right)\in {E}_{R}}{d}_{ij} {y}_{ijk}^{p}\quad \forall k\in K,\forall p\in P,\forall t\in T,$$

(11)

$$\sum_{p\in P}{LT}_{kt}^{p}+\sum_{p\in P}{UT}_{kt}^{p}+\sum_{p\in P}\sum_{(i,j)\in {E}}{t}_{ij} {x}_{ijkt}^{p}\le {T}_{max }\quad \forall k\in K,\forall t\in T,$$

(12)

$$\sum_{\left(j,h\right)\in S}{x}_{jhkt}^{p}\le M\sum_{i\notin V\left[S\right], j\in V\left[S\right]\backslash \{1, n\}}{x}_{ijkt}^{p}\quad \forall S\subseteq E,\forall k\in K,\forall p\in P,\forall t\in T,$$

(13)

$$\sum_{j\in V[E]}{x}_{1jkt}^{1} \ge \sum_{j\in V[E]}{x}_{njkt}^{2}\quad \forall k\in K,\forall t\in T,$$

(14)

$$\sum_{j\in V[E]}{x}_{njkt}^{p} \ge \sum_{j\in V[E]}{x}_{njkt}^{p+1}\quad \forall p\in \left\{2,3,\dots ,|P|-1\right\},\forall k\in K,\forall t\in T,$$

(15)

$$\sum_{\begin{array}{c}\left(1,j\right)\in E\\ j\in V[E]\backslash \{1,n\}\end{array}}{x}_{1jkt}^{p} ={u}_{kt}\quad \forall k\in K,\forall p=1,\forall t\in T,$$

(16)

$$\sum_{\begin{array}{c}\left(j,n\right)\in E\\ j\in V[E]\backslash \{1,n\}\end{array}}{x}_{jnkt}^{p} ={u}_{kt}\quad \forall k\in K,\forall p=1,\forall t\in T,$$

(17)

$$\sum_{\begin{array}{c}\left(n,j\right)\in E\\ j\in V[E]\backslash \{1,n\}\end{array}}{x}_{njkt}^{p} \le {u}_{kt}\quad \forall k\in K,\forall p\in P\backslash \{1\},\forall t\in T,$$

(18)

$$\sum_{\begin{array}{c}\left(j,n\right)\in E\\ j\in V[E]\backslash \{1,n\}\end{array}}{x}_{jnkt}^{p} \le {u}_{kt}\quad \forall k\in K,\forall p\in P\backslash \{1\},\forall t\in T,$$

(19)

$${x}_{ijkt}^{p}\in {Z}^{+}, {y}_{ijkt}^{p}\in \left\{0,1\right\},{u}_{kt}\in \left\{\mathrm{0,1}\right\}\quad \forall \left(i,j\right)\in E , \forall k\in K,\forall p\in P,\forall t\in T,$$

(20)

$${LT}_{kt}^{p}, {UT}_{kt}^{p}\ge 0\quad \forall \left(i,j\right)\in E ,\forall k\in K,\forall p\in P, \forall t\in T.$$

(21)

Objective functions (1), (2), and (3) and (4) stand for the economic, environmental and social aspects, respectively. Objective function (1) minimizes the total waste collection routing cost which is composed of the traversing costs and usage costs of vehicles. Objective function (2) minimizes the total pollution emission released by vehicles. Objective function (3) maximizes the number of hired labor in all periods. Constraint (5) indicates the flow balance in the network for vehicles. Objective function (4) minimizes the workload deviation by balancing the total service time of vehicles in all time periods. In other words, this objective tries to maximize the satisfaction level of workforce. Constraint (6) ensures that required edges are served only by one vehicle. Constraint (7) shows the capacity limitation of vehicles. Constraint (8) states that the required edge can be served by the vehicle traversing it; i.e., vehicles may just traverse an edge without serving it. Constraint (9) indicates that vehicles can be employed if only their usage costs are paid. Constraints (10) and (11) calculate the sum of loading and unloading times for vehicles in each period and trip, respectively. Constraint (12) expresses the maximum available service time for vehicles in each period. Constraint (13) prevents from the potential sub-tours for vehicles in each period and trip. Constraints (14) and (15) form the sequence of vehicles’ trips from 1 to P. Constraints (16) and (17) ensure that the beginning of first trips is the depot and the disposal site is the ending, respectively. Constraints (18) and (19) indicate that if more than one trip is required, the next trips (second, third, fourth, etc.) begin from the disposal site and end at there again. Constraints (20) and (21) display the domain of the variables.

3 Methodology

This section provides the proposed solution techniques to validate, solve and analyze the proposed problem. Accordingly, MOIWOA-MOSA is designed as a hybrid multi-objective optimization algorithm based on MOIWOA and MOSA. Moreover, the ε-constraint method is implemented to solve the model exactly using CPLEX solver/GAMS software.

3.1 Solution representation

In order to design any meta-heuristic algorithm, it is necessary to show how the solutions are represented via a set of vectors or matrices. In this research, a two-part vector is considered for routing of vehicles in each period. The first part includes a permutation of all demand points. In other words, if the number of edges is equal to n, a random order from 2 to n − 1 is created. The second part shows the selected vehicles for each created tour. If a vehicle is used several times, it means that the vehicle has multiple trips.

To interpret this two-part vector, the first vehicle is selected (according to the second part) and based on its capacity, it visits a number of edges (according to the first part) and moves to its final destination. Then the second vehicle is selected from second part and the same process occurs until all the edges are visited. An example with 8 edges and 2 vehicles are depicted in Fig. 2.

Fig. 2

An example of the solution representation

3.2 Heuristic for generating initial solutions

In this section, a heuristic algorithm is developed to yield the initial solutions. The main implementation steps are given as follows:

• Step 1 At the starting point of the algorithm, choose the first planning period. Otherwise, if there is any planning period remaining, choose the next planning period among the existing periods and go to Step 2. Otherwise, go to Step 6.

• Step 2 Select a vehicle randomly among the existing vehicles to begin the first trip from the depot and go to Step 3.

• Step 3 Among the existing required edges, take into account all the required edges that can be added to the trip based on the capacity and available service time limitations of the vehicle, then select one edge with the minimum distance from the depot and go to Step 5. Differently, if there is no eligible required edge fulfilling one of these limitations, move to the disposal site.

• Step 4 If the maximum available service time allows the vehicle to build another trip, begin the next trip from the disposal site and go to Step 3, otherwise, go to Step 5.

• Step 5 If there is a remaining required edge to be served, go to Step 3. Otherwise, select the shortest path to the disposal site and then go to Step 1.

• Step 6 Move to the depot and complete the tour.

• Step 7 Report the obtained solution.

3.3 MOSA

Simulated annealing (SA) is a local search algorithm with a great sufficiency to prevent from being trapped in local optima. It is very effective for solving non-convex or discrete problems. Therefore, SA is employed to solve integer programming problems efficiently (Glover et al., 2006). From its introduction to the present, SA has represented a high efficiency in solving large combinatorial optimization problems (Kubotani & Yoshimura, 2003). Furthermore, the simplicity of implementation, convergence and hill-climbing for eluding local optima are regarded as the main superiority criteria.

Accordingly, MOSA generates non-dominated solutions through a simple probability function that tries to provide the solutions constituting the Pareto optimal front. Based on this probability function, the total space of objective is covered uniformly generating as a wide range of possible well-dispersed and non-dominated solutions. These characteristics have made MOSA a fast reliable algorithm in comparison with the other existing multi-objective algorithms with a variety of applications. The mechanism of the suggested MOSA is adapted from Nam and Park (2000). The superiority of MOSA against other EAs is that there is no need to employ high memory for keeping the population information. The pseudo-code of the suggested MOSA algorithm is given in Fig. 3.

Fig. 3

Pseudo-code of the suggested MOSA (Nam & Park, 2000)

In order to create a neighbor solution in the MOSA, first, one of the solution parts is selected randomly. Next, the current values are replaced with new random values. If the first part is selected, a new random permutation is performed. If the second part is selected, one of the vehicles is randomly selected and replaced in each cell. An example of these changes is shown in Fig. 4.

Fig. 4

An example of creating a neighbor in the proposed MOSA

It should be noted that MOSA is implemented to each initial solution generated by the heuristic algorithm in Sect. 3.2. Finally, we can provide a set of the high-quality initial solutions.

3.4 MOIWO

Invasive Weed Optimization algorithm (IWO) is an evolutionary and intelligent algorithm which was inspired by the processes of propagation, adaptability and survival of weeds. It was first introduced by Mehrabian and Lucas (2006) and then has been employed by different researchers to tackle the complexity of optimization problems. The implementation of IWO has a suitable performance for all test functions. To simulate the colonizing behavior of weeds, the main principal features of the process are regarded as follows:

• i. A certain number of seeds is spread over the search space.

• ii. Each seed grows to a flowering plant and generates childs (seeds) with respect to the fitness value.

• iii. The generated seeds are randomly dispersed over the search space and become new plants.

This process is maintained until the maximum number of plants is provided. Thereafter, only the plants with lower fitness values can survive and generate new seeds and the others are ruined. Accordingly, this process maintained until reaching the maximum number of iterations. Finally, the plant with the best fitness value is regarded as the nearest solution to the optimum.

MOIWOA is a multi-objective variant of the IWO algorithm that was first introduced by Nikoofard et al. (2012) and has been extensively employed by researchers (Goli et al., 2019a; 2019b).

The execution steps of MOIWOA are given as follows:

• Step 1 Generating initial seeds

  The initial solutions are called initial seeds. These initial seeds are the output of the MOSA (Sect. 3.3) applied to the solutions generated by the heuristic algorithm (Sect. 3.2).

• Step 2 Seed reproduction

  In this phase, a member of the population is permissible to generate seeds regarding its own and colony’s lowest and highest fitness values according to Eq. (22):

  $$S={S}_{min}+\left({S}_{max}-{S}_{min}\right)\frac{{f-f}_{worst}}{{f}_{best}-{f}_{worst}}.$$

  (22)

  where S_min and S_max denote the minimum and maximum number of seeds, respectively. Furthermore, ${f}_{worst}$ and ${f}_{best}$ represent the worst and best fitness value, respectively. After specifying the number of seeds, reproduction is performed by executing the following sub-steps:

• Sub-step 1 Swap the routing plan for a certain trip of the vehicle by another vehicle in a period.

• Sub-step 2 Select two various trips of two vehicles randomly. If there is a common edge/edges between these trips, choose one of them randomly, then divide these trips into two parts. The first part of the first trip is combined with the second part of the second trip and the other two parts of the trips are combined. Consequently, two new trips are generated.

• Sub-step 3 Select two trips of a vehicle randomly. If there are two common edges among them, exchange the sequences of the common edges.

• Sub-step 4 Select an edge in a trip of vehicle randomly and reverse its direction.

• Sub-step 5 Select a part of a trip randomly and reverse the direction of the included edges.

An example of seed reputation in the proposed MOIWOA is provided in Fig. 5.

Fig. 5

An example of seed reputation in MOIWOA

• Step 3: Competitive exclusion

  The maximum number of plants can be provided by fast reproduction after several iterations, nevertheless, it is expected that the desired plants are re-generated more than the undesired ones. Moreover, it is critical to reach the maximum number of plants (p_max) to efficiently handle the speed of MOIWOA (Goli et al., 2019a). In this stage, the solutions are sorted based on the non-dominate sorting technique that is described in the following. The weaker solutions are removed to keep p_max solutions for the next step.

• Step 4: Non-dominate sorting

  The following parameters are quantified for each solution in order to rank a set of solutions: n_p as the number of solutions and S_p as the set of solutions dominated by solution p.

According to this technique, the value of crowding distances is the basic factor for sorting the solutions, which is calculated by Eq. (23):

$${d}_{{I}_{j}^{m}}={d}_{{I}_{j}^{m}}+\frac{{f}_{m}^{({I}_{j+1}^{m})}-{f}_{m}^{({I}_{j-1}^{m})}}{{f}_{m}^{max}-{f}_{m}^{min}}\quad \forall m,$$

(23)

where ${f}_{m}^{max}$ and ${f}_{m}^{min}$ stand for the maximum and minimum value of mth objective. Here, ${f}_{m}^{({I}_{j+1}^{m})}$ and ${f}_{m}^{({I}_{j-1}^{m})}$ show the value of mth objective for the next and previous solution of jth solution according to the sorted value of mth objective function, respectively.

Finally, the pseudo-code of the suggested MOIWOA is given by Fig. 6.

Fig. 6

Pseudo-code of the suggested MOIWOA (Goli et al., 2019a)

3.5 MOIWOA-MOSA

The hybrid algorithm of this study works by generating high-quality initial solutions using MOSA to be incorporated into the MOIWOA. Finally, the best possible solution will be found by MOIWOA.

3.6 Taguchi design method

As was clarified in the previous sections, MOSA and MOIWOA have several adjustable parameters that can affect the performance and final results. Hence, providing a suitable combination of these parameters can highly enhance the performance and efficiency of the algorithms. Two approaches for performing the test are known (Taguchi et al., 2005): (i) Standard analysis of variance (ANOVA) and (ii) Signal to noise (S/N) ratio. The value of S/N indicates the amount of dispersion around a certain value, or how the obtained solutions have altered over different experiments. To reduce the dispersion of the objective functions as far as possible, S/N ratio plays an important role. This ratio stands for noise factors accompanying controllable parameters.

By investigating among different Taguchi tables, the relevant table to L₉ is taken into account, which is implemented using Minitab statistical software. The input data is represented in Table 1 and the results of implementing Taguchi design method for MOSA are displayed in Table 2. After applying the test, the optimal values of the parameters of MOSA are given in Table 3 and Fig. 7.

Table 1 Parameter of the MOSA algorithm

Table 2 Results of implementing Taguchi design method for MOSA algorithm

Table 3 Optimal value for the parameter of the MOSA algorithm

Fig. 7

Output results of the Taguchi design for the MOSA algorithm

Next, the Taguchi tests applied with MOIWO algorithm using L₂₇. The selected initial value for MOIWO parameters are shown in Table 4. The results of Taguchi method and the optimal values of MOIWO parameters are represented in Fig. 8 and Tables 5 and 6 and Fig. 8 respectively.

Table 4 Parameter of the MOIWOA

Fig. 8

Output results of the Taguchi design for the MOIWOA algorithm

Table 5 Results of implementing Taguchi design method for MOIWO algorithm

Table 6 Optimal value for the parameter of the MOIWOA

3.7 ε-constraint approach

The ε-constraint approach is known as an effective way to cope with multi-objective optimization problems in the literature, which can generate Pareto solutions. Here, the Pareto front can be achieved by the ε-constraint approach based on Model (24) (Bérubé et al., 2009):

$$ \begin{array}{l} {\text{minimize}}\; f_{{1}} \left( x \right) \hfill \\ {\text{subject to}} \hfill \\ x \in X, \hfill \\ f_{{2}} \left( x \right) \le \varepsilon_{{2}} , \hfill \\ \ldots \hfill \\ f_{n} \left( x \right) \le \varepsilon_{m} . \hfill \\ \end{array} $$

(24)

The execution steps of the ε-constraint approach are given in Fig. 9.

Fig. 9

Execution steps of the proposed ε-constraint

According to the offered ε-constraint method, the first objective function is set to be the main one, and the other objective functions as sub-objective functions. Now, the final formulation relates to the suggested model of the study is displayed in Model (25):

$$ \begin{array}{l} {\text{minimize}}\;f_{{1}} \left( x \right) \hfill \\ {\text{subject to}} \hfill \\ x \in X, \hfill \\ f_{{2}} \left( x \right) \le \varepsilon_{2} , \hfill \\ f_{{2}} \left( x \right) \ge \varepsilon_{3} , \hfill \\ f_{{4}} \left( x \right) \le \varepsilon_{4} . \hfill \\ \end{array} $$

(25)

4 Results

This section provides the model validation and comparisons between the proposed solution techniques through solving 10 problem instances in Tirkolaee et al. (2019b). The parameters took value using a uniform distribution. For example, the demand parameter follows a uniform distribution uniform (1, 4). Moreover, Table 7 shows the input information about the problem instances.

Table 7 Information of the problem instances

The obtained results by the proposed solution techniques are shown in Table 8. With respect to the multi-objective essence of the suggested model, MOSA-MOIWOA and the ε-constraint approach are compared using the mean of ideal distance (MID), number of solutions (NOS), rate of achievement to several objectives (RASO), distancing (D), and CPU time metrics. MID is employed to calculate the mean distance of the Pareto solutions from the ideal solution or the origin of the coordinates. Based on Eq. (25), the lower value of MID indicates more efficiency for the solution methods. Moreover, RAS and D are calculated based on Eqs. (26)–(28).

$$ MID = \frac{1}{NOS}\sum_{s = 1}^{N} {} \sqrt {\sum_{m = 1}^{4} {\,\,f_{s,m}^{2}} } , $$

(26)

$$ RASO = \frac{{\sum\limits_{i = 1}^{n} {\left( {\frac{{f_{i,1} - F_{i} }}{{F_{i} }}} \right)} + \sum\limits_{i = 1}^{n} {\left( {\frac{{f_{i,2} - F_{i} }}{{F_{i} }}} \right)} + \sum\limits_{i = 1}^{n} {\left( {\frac{{f_{i,3} - F_{i} }}{{F_{i} }}} \right) + \sum\limits_{i = 1}^{n} {\left( {\frac{{f_{i,4} - F_{i} }}{{F_{i} }}} \right)} } }}{NOS}, $$

(27)

$$ D = \sqrt {\sum\limits_{m = 1}^{4} {\left( {\mathop {\max}\limits_{s = 1:n} f_{s,m} - \mathop {\min}\limits_{s = 1:n} f_{s,m} } \right)}^{2} } . $$

(28)

where $f_{sol,m}$ represents the mth objective value for the sth solution, and $F_{i} = \min \{ f_{1i} ,\,f_{2i} ,\,f_{3i} ,\,f_{4i} \}$.

Table 8 Comparisons between MOSA-MOIWOA and ε-constraint

Moreover, a run time limitation of 10,000 s is regarded to evaluate the performance of the solution techniques.

Figures 10, 11 and 12 depict the comparisons for the metrics between solution techniques in different problems.

Fig. 10

Comparison of CPU time values

Fig. 11

Comparison of MID values

Fig. 12

Comparison of NOS values

As can be seen in Table 4, the ε-constraint approach couldn’t solve P10 within the run time limitation of 1000 s due to the complexity of the proposed problem in large sizes and the low efficiency of CPLEX solver to solve it. On the other hand, the suggested MOSA-MOIWOA solved the problem instances in less than 85 s. In fact, the comparison of the average CPU time values for these two solution techniques identifies that the MOSA-MOIWOA just spends about 1.2% of the CPU time required for the ε-constraint approach. Figure 5 represents this comparison trends for the CPU time values.

Moreover, by analyzing the values obtained for the two MID and NOS metrics, the quality of the final solutions can be evaluated. The lower MID and the more NOS reflect a more efficient solution technique. Figure 6 displays that these two solution techniques have a close behavior against each other in terms of MID value. Furthermore, the comparison of NOS values in Fig. 12 demonstrates that the MOSA-MOIWOA has a relatively higher efficiency in this index. In order to complete the comparisons, the trend of RASO and D are illustrated in Figs. 13 and 14.

Fig. 13

Comparison of RASO values

Fig. 14

Comparison of D values

As can be seen, the average RASO index in the EPC method and MOSA-MOIWO is 0.94 and 0.74, respectively. Accordingly, it is concluded that the MOSA-MOIWO algorithm has been able to achieve higher quality Pareto solutions. This is while the scattering of Pareto solutions is also of great importance in comparing multi-objective methods (Goli et al., 2020). Accordingly, index D determines the degree of scatter of Pareto solutions. As can be seen in Fig. 14, in all solved problems, the MOSA-MOIWO algorithm has a higher dispersion index than the EPC. Therefore, in general, it is concluded that the designed meta-heuristic algorithm, in addition to the high quality of Pareto solutions, has been able to provide complete coverage of the optimal boundary of the study optimization problem.

The numerical results demonstrate important findings from solving different samples. Due to the high CPU time of the EPC, it can be concluded that the problem is of high complexity. Therefore, MOSA-MOIWO is an efficient method that can solve large-scale problems in a reasonable time. It is also quite obvious that increasing the dimension of the problem leads to a larger value for each objective function. Therefore, it is observed that an upward trend is followed by MID index.

On the other hand, by increasing the dimension of the problem, the solution space also increases, and then, more Pareto solutions can be found. A closer look at RASO index in Fig. 13 and D index in Fig. 14 confirms that the diversity of Pareto solutions is directly related to the problem scale.

4.1 Sensitivity analysis

This section provides a sensitivity analysis on P10 as the largest-sized problem that is the closest one to the real-world scale. For this purpose, the parameter of the maximum available service time (${T}_{max})$ is studied under different change intervals and the behavior of the objective functions is assessed using MOSA-MOIWOA. The obtained results are shown in Table 9 and Figs. 15, 16, 17 and 18.

Table 9 Obtained results for the sensitivity analysis

Fig. 15

Sensitivity analysis of the 1st objective function

Fig. 16

Sensitivity analysis of the 2nd objective function

Fig. 17

Sensitivity analysis of the 3rd objective function

Fig. 18

Sensitivity analysis of the 4th objective function

As can be seen in Table 9 and Figs. 15, 16, 17 and 18, all objective functions reflect remarkable fluctuations against the changes of ${T}_{max}$, but in various behavioral directions. Accordingly, the 1st, 3rd and 4th objective functions decreased by the increase of ${T}_{max}$. This is due to the generated decrease in the usage costs of vehicles and this reduction leads to a significant decrease in the number of hired labor and also deviation of workload. However, the 2nd objective directly increased by the increase of ${T}_{max}$ as it was expected due to more usage of vehicles and less movement from the disposal site to the operational zone.

Table 10 represent different change percentages of the objective functions against the four change intervals considered for ${T}_{max}$. Accordingly, it is revealed that the most increase is related to the 2nd objective function; i.e., total pollution emission, where a 19.50% increase occurs for the 20% increase in ${T}_{max}$. On the other hand, the most decrease is again related to this objective function, where a 18.59% decrease occurs for the 20% decrease in ${T}_{max}$. Overall, increasing the parameter improves the 1st, 3rd and 4th objective functions but it leads to a worse value for the 2nd objective function.

Table 10 Change percentages of the objective functions against parameter’s change intervals

Based on the obtained results, managers can investigate these trade-offs to find the optimal level of the resources to be prepared and allocated in the waste collection system. In fact, sensitivity analysis is a useful tool to identify the optimal policy in a real-world situation.

The results of sensitivity analysis can provide useful managerial insights in the design of SWM systems. In the real world, the parameter T_max is set experimentally and based on historical data. Managers set a specific value for the parameter T_max according to vehicle fuel consumption, driver fatigue, etc. However, with respect to the analyses conducted in this research, this parameter can have several effects in terms of economic, environmental and social aspects to address the sustainable development. It is true that large values for this parameter reduce routing costs, but create unfavorable conditions in terms of pollution and job opportunities. Therefore, it is necessary to find and set the most appropriate value for T_max before designing such a MWM system with a complete and accurate perspective.

5 Conclusions and outlook

This study designed an effective methodology to determine the optimal plans for waste-collection routes and the essential number of vehicles. For this purpose, a novel MILP model was developed to formulate a multi-objective multi-trip sustainable PCARP. The main goals of the problem were to concurrently minimize the total cost, total pollution emission and total workload deviation and also maximize total job opportunity. To validate the proposed mathematical model and cope with its tri-objectiveness, the ε-constraint approach was utilized by CPLEX solver/GAMS software. Moreover, to treat the problem efficiently in large scales, a multi-objective meta-heuristic algorithm, namely, MOSA-MOIWOA was then developed. To enhance the efficiency of the algorithm, the Taguchi design method was also applied to set the parameters. The performance of the MOSA-MOIWOA was then evaluated using various problem instances and three metrics of MID, NOS, D, RASO and CPU time in comparison with the ε-constraint method. The obtained results proved that the MOSA-MOIWOA can provide high-quality solutions within a much shorter CPU time. Moreover, the MID values were close to each other but the proposed algorithm could yield much better NOS value. Eventually, a number of sensitivity analyses were implemented on the maximum available service time of vehicles to assess the behavior of the objective functions and offer useful managerial decision aids. It was obvious that the total cost and total job opportunity have indirect relations with this parameter, but the total pollution emission reflects a direct relation.

According to the main limitations of the study, following suggestions are raised for future research directions:

• (i) Considering multiple depots and disposal sites within the network.

• (ii) Taking into account the locational decisions for depots and disposal sites, accordingly more environmental and social aspects can be incorporated into the model, such as minimizing proximity to inhabited areas and maximizing work safety at disposal sites, respectively (Sagnak et al., 2021).

• (iii) Extending the MSW management network by adding more facilities such as incinerators and recycling/recovery sites.

• (iv) Developing the social sustainability by minimizing the number of garbage pick-up to reduce social discomfort of odor presence time or under pandemics outbreak (e.g., COVID-19 (Tirkolaee et al., 2021)).

• (v) Studying the effect of milk-run strategy in the problem.

• (vi) Applying uncertainty modeling such as fuzzy theory (Tirkolaee et al., 2021), grey systems (Roy et al., 2017), robust optimization (Tirkolaee et al., 2019b), stochastic optimal control (Yılmaz et al., 2020), regime switching (Savku & Weber, 2018).

• (vii) Employing other multi-objective algorithms like non-dominated sorting genetic algorithm III (NSGA-III) (Sangaiah et al., 2020).