1. INTRODUCTION
Location and establishment of facilities are some of the critical issues in the primary stages of industrial system design. The optimal industrial location has always been significant from the perspective of geographers and economists (Shannaq et al., 2013;Golubev et al., 2021;Man et al., 2019). In fact, industrial and service institutions deal with various issues to determine the location of factory or service units, establish equipment and related sections, establish offices across the city and determine distribution centers (ALSoud et al., 2021;Bolano et al., 2020;Mirfani et al., 2021;Heidary et al., 2017). Some modes of discrete location problems and research allocation have been assessed in the literature, including singlefacility location problem, multifacility location problems, and locationallocation problems. In a singlefacility location problem, the goal is to find a new facility location so that there is a minimum total weighted distance between the new and existing facilities (Alwreikat and Rjoub, 2020;Shinkevich et al., 2021). A few simple examples of singlefacility problems include locating a hospital, a fire station, or a library in an urban area and locating a new airport to provide military service (Miandoabchi et al., 2013;Ab Yajid, 2020).
Moreover, multifacility location problems seek to find the optimal location for more than one new facility based on the existing facilities. According to Golubev et al. (2021), this type of problem has several uses, including establishing several warehouses to serve a specific number of regions. In fact, the singlefacility location problem is a certain mode of the multifacility location problem. In location allocation problems, the goal is to find the optimal location for facilities in a way that the cost of transportation from the facilities to customers is minimized. Therefore, a number of optimal facilities must be determined in this type of problem to meet customers’ demand in appropriate locations. In such locationallocation problems, the article deals with the limitation of facility capacity, which prevents the desired facility to meet all demands of a customer point. Accordingly, the entire demand of a customer may be divided between different facilities, each being responsible for meeting only one part of the customer demand. Miandoabchi and Farahani (2011) addressed the problem of designing street directions and lane additions in urban road networks based on reserve capacity. In this problem, the problem was to find the optimum configuration of street directions and twoway street lane allocations and the optimum selection of street lane addition projects to maximize the reserve capacity of the network. They proposed a twolevel mathematical problem, which was solved using a hybrid genetic algorithm and an evolutionary simulated annealing algorithm. These scholars (Shivaprasad, 2021;Astanakulov et al., 2021) then considered the same problem for a threeobjective mode. They presented a mixedinteger programming model and presented three metaheuristic methods (hybrid genetic algorithm, evolutionary simulated annealing algorithm, and artificial bee colony [ABC] algorithm) to solve the model. ZanjiraniFarahani et al. (2013) also reviewed urban transportation network design problems.
Franceschetti et al. (2013) evaluated a timedependent vehicle routing problem, where travel was considered two consecutive periods of traffic congestion and waiting idly. They proposed an integer linear programming formulation. Wang et al. (2017) assessed holiday travel behavior dynamics with integrated multimodal travel information (IMTI) usage. They examined a twoway relationship between holiday travel biography and IMTI usage biography over the life course after controlling for the effects of residential, household structure, employment/education, and car ownership biographies. Based on the webbased life history survey data, statistical characteristics of mobilities in each life biography were first analyzed. Afterward, various randomeffects ordered logistic models were established to investigate the biographical interdependencies from three aspects: intradomain interdependency, interdomain interdependency, and outerdomain interdependency. According to the results, longterm state dependence for different life domains was much more obvious over the life course when explaining holiday travel behavior dynamics and IMTI usage mobilities. To the best of our knowledge, vehicle routing problems are defined on a network in which customers’ arcs and locations are given. According to Huang et al. (2017), arcs somehow represent the distances or expected travel time derived from the underlying road network. They also revealed that when executed, the quality of the solutions obtained from the vehicle routing problem depends largely on the quality of the road network representation. In the mentioned research, path selection in the road network was defined as an integrated decision in the timedependent vehicle routing problem denoted as path flexibility (PF). This means that any arc between two customer nodes has multiple corresponding paths in the road network. Therefore, the decisions to make involve the routing decision and the path selection decision depending upon the departure time at the customers and the congestion levels in the relevant road network. The corresponding routing problem is a timedependent vehicle routing problem with path flexibility (TDVRPPF). Models of this problem are formulated under definite and probable traffic conditions. Considering PF has extremely reduced costs and fuel consumption, having the two concepts of timedependent travel and PF has led to an extensive presentation of possibilities and dynamism in travel time and PF of servers as a natural resource under uncertainty. The route estimation method has been used to generate nearoptimal solutions for TDVRPPF under possible congestion problems. Sun et al. (2018) introduced the timedependent capacitated profitable tour problem with time windows and precedence constraints. This problem concerned determining a tour and its departure time at the depot that maximized the collected profit minus the total travel cost (measured by total travel time). To deal with road congestion, travel times are considered to be timedependent. Shiripour et al. (2017) presented a capacitated locationmultiallocationrouting model for a transportation network with travel times between the nodes represented by links on the network. The concept of multiallocation arises from the possibility of allocating the population in a demand node to more than one server node. In normal conditions, travel time between two nodes is a fixed value. However, since the flow of population in a link can affect the travel time, the impact of the population flow on link time was considered to be simultaneous in the foregoing research. By doing so, the distribution of the population over the network directly impacted the travel link times. The mentioned researchers assumed that all links were twoway and capacities of the server nodes and arcs for accepting the population were limited. The main objective of the aforementioned study was to detect optimal locations of server nodes, optimal allocation of the population in demand nodes to the servers, and optimal allocation of the population of the nodes to different routes to reach the assigned servers so that total transportation time could be minimized. In another study, Shiripour et al. (2016) evaluated a capacitated locationmulti allocationrouting problem with intelligent probabilistic travel times. According to these scholars, the concept of intelligent probabilistic travel times incurred two issues: (1) consideration of some random factors in computing the travel times and (2) impact of the traveling population on these random factors simultaneously.
They considered three random factors: the time spent in traffic, the number of accidents, and road failures. It was assumed that server nodes and arcs have limited capacities for accepting the population. After proposing a function for computing the intelligent probabilistic travel time, we first formulate the problem as a mixedinteger nonlinear programming model and then suitably transform it into a mixedinteger linear programming one. In the mentioned study, the main objective was to determine appropriate server locations among the candidate locations, allocate the existing population in each demand node to server locations, and find the movement path of each member to reach its corresponding server with respect to the simultaneous change of the probabilistic travel times so that the expected total transportation time was minimized.
2. STATEMENT OF THE PROBLEM
Transportation planning is one of the fundamental issues considered in various science fields, including industrial engineering and civil engineering. As mentioned before, the ultimate goal is to achieve a higher service level and better security, save energy, gain economic growth and increase availability. Our problem was a transportation network encompassing several nodes and links. Imagine that G (N, L) is a transportation network, in which N is the set of nodes and L is the set of links. In this network, each node shows a region, department, city, or location that includes a population of individuals, vehicles, or products. The population existing in these locations were considered as travel applicants. In this communication network, all links are twoways and the time required to travel from one node to another node is determined. It is notable that the forward travel time from one node to another may not be equal to its backward time. Since the presence of each person or vehicle in each arc affects the forward and backward travel time from one node to another and given the fact that other people or vehicles are affected by this presence, it seems necessary to evaluate the impact of travels on travel duration along each arc simultaneously. However, since each link has the ability to pass a certain number of people whether in forward or backward travel time in real conditions, a certain capacity was considered for each arc in the proposed model so that, if the population that travels from different demand nodes to one arc exceeds the arc’s capacity, it would block the arc, prevent the presence of additional population in the arc and necessitate the selection of another arc for reaching the destination. Since each link is twoway, the capacity of its forward travel might be different from its backward travel. In addition, the traffic flow on one side of a link might be different from the other side. Therefore, blockage of one side of the link is not interpreted as the blockage of the other side. Moreover, it is logical to assume that the increased travel time of each arc depends on the population passing it. In other words, the more the capacity of the arc is filled, the higher the chance of an increase in the travel time. Furthermore, since each service node is able to serve a limited population, the capacity of all service nodes existing in the network was also capacitated in the present research. In real conditions, the population existing in a node is not necessarily allocated to the same region to receive the service. Therefore, it is possible to allocate different service areas to the population within each area. In this mode, it is possible that a percentage of the population of an applicant node such as A is allocated to a service node such as B, and another percentage is allocated to another service node such as C. Moreover, it should be specified that the percentage of population traveling from node A to node B takes which paths to reach this service node. In addition, the percentage of population traveling from node A to node C takes which paths to reach the C service node. In fact, there a multiple allocation in this mode, and travel applicants existing in a node might be allocated to more than one service location. Having information such as the length of the route and vehicle speed is not enough to predict the travel time (Alwreikat and Rjoub, 2020). This is mainly due to the fact that we may encounter unexpected events such as road accidents and breakdowns on the path, each could affect the travel time. Therefore, focusing on these issues can get us extremely closer to the actual travel time. Here, it is assumed that the travel time of each arc on the network is not certain and the time of each arc depends on several parameters affecting the arc. In each route (arc), parameters such as the number of accidents and vehicle breakdowns can be affected by two factors; one is the damage caused by accident or breakdown that each route (arc) encounters in its normal condition and based on geographical conditions of the region. The other one is the number of accidents or breakdowns that occur based on the population that travels the route (arc). As such, the values of parameters are affected depending on what type of population travels the arc, and this sensitivity of parameters to changes affects the time of the related arc. Therefore, when the time of arcs is not certain, the waiting time of each arc should be considered as the criterion for decision making in the objective function. Therefore, our proposed model should dynamically update the values of the parameters for each person entering each link so that the optimal service nodes could be selected and the optimal allocation of people or vehicles to these nodes and the optimal traffic route could be carried out. Overall, our proposed model sought to find service location(s) and allocate people or vehicles existing in each area of these places and determine the path of transportation of each person or vehicle in each region to reach the related service point while considering the direct effect of random factors and population traffic on random time travels of each arc in a way that all three objectives below are met simultaneously: 1) minimizing the overall transportation time of the network, 2) maximize travel attractiveness for travel applicants, and 3) balance travel applicants’ allocation to each service region to prevent irregular allocation of travel applicants to these regions. In recent studies (Shiripour et al., 2016;Qazani et al., 2021), only the first objective was considered in location, allocation, and routing decisions. Meanwhile, considering two new objectives of maximization of travel attraction and balancing the allocations leads to more real and rational decisions.
Research hypotheses:

1) Travel time in each arc is uncertain.

2) The service providers are finitely capacitated regarding the acceptance of people or vehicles.

3) The arcs are finitely capacitated in terms of the acceptance of people or vehicles.

4) The population existing in each node is not necessarily allocated to one service point to receive services.

5) The maximum number of service locations required is specified.

6) Each node existing in the network can provide services, meaning that there are candidate locations (for choosing service points) per the number of nodes existing in the network.

7) The attractiveness of each arc is determined.
3. PROPOSED MODEL
The main objective of the present study was to select optimal location(s) to provide services, optimally allocate people or vehicles existing in each region to these locations and determine the optimal transportation route of each person or vehicle in each region to reach the related service points while considering the direct effect of random factors and population traffic on possible travel times along each arc in a way that the total transportation time of the system is minimized, the overall attractiveness of different routes is maximized and the allocation of people to various service points is balanced. In fact, we attempted to study the effect of population traffic and some random factors (e.g., road accidents) on the calculation of actual travel time and the decisionmaking process. Evidently, travel time management and making the least increase in travel time will be associated with reduced fuel consumption and decreased travel time, which are among the secondary goals of the present study. Looking around, we come across different networks that pursue goals similar to the proposed issue, such as tourism transportation networks, network of transportation of the injured, freight transportation networks, and student transportation networks. For instance, in the tourism transportation network, we know that deciding about choosing from among the various cities (nodes) that exist at the level of a country (network) as places to visit (tourism) in a particular period of time and choosing a location or route as the place of recreation by a person in each city (node) is extremely important. This is mainly because in a certain period of time, there may be many people on different routes, and usually each person considers the shortest route as their travel path, which increases people’s travel time and blocks some routes. Therefore, it is extremely valuable to provide a model that controls the situation by taking the impact of the current population and random factors on travel time into account and provides each person with the best location and the most appropriate route in terms of time. In this situation, all people can reach their desired destination with the least delay and know about the actual time of their travel from the beginning. Moreover, we considered two other objectives, including balancing the allocations for balanced allocation of people to service locations and maximizing travel attractiveness for those who, in addition to the goal of minimizing time, also care about the attractiveness of the travel route. Indexes, parameters, and variables used in the model are presented below:
Indexes:
Parameters:

Cap^{k} : The capacity of the kth vehicle

d_{i} : The demand of the ith node in

c_{ij} : Distance from the ith node to the jth node

γ_{2}  γ_{1}: Coefficients related to the coordination of goals and their importance

pen : Penalty for each unit of unmet demand

bigM : A large enough number
Variables:

${x}_{ij}^{k}$ : Binary variable, (1) if the ij path is traveled by the kth vehicle; otherwise, (0).

${\alpha}_{i}^{k}$ : Percentage of the ith node demand met by the kth vehicle

D_{i} : Accumulated unmet demand of the ith node

y_{k} : Binary variable. (1) if the kth vehicle is used; otherwise, (0).

${u}_{i}^{k}$ : Auxiliary variable, corresponding to the ith node by the kth vehicle
The objectives are, as follows: (Equation 1)
Constraints (2) show that if a vehicle is used, it can pass the path between two points. Constraints (3) ensure that the vehicle definitely leaves each node it enters. Each vehicle used must exit the depot, which is shown by constraints (4) (assuming that the depot is node 1). Constraints 57 are related to the elimination of subtour and creating rotation. Constraints (8) are related to the capacity of the vehicles. Constraints (9) simultaneously ensure that each satisfaction rate is a value between zero and one and it forces the rate to be valued only in case of the arrival of a vehicle to the desired node. Constraints (10) show that the demand points are met by the vehicles. In addition, constraints (11) determine the level of cumulative accumulated demand while constraints (12) show the binary nature of variable x, and constraints (13) show the nonnegative decision variables.
4. SOLUTION METHOD
In this research, the mathematical model is solved using the cuckoo optimization algorithm (COA), which is a novel algorithm that has less attracted the attention of researchers in the routing field. In order to develop the algorithm and improve its performance, attempts have been made to enhance the mechanisms of this algorithm to improve its performance. One of the weaknesses of the algorithm is the use of a simple heuristic in cuckoo clustering. As observed in the COA steps, the kmeans method is applied, where k number of clusters (members) are first selected randomly from n members as centers of clusters. Afterwards, nk remaining members are allocated to the nearest cluster. In the next stage, all members of cluster centers are recalculated and allocated to clusters based on new centers. This continues as long as the centers of the clusters remain constant. Given the heuristic nature of the kmeans method, the method has unfavorable quality in large scales. Therefore, attempts are made to use metaheuristics for cuckoo clustering. Since the clustering operation is carried out in each COA iteration, the clustering method must be very rapid and present the best clustering in the shortest time possible. In this respect, simulated annealing (SA) is applied for cuckoo clustering. The SA is the thermal process of achieving lowenergy state in a solid. This method moves towards the optimal solution step by step by creating and evaluating consecutive solutions. To move, a new neighborhood is randomly created and assessed. In this method, we evaluate the points near the given point in the search space. The new point will be selected as the new point in the search space in case of being identified as a better point (decreasing the objective function). If it is recognized as worse (increasing the objective function), it is selected based on a probability function. Simply put, search always occurs to reduce the value of the objective function in a cost function minimization process. However, the process sometimes leads to the increase of the price function. Usually, a criterion called the metropolis acceptance criterion is used to accept the next point:
The control parameter in the steel hydration simulation plays the same role as temperature in the physical phenomenon. First, the particle (which represents the current point in the search space) is shown with a very large amount of energy (which indicates the high value of the control parameter C). This high energy allows the particle to escape a local minimization. As the search continues, particle’s energy level decreases (C is reduced). Ultimately, the search will tend to the overall minimum. However, it should be noted that the possibility of algorithm escape from local minimization reduces at low temperature. Therefore, the higher the initial energy, the higher the possibility of reaching a general minimization (15).
This algorithm is primarily selected due to its ability to generate only one solution in each iteration and its extremely simple mechanism. In addition, it can present the best clustering in a very short period. In order to design and implement the proposed algorithm, the strand of solution (the coding system) defined in this study must be introduced. The solution strand of the main problem used in the COA is a vector of the length of the number of customers plus the number of vehicles (N+k). each cell shows a value between zero and one. The first part of this strand shows the solution of customers’ priority to visit, and the second part of the strand shows the percentage of customers for each vehicle. For instance, if N=5 and K=2. An example of this solution strand is shown in Figure 1.
Based on the following equations, the first and second vehicles visit two and three customers, respectively.
In other words, the first vehicle respectively visits customers 3 and 5, whereas the second vehicle visits customers 4, 1 and 2, respectively. Now, if the capacity of the vehicle is more than the customer demand in visiting customers, the unmet amount is reported in the form of the amount of unmet demand and affects the objective function. As mentioned before, clustering is carried out by SA, the solution strand of which is a vector with values between 0 and 1 and the length of the number of cuckoos (M). Accordingly, if the number of desired strands is W, and if the value of each cell is between 0 and $\frac{1}{W}$, it will be allocated to the first strand, and if it is between $\frac{1}{W}$ and $\frac{2}{W}$, it will be allocated to the second strand. This will continue for all clusters.
5. NUMERICAL RESULTS
In this section, the validity and behavior of the proposed model were assessed by presenting an example network, which included 9 nodes and 12 links. In this model, attempts were made to locate a maximum of three service nodes out of nine candidate notes. The population was assumed to be 10 in all applicant nodes and the node capacity to accommodate travel applicants was assumed to be 55. Figure 2 shows the network corresponding to this example.
As observed, when only the first objective in the mathematical model was considered as the goal function, two nodes (nodes 1 and 9) were selected as service nodes, which provided services to 55 and 36 travel applicants, respectively. The entire population in nodes 2, 3, 4, 5, and 40% of the population in node 8 were allocated to the service node 1 and the whole population in nodes 6 and 7 and 60% of the population in node 8 were allocated to the service node 9. Moreover, when the mathematical model considered all three objective functions introduced in the proposed problem simultaneously, 3 nodes (nodes 1, 4, and 9) were selected as service nodes, each of which provided services to 30 travel applicants. According to Figure 2, all populations in nodes 2 and 3 were allocated to service node 1, all populations in nodes 5 and 6 were allocated to service node 4, and all populations in nodes 7 and 8 were allocated to service node 9. In this decisionmaking process, a balance was evident in the allocation of travel applicants to service nodes.
5.1 COA Performance Evaluation
To evaluate the COA, the algorithm is coded in Matlab R2014. Afterwards, 10 problems are generated with different sizes. Notably, demand is generated at each point of the discrete uniform distribution with a boundary below 1000 and a boundary above 2000. Each point is located in a twodimensional space with a length of 150 and a width of 150, and the cost of travel between the two points corresponds to the distance between the two points. In addition, the capacity of vehicles is randomly considered in the range of 800015000. Evidently, vehicle capacity is selected such that the total vehicle capacity is higher than the total demand points. The results obtained from the exact solution of examples generated by GAMS are compared to the COA results. Given the extremely high solution time of GAMS for solving largescale problems, a time limit of 3600 seconds or one hour is considered in this regard. It is notable that in case of a need for more than one hour to solve a problem in GAMS, the software will present a justified (not necessarily optimal) solution, and the program is completed. Table 1 presents a summary of results of comparison of GAMS results with COA results.
As shown in the Table 1, GAMS software cannot find the optimal solution to the last three problems in less than one hour. For these examples, however, the COA is able to find better solution, in relatively less time, compared to GAMS. The solution time of both methods is shown in Figure 3.
As observed, the time of exact solution of this problem with GAMS has an exponential increase, which in creases dramatically with increasing the dimensions of the problem. In the COA, however, the increase rate of the solution time is very small.
Figure 4 shows the diagram of the objective function obtained from the exact solution and the metaheuristic solution of the model.
Figure 4 shows that in different problems, there is no significant difference between COA and exact solution. This algorithm has a 0.05% error in all solved problems with GAMS, which shows the proper performance of the COA in finding the optimal solution to the problem.
6. CONCLUSION
The present study evaluated a locationallocationrouting problem in the transportation network, where travel time along each arc was considered a part of the decision making process. In this research, we introduced a threeobjective model for the problem while considering the travel times that depended on the population existing on the route. The goal of the location of service points was to allocate people or vehicles existing at each point to these centers and determine the route of transportation of each individual or vehicle in each point to reach the related service providers while considering the direct impact of random factors and population traffic on the randomized travel time of each arc, in a way that it would simultaneously minimize the overall transportation time of the network, maximize the travel time for applicants and allocate travel applicants to each service point in a balanced form. In this study, a general form of a threeobjective mathematical programming model was formulated. Moreover, a sample example was assessed to confirm the validity and behavior of the proposed model.