Road Network Design Problem

5.11. Road Network Design Problem#

Attribution

This chapter reuses content Bahman Ahmadi et al. (2024). Find out more here.

During Friday’s session, you will work on the road Network Design Problem (NDP) and explore its implementation using two optimization approaches: Mixed-Integer Linear Programming (MILP) and a Genetic Algorithm (GA). The problem description, mathematical formulation, and Python implementation are provided in this book, should you wish to prepare in advance for the session.

MUDE Exam Information

The road Network Design Problem serves as an example for solving a complex optimization problem. You’re not expected to understand the problem details for the exam.

Problem description#

The road Network Design Problem (NDP) involves determining which road links to build, refurbish, or upgrade to improve the performance of a road network. A specific variation of the problem is the road NDP with capacity expansions, where the objective is to decide which links should have their capacity increased. This is a complex problem that must take into account a variety of factors, including:

  • The current state of the road network

  • The projected future traffic demand

  • The budget available for improvements

  • The impacts of the adjustments (could be environmental, social, etc).

There are a variety of approaches to dealing with the road NDP with capacity expansions. One common approach is to use mathematical optimization models. In this assignment, we use a simplified example to show how optimization can be used to tackle road NDPs. Note that the classical approaches to dealing with road NDPs can be more complicated and will be covered in other courses in the TTE track of the Civil Engineering Master’s program.

In this assignment, the goal is to minimize the total travel time on the network by selecting a predefined number of links for capacity expansion (subject to the available budget). The following (main) assumptions and simplifications are made to make the problem solvable using methods and algorithms that you have learned so far.

2. Route choice behavior#

To assess the quality of solutions for the road capacity expansion problem, it is essential to understand the effect of added capacity on travel time. This requires not only modeling the travel time-flow relationship but also knowing how vehicles distribute themselves across the network.

In congested networks, the route choice behavior of drivers often follows the so-called User Equilibrium (UE) principle. Under UE, each traveler seeks to minimize their own individual generalized travel time. The principle states that, for each origin-destination pair, all routes used have equal and minimal travel times. In this state of equilibrium, no driver can improve their travel time by switching to a different route.

However, calculating UE requires advanced methods, which are beyond the scope of the MUDE course. Instead, we assume that drivers follow the so-called System Optimal (SO) principle. Under SO, route choices are made to minimize the total travel time across the entire network, summed over all drivers. This means some drivers may take longer routes to allow others to save time, a concept also referred to as social equilibrium.

The SO equilibrium is easier to compute than the UE. However, it is important to note that achieving a system-optimal traffic distribution in real-world road networks is nearly impossible, as drivers independently make route choices based on personal preferences and information. In practice, it is not feasible to dictate where drivers should go.

Example of a perfect user equilibrium in a very simple road:

https://files.mude.citg.tudelft.nl/sketchoptimization.png

Fig. 5.47 Route choice behavior according to the user equilibrium (UE) principle#

3. Quadratic terms#

In the formulation of the Road NDP, even after applying the simplifications mentioned earlier, a quadratic term arises. Specifically, the flow (a variable) must be multiplied by the travel time (also a variable), making the problem non-linear.

Fortunately there are different methods to transform quadratic terms to linear variables and constraints with mathematical programming. You do not need to learn these techniques. Most solvers (including Gurobi) can handle this well given some adjustments to the formulation. Gurobi uses the McCormick Envelopes (MCE) to transform each quadratic term into one new variable and four constraints. For more information on MCE and how Gurobi implements them, check these links (but this is not part of your learning goals):

In this assignment, we will proceed with a quadratic term in the objective function. With the earlier simplifications and assumptions, we can formulate and solve the road NDP using the branch-and-bound method, which you have studied before.

We are using the Sioux Falls network which is one of the most used networks in transportation research

https://files.mude.citg.tudelft.nl/SiouxFallsMap_AAA1998.jpg

Fig. 5.48 SiouxFalls map#

Mathematical description#

Decision variables#

We have a set of binary variables \(y_{ij}\), these variables take the value 1 if link \((i,j)\) connecting node \(i\) to node \(j\) is selected for expansion, and 0 otherwise.

We also have two sets of decision variables representing link flows; \(x_{ij}\), representing flow on link \((i,j)\) in cars per hour, and \(x_{ijs}\), representing flow on link \((i,j)\) going to destination \(s\).The first represents the total number of cars passing on that road, while the second represents the number of cars specifically going to destination \(s\). Summing the latter over all \(s\) results in the former for a link \((i,j)\).

Therefore, mathematically we define the domain of the variables as follows:

\[\begin{split}\begin{align} & y_{ij} \in \{0, 1\} \quad \forall (i,j) \in A \\ & x_{ij} \geq 0 \quad \forall (i,j) \in A \\ & x_{ijs} \geq 0 \quad \forall (i,j) \in A, \forall s \in D \\ \end{align}\end{split}\]

As you will see in the code implementation, we have one extra set of variables called x2 (representing \(x^2\)). This is to help Gurobi isolate quadratic terms and perform required transformations based on MCE to keep the problem linear. This is not part of your learning goals.

Test yourself

Objective function#

The objective function of the problem (in its simplest form), is the minimization of the total travel time on the network. This is calculated by multiplying the flow of vehicles on each link by the corresponding travel time and summing over all links (\(A\) represents the collection of all links to simplify the notation):

\[Z = \sum_{(i,j) \in A}{ x_{ij} \cdot t_{ij}} \]

The travel time \(t_{ij}\) is a function of the flow on a link and can be expressed as follows (where \(\beta\) is a parameter):

\[ t_{ij} = t^0_{ij} \cdot \left( 1 + \beta \cfrac{x_{ij}}{c_{ij}}\right) \quad \forall (i,j) \in A \]

Note that the commonly used travel time function based on the Bureau of Public Roads (BPR) is as follows:

\[ t_{ij} = t^0_{ij} \cdot \left( 1 + \alpha \left (\cfrac{x_{ij}}{c_{ij}}\right)^\beta\right) \quad \forall (i,j) \in A \]

Where \(\beta\) usually assumes the value of \(4\) making this function (and the problem) non-linear. Therefore, we use the linear function mentioned before to make the problem manageable using the methods we have learned so far.

The following constraint yields the capacity of each link based on which ones are selected for expansion, when \(y_{ij}\) is 1 there is added capacity as you can see:

\[ c_{ij} = (1 - y_{ij}) \cdot c^0_{ij} + y_{ij} \cdot c^1_{ij} \quad \forall (i,j) \in A \]

This allows us to represent \(t_{ij}\) as:

\[ t_{ij} = t^0_{ij} \cdot \left( 1 + \beta \cdot x_{ij} \cdot \left(\cfrac{1 - y_{ij}}{c^0_{ij}} + \cfrac{y_{ij}}{c^1_{ij}} \right) \right) \quad \forall (i,j) \in A \]

Which leads to the following extended objective function:

\[ Z = \sum_{(i,j) \in A}{ x_{ij} \cdot \left(t^0_{ij} \cdot \left( 1 + \beta \cdot x_{ij} \cdot \left(\cfrac{1 - y_{ij}}{c^0_{ij}} + \cfrac{y_{ij}}{c^1_{ij}} \right) \right) \right)} \]

Now, for Gurobi (and other solvers as well), it is necessary to keep binary variables and quadratic terms clean and separate for the solver to perform the required transformations to linearize the problem. Thus, the following, despite its complexity, is the most solver-friendly formulation of our objective function:

\[ \begin{align} Z = \sum_{(i,j) \in A}{t^0_{ij} \cdot x_{ij}} + \sum_{(i,j) \in A}{\left(t^0_{ij} \cdot \frac{\beta}{c^0_{ij}} \cdot x^2_{ij} \right)} - \sum_{(i,j) \in A}{\left(t^0_{ij} \cdot \frac{\beta}{c^0_{ij}} \cdot x^2_{ij} \cdot y_{ij} \right)} + \sum_{(i,j) \in A}{\left(t^0_{ij} \cdot \frac{\beta}{c^1_{ij}} \cdot x^2_{ij} \cdot y_{ij} \right)} \end{align}\]

Therefore, we use this equation to model our objective function in Gurobi.

Constraints#

We have four sets of constraints for this problem. Let’s go through them one by one and add them to the model.

1. Budget constraint#

We can only extend the capacity of certain number of links based on the available budget. So first, we have to make sure to limit the number of extended links to the max number that can be expanded:

\[ \sum_{(i,j) \in A}{ y_{ij}} = B \]
Test yourself

3. Node flow conservation constraints#

This set of constraints ensures that the total incoming and outgoing flow at each node balances, except at origin and destination nodes where flow is either generated or absorbed. For example, at an intersection, the number of vehicles entering equals the number leaving, except for origins (extra outgoing flow) and destinations (extra incoming flow).

\(d_{is}\) here is the number of travelers from node \(i\) to node \(s\) with the exception of \(d_{ss}\), which is all the demand that arrives at node \(s\).

\[ \sum_{j \in N; (i,j) \in A}{ x_{ijs}} - \sum_{j \in N; (j,i) \in A}{ x_{jis}} = d_{is} \quad \forall i \in N, \forall s \in D \]

The figure gives an example:

https://files.mude.citg.tudelft.nl/equil.png

Fig. 5.49 Constraint example#

\[ \sum_{j \in N; (i,j) \in A}{ x_{ijs}} - \sum_{j \in N; (j,i) \in A}{ x_{jis}} = d_{is} \quad \forall i \in N, \forall s \in D \]
Test yourself

4. Quadratic variable constraints (you do not need to fully understand this)#

These constraints are essentially dummy equations to help Gurobi handle the quadratic terms (introduced earlier as dummy variables). Instead of directly using \(x^2_{ij}\) in the model, we define a new set of decision variables and constraints to represent their values. This informs Gurobi that these are quadratic terms and allows it to replace them with the necessary variables and constraints to maintain linearity. This is not part of your learning goals!

Conclusions#

Here is the complete problem formulation.

\[\begin{split}\begin{align} min \quad & \sum_{(i,j) \in A}{t^0_{ij} \cdot x_{ij}} + \sum_{(i,j) \in A}{\left(\frac{t^0_{ij} \cdot \beta}{c^0_{ij}} \cdot x^2_{ij}\right)} - \sum_{(i,j) \in A}{\left(\frac{t^0_{ij} \cdot \beta}{c^0_{ij}} \cdot x^2_{ij} \cdot y_{ij}\right)} + \sum_{(i,j) \in A}{\left(\frac{t^0_{ij} \cdot \beta}{c^1_{ij}} \cdot x^2_{ij} \cdot y_{ij}\right)} \\ \text{s.t.} \quad \\ & \sum_{(i,j) \in A}{ y_{ij}} = B \\ & \sum_{s \in D}{x_{ijs}} = x_{ij}, \quad \forall (i,j) \in A \\ & \sum_{j \in N; (i,j) \in A}{ x_{ijs}} - \sum_{j \in N; (j,i) \in A}{ x_{jis}} = d_{is}, \quad \forall i \in N, \forall s \in D \\ & y_{ij} \in \{0, 1\}, \quad \forall (i,j) \in A \\ & x_{ij} \geq 0, \quad \forall (i,j) \in A \\ & x_{ijs} \geq 0, \quad \forall (i,j) \in A, \forall s \in D \\ \end{align}\end{split}\]

Adaption for GA#

To be able to use a Genetic Algorithm (GA), we need to define the problem slightly differently. Essentially, we break down the problem into two sub-problems:

  1. The Traffic Assignment (TA) problem: the route choices of the drivers.

  2. The road Network Design Problem (NDP): where we select which links should be upgraded.

We solve the problem by iteratively going between the Traffic assignment and the Design Problem. The idea is for the GA to move to better networks as generations pass which are evaluated by the traffic assignment process that you have learned. We use Gurobi to solve the Traffic Assignment sub-problems, which provide us with the objective function (or fitness function within the context of GA) value of the decision problem (which will be dealt with using GA). This is usually referred to as the iterative-optimization-assignment method since we iteratively improve the objective function value of the NDP using the assignment problem.

Let’s see how this works.

The Network Design sub-problem#

The network design is where we use the genetic algorithm. As explained before, GA uses a population of solutions and iteratively improves this population to evolve to new generations of populations with a better objective function value (being that minimization or maximization). In this problem, the decision variables are links for capacity expansion and the objective function value is the total system travel time that we want to minimize.

\[\begin{split}\begin{align} min \quad & \sum_{(i,j) \in A}{t^0_{ij} \cdot x_{ij}} + \sum_{(i,j) \in A}{\left(t^0_{ij} \cdot \frac{\beta}{c^0_{ij}} \cdot x^2_{ij}\right)} - \sum_{(i,j) \in A}{\left(t^0_{ij} \cdot \frac{\beta}{c^0_{ij}} \cdot x^2_{ij} \cdot y_{ij}\right)} + \sum_{(i,j) \in A}{\left(t^0_{ij} \cdot \frac{\beta}{c^1_{ij}} \cdot x^2_{ij} \cdot y_{ij}\right)} \\ \text{s.t.} \quad \\ & \sum_{(i,j) \in A}{y_{ij}} = B \\ & y_{ij} \in \{0, 1\}, \quad \forall (i,j) \in A \\ \end{align}\end{split}\]

Where the values of \(x_{ij}\) are not decision variables anymore, they will be obtained from solving the Traffic Assignment problem with Gurobi which evaluates the travel times on the network.

The Traffic Assignment sub-problem#

This is just part of the original NDP that assigns traffic to the network based on a set of given capacity values, which are defined based on the values of the DP decision variables (links selected for capacity expansion). The main difference (and the advantage) here is that by separating the binary decision variables, instead of a mixed integer programming problem, which are hard to solve, here we have a quadratic programming problem with continuous decision variables, which will be transformed to a linear problem that Gurobi can solve very fast.

\[\begin{split} \begin{align} min \quad & \sum_{(i,j) \in A}{t^0_{ij} \cdot x_{ij}} + \sum_{(i,j) \in A}{\left(t^0_{ij} \cdot \frac{\beta}{c^0_{ij}} \cdot x^2_{ij}\right)} - \sum_{(i,j) \in A}{\left(t^0_{ij} \cdot \frac{\beta}{c^0_{ij}} \cdot x^2_{ij} \cdot y_{ij}\right)} + \sum_{(i,j) \in A}{\left(t^0_{ij} \cdot \frac{\beta}{c^1_{ij}} \cdot x^2_{ij} \cdot y_{ij}\right)} \\ \text{s.t.} \quad \\ & \sum_{s \in D}{x_{ijs}} = x_{ij}, \quad \forall (i,j) \in A \\ & \sum_{j \in N; (i,j) \in A}{x_{ijs}} - \sum_{j \in N; (j,i) \in A}{x_{jis}} = d_{is}, \quad \forall i \in N, \forall s \in D \\ & x_{ij} \geq 0, \quad \forall (i,j) \in A \\ & x_{ijs} \geq 0, \quad \forall (i,j) \in A, \forall s \in D \\ \end{align} \end{split}\]

Where the values of \(y_{ij}\) are constant and are defined by GA.

Test yourself

Summarizing#

The following is a diagram that shows what you are finally doing to solve the same problem but with a meta-heuristic approach:

https://files.mude.citg.tudelft.nl/GAdiagram.png

Fig. 5.50 Summary meta-heurstic approach#

PyMOO#

PyMOO is a Python library that provides a comprehensive and easy-to-use framework for multi-objective optimization (MOO). For this case, we are going to deal with only one objective; nevertheless, this is a useful tool if you have more objectives. In addition, PyMOO easily allows us to define our optimization problem by specifying the objectives, constraints, and decision variables.

Implementation in Python#

The implementation in Python is shown in the notebooks for MILP and GA. The full notebooks including source data are provided here

We can now start solving our problem!

https://files.mude.citg.tudelft.nl/outputopt.png

Fig. 5.51 A solution of the problem#

Attribution

This chapter is written by Gonçalo Homem de Almeida Correia, Maria Nogal Macho, Jie Gao and Bahman Ahmadi. Find out more here.

Contents