A New Heuristic Method for Transportation Network and Land Use Problem

Our paper deals with the Transportation Network and Land Use (TNLU) problem. It consists in finding, simultaneously, the best location of urban area activities, as well as of the road network design that may minimize the moving cost in the network, and the network costs. We propose a new mixed integer programming formulation of the problem, and a new heuristic method for the resolution of TNLU. Then, we give a methodology to find locations or relocations of some Dakar region amenities (home, shop, work and leisure places), that may reduce travel time or travel distance. The proposed methodology mixes multi-agent simulation with combinatorial optimization techniques; that is individual agent strategies versus global optimization using Geographical Information System. Numerical results which show the effectiveness of the method, and simulations based on the scenario of Dakar city are given.


Introduction
The Transportation Network and Land Use (TNLU) problem is a combination of transportation network optimization problem and land use planning or quadratic assignment problem.This model appear, for example, if we decided to find locations or relocations of some Dakar region amenities (home, shop, work, leisure places), that may reduce travel time or travel distance.By solving this problem we want to have the best layout of the city in the direction of the location of activities and roads linking these activities.A better configuration in the direction of the location of activities and transportation for a city would reduce congestion and pollution.We consider an activity (i) placed on a zone (k) and an important flow (cars) from other activities, in the direction of the activity (i).Suppose that there is traffic lights at intersections, thus flows will be subject to long waiting times before reaching activity (i), hence congestion.During this waiting, passengers of vehicles are subject to breathing harmful gases from cars.It would therefore make sense to ask whether the location of this activity is the best possible (reallocation of this activity or not).Since the 60s many authors are interested in these two problems separately.Indeed we can cite authors such as Scott (1969), Boyce, Farhi and Weischedel (1973), Hoang Hai Hoc (1973) and Billheimer (1970) for network design problems, Lawler(1963), Schlager (1965) and Gordon-McReynolds(1972), Maniezzo(1997), Xia and Yuan(2006), Huizhen Zhang(2010) and many others for the activities's location problems.Among the researchers who are interested in the combined problem we can cite Marc Los and more recently Lin and Feng (2003).They have considered a more complex model in the sense that their formulation allows the location of several activities in the same area.Marc Los gave an exact solving method and heuristic methods.The Los's model assumes the OD flow streams, once the activities have been fixed via the shortest paths in the sense of a fixed cost independent of the flow on the arcs (see Gueye, 2012).
In this paper, the heuristics that we propose are based on Shortest Path (SP), Location of Activity (LA) and Network Design (ND) problems.It gives results corresponding to the best affectation found, by solving a Quadratic Assignment Problem (QAP) which uses shortest path solutions for distance, and the network design (ND) problem associated to these affectations.We use heuristic methods (QAP-heur and ND-heur) for the resolution of QAP and ND problems.
After locating the n activities into the n locations, we develop a new methodology by which a very large set of relocation possibilities can be simulated, analyzed, and the "best" one can be found.The methodology was coded in a prototype software called Simulation-Driven Optimization (SDO), that originated from two projects DAMA (Balac et al., 2014).In this task, SDO uses the multi-agent simulator MATSim (see Balmer et al., 2009).In MATSim, the actors of the modelled system are the agents (i.e the city residents).The agents act according to given "realistic" rules.They try to perform some activities at different places and have learning capabilities.The overall traffic observed in the urban area emerges from the simulation as a consequence of individual agents behaviour, each pursuing his/her individual interests.MATsim basically needs three data to perform a simulation: the transportation network (network.xml),the amenities location (facilities.xml)and the initial agent plans (plans.xml).At the first MATSim iteration, each agent follows one or several possible initial plans contained in the agent plan file.Following a complete MATSim simulation, in SDO we adopt a global (or collective) view which contrasts with the individual behaviour of the agents in the simulation.Our problem is indeed slightly different: it aims at finding some suitable relocations to increase the global accessibility for a set of selected amenities.Let us remark that the MATSim simulations are operating on only one day, such as the global accessibility we seek to improve.So, to be pertinent, the simulated plans should be as representative as possible of what the agent do most frequently.
In the following Section 2, we give a mixed integer programming formulation of the problem.Then, in Section 3 we present one heuristic method, TNLU-heur, to solve it.The Section 4 is devoted to the computation of the lower bound for TNLU.We detail in Section 5 the simulation results of our algorithms.The Section 6 deals with SDO that calculates the optimal locations of activities based on the results of the simulations of MATSim.Finally, we conclude and give some perspectives in Section 7.

Formulation of the Model
We assume to have n locations and n activities to locate in these locations.Consider the graph G(N, Ar) such as N is the set of node consisting different areas (n) and Ar is the set of links joining the n locations in pairs.Let: A = (a ik ) : location cost matrix of the activity i on the area k.Q = (q i j ) : moving flow matrix per unit of times between the activities i and j.C = (c i j ) : construction cost matrix of link or road (i, j).D = (d i j ) : direct distance matrix between two locations i and j or distance of link (i, j).
To determine the locations of activities, we consider the decision variables, called assignment, x i j = 1 if the activity i is affected on the location j and 0 otherwise.To determine the configuration of the network, we need to introduce binary variables y i j = 1 if the link (i, j) is selected and 0 otherwise.Depending on the constructed network, the shortest path within the meaning of the fixed costs of building d i j associated with each arc (i, j) may vary.Then it is necessary to introduce variables of shortest path l kl (y) between k and l and variables z kl a which take value 1 if the shortest path between k and l goes by the link a = (i, j) and take 0 otherwise.It is also estimated that if the activity i is placed on the node k and activity j on the node l, then the cost of transportation of a vehicle from k to l is equal to the distance l kl (y).On the other hand, the locations of activities induce fixed costs a ik on the back of construction costs of roads c i j for the network construction (for more details see Gueye, 2012).
The sum of the location, transportation and construction costs is expressed mathematically by the following function: This is a cubic integer model.The "Optimal" network and activity locations are obtained by minimizing the above function under a set of constraints described below.
where m is a intermediate node different to k and l taken respectively as origin and destination of the considered path, δ The constraints (1) and (2) represent assignment constraints ensuring that only one activity is localized in one area only.
For the shortest paths, are used the flow variables z kl a where are associated the constraints of flow conservation.Constraint (5) can give to variable l kl value of the shortest path from k to l.The constraint (6) ensures that only are taken into account in the calculation of shortest paths, the link actually constructed.
Thus, in our model the variables are mainly x ik , y i j and l kl , involved in the objective function, and z kl a which occurs only in the constraints.This model is an extension of Koopmans-Beckmann's problem.In addition, our problem, which generalized quadratic assignment problem1 is NP-complete (Garey et al., 1979;Sahni et al., 1976).Indeed, if we fix all the variables y i j then the problem becomes of finding activities assignments minimizing the sum of travel costs.This is equivalent to solve a quadratic assignment problem.In sequel, we can reformulate the model by replacing the arc a by (i, j).

A Decomposition Method
We propose in this section to solve three sub problems for TNLU, which are: • Shortest path problems; • A location of activity problem; • A network design problem.
It starts with Shortest Path problems (SP) which solutions are the shortest path between the different locations and denoted γ kl .Then the solutions of shortest path problems are introduced into the Location of Activity problem (LA), which returns as solution (x * ik ) the allocation of different activities in different locations.Finally, x * ik are used in the Network Design problem (ND), which returns as solution the shortest path and roads built, respectively denoted l * kl and y * i j .So (x * ik ), (l * kl ) and (y * i j ) are solution of the decomposition method of TNLU.The TNLU-heur (Figure 1) is a decomposition method of TNLU which uses heuristic methods (QAP-heur and ND-heur in section 3.2 and 3.3) for QAP and ND problem.

Shortest Path Problems (SP)
This section contains the resolution of the SP problems which is formulated as follows: Thus we solve n 2 − n SP problems using IBM CPLEX.Finally, at the end of the iterations (of these SP problems), we

Activity Location Problem (AL)
For the LA problem we use solutions of previous problems(γ kl ), and we formulate the problem as following : If we solve the Quadratic Assignment Problem (QAP), we find the solutions x * ik .Due to the fact that QAP is NP-hard, we propose a heuristic to solve it.
• A heuristic for QAP (QAP-heur) We are given n initial assignments in a random choice.The value of the objective function is evaluated, then two assignments are chosen according to a specific criterion.Then the value of the objective function is evaluated in the case of swapping these two assignments.If there is a reduction of the cost then this permutation is preserved if not the previous location are kept.The process is repeated until a given number of iterations.Thus for an initial allocation, we choose according to a criterion the best neighbour of the 2-opt permutation.The selection criterion is the following (reducing the risk that the algorithm does not stop on a local solution during search): } which maximizes the function : • Otherwise, if it maximizes the value: We put a counter that gives the number of iterations since the last best solution, and if it exceeds a given threshold, then we consider that the algorithm is stuck on a local solution.The threshold is empirically and randomly chosen in the interval [n 2 /2, 2n 2 ].The algorithm is a 2-opt neighbourhood search method associated with local search as tabou-search and stochastic search.
Data: n number of locations and activities, m number of iterations Random initialization of x; Algorithm 1: QAP-heur

Network Design Problem (ND)
We solve the following problem : The resolution of ND give us the solutions l * kl , y * i j and z kl * i j .We can simplify the objective function by c i j y i j without changing the problem as the sum a ik x * ik is constant.
• A heuristic for Network Design (ND-heur) In this sub-section we propose a heuristic for the Network Design problem (ND).This heuristic is known in the literature as deletion method (see Gamvros).The ND-heur consists of: 1.The method start with all edges in accordance with the instance solved.This mean that ∀i j ∈ {1, . . ., n}, if y * i j indicates the constracted roads then : 2. After that, choose y in the neighbourhood of y * such that : Set y k 0 l 0 = 0 and ∀(k, l) (k 0 , l 0 ) y kl = y * kl 3.Then, solve the following problem : where x * is a fixed activities allocation, to compute the shortest path associated to ȳ.So this problem is the network design problem associated with the activities assignment x * .We solve this problem by using the Flyod-Warshall algorithm (Floyd, 1962).
4. After resolution of the problem, if there is a feasible solution and a cost reduction, update y * = y.
5. Finally, repeat all the process from phase 2 up to n 2 − n times.
The step 2 specifies that a link which has the most contribution in the cost n ∑ i, j=1,i j c i j y i j is a candidate for deletion and is effectively deleted in the step 3 ( if there is a reduction of the global cost : ) .

Lower Bound Based on Distances
In this section, the proposed bound is based on a linearization of the cubic term of TNLU problem formulation (see section 2).We replace n ∑ k=1 n ∑ l=1 l kl x ik x jl by a new variable D i j ∀i, j = 1, ..., n, i j, with additional constraints to keep the problem equivalent.So D i j have a structure of distance, which allow us to add triangular inequality to the constraints (equation ( 3)).We obtain the following linear problem: with X the domain of affectation variables presented in subsection 3.2, and γ kl = min The relaxation of TNLU-LIN is a lower bound of TNLU, and we solve the problem using IBM CPLEX.Depending on the structure of the instances, to have a tighter bound, we can replace the inequalities of the constraints (6), ( 8) and ( 9) of (D ′ γ ) by equalities.

Numerical Simulations on Test Problems
The behavior of TNLU-heur is evaluated using test problems on academics instances and Dakar city.For academics instances we use data from Los (1978) and QAPLIB (Burkard et al., 1996) benchmark.The QAPLIB problems are instances of Nugent, Vollman, and Ruml (denoted by Nug), Elshafei (Els), Wolkowicz (Had), Krarup and Pruznan (Kra).More precisely we construct the instances from flows and distances of "Nug", "Els", "Had" and "Kra", where there is no localization and construction costs.In addition, the costs of localization are null and construction costs are obtained by multiplying the distances by a mean cost 10 ie c i j = 10.0 × d i j ∀ i, j = 1, . . ., n for instances "Nug" and " Had", by 100 for instances "Els", and by 1000 for instances "Kra".The mean construction costs are fixed according to the order of size of the values of distances.
Recall the considered data: n : number of locations and activities.
A = (a ik ) : location cost matrix of the activity i on the area k.Q = (q i j ) : moving flow matrix per unit of times between the activities i and j.C = (c i j ) : construction cost matrix of link or road (i, j).D = (d i j ) : direct distance matrix between two locations i and j or distance of link (i, j).
The lower bound has been solved with IBM Cplex Optimizer 12.5 solver interfaced with C++.The tests have been performed on a HP ProBook 450 Intel (R) Core (TM) i5 3230M CPU @2.40GHz.
For simulations on academic instances and Dakar city, we have executed ten times and have kept the best solution.This is due to the stochastic nature of the algorithm.

Tests on Academic Instances
In all tables, Los x means data of Los for n = x, Nug x means data of Nug for n = x and so on.In order to evaluate the result obtained by our heuristic we calculate the gap as follows: where opt is the best value found (upper bound) and LB is the lower bound based on distances.The following table 1 summarizes results for different step sizes.We see that the heuristic TNLU-heur provides good upper bounds for the TNLU problem.Indeed, for some instances of the academic test the gap is quite small.For instances Had and Nug, the CPU time grows at a linear rate with n.However, for instances Los the CPU time remains nearly constant, firstly for n=8, 10,12; and secondly for n=14,15.

Tests on Dakar City
For real life simulation of Dakar we could get from the DGTC (Direction des Travaux Géographiques et Cartographiques), a zonal division of the region into zones corresponding to the 2000 Census Report information with projection to 2015.And from this Census Report we could have information on the flow of movement (q i j ) between areas.For data of cost construction of roads (c i j ), we obtained information about mean costs updated in January 2008, and our choice fell on the road type monolayer whose width is 7m and the cost per km (all taxes included) is 33 246 500 FCFA (≃ 21.7764575e + 09 euro).So if we denote by C M the mean cost of road construction, then the cost of construction of the road is calculated by c i j = C M * d i j .
To calculate the distance we used OSRM project (Open Source Routing Machine) which is a C++ implementation of a routing engine high performance for shortest paths in road networks.It combines sophisticated routing algorithms with the data of the road network and free from OpenStreetMap (OSM).OSRM is able to calculate and output a shortest path between any origin and destination in milliseconds.
For the cost of localization activities (a ik ), we did not find an institution that can provide us with this information.Therefore we decided to build ourselves these costs realistically.We start with a choice of activities to locate in zones.The tests were done for a maximum number of zones equal to n = 20.Table 2 illustrates the concerned activities and zones.
In addition, to give the cost of locating an activity on a zone, we evaluated the mean cost of land for this activity in this zone and the cost of construction of the building for the same activity.The sum of these two costs gives the cost of location.So we look for the mean cost of 150m 2 and for a considered activity we estimate approximately a mean size of the area, then we can obtain a mean cost of land for the activity in this zone.
Table 7, where DKx means data on Dakar city intances with size n=x.With this data (information) the simulations has been performed.Figure 3 shows the logical connections between the different areas, given by the optimal solution for We can see that Diamniadio must be connected to Rufisque-Traditionnel and itself joined to Medina and Plateau-Sud, through direct links or roads.That is, it would be interesting to join the suburbs to downtown of Dakar, through direct connections to make the road network more fluid.In order not to overload the road network, we can develop a network of maritime shuttle linking these suburb areas and the downtown.And this is in perfect harmony with the geometrical shape of the entire south-eastern part of the region of Dakar.We see that in Table 7 the good result for simulations, where the gap is quite small for all sizes of tested problem.
In addition to the analysis of the SP, AL and ND problems, a new analysis of the set of relocation possibilities is proposed in the next section.

Multi-agent Simulations
In this section, the behavior of SDO, only limited to the test of localisation, is evaluated using test problems on Dakar city.We show how the simulator has been used to give the SDO software that can be seen as an over-layer of MATSim.
We use MATSim as a multi-agent simulator system, and need for that to generate agent plans.The necessary adaptation elements of the simulations for Dakar city are presented in sub-section 6.3.
Figure 3. Logical links between nodes

MATSim
A " multi-agent system (MAS)" is a system composed of a set of agents located in a certain environment and interacting according to certain relations.An agent is an entity characterized by the fact that it is, at least partially, autonomously.
The implementation of this system on a computer by computer programs, translating these interactions is the multi-agent simulation.In our paper, the agents are individuals (flows of activities) of a city that interact during the displacement.MATSim is a simulation tool multi-agent especially in transport system simulation.The tool is designed since a decade thanks to the collaboration of volunteers and researchers.This is one of the biggest free systems in that category.With adequate computer resources, MATSim is able to simulate the behavior of millions of individuals on very large transportation networks.It is compatible with the data sets of various geographic information systems as OpenStreetMap.And being free, it is designed to be extensible and customizable.
Very schematically, the MATSim simulations are conducted in three steps (see Figure 4).Each agent (or individual) has at the beginning of the simulation a plan of activities that he wants to perform on the day.This can be for example: leave home at 8:00 am to bring his children to school and then go to work, take back their children at 4:30 pm to lead them in a sporting activity, retrieve them to 6:00 pm to do some shopping and then come back to home.Initially, the plan thus specifies the sequence of tasks to be done by individuals but does not mention where exactly and how these activities will be made.The execution of the plan will specifically consists in making choices on these specific places and displacement routes to perform these activities.
After the first iteration of the simulation, a phase called "scoring" is performed.It will be to evaluate, to put a note (or a score) on the selected locations and routes.This score is a function called utility in the mathematical sense (function) and economic (utility) of the term that the alleged agents thinking and all acting in the same way (which is restrictive,...) will seek to maximize.
In this objective of maximizing that a replanning phase follows that of "scoring".Replanning allows agents to change eventually their plans by deciding other execution venues of their activities, other routes, or of elimination of activities, helping to increase the utility of their displacements.
The system iterate through these three steps until a maximum number of iterations, or until observation of minimal changes of utility scores.Note, however, that from a mathematical point of view, there is no evidence, in this system, of convergence to a stable state.MATSim simulation requires data on the transport network, the plans of the agents and the location of all the activities they will perform.

Simulation-Driven Optimization (SDO)
The SDO is an overlay of MATSim as part of a similar project called DAMA (Balac et al., 2014).The tool calculates the optimal locations of activities based on the results of the simulations of MATSim.The general operation is described in the Figure 5  simulation and propose a relocation of some activities (possibly all) to increase the utility scores.The problem that solves SDO is whether where should be located the activities so that the sum of utility scores of the agents is the highest possible.Due to the complexity of the function used by MATSim, we substitute in SDO, to that function the sum of displacement time or traveled distances.Intuitively, more displacement will be shorter in time and distance, and more one can hoped that the utility for the agent will be strong.
Thus we have the same problematic of combinatorial optimization designated by the acronym LU.For the time, SDO only limited to the aspect of localization.We do not treat the problem TNLU discussed in Section 2 but only a part (AL) thereof.The full integration of methods developed for TNLU will be the subject of further developments.
The large quantity of possible activities in a city makes necessary the execution of a previous phase called "clustering" where activities will be grouped into groups (clusters) of activities.Each group will be considered as representative of a single entity.This has the effect of reducing the number of activities to be taken into account.The total amount of displacement between these groups, as well as the displacement times are calculated in the phase known as "aggregation".These data represent matrices of flow and distance of the problem AL.Due to them, the heuristics developed for solving the problem AL are executed, after which an "optimal" relocation is proposed.This relocation changes the geographical coordinates of the proposed activities.Then, a new MATsim simulation is restarted, with these new localizations, and so on until a fixed number of iterations.Analogously to the MATSim simulation, there is no mathematical proof that the loop which concatenates simulation and optimization of localizations converges to a stable state of the displacements.

SDO computations on Dakar City
The use of SDO tool to the cases of Dakar requires to acquire the necessary files for MATSim simulations.Three types of files are generally used: the network file, the file containing the original plans of the agents, and one containing the list of activities and localizations.The variations (in reverse directions) of curves of distance and time appear to show that it is this phenomenon that has occurred for the case study in Dakar.

Conclusion
In this paper, we have proposed and analysed for the first time a new class of heuristic (TNLU-heur) for the TNLU problem.The new mixed integer formulation is a reformulated version of the Los's one in a more general sense.Unlike the model of Los we do not impose symmetry on the roads.Our new algorithm for TNLU was presented and applied to problems from the literature designed to highlight some intrinsic difficulties of TNLU, and from Dakar city.In addition to the convergence analysis, a new methodology of relocation possibilities is proposed.The SDO perfomance is evaluated by studying the quality of solution set in terms of optimal locations of activities.We investigated the average of displacement times, distances and scores; and established an interesting result characterizing the Utility, Travel Time and Travel Distance.In future work, we plan to study the combination between TNLU-heur and genetic or greedy algorithm.Another aspect that we wish to study will comprise the full integration of SP, AL and ND problems for further development.

Figure 7 .
Figure 7. Average of the utilities

Table 1 .
Results of TNLU-heur For each instance, the Time CPU is the minimal execution time in the ten experiments.