A hybrid GA-BFO algorithm for the profit-maximizing capacitated vehicle routing problem under uncertain paradigm

Abstract

This study designs a new variant of the capacitated vehicle routing problem (CVRP) under a fuzzy environment. In CVRP, several vehicles start their journey from a central depot to provide services to different cities and finally return to the depot. This paper introduces an additional time beyond the service time at each city to fulfill the pre-ordered demands. The need for this excess service time is to provide the services to new customers who are not enlisted at the start of the process. It is a market enhancement step. The proposed model’s main objective is to find the maximum time-dependent profit by using the optimum number of vehicles in an appropriate route and spending optimum excess service time in each city. The model considers travel time and travel cost as fuzzy numbers. An expected value model (EVM) is formulated using the credibility approach on fuzzy variables. A hybrid meta-heuristic method combining a genetic algorithm (GA) and bacteria foraging optimization algorithm (BFOA) is designed to solve the proposed model. The proposed model is explained with the help of some numerical examples. Sensitivity analyses based on different independent parameters of the algorithms are also conducted.

Keywords

Capacitated vehicle routing problem profit maximization fuzzy credibility theory hybrid algorithm genetic algorithm bacteria foraging optimization algorithm

1 Introduction

In the competitive world of business, one of the most critical issues is to deliver the goods to the customers within a reasonable time and preferably as fast as possible once the order is placed. This issue is a significant operational process, and it has become a natural area of research in the field of supply chain and logistics. An efficient supply and distribution management can reduce a significant drop in the running cost for a company. In this regard, the vehicle routing problem (VRP) is a very well-known problem. In VRP, the location and demands are known for all the cities. In this problem, some vehicles start from a single fixed depot after serving all the customers in an appropriate route by maintaining the vehicle’s capacity constraints and returns to the starting point. It focuses on the minimization of the cost or distance traveled by those sets of vehicles.

This paper proposes a new variant in the domain of capacitated vehicle routing problems. It introduces an additional service time beyond its stay time to fulfill the previously ordered demands. In the real-world scenario, several times, while supplying the products ordered, customers may ask for more quantities of the same or different products to fulfill their recent excess demands. In this regard, we consider the above mentioned additional service time. The proposed model also considers that each city may have several customers. Some are enlisted in the system from the beginning of the process, and some are not. The basic need for this excess service time is to provide the services to new customers who are not enlisted at the start of the process. This is a step to increase the market of any product(s). However, there is no guarantee of fulfilling the excess demand as each vehicle has some fixed capacity, and for each vehicle, there is a maximum time limit. This excess service time in each city will return some extra earning and incur some expenditure. This paper considers a closed CVRP model having a single depot. Now the proposed model is illustrated below.

There is a depot or a central position from where all the vehicles start their journey in the capacitated vehicle routing problem. The proposed model is applicable in a salesman or medical representative problem where the depot denotes the medical house or warehouse locations. Every vehicle or salesperson spends some predetermined period in each city to deliver or to market their products. In this period, they are supposed to have some earnings or returns, as well as some expenses. So the revenues minus the costs to serve a particular city yields the profit in that city. This amount of profit in a city depends on the time a salesman spends there. All the cities must be served with base demands. This quantity is a deterministic term, and it is known earlier before the solution procedure is going to start. All the salesmen or vehicles begin at the same time, and each of them has different routes. All the routes are mutually exclusive, i.e., each city must be served by only one vehicle. All the vehicles return to the depot after the completion of their delivery. After serving the base demand of a city, the vehicle or salesman spends extra time to market or supply more products. This additional quantity of supply is called “excess demand.” The service time or unloading time for the base demand is known earlier, but the excess time is not concretely predictable. The total need, consisting of the base demands and extra demands for all the cities for a particular route, must be less than the vehicle’s capacity. There is an upper limit on the time before which all the conveyances have to return to the depot. The traveling time or cost from one city to another is equal for all the vehicles. Also, all the carriages have the same capacity. A specific part of each conveyance’s capacity is filled up for the base demand, and the remaining portion is used for the excess demand.

The solution methodology of any vehicle routing problem can be categorized in two different ways. One is the exact method, and the other is the meta-heuristic method. The exact methods are suitable for small-scale problems, but they used to take a long time to find the solution. Whereas the meta-heuristic processes can deal with large-scale problems that may not find the best solution, it may lead to a very close solution within very little time compared to the exact method. Here this paper uses a hybrid heuristic method combining the strength of genetic algorithm (GA) and bacteria foraging algorithm (BFO).

The significant contribution of this paper is listed below.

A new variant of the capacitated vehicle routing problem is presented here.

The system can work fine in an uncertain environment.

A fuzzy model is also constructed for uncertain travel time and cost.

A hybrid metaheuristic algorithm based on GA and BFO algorithm is designed to solve the proposed model.

Sensitivity analyses are performed for all the independent parameters in the model.

The paper is organized using several sections. The introduction of the work is presented in the first section. The literature review is provided in Section 2. Section 3 deals with the motivation of this work. Some preliminaries required to understand the model are described in section 4. Section 5 presents the mathematical model formulation. The proposed GA-BFO based hybrid algorithm is described in Section 6. Section 7 consists of some numerical illustrations of different models. Section 8 deals with the result analysis. Finally, Section 9 concludes the work.

2 Literature review

The basic model of VRP is also called the capacitated vehicle routing problem (CVRP) [1]. There are several variations of CVRP problems. Based on the dynamicity of introducing new customers and their demand, CVRP problems can be classified into two groups. One of these is the static VRP [2], where the customers, as well as their needs, are fixed. Also, the demands are known before the vehicles start. The dynamic VRP [3] considers a part of the total customers’ demands to be known when the vehicles have already started. Based on a return to the depot after the completion of service, the CVRP problems can be classified into two groups. The closed VRP [4] allows all the vehicles to return to the depot after completing the service. Simultaneously, the open VRP (OVRP) [5] restricts the vehicles to return to the central depot after visiting the customers. This model is generally used when an organization uses the third party logistics for its supply and distribution of goods. In the VRP with time window (VRPTW) [6], the products have to be delivered to each customer within a predefined time window. In literature, some of the articles related to VRP consider the crisp time window; also, few studies in this domain consider fuzzy time window [7].

Based on the number of depots in the problem, VRP problems are classified as either single-depot VRP [8] or multi-depot VRP (MDVRP) [9]. There is another variation of VRP called green VRP [10], where the main area of concern is the reduction in carbon emission from vehicles. In 2014, a research work [11] was published on environmental VRP, focusing on the minimization of CO2 emission due to transportation. Here a hybrid artificial bee colony (ABC) algorithm was used and it performed slightly better than the existing method in the literature. Schyns [12] has presented a new variety of dynamic VRP that considered time windows, partial split delivery, and heterogeneous fleets. This problem was solved using an ant colony optimization (ACO). Hu et al. [13] proposed the advancement of simultaneous pickup and delivery VRP under closed-loop logistics by imposing the dynamic nature, uncertain pickup, and delivery of incompatible goods. In 2016, a heterogeneous VRP [14] with mixed backhaul was solved using label-based ACO by Wu et al. This problem also considers the time window. Another work on minimizing the time-dependent CO₂ emission in VRP [15], especially in urban areas, was published in 2016. For solving the problem, a tabu search method was used. Du et al. [16] proposed a problem on multi-depot VRP in 2017. In this study, the authors considered hazardous materials. It was solved using fuzzy bi-level programming. Huang et al. [17] investigated a particular type of VRP that supports path flexibility, which allows selecting any of the direct paths from one city to its neighboring city depending upon several conditions like congestion in the network and the journey. In 2017, a new variety of VRP for electric vehicles was first introduced to overcome the limited driving range. It was solved by proper planning of battery charging stations having nonlinear charging functions [18]. Ng et al. [19] proposed a CVRP having strategies under time-dependent traffic congestion. It was solved using the multiple colonies ABC algorithm. In 2018, a new concept of time window assignment VRP having two-stage stochastic optimization [20] was published. Here, the product dependent time window and the delivery schedule were both considered. The work also includes sensitivity analysis on three different models. In the same year, another study [21] was investigated on VRP having a heterogeneous fixed fleet vehicle. This problem also considers fuel and carbon emissions.

3 Motivation

There are several works on different types of VRP that have already been published. However, none of the works has dealt with this type of extra service time that can fulfill some other customers’ requirements that may not be known a priory but belong to the same city where the vehicle is waiting for the service. In this paper, we consider extra service time to handle such a scenario. Sometimes, the enlisted customers may demand some more requirements beyond the pre-ordered amount because of several issues like weather effect, sessional, fashion trends, etc. This incremented demand produces more profit for the organization without investing in fuel costs to move vehicles to other cities. Most of the papers on VRP does not consider any time budget for the vehicles. Here we also consider maximum time budget constraints in our proposed model. It motivates us to develop a new model that can implement the above facts. Besides, most of the papers published on VRP are in a crisp environment. However, real-life problems are uncertain in nature. The fuzzy set theory has the flexibility to deal with this kind of uncertainty based on expert judgment, past data, and observations. Here, travel time, travel cost, and capacity constraint parameters are considered fuzzy numbers to implement the above model.

4 Mathematical preliminaries

This section presents some basic definitions and properties related to fuzzy numbers.

Fuzzy number

A fuzzy number [22] is a normal and convex fuzzy subset of $R$ . Here, the normality indicates that there is a point x₀ such that $μ_{\tilde{A}} (x_{0}) = 1$ and convexity means that $μ_{\tilde{A}} ({α x}_{1} + (1 - α) x_{2}) ⩾ \min (μ_{\tilde{A}} (x_{1}), μ_{\tilde{A}} (x_{2}))$ (1) for any α ∈ [0, 1] and $x_{1}, x_{2} \in R$ .

Triangular fuzzy number

It is a fuzzy number represented as $\tilde{A} = (a_{1}, a_{2}, a_{3})$ where a₁ ⩽ a₂ ⩽ a₃ and the membership function of the number is given by $μ_{\tilde{A}} (x) = {\begin{matrix} 0 & if & x 〈 a_{1} and x 〉 a_{3} \\ \frac{x - a_{1}}{a_{2} - a_{1}} & if & a_{1} ⩽ x ⩽ a_{2} \\ \frac{a_{3} - x}{a_{3} - a_{2}} & if & a_{2} ⩽ x ⩽ a_{3} \end{matrix}$ (2)

Possibility, Necessity, and Credibility measure

If $\tilde{A}$ and $\tilde{B}$ be two fuzzy variables, then the possibility of $\tilde{A} * \tilde{B}$ is represented by $Pos (\tilde{A} * \tilde{B})$ where * represents any relation 〈, 〉, = , ⩾ , ⩽ and this is defined as $\begin{matrix} Pos (\tilde{A} * \tilde{B}) = \\ {sup {min (μ_{\tilde{A}} (x), μ_{\tilde{B}} (y))}, x, y ɛ R, x * y} \end{matrix}$ (3)

The necessity of $\tilde{A} * \tilde{B}$ is represented by $Nes (\tilde{A} * \tilde{B})$ where * represents any relation 〈, 〉, = , ⩾ , ⩽ and this is defined as $\begin{matrix} Nes (\tilde{A} * \tilde{B}) = \\ {inf {max (1 - μ_{\tilde{A}} (x), μ_{\tilde{B}} (y))}, x, y ɛ R, x * y \end{matrix}$ (4)

The relationship between possibility and necessity is given by $Nes (\tilde{A} * \tilde{B}) = 1 - Pos (\bar{\tilde{A} * \tilde{B}})$ (5)

The credibility of the fuzzy event $\tilde{A}$ is defined as the average of possibility and necessity of that event [23].

The credibility measure is defined as $Cr (\tilde{A}) = \frac{1}{2} [Pos (\tilde{A}) + Nes (\tilde{A})]$ (6)

The credibility of a triangular fuzzy variable $\tilde{A} = (a_{1}, a_{2}, a_{3})$ can be defined as $Cr {\tilde{A} ⩾ r} = {\begin{matrix} 1 ifr ⩽ a_{1} \\ \frac{2 a_{2} - a_{1} - r}{2 (a_{2} - a_{1})} & if & a_{1} ⩽ r ⩽ a_{2} \\ \frac{a_{3} - r}{2 (a_{3} - a_{2})} & if & a_{2} ⩽ r ⩽ a_{3} \\ 0 & if & r ⩾ a_{3} \end{matrix}$ (7) $Cr {\tilde{A} ⩽ r} = {\begin{matrix} 0 & if & r ⩽ a_{1} \\ \frac{r - a_{1}}{2 (a_{2} - a_{1})} & if & a_{1} ⩽ r ⩽ a_{2} \\ \frac{2 a_{3} - a_{2} - r}{2 (a_{3} - a_{2})} & if & a_{2} ⩽ r ⩽ a_{3} \\ 1 & if & r ⩾ a_{3} \end{matrix}$ (8)

The expected value of a fuzzy variable

The expected value of a fuzzy variable $\tilde{A}$ is represented as $E (\tilde{A})$ and defined as $E (\tilde{A}) = \int_{0}^{\infty} Cr {\tilde{A} ⩾ r} dr - \int_{- \infty}^{0} Cr {\tilde{A} ⩽ r} dr$ (9)

Now for any two bounded fuzzy variable $\tilde{x}$ , $\tilde{y}$ and for any real number a, b, $E (a \tilde{x} + b \tilde{y}) = aE (\tilde{x}) + bE (\tilde{y})$ (10)

If $\tilde{A} = (a_{1}, a_{2}, a_{3})$ be a triangular fuzzy number. Then the expected value of $\tilde{A}$ i.e., $E (\tilde{A})$ will be $\frac{a_{1} + 2 a_{2} + a_{3}}{4}$ .

Theorem 1. [24]: If $\tilde{A} = (a_{1}, a_{2}, a_{3}, a_{4})$ be a trapezoidal fuzzy variable with a₁ ⩽ a₂ ⩽ a₃ ⩽ a₄

Then for a predefined credibility level β with 0 ⩽ β ⩽ 1 $Cr {h (x) ⩽ \tilde{A}} ⩾ β$ is equivalent to $h (x) ⩽ {\begin{matrix} (1 - 2 β) a_{4} + 2 β a_{3} & if & β ⩽ 0.5 \\ 2 (1 - β) a_{2} + (2 β - 1) a_{1} & if & β > 0.5 \end{matrix}$ (11)

Now for a Trapezoidal fuzzy variable $\tilde{A} = (a_{1}, a_{2}, a_{3}, a_{4}$ ) if a₂ = a₃. Then the fuzzy variable is converted to a triangular fuzzy variable. In that case, for a fuzzy variable $\tilde{A} = (a_{1}, a_{2}, a_{3})$ with a₁ ⩽ a₂ ⩽ a₃ for a predefined credibility level β with 0 ⩽ β ⩽ 1 $Cr {h (x) ⩽ \tilde{A}} ⩾ β$ is equivalent to $h (x) ⩽ {\begin{matrix} (1 - 2 β) a_{3} + 2 β a_{2} & if & β ⩽ 0.5 \\ 2 (1 - β) a_{2} + (2 β - 1) a_{1} & if & β > 0.5 \end{matrix}$ (12)

Fuzzy Linear Programming problem and Expected value model: $Let MaxZ = f (u, \tilde{A})$ (13)

Subject to $g_{j} (u, \tilde{A}) ⩽ 0 for j = 1, 2, \dots . k .$ (14) be a fuzzy linear programming problem with u and $\tilde{A}$ as decision variables and the fuzzy variable, respectively.

Then by the expected value model [23] the problem can be written as $MaxZ = E [f (u, \tilde{A})]$ (15)

Subject to $E [g_{j} (u, \tilde{A})] ⩽ 0] for j = 1, 2, \dots, k .$ (16)

5 Model formulation

The concept of CVRP was first described by Dantzig et al. [25] in their work on the truck dispatching problem. In general, the CVRP can be stated as follows: Let G = (V, E) be a graph, V ={ j₀, j₁, j₂, …… , j_n } being the vertex set, where j₀ refers to the depot, the customers are indexed j₀, j₁, j₂, …… , j_n, and E ={ (j_m, j_n) : j_m, j_nV } is the edge set. Let n be the number of customers, N ={ 1, 2, . . . , n } be the set of customers, and {0, 1, . . . , n} be the set of all customers and the depot, which is identified by 0. Let K = {1, 2, . . . , k} be the set of vehicles. Let d_i denote the demand of customer iN, L_h indicate the load of vehicle hK. Furthermore, $\underset{i}{maxd} ⩽ \underset{h}{maxL}$

Let us define $y_{ih} = {\begin{matrix} 1, if point i is visited by vehicle h \\ 0, otherwise \end{matrix}$ (17) $x_{ijh} = {\begin{matrix} 1, if vehicle h travels directly \\ from node i to node j \\ 0, otherwise \end{matrix}$ (18)

Let c_ij be the transport cost or distance between customer i and j. The mathematical model for the capacitated vehicle routing problem is as follows: $Min Z = \sum_{i = 0}^{n} \sum_{j = 0}^{n} \sum_{h = 1}^{k} c_{ij} x_{ijh}$ (19)

Subject to $\sum_{i = 0}^{n} d_{i} y_{ih} ⩽ L_{h}, \forall h$ (20) $\sum_{h = 1}^{k} y_{ih} = 1, i = 1, 2, \dots \dots . ., n$ (21) $\sum_{i = 0, i \neq j}^{n} x_{ijh} = y_{jh}, j = 1, 2, \dots \dots . ., n \forall h$ (22) $\sum_{j = 0, j \neq i}^{n} x_{ijh} = y_{ih}, i = 1, 2, \dots \dots . ., n \forall h$ (23)

In the above model, expression (19) represents the objective function, which minimizes the total distance traveled. The inequality (20) refers to the vehicle capacity constraints; that is, certain vehicles must deliver goods to customers within vehicle capacity; (21) represents that every demand node must be served. Equations (22) and (23) ensure that each demand node is served by the same vehicle.

The above model is used for minimization of cost based on distance. A little modification in the model mentioned above can be used for the minimization of time. Let t_ij be the time needed to reach from customer i to customer j. The objective function, in this case, will be as follows: $Min Z = \sum_{i = 0}^{n} \sum_{j = 0}^{n} \sum_{h = 1}^{k} t_{ij} x_{ijh}$ (24) subject to the constraints (20) –(23).

But, it has been observed that, though in the cost minimization problem, the best solution incurs at the least amount of cost, it takes more time. Whereas, in the case of the time minimization problem, the best solution occurs at least time but assumes more cost. This multi-objective trade-off is optimized using the weighted average method. The objective function will be as follows: $Min Z = \sum_{i = 0}^{n} \sum_{j = 0}^{n} \sum_{h = 1}^{k} c_{ij} x_{ijh} + c_{1} \sum_{i = 0}^{n} \sum_{j = 0}^{n} \sum_{h = 1}^{k} t_{ij} x_{ijh}$ (25) subject to the constraints (2) - (7).

Here c₁ is the cost of traveling for one unit of time. This model yields optimized results in terms of cost awell as time. That is why this model is better to use in the case of profit maximization.

5.1 Proposed Profit Maximization CVRP crisp model

The list of parameters and variables referred to in the model with their description is presented below:

Parameters

C_ij Cost of traveling from i^th city to j^th city

R_i (t) Return or earning for spending t time in the i^th city

E_i (t) Expenses for spending t time in the i^th city

P_i (t) Profit on the i^th city after spending t time that is equal to the minus of E_i (t) from R_i (t)

u_i Base demand in city i.

v_i Excess demand in city i per unit of time

d_i Total demand in city i, so d_i = (u_i + v_it_i)

s_i Service time to supply the base demand in i^th city

t_i Extra service time to supply the excess demand in i^th city

T_i Total service time in i^th city, i.e., the addition of S_i and t_i

t_ij Travel time from i^th city to j^th city

M Total tour time

N Count of all the cities

K Count of all the vehicles or salespeople

C Capacity of each vehicle

α Percentage of the capacity of the vehicle that is filled up for the base demands of the route of each vehicle

C₁ Cost of traveling for one unit of time

Decision variables $x_{ij}^{k} = {\begin{matrix} 1 if the vehicle k travels \\ from city i to city j \\ 0 otherwise \end{matrix}$ (26) $y_{j}^{k} = {\begin{matrix} 1 if the city j is visited by the vehicle k \\ 0 otherwise \end{matrix}$ (27)

Now the proposed model for a crisp variable can be formulated as follows, $\begin{matrix} Maximize Z = \sum_{i = 1}^{N} P_{i} (t) - \\ (\sum_{i = 1}^{N} \sum_{j = 1}^{N} \sum_{k = 1}^{K} C_{ij} x_{ij}^{k} + C_{1} \sum_{i = 1}^{N} \sum_{j = 1}^{N} \sum_{k = 1}^{K} t_{ij} x_{ij}^{k}) \end{matrix}$ (28)

Equation (28) represents the objective function, i.e., maximization of profit, which is the difference between the profit and expenditure.

Where $P_{i} (t_{i}) = r_{i} (t_{i}) - e_{i} (t_{i}), i = 1, 2, \dots \dots ., N$ (29)

Equation (29) gives the profit of visiting the i^th city by fulfilling some extra demand. t_i is the extra service time needed to fulfill the excess demand. r_i (t_i) and e_i (t_i) are the return and expenditure respectively during the extra service time t_i.

Subject to $\begin{matrix} \sum_{i = 1}^{N} x_{ij}^{k} = y_{j}^{k}, j = 1, 2, \dots \dots . ., N, \\ k = 1, 2, \dots \dots, K \end{matrix}$ (30) $\begin{matrix} \sum_{j = 1}^{N} x_{ij}^{k} = y_{i}^{k}, i = 1, 2, \dots \dots . ., N \forall k, \\ k = 1, 2, \dots \dots, K \end{matrix}$ (31)

Equations (30) and (31) indicate every city must be served by exactly one vehicle. $\sum_{k = 1}^{K} y_{i}^{k} = 1, i = 1, 2, \dots \dots . ., N$ (32)

As $y_{i}^{k}$ is a binary decision variable, at least one $y_{i}^{k}$ must be 1. Therefore, Equation (32) refers to every city must be served $\begin{matrix} (\sum_{i = 1}^{N} \sum_{j = 1}^{N} t_{ij} x_{ij}^{k} + \sum_{i = 1}^{N} (s_{i} + t_{i}) y_{i}^{k}) \\ ⩽ M, \forall k = 1, 2, \dots \dots, K \end{matrix}$ (33)

M is the maximum tour limit for each vehicle. The left-hand side of in Equation (33) represents each vehicle’s total tour time, which includes travel time from one city to another in its route, service time, and excess service time at each city. Therefore, in Equation (33) ensures the time budget constraints. $\sum_{i = 1}^{N} u_{i} y_{i}^{k} ⩽ α C, k = 1, 2, \dots \dots, K$ (34) $\sum_{i = 1}^{N} v_{i} t_{i} y_{i}^{k} ⩽ (1 - α) C, k = 1, 2, \dots \dots, K$ (35) $\sum_{i = 1}^{N} d_{i} y_{i}^{k} ⩽ C, k = 1, 2, \dots \dots, K (36)$ (36)

The inequality in (36) shows the excess demand constraint, which ensures that the total demand served by one vehicle must be less than the vehicle capacity. Condition (36) derived from the two restrictions (34) and (35). where (35) refers to excess demand constraint, and the in Equation (16) indicates base demand constraint. $\begin{matrix} x_{ij}^{k} = 0 or 1, i, j = 0, 1, \dots, N, i \neq j, \\ k = 1, 2, \dots . ., K \end{matrix}$ (37) $\begin{matrix} y_{i}^{k} = 0 or 1, i = 1, 2, . . . ., N, k = 1, 2, \dots ., K \end{matrix}$ (38)

Equations (37) and (38) shows that $x_{ij}^{k}$ and $y_{i}^{k}$ are binary decision variables. $\sum_{i = 1}^{N} x_{i 0} ⩽ K$ (39) $\sum_{j = 1}^{N} x_{0 j} ⩽ K$ (40)

Constraints (39) and (40) are used to maintain the maximum vehicle number constraints, i.e., the total number of vehicles that start from the depot and return to the same must not exceed K. $\sum_{k = 1}^{K} \sum_{i \in S} \sum_{j \in S, i \neq j} x_{ij}^{k} ⩽ | S | - 1 \forall S \subseteq V \ {0}$ (41)

Sub-tour elimination constraint is represented by constraint (41).

5.2 Fuzzy model for profit maximization CVRP

Generally, in a real-life scenario, the profit-maximization mentioned above the CVRP problem is not deterministic because there are different uncertain parameters like the traveling cost and traveling time from one city to another city. In this paper, t_ij, c_ij and α are taken as TFN (Triangular Fuzzy Number) represented by ( $t_{ij}^{l}, t_{ij}^{m}, t_{ij}^{n}$ ), ( $c_{ij}^{l}, c_{ij}^{m}, c_{ij}^{n}$ ) and (α^l, α^m, αⁿ) respectively. Let the membership function of the above mentioned fuzzy numbers are μ_t(x), μ_c(x), and μ_α(x) respectively and given by $μ_{t} (x) = {\begin{matrix} 0 & if & x 〈 t_{ij}^{l} and x 〉 t_{ij}^{n} \\ \frac{x - t_{ij}^{l}}{t_{ij}^{m} - t_{ij}^{l}} & if & t_{ij}^{l} ⩽ x < t_{ij}^{m} \\ \frac{t_{ij}^{n} - x}{t_{ij}^{n} - t_{ij}^{m}} & if & t_{ij}^{m} ⩽ x ⩽ t_{ij}^{n} \end{matrix}$ (42) and similarly, μ_c(x) and μ_α(x) can be defined.

Therefore, the fuzzy model can be written as $\begin{matrix} Maximize Z = \sum_{i = 1}^{N} P_{i} (t) - \\ (\sum_{i = 1}^{N} \sum_{j = 1}^{N} \sum_{k = 1}^{K} {\tilde{c}}_{ij} x_{ij}^{k} + C_{1} \sum_{i = 1}^{N} \sum_{j = 1}^{N} \sum_{k = 1}^{K} {\tilde{t}}_{ij} x_{ij}^{k}) \end{matrix}$ (43)

Subject to constraints (30) to (33) and (35) to (41) and $\sum_{i = 1}^{N} u_{i} y_{i}^{k} ⩽ \tilde{α} C, k = 1, 2, \dots \dots, K$ (44)

Now using the fuzzy expected value of the objective the above fuzzy model can be written as $\begin{matrix} Maximize Z = \\ E [\sum_{i = 1}^{N} P_{i} (t) - (\sum_{i = 1}^{N} \sum_{j = 1}^{N} \sum_{k = 1}^{K} {\tilde{c}}_{ij} x_{ij}^{k} \\ + C_{1} \sum_{i = 1}^{N} \sum_{j = 1}^{N} \sum_{k = 1}^{K} {\tilde{t}}_{ij} x_{ij}^{k})] \end{matrix}$ (45)

Subject to constraints (30) to (33) and (35) to (41) and $Cr (\sum_{i = 1}^{N} u_{i} y_{i}^{k} - \tilde{α} C) ⩽ 0, k = 1, 2, \dots \dots, K$ (46)

As the fuzzy variable ${\tilde{c}}_{ij}$ and ${\tilde{t}}_{ij}$ are independent fuzzy variable, using chance constraints method [26] the objective function (43) can be written as $\begin{matrix} Maximize Z = \sum_{i = 1}^{N} P_{i} (t) - \\ (\sum_{i = 1}^{N} \sum_{j = 1}^{N} \sum_{k = 1}^{K} Cr [{\tilde{c}}_{ij}] x_{ij}^{k} + \\ C_{1} \sum_{i = 1}^{N} \sum_{j = 1}^{N} \sum_{k = 1}^{K} Cr [{\tilde{t}}_{ij}] x_{ij}^{k}) \end{matrix}$ (47)

For predefined credibility β, the fuzzy constraint can be written as $\begin{matrix} \sum_{i = 1}^{N} u_{i} y_{i}^{k} ⩽ \\ {\begin{matrix} ((1 - 2 β) α_{3} if β ⩽ 0.5 \\ + 2 β α_{2}) C \\ (2 (1 - β) α_{2} if β > 0.5 \\ + (2 β - 1) α_{1}) C \end{matrix}, k = 1, 2, \dots, K \end{matrix}$ (48)

When we apply the expected value model, the objective (43) changes to (45), and the constraint (44) will be as follows $\sum_{i = 1}^{N} u_{i} y_{i}^{k} ⩽ CE (\tilde{α)}, k = 1, 2, \dots \dots, K$ (49)

6 The proposed GA-BFO based hybrid algorithm for the proposed model

The purpose of the proposed GA-BFO hybrid algorithm is to maximize the profit by routing the vehicles in the optimized paths concerning the time, and the distance traveled and finding the excess service time to invest in each city. It combines the capability of both GA and BFO algorithms. The proposed model’s solution involves two decisions; first, each vehicle’s route serves all the cities, and second is the extra time to be spent on each city to gain more profit. Therefore, we use GA for finding the route of each vehicle and BFO to find the extra time to stay in each city. This paper uses natural number encoding for the genetic algorithm. Each chromosome is broken into some sub-routes based on the different destination cities’ demands and the vehicles’ capacity. Figure 1 represents the natural number encoding scheme, i.e., the depot is represented by 0, city 1 is represented by 1, and so on.

Fig. 1

Encoding rule.

Here, GA finds the best-optimized paths concerning the time and the distance traveled. Each sub-route is passed through the BFO algorithm to get the excess service time spent in each city to get the highest return, excluding the cost. Finally, it concatenates all the extra service times to generate a complete list of excess time for the cities in different routes. This algorithm finds the best chromosome in each generation from the GA part, and then it passes all the sub-routes to the BFO part to get the best optimized excess service time to wait in each city. In this way, for each generation of GA, the BFO executes the same number of times as the number of vehicles used in the best path found through the GA. A flowchart of this proposed algorithm is depicted in Figs. 2 and 3.

Fig. 2

Flowchart of the GA part.

Fig. 3

Flowchart of the BFO part.

Genetic Algorithm:

In the literature, several algorithms originate from the concept of biology. The Genetic algorithm is a state-of-the-art algorithm based on evolutionary strategy. The proper implementation of the algorithm helps to get better solutions in each successive iteration and these iterations continue until it converges or reaches the desired accuracy level. In this paper, we use a combination of natural numbers to encode the solution, where each natural number represents a city. Here we use tournament selection with selection pressure four, and for crossover operation, we choose partially mapped crossover (PMX) [27]. The mutation prevents the algorithm from being trapped in local minima. To implement GA, we use inversion mutation [28]. Two essential parameters of GA are crossover and mutation probability. Here we use adaptive crossover and mutation probability [29]. The first part of Fig. 2 shows the flowchart of GA.

Fitness function:

The fundamental objective of VRP is to minimize the cost of the distance traveled and the total time of the traversal. Based on this concept, the fitness function used in GA is designed as follows. $Fitness = \sum_{i = 0}^{n} \sum_{j = 0}^{n} \sum_{h = 1}^{k} c_{ij} x_{ijk} + c_{1} \sum_{i = 0}^{n} \sum_{j = 0}^{n} \sum_{h = 1}^{k} t_{ij} x_{ijk}$ where C₁ is the traveling cost for one unit of time.

Bacterial Foraging Optimization Algorithm:

The social behaviors of animals like food searching by ants, bird-flocking, and fish-schooling inspire scientists to develop so many meta-heuristic algorithms like ant colony optimization, particle swarm optimization, etc. The BFOA is also an example of such an algorithm. Passino proposed it [30], and it was a new nature-inspired algorithm in this series. This algorithm mimics the foraging behavior of E. Coli bacteria. At the time of foraging, the bacteria used to move with the help of a set of flagella. This movement of bacteria is generally of two types, tumble or swim. Depending on the food value or the health value, the new position the bacteria used to move straight (i.e., called swim) or in a different direction (tumble). These two steps are called chemotaxis movement and the purpose of these movements is to move towards a nutrient gradient and avoid a noxious environment. The problem with this algorithm is that it suffers from local minima, and its convergence rate is prolonged. We use an adaptive search mechanism and comprehensive learning strategy to increase the performance of the original BFO, called adaptive-comprehensive learning bacterial foraging optimization (ACLBFO) [31]. Here the extra service time spent on the cities of a vehicle’s route represents a bacteria, and the total profit in the corresponding cities of a route forms the health value. In the original BFO, a constant chemotaxis step length is used. This paper used an adaptive chemotaxis step length based on a non-linearly decreasing modulation model. The BFO algorithm is subdivided into four parts: chemotaxis, swarming, reproduction, and elimination and dispersal.

The steps of the BFO algorithm are described as follows:

First, initialize all the parameters and positions: population (S), Number of chemotactic steps in life (N_c), Maximum number of continuous swim (N_s), Number of reproduction steps to be taken (N_re), Number of elimination dispersal steps to be considered (N_ed), Probability of elimination dispersal (P_ed), Chemotaxis step length (C), etc.

Evaluate the health value of the bacteria of the initial population.

Find the pbest of each bacteria with the present values.

For (elimination dispersal loop- N_ed iterations)

For (reproduction loop –N_re iterations)

For (chemotaxis loop- N_c iterations)

Calculate the chemotaxis step size using the Equation (33).

Update the position of all the bacteria either by swim or tumble based on the health value and the value of N_s.

Evaluate the health value of all the new bacteria.

Update the gbest and the pbest.

End For (Chemotaxis loop).

Sort bacteria in descending order based on the health values.

Replace the last S/2 bacteria with the first S/2 bacteria.

End For (Reproduction loop).

Do Eliminate and dispersal step with probability P_ed.

End For (elimination dispersal loop).

End

6.1 Proposed GA-BFO algorithm

The steps of this proposed algorithm are described as follows:

Initialize parameters: Population size P_SZ; Maximum Generation MX_GN ; Adaptive crossover probability and mutation probability parameters; Number of destinations NOD.

Generate a random population using the natural number encoding method.

Let Gn_count = 1.

If Gn_count ⩽ MX_GN, then go to step5, else go to step 24.

Decode chromosomes and calculate the corresponding fitness values.

Let C = 2.

Selection of parents using the Tournament selection method.

Calculate p_c and apply crossover.

Calculate p_m and apply the mutation.

Store the offspring to form a new population.

C = C + 2.

If C ⩽ P_SZ, then go to step 7, else go to step 13.

Search the new population and store the best chromosome.

Calculate the number of vehicles TNV required to fulfill the demands of the destination cities present in the best chromosome.

Let VC = 1.

If VC ⩽ TNV then go to step 17, else go to step 21.

Apply the BFO algorithm on the VC^th sub route.

Store the extra service time for each city for VC^th sub route.

VC = VC + 1.

Go to step 16.

Store the concatenated list of extra service time and profit in BST_SL.

Gn_count = Gn_count + 1

go to step 4.

Find the best solution from BST_SL in respect of profit.

6.2 Constraint handling

Constraint handling mechanism is a critical issue of solving constrained optimization problems while using the meta-heuristic algorithms. Most of the constraints, except the restrictions (30) to (36) and (39) to (41), are already explained in the mathematical model section. The others are defined in the literature. Therefore, only the above constraints are needed to be explained about its handling method.

In this paper, we use natural number encoding (a combination of 1 to n, where n is the number of the city in the problem) to represent chromosomes for finding the routes of each vehicle using GA. Constraints (30), (31), and (32) are automatically satisfied due to the encoding scheme used in this paper. We encode a bacterium in BFO with real numbers to optimize the excess service time. The total extra service time for a vehicle is equal to the difference between the predefined maximum allowable time and that vehicle’s tour time. We randomly split the whole excess service time among the cities in the route. Therefore, the constraint (33) is also preserved using the encoding scheme of the bacterium. Restrictions (34), (35), and (36) are the constraints related to vehicle capacity. When we evaluate each chromosome’s fitness, we take care of vehicle capacity constraints by varying the vehicle number. Conditions (39) and (40) are related to the number of maximum allowable vehicles. We make a significant penalty on fitness value if these two constraints are violated. Restraint (41) is used for sub-tour elimination; this is also maintained through the encoding scheme and binary decision variables.

7 Numerical Illustration

7.1 Test data for the crisp model

The proposed model is tested for a different set of cities (n = 10, 19, 28, 37, 46, 55, 64, 73, 82, 91). The traveling cost and traveling time between different cities are shown in Tables 1 and 2. The return and expenditure at the i^th city for the extra service-time t is calculated using the following expressions. $R_{i} (t) = a_{i} + b_{i} t - c_{i} t^{2}$ (50) $E_{i} (t) = e_{i} t$ (51) where a_i is the fixed return from the visit of i^th city, b_i and c_i are the constants used to determine the return from the excess stay of time t on i^th city, e_i is the rate of expenses per unit time of excess stay at the i^th city. The values of the constants mentioned above are given in Table 3.

Table 1

The cost to travel between different cities (D for depot)

City	D	1	2	3	4	5	6	7	8	9
D	0	25	28	32	20	26	37	8	29	20
1	37	0	20	28	35	40	30	42	28	4
2	42	28	0	30	25	35	9	32	40	30
3	28	30	7	0	20	25	30	35	22	37
4	37	22	35	30	0	20	25	30	9	28
5	25	30	25	8	28	0	32	40	32	30
6	28	25	30	22	37	40	0	10	32	20
7	20	5	32	40	35	25	40	0	22	37
8	30	40	35	25	20	22	37	32	0	28
9	28	30	28	20	11	32	37	40	30	0

Table 2

The time to travel between different cities (D for depot)

City	D	1	2	3	4	5	6	7	8	9
D	0	5	5	4	6	5	3	7	5	6
1	3	0	6	5	3	2	5	2	6	7
2	2	6	0	5	7	4	9	5	2	5
3	6	5	9	0	8	7	5	4	8	3
4	3	8	4	5	0	8	7	5	9	6
5	7	5	7	9	6	0	5	2	5	5
6	6	7	5	8	3	2	0	9	5	8
7	8	9	5	2	4	7	2	0	8	3
8	5	2	4	7	8	8	3	5	0	6
9	6	5	6	8	7	5	3	2	5	0

Table 3

Values of different constants and BD and EDR

	D	1	2	3	4	5	6	7	8	9
a_i	0	60	140	70	80	110	90	100	60	140
b_i	0	33	60	47	65	64	31	18	110	116
c_i	0	1	3	4	2	3	2	1	4	5
e_i	0	9	12	7	13	10	7	8	14	6
BD	0	4	2	3	4	1	5	4	2	3
EDR	0	0.6	0.3	0.1	0.7	0.2	0.5	0.4	0.6	0.4

(Note: D for the depot, BD for Base Demand, EDR for Excess Demand Rate).

7.2 Test data of the fuzzy model

The proposed fuzzy model is tested for the 10 city problem, including depot. The traveling cost and traveling time between different cities are shown in Tables 4 and 5. The return and expenditure at the i^th city for the extra service-time t is the same as the crisp problem. The fuzzy distance presented in Table 4 are generated based on crisp distance using the following formulae d_l = d - random number between (1, 2, 3), d_m = d and d_u = d + random number between (1, 2, 3)

Table 4
Fuzzy distance between different cities (D for depot)

D 1 2 3 4 5 6 7 8 9

D 0 22, 25, 27 27, 28, 29 30,32, 33 19, 20, 22 24, 26, 28 35, 37, 38 7, 8, 9 28, 29, 31 18, 20, 21

1 36, 37, 39 0 18, 20,21 27, 28,30 34, 35, 36 38, 40, 41 29, 30, 31 41, 42, 43 27, 28, 30 3, 4, 5

2 40, 42,43 27, 28, 30 0 28, 30, 31 24, 25, 26 34, 35, 37 8, 9, 10 31, 32, 33 39, 40, 41 29, 30,31

3 27, 28, 30 29,30, 31 6, 7, 8 0 19, 20, 21 23, 25, 26 28, 30, 32 34, 35, 36 21, 22, 24 35, 37, 39

4 27, 28, 30 21, 22, 23 33, 35, 37 29, 30, 31 0 19, 20, 22 24, 25, 26 29, 30, 31 8, 9, 10 27, 28, 30

5 23, 25, 26 29, 30, 31 23, 25, 26 7, 8, 10 26, 28, 30 0 31, 32, 33 38, 40, 41 31, 32, 33 28, 30, 31

6 27, 28,29 23, 25, 26 29, 30, 31 21, 22, 23 35, 37, 39 39, 40, 41 0 9, 10, 11 31, 32, 34 19, 20, 22

7 19, 20, 22 4, 5, 6 31, 32, 33 39, 40, 42 34, 35, 36 24, 25, 26 39, 40, 42 0 21, 22, 23 36, 37, 39

8 28, 30, 31 39, 40, 41 34, 35, 36 24, 25, 26 18, 20, 21 20, 22, 24 35, 37, 40 31, 32, 33 0 27, 28, 29

9 27, 28, 30 29, 30, 31 27, 28, 30 18, 20, 22 10, 11, 12 31, 32, 33 36, 37, 38 39, 40, 41 29, 30, 31 0

	D	1	2	3	4	5	6	7	8	9
D	0	22, 25, 27	27, 28, 29	30,32, 33	19, 20, 22	24, 26, 28	35, 37, 38	7, 8, 9	28, 29, 31	18, 20, 21
1	36, 37, 39	0	18, 20,21	27, 28,30	34, 35, 36	38, 40, 41	29, 30, 31	41, 42, 43	27, 28, 30	3, 4, 5
2	40, 42,43	27, 28, 30	0	28, 30, 31	24, 25, 26	34, 35, 37	8, 9, 10	31, 32, 33	39, 40, 41	29, 30,31
3	27, 28, 30	29,30, 31	6, 7, 8	0	19, 20, 21	23, 25, 26	28, 30, 32	34, 35, 36	21, 22, 24	35, 37, 39
4	27, 28, 30	21, 22, 23	33, 35, 37	29, 30, 31	0	19, 20, 22	24, 25, 26	29, 30, 31	8, 9, 10	27, 28, 30
5	23, 25, 26	29, 30, 31	23, 25, 26	7, 8, 10	26, 28, 30	0	31, 32, 33	38, 40, 41	31, 32, 33	28, 30, 31
6	27, 28,29	23, 25, 26	29, 30, 31	21, 22, 23	35, 37, 39	39, 40, 41	0	9, 10, 11	31, 32, 34	19, 20, 22
7	19, 20, 22	4, 5, 6	31, 32, 33	39, 40, 42	34, 35, 36	24, 25, 26	39, 40, 42	0	21, 22, 23	36, 37, 39
8	28, 30, 31	39, 40, 41	34, 35, 36	24, 25, 26	18, 20, 21	20, 22, 24	35, 37, 40	31, 32, 33	0	27, 28, 29
9	27, 28, 30	29, 30, 31	27, 28, 30	18, 20, 22	10, 11, 12	31, 32, 33	36, 37, 38	39, 40, 41	29, 30, 31	0

Table 5

Fuzzy Time to travel between different cities (D for depot)

City	D	1	2	3	4	5	6	7	8	9
D	0	4,5,6	4,5,6	3,4,5	5,6,7	4,5,6	2,3,4	6,7,8	4,5,6	5,6,7
1	2,3,4	0	5,6,7	4,5,6	2,3,4	1,2,3	4,5,6	1,2,3	5,6,7	6,7,8
2	1,2,3	5,6,7	0	4,5,6	6,7,8	3,4,5	8,9,10	4,5,6	1,2,3	4,5,6
3	5,6,7	4,5,6	8,9,10	0	7,8,9	6,7,8	4,5,6	3,4,5	7,8,9	2,3,4
4	2,3,4	7,8,9	3,4,5	4,5,6	0	7,8,9	6,7,8	4,5,6	8,9,10	5,6,7
5	6,7,8	4,5,6	6,7,8	8,9,10	5,6,7	0	4,5,6	1,2,3	4,5,6	4,5,6
6	5,6,7	6,7,8	4,5,6	7,8,9	2,3,4	1,2,3	0	8,9,10	4,5,6	7,8,9
7	7,8,9	8,9,10	4,5,6	1,2,3	3,4,5	6,7,8	1,2,3	0	7,8,9	2,3,4
8	4,5,6	1,2,3	3,4,5	6,7,8	7,8,9	7,8,9	2,3,4	4,5,6	0	5,6,7
9	5,6,7	4,5,6	5,6,7	7,8,9	6,7,8	4,5,6	2,3,4	1,2,3	4,5,6	0

Where (d_l, d_m, d_u) as fuzzy distance and d as crisp distance.

The fuzzy demand in Table 5 are generated based on crisp demand using the following formulae

$\begin{matrix} {Demand}_{l} = demand - 1, {Demand}_{m} = demand \\ and {Demand}_{u} = demand + 1 \end{matrix}$ (52)

Where (Demand_l, Demand_m, Demand_u) is the fuzzy demand and demand is the crisp demand.

8 Results and discussion

The proposed profit maximization VRP model is solved using four different metaheuristic algorithms, viz. GA, ALO, BFO, and proposed hybrid GA-BFO algorithm. All the algorithms are implemented using C language in the Linux environment running on an Intel Core i5 CPU T6500 of 2.44 GHz. Based on the several trial and error, we choose the algorithms’ optimal parameter settings and are presented in Table 6.

Table 6
Parameter Settings of the algorithms

Algorithm Parameter Name Value

Genetic Algorithm Population Size 120

Generation 250

Selection Pressure 4

Crossover Probability 0.95

Mutation Probability 0.05

Bacterial Foraging Optimization Algorithm Population Size 40

Number of chemotactic steps in life 250

Maximum number of continuous swim 6

Number of reproduction steps to be taken 5

Number of elimination dispersal steps to be considered 2

Probability of elimination dispersal 0.05

Chemotaxis step length Min = 0.01, Max = 0.1

Algorithm	Parameter Name	Value
Genetic Algorithm	Population Size	120
	Generation	250
	Selection Pressure	4
	Crossover Probability	0.95
	Mutation Probability	0.05
Bacterial Foraging Optimization Algorithm	Population Size	40
	Number of chemotactic steps in life	250
	Maximum number of continuous swim	6
	Number of reproduction steps to be taken	5
	Number of elimination dispersal steps to be considered	2
	Probability of elimination dispersal	0.05
	Chemotaxis step length	Min = 0.01, Max = 0.1

8.1 Analysis for crisp model

In the case of classical CVRP, it is observed that if we focus on the cost minimization, it takes a long time to complete the tour, whereas, if the focus is on the travel time minimization, then the cost due to distance increases. Table 7 and Fig. 4 show this trade-off of cost and time for the 10 city problem.

Table 7
The cost-time trade-off for ten city problem

Time 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29

Cost 420 410 399 384 358 345 304 297 289 279 283 277 262 265 270 250

Fig. 4

The cost - time trade-off for ten city problem.

The proposed model is solved for the above-mentioned input data set using the proposed GA-BFO hybrid algorithm. The results for the ten cities, including the depot, are presented in Table 8, and here the total tour time M is 30. The second column is used to present the derived path for the corresponding run. Here the 0′s in the path are used to indicate the depot, and the number of vehicles used is presented at the end of the path within parenthesis. The third column presents the excess time to wait in each city, corresponding to the path. For example, consider the first run, the path is 0962015730840 (3) and excess time in the cities are 7.80 0.10 6.10 2.20 6.40 0.20 3.20 8.30 5.70. It means the above path includes three vehicles and the paths are (Depot, 9, 6, 2, Depot),) and (Depot, 8, 4, Depot). The excess waiting time in city9 is 7.80, and the excess waiting time in city6 is 0.10. Similarly, the excess waiting time in city4 is 5.70. The fourth column indicates the total earnings for the corresponding path. The fifth column represents the cost involved in the corresponding path. Finally, the sixth column shows the net profit, i.e., the [return (fourthcolumn) - cost (fifthcolumn)]. The average value of the profits of all runs is 2306.41, the standard deviation is 14.4851, and the best profit obtained is 2323.49.

Table 8

Results for the ten cities including depot when M (total tour time) = 30

Run No	Path (Number of Vehicles)	Excess Time in City	Return	Cost	Profit	Mean/SD
1	0962015730840 (3)	7.800.106.102.206.400.203.208.305.70	2699.69	388	2312.69	2306.41 / 14.4851
2	0742039106580 (3)	0.304.904.814.209.501.300.505.309.20	2699.90	395	2305.90
3	0973062801540 (3)	10.040.173.800.104.207.810.507.506.00	2691.19	395	2296.19
4	0973015406280 (3)	10.040.173.790.707.505.800.104.507.66	2691.38	395	2296.38.
5	0396015720840 (3).	4.408.900.001.606.900.105.408.305.70.	2722.49.	399.	2323.49
6	0514062807390 (3)	7.500.506.000.004.207.900.002.789.22	2664.88	378	2287.88
		1.707.600.304.404.208.900.028.305.70	2714.59	399	2316.59
8	0812065403970 (3)	8.200.106.700.107.005.053.339.670.00	2687.66	402	2286.66
9	0157203960840 (3)	1.406.600.105.904.308.900.008.305.70	2720.80	399	2322.80
10	0962015730840 (3)	7.800.106.101.807.400.103.708.305.70	2709.52	388	2322.52

After the value 30 for M (total tour time), it is increased step by step up to 80 with the gap of 10, resulting the profit is increasing up to a particular value of M (60), and then it is decreasing as the excess waiting after a specific time is not beneficial anymore. This fact is presented in Table 9, and it is clearly observed in Fig. 5.

Table 9

Results for the ten cities including depot with varying M (total tour time)

M	Path	Time	return	cost	Profit
30	0 3 9 6 0 1 5 7 2 0 8 4 0	4.40 8.90 0.00 1.60 6.90 0.10 5.40 8.30 5.70	2722.49	399	2323.49
40	0 5 4 6 0 3 1 9 0 7 8 2 0	4.90 0.02 3.60 1.50 9.35 0.30 9.10 4.70	2717.68	308	2409.68
50	0 9 5 6 0 7 1 4 0 3 2 8 0	9.50 7.70 1.32 0.00 0.00 4.27 4.60 6.30 9.40	2727.09	306	2421.09
60	0 9 2 6 0 7 1 8 0 4 5 3 0	8.71 4.17 0.25 0.00 0.00 8.29 6.80 9.30 3.51	2698.51	232	2466.51
70	0 7 8 6 0 2 1 9 0 4 3 5 0	0.00 6.62 0.00 5.80 0.63 9.70 6.90 4.80 8.32	2704.01	283	2421.01
80	0 9 5 6 0 7 8 1 0 4 3 2 0	9.60 7.60 1.27 0.10 8.17 0.03 6.30 2.40 4.49	2704.73	318	2386.73

Fig. 5

Total tour time vs. Profit for ten city problem (crisp model).

The proposed model and the algorithm are also tested on a different virtual dataset with 19, 28, 37, 46, 55, 64, 73, 82, and 91 cities. The values of cost c_ij and time t_ij are chosen randomly within a range [t_l, t_u] and [c_l, c_u] using the time-cost relationship calculated from the data set used in the 10 city problem. In this paper, the value of t_l = 0 and t_u = 0 are taken into consideration. Now, for the values of a_i, b_i, c_i, e_i, BD and EDR are used repetitively, i.e., The a_i values for city1, city11, city20, and city 29 are the same; similarly, all the other values are generated. The detailed results for 19 cities and the 37 cities problem with M = 80 are presented in Table 10 as an example.

Table 10

Results for a different number of cities when M (total tour time) = 80

No of City	Path	Time	return	cost	profit
19	0 2 7 9 0 10 18 12 0 15 8 11 0 1 17 16 0 3 6 13 0 5 14 4 0	6.70 0.30 9.67 0.97 8.79 2.24 0.90 7.92 2.50 0.00 8.01 0.00 4.40 0.11 3.57 8.20 8.30 8.12	5335.93	707	4628.93
37	0 16 23 6 0 9 22 25 0 1 14 35 30 0 10 3 20 21 0 36 13 7 0 12 18 24 0 19 4 8 0 27 15 34 0 11 29 28 2 0 31 5 32 33 0 17 26 0	2.70 8.60 4.40 8.50 0.82 0.00 0.10 5.20 6.11 2.10 1.30 3.90 5.03 2.74 8.20 1.01 0.01 4.20 8.90 0.03 0.00 1.72 6.30 7.03 0.00 0.00 4.70 4.80 1.22 4.71 2.70 5.30 5.21 0.01 9.00 9.01	9888.81	1256	8632.81

The result of the proposed hybrid method is compared with the GA, BFO, and ant lion optimization (ALO) algorithm. The comparisons are presented in Table 11.

Table 11

Comparison of the results for different algorithms

Instances	The proposed hybrid algorithm	Genetic algorithm	Bacteria foraging algorithm	Ant Lion optimization
Prob-1 10 city	2704.73	2645.92	2655.46	2644.65
Prob-2 19 city	5335.93	5035.88	5139.82	5087.72
Prob-3 28 city	7817.89	7308.55	7286.03	7261.29
Prob-4 37 city	9888.81	9297.28	9299.3	9250.65
Prob-5 46 city	10973.53	10921.41	10893.45	10889.58
Prob-6 55 city	12394.74	12388.05	12391.16	12390.55
Prob-7 64 city	13605.21	13600.11	13595.86	13606.8
Prob-8 73 city	14617.67	14610.55	14617.34	14604.4
Prob-9 82 city	15420.57	15421.6	15417.29	15414.33
Prob-10 91 city	17647.82	17643.61	17634.2	17644.41

It is clear from Table 11 that the proposed hybrid method yields better results compared to other methods we have considered. Out of 25 independent runs of each algorithm, the best are found for the comparison. The descriptive statistics of the results of 37 cities and 28 city problems are presented in Tables 12 and 13. A non-parametric statistical test, namely the Wilcoxon rank-sum test, is also conducted to compare the proposed algorithm with the other algorithms at a 5% level of significance. The null hypothesis is as follows:

$\begin{matrix} H_{0} : The proposed Algorithm is statistically \\ equivalent to the other algorithms . \end{matrix}$

Table 12

Descriptive statistics of the results of 37 cities problems

Method	Mean	Standard Deviation	MSE	Coefficient of variance
GA-BFO Hybrid	9694.421	53.34324209	16.86861	0.550246808
GA	9540.778	98.31236011	31.0891	1.03044385
BFO	9606.142	68.27450863	21.5903	0.710738074
ALO	9585.956	70.89016437	22.41744	0.73952107

Table 13

Descriptive statistics of the results of 28 cities problems

Method	Mean	Standard Deviation	MSE	Coefficient of variance
GA-BFO Hybrid	7531.525	95.47358984	30.1914	1.26765283
GA	7452.05	167.4103795	52.9398103	2.246501023
BFO	7395.926	148.6144672	46.9960209	2.009409872
ALO	7450.457	135.8193314	42.9498438	1.822966449

The detailed results (p-values) are presented in Table 14.

Table 14

The p- values of the Wilcoxon test for the different instances on a comparison of the proposed algorithm with the other algorithms

Instance	p-values of Wilcoxon test
	GA	BFO	ALO
10 City	0.0086	0.0137	0.0194
19 City	0.0054	0.0110	0.1040
28 City	0.0122	0.0136	0.2984
37 City	0.0074	0.0054	0.0136
46 City	0.0033	0.0130	0.0039
55 City	0.0012	0.0040	0.0159
64 City	0.0036	0.0200	0.0091
73 City	0.0028	0.0159	0.0115
82 City	0.0047	0.0167	0.0084
91 City	0.0010	0.0166	0.0094

The average fitness of each generation using the four different algorithms for “Prob-10 91 city” is depicted in Fig. 6.

Fig. 6

Generation vs. average fitness plot for Prob-10 91 city.

Finally, sensitivity analysis is performed on independent parameters, which shows the changes in the main output concerning the changes in different important parameters. It also helps to identify the most critical parameters. Here the sensitivity analysis is conducted with the parameters a_i, b_i, c_i and e_i. Here each of the parameters a_i, b_i, c_i and d_i are varied from -5% to 5% with a step-length of 2.5%. The sensitivity is presented in Table 15 and depicted in Fig. 7.

Table 15

Sensitivity Table for different parameters

% of the change in Parameter	% of the change in profit with the change in
	a _i	b _i	c _i	e_i.
–5	–1.305925659	–7.538122016	1.132495957	0.933051799
–2.5	0.414930563	3.520622959	–1.879110826	–0.062868267
2.5	0.24453588	3.104391674	–1.726492688	–0.701523146
5	1.893418776	6.121634922	–2.792218209	–2.683391071

Fig. 7

Sensitivity Graph.

8.2 Alysis for fuzzy model

In the case of fuzzy profit maximization CVRP, we observe that it shows a better profit compared to the results of a crisp model. The fuzzy model is solved using the expected value model (EVM), as mentioned earlier in this paper. The results of fuzzy profit maximization CVRP are given below for the different values of credibility level β.

The results show that the expected profit is directly proportional to the credibility level. This is presented in Table 16, and this fact is observed in Fig. 8.

Table 16
Results for the ten cities including depot with varying credibility level (Value)

β Value Expected profit

0.35 2296.65

0.45 2297.21

0.55 2392.67

0.65 2480.35

0.75 2623.07

0.85 2608.34

0.95 2622.67

β Value	Expected profit
0.35	2296.65
0.45	2297.21
0.55	2392.67
0.65	2480.35
0.75	2623.07
0.85	2608.34
0.95	2622.67

Fig. 8

Credibility level vs. Expected Profit for ten city problem (Fuzzy model).

8.3 Discussion

This paper proposes a VRP with a tour time constraint to maximize profit in both crisp and fuzzy environments and solved using a GA-BFO hybrid metaheuristic algorithm. Here in Table 7 we presented a trade-off between cost and time while serving the customers, and we observe that time and cost are inversely proportional. From Table 8 we found that the algorithm is producing stable results with less variance in different executions. Table 9 presents the solution of 10 city problem for different tour time. We found the profit initially increasing with the increase of tour time but after a certain limit profit started decreasing. Tables 10 and 11 present the results of the proposed model for a different number of cities. From Tables 12 and 13, we observe that the coefficient of variance for the proposed hybrid GA-BFO algorithm is the smallest among the four algorithms, i.e., the proposed algorithm is more consistent than the other algorithms used here. Figure 6 presents a convergence graph of the algorithms used here. We visualize the steady and faster convergence of the proposed algorithm among the used algorithms from this figure. Table 14 gives the p-value matrix of the Wilcoxon’s Signed-rank test. From this table, we observe that in most cases, the proposed hybrid GA-BFO algorithm is statistically significant than the other algorithms. As we have used some independent parameters to construct the proposed model, we have performed the sensitivity analysis to find the impact of the independent parameters. From Table 15 and Fig. 7 we observe that the parameter b_i has the most significant impact on the proposed model. The results of the proposed fuzzy model using chance constraints are presented in Table 16. The result of Table 16 is also visualized in Fig. 8; from this result, we observe that the profit is approximately proportional to the credibility level.

9 Conclusions

This paper has proposed a new variant of vehicle routing problem applied in several realistic scenarios. Here the VRP problem is considered as a profit maximization problem. This problem has also been illustrated with numerical data in a crisp and fuzzy environment. It shows that this problem is more suitable even in an imprecise paradigm, which fits better in the real-life implications since it produces better results in such cases. For the first time, GA-BFO based hybrid meta-heuristic algorithm has been successfully applied to the VRP problem with the extra service-time concept. The model formulation and algorithm is quite simple and general. The proposed model can be used to implement the sales and marketing of different product-based and service-based companies that want to spend more time in different cities to increase the profit. This profit-maximization VRP can also be extended to rough, fuzzy-rough, random, etc. uncertain environments.

This study can establish a perfect direction for future research to handle many profit maximization problems in a real-world scenario like the distribution of bakery products, dairy products, and perishable food products.

References

Chiang

T.C.

and Hsu

W.H.

, A knowledge-based evolutionary algorithm for the multi-objective vehicle routing problem with time windows, Computers & Operations Research 45 (2014), 25–37.

Laporte

, The Vehicle Routing Problem: An Overview of Exact and Approximate Algorithms, European Journal of Operational Research 59(3) (1992), 345–358.

Psaraftis

, Dynamic vehicle routing: status and prospects, Annals of Operations Research 61 (1995), 143–164.

Cordeau

J.F.

, Gendreau

, Hertz

, Laporte

and Sormany

J.S.

, New Heuristics for the Vehicle Routing Problem, In: A. Langevi and D. Riopel, Eds., Logistics Systems: Design and Optimization, Springer, New York, (2005), 279–297.

Brandao

, A tabu search heuristic algorithm for open vehicle routing problem, European Journal of Operational Research 157 (2004), 5522–5641.

Oliveira

H.C.B.de.

and Vasconcelos

G.C.

, A hybrid search method for the vehicle routing problem with time windows, Annals of Operations Research 180 (2008), 125–144.

Tang

, Pan

, Fung

R.Y.K.

and Lau

, Vehicle routing problem with fuzzy time windows, Fuzzy Sets and Systems 160(5) (2009), 683–695.

Laporte

, Fifty Years of Vehicle Routing, Transportation Science 43(4) (2009), 408–416.

Renaud

, Laporte

and Boctor

F.F.

, A tabu search heuristic for the multi-depot vehicle routing problem, Journal of Computers & Operations Research 23(3) (1996), 229–235.

10.

Lin

, Choy

K.L.

, Ho

G.T.S.

, Chung

S.H.

and Lam

H.Y.

, Survey of green vehicle routing problem: past and future trends, Expert Systems with Applications 41(4) (2013), 1118–1138.

11.

Zhang

, Lee

C.K.M.

, Choy

K.L.

, Ho

and Ip

W.H.

, Design and development of a hybrid artificial bee colony algorithm for the environmental vehicle routing problem, Transportation Research Part D: Transport and Environment 31 (2014), 85–99. DOI: 10.1016/j.trd.2014.05.015

12.

Schyns

, An ant colony system for responsive dynamic vehicle routing, European Journal of Operational Research 245(3) (2015), 704–718. DOI: 10.1016/j.ejor.2015.04.009

13.

, Sheu

, Zhao

and Lu

, A dynamic closed-loop vehicle routing problem with uncertainty and incompatible goods, Transportation Research Part C: Emerging Technologies 55 (2015), 273–297. DOI: 10.1016/j.trc.2015.01.010

14.

, Tian

and Jin

, A label based ant colony algorithm for heterogeneous vehicle routing with mixed backhaul, Applied Soft Computing 47 (2016), 224–234. DOI: 10.1016/j.asoc.2016.05.011

15.

Ehmke

J.F.

, Campbell

A.M.

and Thomas

B.W.

, Vehicle routing to minimize time-dependent emissions in urban areas, European Journal of Operational Research 251(2) (2016), 478–494. DOI: 10.1016/j.ejor.2015.11.034

16.

, Li

, Yu

, Dan

and Zhou

, Multi-depot vehicle routing problem for hazardous materials transportation: A fuzzy bilevel programming, Information Sciences 399 (2017), 201–218. DOI: 10.1016/j.ins.2017.02.011

17.

Huang

, Zhao

, Woensel

T.V.

and Gross

, Time-dependent vehicle routing problem with path flexibility, Transportation Research Part B: Methodological 95 (2017), 169–195. DOI: 10.1016/j.trb.2016.10.013

18.

Montoya

, Guéret

, Mendoza

J.E.

and Villegas

J.G.

, The electric vehicle routing problem with nonlinear charging function, Transportation Research Part B: Methodological 103 (2017), 87–110. DOI: 10.1016/j.trb.2017.02.004

19.

K.K.H.

, Lee

C.K.M.

, Zhang

S.Z.

, Wu

and Ho

, A multiple colonies artificial bee colony algorithm for a capacitated vehicle routing problem and re-routing strategies under time-dependent traffic congestion, Computers & Industrial Engineering 109 (2017), 151–168. DOI: 10.1016/j.cie.2017.05.004

20.

Neves-Moreira

, Silva

D.P.

, Guimarães

, Amorim

and Almada-Lob

, The time window assignment vehicle routing problem with product dependent deliveries, Transportation Research Part E: Logistics and Transportation Review 116 (2018), 163–183. DOI: 10.1016/j.tre.2018.03.004

21.

, Wang

and Zhang

, Heterogeneous fixed fleet vehicle routing problem based on fuel and carbon emission, Journal of Cleaner Production 201 (2018), 896–908. DOI: 10.1016/j.jclepro.2018.08.075

22.

Grzegorzewski

, Nearest interval approximation of a fuzzy number, Fuzzy Sets and Systems 130(3) (2002), 321–330.

23.

Liu

and Liu

Y.K.

, Expected value of fuzzy variable and fuzzy expected value models, IEEE Transactions on Fuzzy Systems 10 (2002), 445–450.

24.

Kundu

, Kar

and Maiti

, Multi-objective solid transportation problems with budget constraint in uncertain environment, International Journal of Systems Science 45(8), 1668–1682.

25.

Dantzig

and Ramser

, The truck dispatching problem, Management Science 6 (1959), 80–91.

26.

Pishvaee

M.S.

, Torabi

S.A.

and Razmi

, Credibility-based fuzzy mathematical programming model for green logistics design under uncertainty, Comput Ind Eng 62 (2012), 624–632.

27.

Goldberg

and Lingle

, Alleles, loci, and the traveling salesman problem, In Proceedings of the 1st international conference on genetic algorithms and their applications (1985), 154–159. Lawrence Erlbaum Associates, Publishers.

28.

Deep

and Mebrahtu

, Variant of partially mapped crossover for the Travelling Salesman problems, International Journal of Combinatorial Optimization Problems and Informatics 3(1) (2012), 47–69. ISSN: 2007-1558.

29.

Srinivas

and Patnaik

L.M.

, Adaptive probabilities of crossover and mutation in genetic algorithms, IEEE Transactions on Systems Man and Cybernetics 24(4) (1994), 656–667.

30.

Passino

K.M.

, Biomimicry of Bacterial Foraging for Distributed Optimization and Control, IEEE Control Systems Magazine (2002), 52–67.

31.

Niu