Methodology to Solve the Combination of the Generalized Assignment Problem and the Vehicle Routing Problem : A Case Study in Drug and Medical Instrument Sales and Service

This article presents algorithms for solving a special case of the vehicle routing problem (VRP). We define our proposed problem of a special VRP case as a combination of two hard problems: the generalized assignment and the vehicle routing problem. The different evolution (DE) algorithm is used to solve the problem. The recombination process of the original DE is modified by adding two more sets of vectors—best vector and random vector—and using two other sets—target vector and trial vector. The linear probability formula is proposed to potentially use one out of the four sets of vectors. This is called the modified DE (MDE) algorithm. Two local searches are integrated into the MDE algorithm: exchange and insert. These procedures create a DE and MDE that use (1) no local search techniques, (2) two local search techniques, (3) only the exchange procedure, and (4) only the insert procedure. This generates four DE algorithms and four MDE algorithms. The proposed methods are tested with 15 tested instances and one case study. The current procedure is compared with all proposed heuristics. The computational result shows that, in the case study, the best DE algorithm (DE-4) has a 1.6% better solution than that of the current practice, whereas the MDE algorithm is 8.2% better. The MDE algorithm that uses the same local search as the DE algorithms generates a maximum 5.814% better solution than that of the DE algorithms.


Introduction
Marketing techniques are a major factor that can affect sale volume due to the high competitiveness in the real global market.Many techniques have been used in various types of markets to sell different products.Drug and medical instrument sales and service (DMIS) are products that require special marketing techniques.The marketing technique used in DMIS is one where the company sends a seller to meet the customers in their working place to entice, convince, and provide new information to the customers.This can be interpreted as direct sales as the seller attends the customer's location to sell the product.The DMIS customers include hospitals, medical clinics, beauty clinics, drug stores, etc.The profit from DMIS is the number of products the company can sell and minus it the cost incurred.The cost of DMIS comprises production, labor, travel, administrative, accommodation, and managerial costs etc.In this study, the production, administrative, overhead, and managerial costs are approximated as the operation cost which is shown in cost per unit of the product and we assume that other cost terms are not significant and are included in the overhead cost.A certain number of sellers will be sent to sell and service a group of customers.In assigning the different customers to the seller, the sales may change as sellers have different skill levels (experience) and the customers have different stimulation probability (probability of the customer increasing the volume of product when the sellers convince them).When assigning a number of customers to be served by one seller, there is the potential for a generalized assignment problem (GAP).This may occur when the seller constructs their travel plan.The road that connects the customers has different conditions, which affect the driving speed.The driving speed affects both the fuel consumption and the traveling time.This problem can be interpreted as a GAP because we assign the customers to a group of sellers and the traveling times are affected.In one day, the sellers have only eight hours of traveling.Therefore, this problem can also be categorized as a GAP.
The GAP problem is one of how to assign the customers to sellers.When assigning the different customers to the sellers, the traveling route, which has different traveling times (consuming different resources), is affected.Assigning a customer to be served by different sellers may generate different sales volumes, which affects the total profit obtained by the firm.Generally, in a GAP, only the resource consumption difference, when assigning different customers to a seller, is determined but in the proposed problem, the objective function is also affected.Therefore, the proposed heuristic is a special case of GAP.
Many studies have attempted to find a good solution to the traditional GAP.Several exact methods have been used, such as branch and bound (Ross and Soland 1975), the multiplier adjustment method (Fisher and Jaikumar 1981), and the branch and price method (Savelsbergh 1997).Heuristic and metaheuristic methods have been widely used to solve GAP due to requiring much less computational time than the exact methods, such as the differential evolution algorithm (Sethanan and Pitakaso 2016a), Simulated Annealing (Osman 1995), Tabu Search (Laguna et al. 1995;Dıaz and Fernández 2001), genetic algorithm (Chu and Beasley 1997), and bee algorithm (Özbakir et al. 2010).
After the customers have been assigned to a group of sellers, the sellers in that group create their traveling plan.They have to visit all assigned customers within five days, and some customers may require more than one visit during the trip.The sellers may sleep in the city where the last customer in their daily plan is located.This is interpreted as a vehicle routing problem (VRP).The goal of the VRP is to construct a route to minimize traveling distance.The road conditions that affect the driving speed are not considered traditionally.The driving speed affects the fuel consumption rate, which affects the total cost to visit all customers.In this study, the objective function was changed from minimizing the total distance to minimizing the cost to the firm to visit all customers.
VRP involves a defined number of customers and sellers.The solution to the VRP is the optimal traveling plan for each seller.After the seller finishes their route, they return to the depot or the head office.VRP was introduced by Dantzig and Ramser (Dantzig and Ramser 1959), and Lenstra and Kan (Lenstra and Kan 1981) proved that VRP is a Non Polynomial-hard problem (NP-hard problem).Various types of VRP have been consecutively proposed after Dantzig and Ramser (Dantzig and Ramser 1959), such as VRP with time windows, heterogeneous fleet, or Multi-Depot Heterogeneous Vehicle Routing Problem with Time Windows (Bettinelli et al. 2011;Xu et al. 2012).Braekers et al. (2015) provided an overview of the scenario and the problem's physical characteristics, extending the basic incapacitated VRP, which was most often considered in the reviewed articles.The special case of the VRP problem was proposed by Kaewman and Akararungruangkul (2018) and Akararungruangkul and Kaewman (2018), where one pickup point can be visited more than once and the road condition is considered as the objective function.
The proposed problem is one where the GAP and VRP are integrated.To the best of our knowledge, no study has addressed this problem.The VRP is usually integrated with the scheduling problem (Han et al. 2017;Zhan et al. 2015;Nair et al. 2018).The most relevant work was proposed by Hervert-Escobar et al. (2016), who presented an integrated approach to the assignment, scheduling, Adm. Sci.2019, 9, 3 3 of 21 and routing problems in a sales territory business plan.This paper determined the minimum number of sellers required to attend a set of customers located in a certain region, considering the weekly schedule plan of the visits and forming the optimal route.In our study, the customers are not yet grouped, and the sellers work in a group in which another seller can visit a customer instead of the usual seller when the usual seller is not free.The customers may be visited more than once; they can ask to be visited more than once in a period by a seller.Our routing phase of the problem is also different from that proposed by Hervert-Escobar et al. (2016).The proposed problem integrated the accommodation cost.Therefore, the decision maker must return to the head office as soon as possible so that the accommodation cost is minimized.Furthermore, the sellers' skill level, customers' stimulation level, road condition, and average driving speed, which affects the fuel consumption rate (Fuel Economy in Automobiles (Wikipedia 2018)), were introduced in our problem.Finally, the total profit was used instead of the total cost, as in other studies, as the introduced terms affect the sales volume or the costs.Overall, the problem is a very hard combinatorial optimization problem.
In developing the methodology to solve the VRP and GAP problem, many methods have been proposed, such as Simulated Annealing (Osman 1995), Tabu search (Laguna et al. 1995;Dıaz and Fernández 2001), BEE algorithm (Özbakir et al. 2010), the evolutionary algorithm (Weise et al. 2010), and the genetics algorithm (Liu et al. 2012).Recently, the differential evolution (DE) algorithms have been successfully applied to solve VRP, GAP, and other combinatorial optimization problems.
The DE algorithm applies population-based metaheuristics, which is one of the most powerful and interesting evolutionary algorithms.Generally, it has four common steps: (1) generate the initial solution; (2) the mutation process; (3) the recombination process, and ( 4) the selection process.DE has been successfully applied in several fields, such as production scheduling (Pitakaso 2015;Pitakaso and Sethanan 2015) and manufacturing problems (López Cruz et al. 2003).Dechampai et al. (2015) proposed a DE to solve the capacitated VRP with the flexibility of mixing pickup and delivery services and maximum duration of a route in the poultry industry.Akararungruangkul and Kaewman (2018) proposed a DE to solve the special case of VRP in a student pickup system.They also modify the mutation process by using a different formulae to increase the search capacity of DE.The computational result showed that the new formula outperforms the traditional DE formula.Sethanan and Pitakaso (2016b) improved DE by adding two more steps, the reincarnation and survival processes, to improve the intensification search of the DE.These steps were added after the selection process.Sethanan and Pitakaso (2016a) improved the DE algorithms by adding effective local searches to increase the search mechanism of DE for solving the generalized assignment problem.The use of different mutation and recombination processes and matching to improve the solution quality has shown that using different pairs of mutation and recombination processes produced different quality solutions (Sethanan and Pitakaso 2016b;Teoh et al. 2013).From the literature, DE is more effective if the local search is included, but this increases the computational time.DE has been modified by adding some more processes to the original version, such as recombination, reborn, and reincarnation processes, to obtain a better solution.These processes can improve the solution quality of the original DE by enhancing the search capability of the original system.
In this study, we attempted to enhance the search capacity of the DE by adding both diversification and/or exploration and intensification.We modified the recombination process by adding more choice to obtain new trial vectors.Originally, two types of vectors were used: target and mutant vectors.We added two more vectors: (1) best vector, which is used to increase the intensive search capability, and (2) random vector, which is used to increase the exploration capability.The chance of using these four types of vectors is controlled by the newly designed equation.
This article is organized as follows.In Section 2, we explain the case study of the proposed problem and in Section 3 outline our procedure.Section 4 presents the proposed heuristics.Section 5 provides the computational results, and Section 6 outlines our conclusions.

The Case Study and Problem Definition
The company we identified is one which sells drugs to hospitals and drug stores in Thailand.The company is located in Bangkok.There are 133 hospitals and 227 drug stores to be serviced.The customers are located in every province of Thailand.The number of sellers that the company hires is 50.The maximum number of sellers in each group cannot exceed 10.Thus, the minimum number of groups has to be at least equal to the number of regions in Thailand.In Thailand, there are five regions which are North, North-east, South, West, and Center-East.Therefore, 50 sellers have to be divided into, at least, five groups.The problem of the company is how to find the daily optimal route that generates the maximum profit for the company.The actors who are involved are customers, sellers, and the cities in which the customer is located.To make the case study are more understandable for the reader, we will use a small example to explain the case study as follows.As shown in Figures 1  and 2, there are 12 customers and three sellers, as an example of the DMIS problem.The minimum demand, the probability of the customer buying more than the minimum demand (stimulation level), the number of visits required, and the order in a particular week are shown in Table 1.Notably, it is possible that a customer needs more than one visit during one planning period.If so, this customer cannot be visited on two consecutive days; at least one day is needed to allow the customer time to consider before meeting the seller again.The unit's selling price is 120 baht, and the total production cost per unit is 40 baht.Therefore, the profit margin of this product is 80 baht per unit.

The Case Study and Problem Definition
The company we identified is one which sells drugs to hospitals and drug stores in Thailand.The company is located in Bangkok.There are 133 hospitals and 227 drug stores to be serviced.The customers are located in every province of Thailand.The number of sellers that the company hires is 50.The maximum number of sellers in each group cannot exceed 10.Thus, the minimum number of groups has to be at least equal to the number of regions in Thailand.In Thailand, there are five regions which are North, North-east, South, West, and Center-East.Therefore, 50 sellers have to be divided into, at least, five groups.The problem of the company is how to find the daily optimal route that generates the maximum profit for the company.The actors who are involved are customers, sellers, and the cities in which the customer is located.To make the case study are more understandable for the reader, we will use a small example to explain the case study as follows.As shown in Figures 1  and 2, there are 12 customers and three sellers, as an example of the DMIS problem.The minimum demand, the probability of the customer buying more than the minimum demand (stimulation level), the number of visits required, and the order in a particular week are shown in Table 1.Notably, it is possible that a customer needs more than one visit during one planning period.If so, this customer cannot be visited on two consecutive days; at least one day is needed to allow the customer time to consider before meeting the seller again.The unit's selling price is 120 baht, and the total production cost per unit is 40 baht.Therefore, the profit margin of this product is 80 baht per unit.

The Case Study and Problem Definition
The company we identified is one which sells drugs to hospitals and drug stores in Thailand.The company is located in Bangkok.There are 133 hospitals and 227 drug stores to be serviced.The customers are located in every province of Thailand.The number of sellers that the company hires is 50.The maximum number of sellers in each group cannot exceed 10.Thus, the minimum number of groups has to be at least equal to the number of regions in Thailand.In Thailand, there are five regions which are North, North-east, South, West, and Center-East.Therefore, 50 sellers have to be divided into, at least, five groups.The problem of the company is how to find the daily optimal route that generates the maximum profit for the company.The actors who are involved are customers, sellers, and the cities in which the customer is located.To make the case study are more understandable for the reader, we will use a small example to explain the case study as follows.As shown in Figures 1  and 2, there are 12 customers and three sellers, as an example of the DMIS problem.The minimum demand, the probability of the customer buying more than the minimum demand (stimulation level), the number of visits required, and the order in a particular week are shown in Table 1.Notably, it is possible that a customer needs more than one visit during one planning period.If so, this customer cannot be visited on two consecutive days; at least one day is needed to allow the customer time to consider before meeting the seller again.The unit's selling price is 120 baht, and the total production cost per unit is 40 baht.Therefore, the profit margin of this product is 80 baht per unit.

Customer Head Office
Customer Head Office Hotels used   The distance between customers (road) and the average speed on that road, which depends on the condition of the road, are used to calculate the fuel consumption rate of the road.An example of the proposed problem is as follows.There are three sellers, sellers A, B and C, and the standard salary of a seller is 24,000 Baht; therefore, weekly wages are 6000 baht.The skill levels of A, B, and C are 1.2, 1.0, and 1.3, respectively.The weekly wages of the seller are calculated from the base wages multiplied by the skill level.Therefore, the weekly wage of A, B, and C are 7200, 6000, and 7800 baht, respectively.The required solution of the problem is shown in Table 2. Table 2 provides an example of the expected solution of the DMIS.The 12 customers are divided into 2 groups: the first group includes customers 1, 2, 3, 4, 5, 7, 9, 10, and 11, while the second group includes customers 6, 8, and 12.The first group is serviced by sellers A and B, while the second group of customer is serviced by seller C. The seller starts from customer one's area and returns to customer one as well.In a day, the city in which the last customer in that route is located is used as the rest area (accommodation cost incurred) of the seller before they start the next route the next day.
Table 2 shows that sellers A and B work five days while seller C works four days.The order collected is the minimum order for the customer.The expected order is calculated by the probability of increasing demand of that customer multiplied by the experience level of the seller and the minimum demand of that customer.For example, on day two, seller A has route seven to three, the skill level of seller A is 1.2, and the probability of increasing the demand of customer 3 is 0.12, and demand of customer 3 is 210; therefore, the expected order is 30.24units.The total product order is the order collected plus the expected order collected.Therefore, on day two, the total product order is 240.24 units.The result of time used, minimum demand collected, and expected increasing demand is shown in Table 2.
The total orders collected (including expected order) of sellers A, B, and C are 2080.76,1517.02, and 1273.25 units, respectively.Rounded up, these numbers are 2081, 1518, and 1274 units, respectively.Thus, the profit margin of this company will be 4872 × 80 = 389,760 baht.Regarding road travel, the fuel consumed for all trips is 201.474L. This was determined by the distance between each connection multiplied by the fuel consumption rate for that connection (road).For example, the distance between customers 7 and 3 is 150 km (given information), and the fuel consumption rate of this connection is 0.09.The fuel used to travel this connection is, therefore, 13.5 L. The fuel cost is 33 baht per liter; therefore, the fuel cost is 201.474 × 33 = 6648.642baht.
The base weekly wage of the seller is 6000 baht, and the income is this base multiplied by the skill level.Thus, sellers A, B, and C have weekly wages of 7200, 6000, and 7800 baht, respectively.Total weekly seller wages are 21,000 Baht.
The next expense of the company is the sellers' commission, which can be calculated from the total sales volume multiplied by the commission rate.If the commission rate is 5% of total sales volume and the total sales are 4872 units, the total income is 4872 × 120 × 0.05 = 29,232 Baht.
The last expense is the hotel cost.The number of nights the sellers have to sleep in a hotel occurs when the last customer that a particular seller visits in a day is not customer 1.From Table 2, there are 11 nights.The hotel cost is 550 Baht per night.Thus, the total hotel cost is 6050 baht.
The total expense comprises fuel cost, weekly wage of the seller, accommodation cost, and commission, and is 6648.642+ 21,000 + 29,232 + 6050 = 62,930.642Baht.The total profit of the company for this week is 389,760 − 62,930.642= 326,829.358Baht.
Figure 3 represents the income component of the company, which is the total sale units (TS) multiplied by the selling price (Pr).The total sale units comprises the minimum sale units, which is the minimum guaranteed order of the customer plus the expected sales units, which is the expected increased order quantity, dependent on the chance to increase the order of the customer and the seller experience.The total cost is subtracted from the income, which includes fuel cost, unit cost, hotel cost, weekly income, and the commission.The total cost component is shown in Figure 4, and the total cost calculation is shown in Figure 5.
of customer is serviced by seller C. The seller starts from customer one's area and returns to customer one as well.In a day, the city in which the last customer in that route is located is used as the rest area (accommodation cost incurred) of the seller before they start the next route the next day.
Table 2 shows that sellers A and B work five days while seller C works four days.The order collected is the minimum order for the customer.The expected order is calculated by the probability of increasing demand of that customer multiplied by the experience level of the seller and the minimum demand of that customer.For example, on day two, seller A has route seven to three, the skill level of seller A is 1.2, and the probability of increasing the demand of customer 3 is 0.12, and demand of customer 3 is 210; therefore, the expected order is 30.24units.The total product order is the order collected plus the expected order collected.Therefore, on day two, the total product order is 240.24 units.The result of time used, minimum demand collected, and expected increasing demand is shown in Table 2.
The total orders collected (including expected order) of sellers A, B, and C are 2080.76,1517.02, and 1273.25 units, respectively.Rounded up, these numbers are 2081, 1518, and 1274 units, respectively.Thus, the profit margin of this company will be 4872 × 80 = 389,760 baht.Regarding road travel, the fuel consumed for all trips is 201.474L. This was determined by the distance between each connection multiplied by the fuel consumption rate for that connection (road).For example, the distance between customers 7 and 3 is 150 km (given information), and the fuel consumption rate of this connection is 0.09.The fuel used to travel this connection is, therefore, 13.5 L. The fuel cost is 33 baht per liter; therefore, the fuel cost is 201.474 × 33 = 6648.642baht.
The base weekly wage of the seller is 6000 baht, and the income is this base multiplied by the skill level.Thus, sellers A, B, and C have weekly wages of 7200, 6000, and 7800 baht, respectively.Total weekly seller wages are 21,000 Baht.
The next expense of the company is the sellers' commission, which can be calculated from the total sales volume multiplied by the commission rate.If the commission rate is 5% of total sales volume and the total sales are 4872 units, the total income is 4872 × 120 × 0.05 = 29,232 Baht.
The last expense is the hotel cost.The number of nights the sellers have to sleep in a hotel occurs when the last customer that a particular seller visits in a day is not customer 1.From Table 2, there are 11 nights.The hotel cost is 550 Baht per night.Thus, the total hotel cost is 6050 baht.
The total expense comprises fuel cost, weekly wage of the seller, accommodation cost, and commission, and is 6648.642+ 21,000 +29,232 + 6050 = 62,930.642Baht.The total profit of the company for this week is 389,760 − 62,930.642= 326,829.358Baht.
Figure 3 represents the income component of the company, which is the total sale units (TS) multiplied by the selling price (Pr).The total sale units comprises the minimum sale units, which is the minimum guaranteed order of the customer plus the expected sales units, which is the expected increased order quantity, dependent on the chance to increase the order of the customer and the seller experience.The total cost is subtracted from the income, which includes fuel cost, unit cost, hotel cost, weekly income, and the commission.The total cost component is shown in Figure 4, and the total cost calculation is shown in Figure 5.

Current Practice Procedure (CPP)
The current practice procedure (CPP) is the method that the company is currently using to solve the seller scheduling problem (Figure 1).The CPP procedure is shown in Figure 6.

Day Seller
After the daily travel plan of each seller has been determined in Step 5, Step 6 is applied using the criteria of grouping the sellers.Seller B will pick seller C to join their group.The next step is rerouting the grouped sellers using nearest neighbor heuristics.The nearest neighbor is the heuristics where we select the next customer to visit after the current customer by selecting the closest customer to the current customer.The result is shown in Table 4.

Current Practice Procedure (CPP)
The current practice procedure (CPP) is the method that the company is currently using to solve the seller scheduling problem (Figure 1).The CPP procedure is shown in Figure 6.

Day Seller
After the daily travel plan of each seller has been determined in Step 5, Step 6 is applied using the criteria of grouping the sellers.Seller B will pick seller C to join their group.The next step is rerouting the grouped sellers using nearest neighbor heuristics.The nearest neighbor is the heuristics where we select the next customer to visit after the current customer by selecting the closest customer to the current customer.The result is shown in Table 4.

Current Practice Procedure (CPP)
The current practice procedure (CPP) is the method that the company is currently using to solve the seller scheduling problem (Figure 1).The CPP procedure is shown in Figure 6.
From the example given in Section 2, the result of applying the CPP method is shown in the following details.From Table 1, there are 19 visits, therefore, B = 19.Then, using Equation (1), the maximum number of visits for each seller is P = 19 3 = 6.33 ≈ 7. When B and P are obtained, the next step of the CPP method is to find O S and O S .Therefore, O C = {1, 7, 6, 3, 2, 8, 9, 4, 12, 5, 10, 11}, due to having the following demand, {200, 930, 450, 450, 400, 330, 300, 300, 290, 240, 200, 190}, and O S = {C, A, B} due to the following skill level {1.3, 1.2, 1.0}.When Step 5 is executed while maintaining conditions 1-5, the result is shown in Table 3.After the daily travel plan of each seller has been determined in Step 5, Step 6 is applied using the criteria of grouping the sellers.Seller B will pick seller C to join their group.The next step is re-routing the grouped sellers using nearest neighbor heuristics.The nearest neighbor is the heuristics where we select the next customer to visit after the current customer by selecting the closest customer to the current customer.The result is shown in Table 4.The total profit was calculated as shown in Section 2. The result is shown in Table 5.The total profit of this planning is 368,484.523Baht.We programmed the current practice procedure in C++.In the first iteration of the current practice that we programmed was exactly the same as Steps 1-8.In the next iterations, a random number [0, 1] was picked to multiply the demand of the customer and the skill level of the seller.Therefore, in each iteration, the solution changed because the order of the customer and seller have changed.This enables the current practice to be compared with the proposed heuristics, which will be explained in the next session.
Step 1: Calculate B, total number of times that all customers need to be visited.
Step 2: Calculate P, maximum visiting times for each seller using following Equation  =   when R is total number of seller Step 3: Find the order of customers O C which is the order of customers sorted from maximum demand to minimum demand.
Step 4: Find the order of sellers O S in regards to their skill level from maximum to minimum demand Step 5: Find the traveling plan for each seller as following way; 5.1 Make free schedule table of 5 working days of each seller 5.2 Assigned the customer into all working days by assign the customer which is in the list O C to the seller which is in order O S step by step.In the assigning of the customers to the sellers the following condition must be kept.
Condition 1: The customer needs more than 1 times to be visited, the customer must not be visit on two consecutive days.
Condition 2: Total working time must less than 10 h or 600 min per day.
Condition 3: For each day, the seller will sleep in the last city where the last customer is located.
Condition 4: When we assign one customer to a seller on a particular day, if the next assignment (the closest customer to the last customer in list O C ) means that seller has not enough time, we can pick the next customer in the order O C to replace it.
Condition 5: When the last customer that a particular seller has to service is visited, the seller must return to depot.The time including the travel back to depot must be less than 600 min.
Step 6: Group the sellers to reduce the working day.The criteria for choosing the sellers to group together is that the seller that has least time used will pick the seller that has greatest time used to join.
Step 7: Re-routing the travelling plan of all available groups is obtained in step (6) using nearest neighborhood heuristics Step 8: Calculate the total profit (using income and cost component as shown in Figure 3 and 4. The total profit was calculated as shown in Section 2. The result is shown in Table 5.The total profit of this planning is 368,484.523Baht.We programmed the current practice procedure in C++.In the first iteration of the current practice that we programmed was exactly the same as Steps 1-8.In the next iterations, a random number [0, 1] was picked to multiply the demand of the customer and the skill level of the seller.Therefore, in each iteration, the solution changed because the order of the customer and seller have changed.This enables the current practice to be compared with the proposed heuristics, which will be explained in the next session.

The Proposed Heuristic
Generally, a differential evolution algorithm (DE) comprises four general steps: (1) generate the initial solution; (2) implement the mutation process; (3) execute the recombination process, and ( 4) complete the selection process.The DE was developed in an attempt to increase the search performance of the original DE algorithm.
The proposed heuristics is called the modified differential evolution algorithm (MDE).Two sets of vectors were introduced to the recombination process: a set of best vectors (BV), and a set of random vectors (RV).Therefore, four sets of vectors were used to process the recombination of the vectors: (1) target vector (TV), ( 2) mutant vector (MV), (3) best vector, and ( 4) random vector.To execute the recombination process at each iteration, two values of crossover rate (CR) (CR1 and CR2) were used to determine which set of vectors would be used.The proposed algorithm is explained step-by-step in the following section.

Generate Initial Solution
In this section, two main procedures are explained: (1) the encoding and ( 2) the decoding methods.The encoding method is where we design and present the DE's solution, which is used to execute each step of the DE.As the DE mechanism was first designed to use with continuous optimization, but the proposed problem is a combinatorial optimization, a decoding method was needed to transform the DE solution into the proposed problem's solution.

Encoding Method
A vector that represents the problem comprises two entities: customer's and seller's vector, called target vectors.These two vectors are used to decode to obtain one solution of the proposed problem.The customer vectors are shown in Table 6.There are five customer vectors or number of populations of the customer's vector (NPc).Each vector has the size of 1 × C, where C is number of customers.A customer's vector includes C (in Table 6, C = 12) positions.The value in each position is randomly generated in the first iteration, whereas other iterations are obtained from the DE process.Another type of vector that was used to construct the solution for the proposed problem is the seller's vector, which is shown in Table 7.A seller's vector has the size of 1 × S and has NP vectors (NP is the number of populations of the seller's vector).This vector is also randomly generated.In Table 7, NPs equals 4 while S = 3.From Table 7, there are three sellers and five vectors.For example, vector 5 includes three positions: 0.65, 0.14, and 0.94.The solution of the proposed problem is determined by the decoding method using these two sets of vectors, which will be explained in Section 4.1.2.

Decoding Method
The decoding method was used to transform the vectors shown in Tables 6 and 7 into the proposed problem's solution, as shown in Figure 7. From the example in Section 2, using the procedure given in Figure 7, G is calculated using Equation ( 2) when TS = 3 and if MS = 2 when TS = total number of seller, and MS = maximum number of the sellers that can be assigned to each group.Therefore, G equals to 2 groups.
Step 2 provides the order of the seller (O S ).This order is obtained by sorting the sellers according to their value in the position of that vector.An example of finding the O S of vector 2 is shown in Table 8.In Step 3, using the information shown in Table 8, the seller in order O S is assigned to 2 groups of sellers.When MS = 2, group 1 comprises two sellers (sellers B and C) and group 2 has one member, which is seller A. In Step 4, order O C is obtained.The sorting is conducted by the value in the position of the vector.An example of the O C of vector 1 is shown in Table 9.In Step 5, use Equation (3) to calculate Nv.From the example, Nv = 19/3 = 6.33 ≈ 7. The next step (step 6) constructs the route for each seller.An example of the route construction is shown in Table 10.Step 2: Find   , the order of the seller obtained by sorting the seller according to their value in position of a seller's vector.
Step 3: Assign the sellers into the each group according to   while the maximum number of sellers in each group should not exceed MS.
Step 4: Find the order of customers O C with regard to their value in position of a customer's vector.
Step 5: Calculate Nv which is the average number of visits per seller using the following equation When TV = total number of visits required by all customers and NS is number of seller.
Step 6: Construct the route of each seller in each day during the planning period by using priority.
The priority rules are: 6.1.The customer that has the maximum number of visit will be assigned to the first seller in each group first.

The remaining traveling plan (routing) using the order of O C
When we construct the route the following conditions need to be satisfied:  From Table 10, we can observe that we expected to sell a total of 5102.28 ≈ 5103 units.The total fuel used was 206.638 (calculated from the total distance multiplied by the type of road).The cost calculation of the routes constructed from the previous steps is shown in Table 11.
Table 11.Cost calculation of the visiting plan shown in Table 10.

Income Expense Profit
Total income (1) 5103 × 120 612,360 Total production cost (2) 5103 × 40 = 204,120 Income (1)-( 2 The last step of the proposed procedure is applying a local search to the solution shown in Table 10.The first local search is the exchange algorithm (EA).EA involves picking two customers and exchanging their position in the route.For example, on day 3, group 1 follows the route 4-2-12 and the next day, 12-11-6.This route consumes 36.413L of fuel.If customers 2 and 12 are selected to be exchanged, the new route will be 4-12-2 and 2-11-6, and this route consumes 41.9 L of fuel.The fuel consumption was not reduced.Thus, this route is declined.We continued to perform EA until all customers were exchanged.The second local search is the insertion algorithm (IA).This algorithm was used to move a customer from a route into another route such as 4-2-12 and 12-11-6.If customer 2 is moved out from their route to another, for example, between 12-11, then the new routes will be 4-2 and 12-2-11-6.This insertion is needed to check the time limit, and all conditions of the assignment must be satisfied before the move will be permanently inserted.Finally, the insert is only carried out when the fuel usage is reduced.

Execute Mutation Process
Some target vectors, which are provided in Tables 6 and 7, were randomly selected to execute the mutation process.The result of the mutation process is called the mutant vector.Mutant vectors were produced using Equation ( 4).We first randomly selected three vectors out of the NPc/NPs target vectors.Customer vectors and seller vectors were treated separately.Let r 1 , r 2 , and r 3 denote the randomly selected vectors from the target vector (both sets of vectors).F is the scaling factor.In the proposed heuristics all F were set to 0.8; i is the vector number; j is the vector's position; G is the current iteration of the simulation; X r 1 ,j,G , X r 2 ,j,G , and X r 3 ,j,G are the values in the position of the target vector r 1 , r 2 , and r 3 , respectively; and V i,j,G+1 is the mutant vector generated for vector i, position j in iteration G + 1.
Each group of vectors (customer or seller) use the same controlled parameters, such as F and CR (recombination rate).
An example of how to use Equation ( 1) is explained as follows.We assume that F = 2.To obtain the first mutant vector we randomly selected three vectors (r 1 , r 2 , and r 3 ) which are given in Table 7.The selected target vectors are 4, 2, and 3.The first position of mutant vector #1 (V 1,1,G+1 ) is 0.79 + 2 × (0.61 − 0.95) = 0.11.This mechanism will be applied to all position (j = 1, 2, and 3).Although, the mutant vector 1 is 0.11, 1.37, and 0.48.All four mutant vectors (equals to number of NPs) will be executed in the same manner.

Execute Recombination Process
Originally, Equation ( 2) is used to construct the trial vector (U i,j,G ) using the information from the mutant vector (V i,j,G ) and target vector (X i,j,G ).
where rand i,j is the random number of position j of vector i.
The MDE in Equation ( 2) is rewritten as in Equaltion (3) where V BV i,j,G is the set of best vectors, V RV i,j,G is the set of random vectors, and CR1, CR2, and CR3 are predefined parameters.CR1, CR2, and CR3 control the moving of the search area of the algorithm.When CR1 is high, the solution is guided by the current best solution.If the space between CR1 and CR2 is high, the new exploration behavior is enhanced.CR3 level controls the use of the old target vector, or mutant vector.
The self-adaptive CR1, CR2, and CR3 are used to adjust the parameter values.The adaptive CR1, CR2, and CR3 are calculated as follows.Let CR2 be the center probability of selecting all types of vectors.Set minimum CR2 (B CR ) and maximum CR2 (M CR ).From the preliminary test, the best B CR is 0.2, and the best M CR is 0.8.CR2 is iteratively calculated using Equation ( 4), where M T is the maximum iteration, and G is the current iterations.
If B CR is 0.2 and M CR is 0.8, CR3 is calculated using Equation ( 5): M CR3 was set to the maximum value of CR3, which was 0.9.
Adm. Sci.2019, 9, 3 14 of 21 The using of best vector (BV) and random vector (RV) is controlled by CR1.CR1 is calculated using Equation ( 6): where G NB is the number of iterations when the new best solution is not found and Min CR1 is the minimum value of CR1 that is allowed.CR1 starts from CR2 2 and resets itself when the new best solution is found.

Perform Selection Process
The result of the selection process is the new target vector (X i,j,G+1 ), which uses Equation ( 7) to select the new target from the previous target vector or the trial vector: Equation ( 7) is used to select the best vector between the trial vector and the current target vector to be the new target vectors.
Steps 4.2 to 4.4 are iteratively executed.We concluded the framework of the MDE as a picture as shown in Figure 8.

Perform Selection Process
The result of the selection process is the new target vector ( , , ), which uses Equation ( 7) to select the new target from the previous target vector or the trial vector: Equation ( 7) is used to select the best vector between the trial vector and the current target vector to be the new target vectors.
Steps 4.2 to 4.4 are iteratively executed.We concluded the framework of the MDE as a picture as shown in Figure 8.

Computational Framework and Result
The proposed heuristics was coded in Dev C++ using PC Intel Core i3 CPU 3.70 GHz Ram DDR4 8 GB.Fifteen randomly generated data were tested to compare the proposed heuristics with the current practice procedure that the company was currently using.A real-world case study was solved as well (test instance 16).The details of each test instances are given in Table 12.The detail of the End Set NP, CR, F, NP (size of vector)

Begin
For G = 1 to Gmax when G = iterations and Gmax = Maximum iteration For N = 1 to NP Random Generate Target vector   1 , , and RV and update BV Produce Mutant Vector N (Mutation Process) using the following equaltion.


Using the following equaltions.

Computational Framework and Result
The proposed heuristics was coded in Dev C++ using PC Intel Core i3 CPU 3.70 GHz Ram DDR4 8 GB.Fifteen randomly generated data were tested to compare the proposed heuristics with the current practice procedure that the company was currently using.A real-world case study was solved as well (test instance 16).The details of each test instances are given in Table 12.The detail of the proposed heuristics is shown in Table 13.The simulation was executed five times with different stopping criteria depending on the aim of testing.The best solution out of the five runs was recorded and is reported in Tables 14-18.Note: # customer is number of customers, # seller is number of sellers, and # maxSeller is the maximum number of sellers in each group.
From Table 12, the number of customers started at 30 and increased to 230 customers.The case study had 350 customers.The number of sellers ranged from 4 to 50 in the case study, whereas the number of sellers in each group ranged from 2 to 10.We have eight proposed heuristics: DE-1, DE-2, DE-3, DE-4, MDE-1, MDE-2, MDE-3, and MDE-4, the detail of which is shown in Table 13.same local search , MDE outperformed DE by using a new recombination process to improve the quality of the original DE.
Table 14 provides the result of simulation when all proposed heuristics use the same number of iterations (1000 iterations) as the stopping criteria.Table 16 provides the result of the simulation when we set the computational time to 90 min.We aimed to determine the performance of all heuristics when using the same computational time.The statistical test has been executed using Wilcoxon Sign Rank test with 95% confident interval, and the computational result is shown in Table 17.17, the results correspond to the results shown in Table 15.All MDE algorithms outperformed the DE algorithms.Two local search procedures, exchange and insertion, were not significantly different in finding a solution, but when used together, they performed well.
Table 18 shows the percent difference in the heuristics, and the best solution found compared to all heuristics.From the results in Table 19, we plotted a graph to visualize the progression of the simulation result when using different computational times (Figure 9).From Figure 9, all MDE algorithms found a good solution faster than the original DE algorithms.The algorithms that used local search started with a better solution than the algorithms without local search.The other value that is interesting is the difference between the solution generated from DE and MDE algorithms when using 10 min and 120 min.The large gap in the solution given different computational time means the capability of searching for a new solution is better than the algorithms with a smaller gap.The gap of each method is shown in Figure 10.
result when using different computational times (Figure 9).From Figure 9, all MDE algorithms found a good solution faster than the original DE algorithms.The algorithms that used local search started with a better solution than the algorithms without local search.The other value that is interesting is the difference between the solution generated from DE and MDE algorithms when using 10 min and 120 min.The large gap in the solution given different computational time means the capability of searching for a new solution is better than the algorithms with a smaller gap.The gap of each method is shown in Figure 10.x 10000 Best solution found using different computational time Simulation time (minutes) Profit (Baht) From Figure 10, the algorithms with local search had a larger gap between the first and last solution.This means that the local search enhances the searching capability of the algorithm.The MDE algorithms had larger gaps than the DE algorithms.This means the new recombination formula improves the performance of the original DE for finding a better solution.

Conclusions and Outlook
This article presented algorithms to solve a special case of the vehicle routing problem.The special case of the VRP had the following attributes, which are different from the original and other variants of VRP: (1) One customer can be visited more than once and if more than one visit is required, the seller cannot visit that customer on two consecutive days.(2) The customers are grouped, and each group includes some sellers.The sellers in the group can share the customers.(3) On each day, the seller sleeps in the last city they visit.( 4) All sellers must return the head office when they finish their travel plan.( 5) When traveling from one customer to another customer, the road condition controls the average speed on that route, which affects the fuel consumption, which may affect the travel plan.
A different evolution (DE) algorithm was used to solve the problem.The recombination process was modified by adding two more vector sets: best vector (BV) and random vector (RV), to additionally use the set of mutant or target vectors.Linear probability was proposed to enable the use of four sets of vectors.Two local search techniques were added to the modified differential evolution algorithm (MDE).
The computational result showed that the MDE algorithms outperform the DE algorithms due From Figure 10, the algorithms with local search had a larger gap between the first and last solution.This means that the local search enhances the searching capability of the algorithm.The MDE algorithms had larger gaps than the DE algorithms.This means the new recombination formula improves the performance of the original DE for finding a better solution.

Conclusions and Outlook
This article presented algorithms to solve a special case of the vehicle routing problem.The special case of the VRP had the following attributes, which are different from the original and other variants of VRP: (1) One customer can be visited more than once and if more than one visit is required, the seller cannot visit that customer on two consecutive days.(2) The customers are grouped, and each group includes some sellers.The sellers in the group can share the customers.(3) On each day, the seller sleeps in the last city they visit.( 4) All sellers must return the head office when they finish their travel plan.( 5) When traveling from one customer to another customer, the road condition controls the average speed on that route, which affects the fuel consumption, which may affect the travel plan.
A different evolution (DE) algorithm was used to solve the problem.The recombination process was modified by adding two more vector sets: best vector (BV) and random vector (RV), to additionally Adm. Sci.2019, 9, 3 20 of 21 use the set of mutant or target vectors.Linear probability was proposed to enable the use of four sets of vectors.Two local search techniques were added to the modified differential evolution algorithm (MDE).
The computational result showed that the MDE algorithms outperform the DE algorithms due to different recombination process.The new recombination process enhances the intensification capability of the original DE.A good solution was guided by BV, and when it was stuck on the local optimal, the RV helped the heuristics to escape the local optimum.New mechanisms were used while we maintained the excellent search ability of the original DE because we still used the TV and MV in the MDE algorithm.
Two local searches were presented.The computational result shows that DE and MDE algorithms that use local search produce better solutions than algorithms without local search.Local searches that are exchange or insertion algorithms performed no differently in improving the solution quality, but when combined, the search capability was better than using only one local search.
In the case study, the best DE algorithm (DE-4) produced a solution 1.6% better than the current practice procedure, whereas the MDE algorithm solution was 8.2% better from the current practice procedure.Pairwise comparison of the different DE and MDE algorithms that did not use local search, used exchange, used insertion, and used both exchange and insertion algorithms produced values of 2.244%, 5.796%, 5.814%, and 6.446%, respectively, which means the recombination process used in MDE algorithms enhances the capability of the DE algorithms for all local search types.
Overall, we conclude that both local search and the new scheme of recombination work together well and increase the search ability of the original DE.This method may enable the invention of new recombination schemes for the original DE, which can increase the search capacity of the DE.An effective local search should increase the searching ability of the DE as well.These two themes are an interesting area of future research.
Notes: #C = Customer name, #O = Number of orders placed by a particular customer, #D = Demand of the order, and Prob = probability to increase the demand assumed to be linearly dependent on the seller skill.#S = Service time, which is the direct sell time (minutes).

Figure 9 .
Figure 9. Best solution found using various computational time.

Figure 10 .
Figure 10.Percent difference in the best solution found using 10 and 120 min of computational time.

Table 1 .
Details of the customers.

Table 2 .
The schedule and route planning of 12 customers and three sellers.

Table 3 .
Route results of Step 5.

Table 4 .
The new routing for each group of sellers.

Table 3 .
Route results of Step 5.

Table 4 .
The new routing for each group of sellers.

Table 3 .
Route results of Step 5.

Table 4 .
The new routing for each group of sellers.

Table 5 .
Total profit calculation in Thai Baht.
Note: # NPc is the vector label, which ranges from 1 to 5.

Table 8 .
Results of sorting the seller's vector.

Table 9 .
Results of the calculation of cumulative visiting times.

Table 10 .
Results of Step 4. Calculate G which is the number of group of the sellers/customers using the following Equation

Table 12 .
Characteristics of the test instances and the case study.

Table 16 .
Profit generated in 90 min using different methods."Best" Best is the best solution amongst all compared heuristics.

Table 17 .
Statistical test of results shown in Table16.

Table 19 .
Profit of the case study using different controlled simulation times in Thai Baht.