UC Berkeley Math 228B, Spring 2024, Problem Set 6
PerOlof Persson
persson@berkeley.edu
Due April 19
This content is protected and may not be shared, uploaded, or distributed.
Consider the traffic flow problem, described by the nonlinear hyperbolic equation: \[\begin{aligned} \frac{\partial\rho}{\partial t} + \frac{\partial\rho u}{\partial x} = 0 \end{aligned}\] with \(\rho=\rho(x,t)\) the density of cars (vehicles/km), and \(u=u(x,t)\) the velocity. Assume that the velocity \(u\) is given as a function of \(\rho\): \[\begin{aligned} u=u_\mathrm{max}\left( 1\frac{\rho}{\rho_\mathrm{max}} \right), \end{aligned}\] with \(u_\mathrm{max}\) the maximum speed and \(0\le\rho\le\rho_\mathrm{max}\). The flux of cars is therefore given by: \[\begin{aligned} f(\rho) = \rho u_\mathrm{max}\left( 1  \frac{\rho}{\rho_\mathrm{max}} \right). \end{aligned}\] We will solve this problem using a first order finite volume scheme: \[\begin{aligned} \rho_i^{n+1} = \rho_i^n\frac{\Delta t}{\Delta x} \left( F_{i+\frac{1}{2}}^n  F_{i\frac{1}{2}}^n \right). \end{aligned}\] For the numerical flux function, we will consider two different schemes:
 Roe’s Scheme

The expression of the numerical flux is given by: \[\begin{aligned} F_{i+\frac{1}{2}}^R = \frac{1}{2}\left[ f(\rho_i)+f(\rho_{i+1}) \right] \frac{1}{2}\left a_{i+\frac{1}{2}} \right (\rho_{i+1}\rho_i) \end{aligned}\] with \[\begin{aligned} a_{i+\frac{1}{2}} = u_\mathrm{max}\left( 1\frac{\rho_i+\rho_{i+1}}{\rho_\mathrm{max}} \right). \end{aligned}\] Note that \(a_{i+\frac{1}{2}}\) satisfies \[\begin{aligned} f(\rho_{i+1})f(\rho_i) = a_{i+\frac{1}{2}} (\rho_{i+1}\rho_i). \end{aligned}\]  Godunov’s Scheme

In this case the numerical flux is given by: \[\begin{aligned} F_{i+\frac{1}{2}}^G = f\left(\rho\left(x_{i+\frac{1}{2}},t^{n+}\right)\right) = \begin{cases} \min_{\rho\in[\rho_i,\rho_{i+1}]} f(\rho), & \rho_i<\rho_{i+1} \\ \max_{\rho\in[\rho_i,\rho_{i+1}]} f(\rho), & \rho_i>\rho_{i+1}. \end{cases} \end{aligned}\]
For both Roe’s Scheme and Godunov’s Scheme, look at the problem of a traffic light turning green at time \(t=0\). We are interested in the solution at \(t=2\) using both schemes. Use the following problem parameters: \[\begin{aligned} u_\mathrm{max}= 1.0, \quad \rho_\mathrm{max}= 1.0, \quad \rho_L = 0.8 \end{aligned}\] Solve on the domain \(x\in[2,2]\) using \(401\) cells, that is, the cell midpoints are located at \(x_j = 2+j\Delta x\) for \(j=0,\ldots,N\), with \(N=400\) and \(\Delta x = 4/N\). Note that the cell centers are located between the cell centers, so e.g. the left endpoint of the middle cell is located at \(x_{1/2} = \Delta x / 2\). To satisfy the CFL condition, use a timestep of \(\Delta t = \frac{0.8 \Delta x}{u_\mathrm{max}}\).
The initial condition at the instant when the traffic light turns green is \[\begin{aligned} \rho(x,t=0) = \begin{cases} \rho_L, & x<0 \\ 0, & x\ge 0 \end{cases} \end{aligned}\]
Finally, you need to specify boundary conditions, that is, the numerical fluxes at the endpoints of the domain. For simplicity, simply freeze \(\rho\) in the first and the last cells, which means no boundary fluxes are needed (this will implement a socalled characteristic boundary condition, where your numerical flux function will decide if the boundary is an inflow or an outflow).
What do you observe for each of the schemes? Explain briefly why the behavior you get arises.
For problems 2  3, use only the scheme(s) which are valid (that is, converging to the weak entropy solution).Simulate the effect of a traffic light at \(x=\frac{\Delta x}{2}\) which has a period of \(T=T_1+T_2=2\) time units. Assume that the traffic light is \(T_1=1\) units on red and \(T_2=1\) units on green. Set \(\rho_L=\frac{\rho_\mathrm{max}}{2}\) on the left boundary, giving a maximum flux and a high density of incoming cars. Determine the average flow, or the capacity of cars over a time period \(T\).
The average flow can be approximated as \[\begin{aligned} \dot{q}=\frac{1}{N_T} \sum_{n=1}^{N_T} f^n = \frac{1}{N^T}\sum_{n=1}^{N_T} \rho^n u^n, \end{aligned}\] where superscripts indicate the timestep and \(N_T\) is the number of time steps for each period \(T\). You should run your computation until \(\dot{q}\) over a time period does not change. Note that because of continuity, this can be evaluated at any cell of the interior of the domain (not at the boundaries since they are frozen).
Note: A red traffic light can be modeled by simply setting \(F_{i1/2}=0\) at the position where the traffic light is located.
Assume now that we simulate two traffic lights, one located at \(x=0\), and the other at \(x=0.15\), both with a period \(T\). Calculate the road capacity (= average flow) for different delay factors. That is if the first light turns green at time \(t\), then the second light will turn green at \(t+\tau\). Solve for \(\tau=k\frac{T}{10}\), \(k=0, \ldots, 9\). Plot your results of capacity vs \(\tau\) and determine the optimal delay \(\tau\).
Optional, not graded. Consider the Eikonal equation for first arrivals/optimal path planning problems. The equation can be written \[\begin{aligned} F(x,y)  \nabla \phi(x,y)  = 1, \end{aligned}\] where \(F(x,y)\) is the speed function. We will use the solution \(\phi(x,y)\) to determine the optimal path from the departure point \((x_D,y_D)\) to the arrival point \((x_A, y_A)\), where \(\phi(x_A,y_A)=t\). The level set with value \(t\) of the function \(\phi(x,y)\) gives the maximum distance away from our departure point that can be traveled in a time \(t\). In addition, the optimal path between the departure point and the arrival point is determined by traveling in the normal direction of the level sets.
We will solve the Eikonal equation using a level set method and a timestepping approach. The equation \[\begin{aligned} \phi_t + F  \nabla \phi  = 1,\qquad \phi(x_D,y_D)=0 \end{aligned}\] is integrated in time until a steadystate is reached. This is not an efficient method for solving the Eikonal equation, but it will be sufficient for this problem set (and it illustrates how to solve more general timedependent problems). If you are interested in more sophisticated solvers, feel free to implement the more efficient fast marching method instead.
Write a computer code to solve the Eikonal equation on the unit square \(x,y\in [0,1]\). Use the firstorder upwinded scheme in space, and appropriate treatment of the boundaries.
Write a computer code to find the optimal path between the departure/arrival point, by solving the ODE \(d\boldsymbol{r}/dt = \boldsymbol{n}\), where \(\boldsymbol{r}\) is the current position on the path and \(\boldsymbol{n}\) is the normal vector.
Run your codes with grid spacing \(h=1/100\) for the following cases and plot both the solutions (e.g. as contour curves of \(\phi(x,y)\)) and the optimal paths:
\((x_D,y_D)=(0.2, 0.2)\), \((x_A,y_A)=(0.8, 0.8)\), \(F(x,y)=1\) (for testing).
\((x_D,y_D)=(0.2, 0.2)\), \((x_A,y_A)=(0.8, 0.45)\), and \[\begin{aligned} F(x,y)= \begin{cases} 1.0 & \text{if }y\ge 0.5, \\ 0.5 & \text{if }y< 0.5. \end{cases} \end{aligned}\]
\((x_D,y_D)=(0.2, 0.2)\), \((x_A,y_A)=(0.8, 0.8)\), and \[\begin{aligned} F(x,y)=10.9\cdot\cos(4\pi x)\cdot e^{10\,((x.5)^2+(y0.5)^2)} \end{aligned}\]
Make up your own speed function \(F\), only returning the values 0.01 and 1 but with a nontrivial optimal path.