# D-OccNet: Detailed 3D Reconstruction Using Cross-Domain Learning

Minhaj Uddin Ansari1, Talha Bilal2, Naeem Akhter3  1School of Computing
Queen’s University, Kingston, Ontario K7M2W8
Email<EMAIL_ADDRESS>2Escuela Politécnica Superior (EPS)
Universidad Autónoma de Madrid, Madrid, Spain 28049
Email<EMAIL_ADDRESS>3Department of Computer and Information
Sciences
Pakistan Institute of Engineering and Applied Sciences, Islamabad, Pakistan
45650
Email<EMAIL_ADDRESS>
###### Abstract

Deep learning based 3D reconstruction of single view 2D image is becoming
increasingly popular due to their wide range of real-world applications, but
this task is inherently challenging because of the partial observability of an
object from a single perspective. Recently, state of the art probability-based
Occupancy Networks reconstructed 3D surfaces from three different types of
input domains: single view 2D image, point cloud and voxel. In this study, we
extend the work on Occupancy Networks by exploiting cross-domain learning of
image and point cloud domains. Specifically, we first convert the single view
2D image into a simpler point cloud representation, and then reconstruct a 3D
surface from it. Our network, the Double Occupancy Network (D-OccNet)
outperforms Occupancy Networks in terms of visual quality and details captured
in the 3D reconstruction.

###### Index Terms:

Mesh; Occupancy Network; Point Cloud; Reconstruction;

##  I Introduction

3D image reconstruction from single view 2D images is a fundamental challenge
in computer vision. In recent days, demand for 2D to 3D reconstruction has
increased due to its wide range of applications in autonomous vehicles,
virtual reality, medicine, augmented reality and robotics. Although this task
is the basic functionality of the human brain, it is a very complex problem
for a computer. Fortunately, with the availability of large open source 3D
datasets like ShapeNet [16], deep learning methods have gained rapid
advancement in this area. Deep learning algorithms are now capable of
developing 3D images in a variety of domains due to supervised learning on 3D
datasets. Presently, deep learning algorithms for reconstructing 3D images are
divided into three major categories: voxel, point cloud and mesh based on
output domain representation of the reconstructed image. We will briefly go
through the research contributions in all three groups.

Voxels are small cubic units depicting pixels in 3D domain, where each
individual pixel has a length, width, and height. In recent years, there has
been a lot of work done in voxel domain due to its simplicity. However, since
each voxel must occupy volume, 3D images in the voxel domain take up a lot of
memory and have coarse boundaries and surfaces, inevitably resulting in lower
resolutions. C.B Choy et al [1] introduces 3D Reconstruction Neural Network
(3D-R2N2) which is an extension of the standard Long Short Term Memory (LSTM)
network. 3D-R2N2 unifies single-view and multi-view 3D reconstruction in a
single framework and takes the advantage of the power of the LSTM to retain
previous observations and incrementally refine the output reconstruction as
more observations become available. H. Xie et al [17] uses context aware
fusion module in its encoder-decoder structure for high quality reconstruction
in voxel domain. X. Yan et al [18] also uses encoder-decoder architecture and
introduces a 2D silhouette loss function based on perspective transformations.
J. Wu [15] leverages the power of 3D convolutional networks and generative
adversarial nets to accomplish this task of 3D reconstruction in an
adversarial fashion. MarrNet [14] is an end to end trainable model consisting
of a 2.5D sketch estimator, a 3D shape estimator and a reprojection
consistency function to enforce alignment between estimated 3D shape and
inferred 2.5D sketches.

Figure 1:  Overview of the D-OccNet

Point clouds are points distributed throughout a bounded 3D space. Since the
points are distributed with a high degree of freedom, and there aren’t any
connections between them to ensure they are in the proper position relative to
one another, it is difficult to control the 3D images in the point cloud
domain generated by deep learning algorithms. H. Fan et al [2] introduces a
novel conditional shape sampler, capable of predicting multiple plausible 3D
point clouds from an input image. P. Mandikal et al [10] proposes an latent
embedding matching approach for 3D reconstruction. They first train a 3D point
cloud auto-encoder and then learn a mapping from the 2D image to the
corresponding learnt embedding. They introduced a probabilistic training
scheme in the image encoder to obtain multiple plausible 3D mappings from a
single input image. C. Lin et al [8] uses 2D convolutional operation based
structure generator to predict 3D structure and introduces a pseudo-rendering
pipeline to serve as a differentiable approximation of true rendering. They
also use the pseudo-rendered depth images for 2D projection optimization for
learning dense 3D shapes. Zhao et al [19] introduced Variational Auto Encoder
(VAE) based 3D Point Capsule Network architecture to process sparse 3D point
cloud.

A mesh a set of vertices and edges connected to form surfaces. Meshes are the
best choice for 3D images because they aren’t as memory intensive as voxel
images, they have smoother surfaces compared to voxel images, and the vertices
are connected by edges to one another, so the structure and texture rendering
is more stable compared to point cloud images. C. Kong et al [6] uses
dictionary of dense CAD models and propose a two-step strategy. First it
employs orthogonal matching pursuit to rapidly choose the closest single CAD
model in the dictionary to the projected image. Second, it employs a novel
graph embedding based on local dense correspondence to allow for sparse linear
combinations of CAD models. T. Groueix et al [3] introduces a novel 3D surface
generation network AtlasNet, which is composed of a union of learnable
parametrizations. These learnable parametrizations transform a set of 2D
squares to the surface. A recent trend is to use deformation network for 3D
reconstruction. D. Jack et al [5] introduces a network that learns to deform
points sampled from a high-quality mesh and produces meshes with fine-grained
geometry. Pixel2Mesh [13] uses graph convolutions to deform an initial
ellipsoid according to features extracted from the images. This technique is
really limited to reconstruct a few categories because it highly depends upon
shape of initial 3D shape (ellipsoid).

Many techniques in Computer Vision focus on 3D reconstruction from more than
one perspective (multiple 2D views) of an object, but in many cases we are
only available with one perspective (single 2D view) of the object. 3D
reconstruction using single view is more difficult task because you are only
provided with partial observation of object which can lead to infinite
representations. However it is also the most natural considering humans can
easily visualize the 3D representation of an object by looking at it from a
single perspective. In this paper we focus on 3D reconstruction from only a
single view of the object.

Our Double Occupancy Network (D-OccNet) is inspired by the state of the art
Occupancy Networks (OccNet) [11] that learns boundary of object in 3D space
using occupancy values. The difference between OccNet and D-OccNet is that
OccNet only uses a single domain for image to mesh conversion, whereas
D-OccNet is based on cross domain model learning. In D-OccNet, the input 2D
single view image is first transformed into a point cloud and then post-
processor transforms the point cloud into the mesh domain. This two-step
transformation helps better reconstruct an image in the 3D domain. While
quantitative results of D-OccNet are comparable to OccNet, it outperforms
OccNet in qualitative comparison.

The rest of this paper is organized as follows: Section 2 describes the
working of Occupancy Network; Section 3 explains the architecture of learning
based D-OccNet model; Section 4 discusses the training setup and experimental
results; and, finally, Section 5 concludes this work.

##  II Occupancy Networks

The Occupancy Network (OccNet) is an encoder decoder convolution neural
network architecture. The encoder encodes the input into a vector and the
decoder generates a 3D mesh from that vector. Occupancy Networks can convert
2D single view images, point clouds and voxels to probability occupancy maps
and these maps are used to extract corresponding 3D meshes. For different
types of inputs, different types encoders are used in the Occupancy Network;
however, in all cases, the decoder remains the same.

###  II-A Encoder

As discussed in the previous section, an occupancy network can encode either
an image, or a point cloud or a voxel. A 2D classification networks encoder,
such as ResNet18 [4] can be used if the input is a single view 2D image. A
PointNet [12] encoder is used if the input is a point cloud and a 3D
convolution encoder network is used if the input is a voxel. When a 2D single
view image is passed to ResNet18; a point cloud is passed to a PointNet; or a
voxel is passed to a 3D convolution encoder network, the encoder outputs a
vector. The dimensions of the vector can vary, but the decoder can be adjusted
to learn the features encoded in the vector regardless of its dimensionality.

###  II-B Decoder

The decoder uses CBN (Conditional Batch Normalization), fully connected layers
and Resnet blocks to generate an occupancy map from the input encoded vector.
Nx3 3D-coordinate points are passed into a fully connected layer to generate a
feature vector of Nx256 dimensions where each point is associated with a
256-dimensional vector. This feature vector is then passed into a Resnet
block. The Resnet block starts with applying CBN layer, followed by ReLU
activation and then a fully connected layer. The output is sent into a second
CBN layer, followed by ReLU activation and a second fully connected layer. The
output is then added to the Resnet input. This process repeats five times.
After five iterations, the output is again passed through a CBN and ReLU
layer, and finally a fully connected layer reduces the output to a Nx1
1-dimensional vector. This vector is passed through a sigmoid function to
produce the occupancy probabilities corresponding to each point in the input
Nx3 vector. The CBN layer receives the encoded vector c from the encoder and
passes it through two fully connected layers $\beta$(c) and $\gamma$(c). The
input feature vector is normalized using first and second order moments,
multiplied with the output of fully connected layer $\gamma$(c) and then added
to the output of fully connected layer $\beta$(c). In order to calculate
occupancy value, cross entropy loss is used which is calculated on each
minibatch.

$f_{out}=\gamma(c){(f}_{in}-\mu)/\sqrt{\sigma^{2}+\epsilon}+\beta(c)$

###  II-C Mesh Extraction

The decoder outputs the occupancy values associated with each point in the Nx3
3D-coordinate points matrix initially sent as input. Points that have high
occupancy value, or high probability, are considered lying inside the object
and points that have lower occupancy values are considered lying outside the
object. Multiresolution Iso-Surface Extraction (MISE) [7][9] is used to
develop cube around these points. The cubes can then be broken down to the
desired resolution. After the desired resolution is obtained, the marching
cube algorithm [20] is used to convert the voxel to a mesh.

##  III Method

We found that 2D image to 3D mesh conversion by the OccNets resulted in meshes
that weren’t as perfect compared to their corresponding ground truths. This is
because the Occupancy Network encoder CNN could only generate latent vectors
containing information of the 2D image fed to it. There wasn’t any 3D
information embedded into the latent vector. We hypothesized that if we used a
point cloud representation of the input image, thereby integrating some 3D
information, and fed it into the Occupancy Network encoder, the latent vector
produced by the encoder would contain some 3D information of the input,
resulting in better quality 3D meshes. Using this hypothesis, we introduced
Double Occupancy Network (D-OccNet) as an extension of OccNet. The ”D” in
D-OccNet refers to the architectural design of this network. It contains two
OccNets connected together by a small processing unit in between the output of
the first OccNet and the input of the second OccNet. The main job of this
small processing unit is to convert the mesh output of the first OccNet into a
suitable point cloud input for the second OccNet. We will discuss the
architecture in detail in the upcoming section.

Figure 2:  Internal architecture of D-OccNet

D-OccNet contains three steps in the following sequence (1) image to mesh
reconstruction, (2) mesh to point cloud conversion and (3) point cloud to mesh
reconstruction.

###  III-A Image to Mesh Reconstruction

A single view 2D image is fed into a pre-trained Resnet18 model. Resnet18
encodes the image into a latent vector of size 256x1. Resnet18 can be replaced
with any 2D CNN network to generate a latent vector. This latent vector is
sent to the OccNet’s decoder to produce a 3D mesh.

###  III-B Mesh to Point Cloud Conversion

The 3D mesh generated in the first step contains approximately 8000-12000
vertices connected by edges. This is as far as the original OccNet goes. The
PointNet encoder has been trained on 300 points point cloud as input. To make
the input suitable for PointNet, we use MeshLab, an open source 3D mesh
processing software to reduce the number of points in the 3D mesh to 300. We
save the 300 points of the generated mesh as a point cloud, thereby removing
the edge information.

###  III-C Point Cloud to Mesh Reconstruction

The point cloud from the previous step is fed into the PointNet encoder.
PointNet generates a latent vector of size 512x1, which is sent into the
OccNet’s decoder. The decoder generates a 3D mesh from the latent vector. This
3D mesh is the final output of our D-OccNet.

##  IV Training and Results

Figure 3:  Comparison of IoU, Chamfer and Normal Consistency between OccNet
and D-OccNet

In our implementation of D-OccNet, the training of two OccNets are involved.
We train the first Occnet to generate meshes from single view 2D images. The
first OccNet contains a ResNet18 encoder and OccNet’s decoder. We then train
our second OccNet to generate meshes from point cloud.The second OccNet
contains a PointNet encoder and the OccNet’s decoder.

We use the same dataset as as used by the authors in [1]. This dataset
contains 13 object categories of 3D CAD models taken from the ShapeNet
dataset. ShapeNet is a large annotated dataset that contains 12,000 CAD models
belonging to 270 different object categories. It is one of the most widely
used dataset for training supervised deep learning models.

We augment the input data by scaling it by a random factor between 0.75 and 1
and then taking a random crop. Train/Test split ratio is set the same as in
[2]. Batch size is set to 64 and the Adam optimizer is set with a learning
rate of 1e-4, $\beta$1 = 0.9 and $\beta$2 = 0.999. Validation accuracy is set
as the stopping criteria for training. The code has been written in Python2.
For the training, we use NVIDIA GTX 1070 GPU in a 64 GB RAM Core i7 computer
running on Ubuntu OS.

We compare the output meshes from D-OccNet with the original OccNet. While
surface generation methods struggle to capture shape details such as table
surface, airplane wings and car wheels, our method preserves these details
that reside in the source mesh. We evaluate D-OccNet’s output both
quantitatively and qualitatively and draw a comparison with the original
OccNet. For Quantitative comparison we use IOU, Chamfer-L1 and Normal
Consistency metrics. IoU is the intersection over union metric which computes
the overlapping of two meshes, in our case the ground truth and predicted
mesh, with the union. If the overlapping is same, it will give a score closer
to 1, meaning the final output is closer to the ground truth. As evident by
the Fig. 4, our D-OccNet produces meshes having IOU comparable to the original
OccNet. The Chamfer loss compares the distance between points in the output
and predicted meshes. The smaller the distance, the lesser the dissimilarity
is between the meshes. Again, our meshes outperform OccNet’s meshes when
Chamfer L1 losses are compared. Finally, Normal Consistency metric shows that
the normal of our mesh are closer to the normal of the ground truth mesh.

Quantitatively D-OccNet is comparable to OccNet but our meshes have visual
superiority. Figure 5 shows a visual comparison between the outputs of both
networks. In case of the airplane mesh, our D-OccNet reconstructed the front
wheel of the plane, which wasn’t reconstructed by OccNet. The engine fans are
also much smoother and sharper. The table mesh has sharper edges and finer
corners, and the car mesh has sharper passenger seats and includes rearview
mirrors.

The main reason for similar quantitative result was due to choice of metrics.
There is no available metric that can differentiate the little details between
two meshes. For example, if OccNet generates a 3D mesh of a car that has no
rearview mirrors but is of the right volume, and if D-OccNet produces the same
car mesh with rearview mirrors, but the car is slightly larger in volume than
the ground truth, the IoU metric values obtained by comparing the ground truth
of both networks will be similar because the volume overlapping with the
ground truth is almost identical in both cases. In other words, the IoU of
OccNets mesh output missing the rearview mirror is same is the IoU of
D-OccNets mesh output having the rearview mirror but the body being larger in
size. Although this might be critical in medical applications where small
changes of volume can hold great significance, other applications such as
recreating an environment based on the images captured by a self driving car
camera or creating a hologram of a product from an online shopping store might
not need this level of precision. Instead, the greater focus in these type of
applications may be to correctly reproduce the image visually. A novel metric
is needed to capture this difference.

OccNet’s performance was not good in reconstructing 3D objects from single
view 2D images because it did not get any 3D information about the input. We
exploited this draw back by first generating 3D information about the image in
the form of point clouds and then using this 3D information for reconstruction
purposes which significantly improved our output quality.

Figure 4:  Qualitative comparison between OccNet and D-OccNet

##  V Conclusion

In this work, we designed a method for a detailed and precise 3D
reconstruction from a single 2D perspective of an object. Our three-step
network converted a 2D image into a point cloud, and then created a 3D object
using the point cloud as input. We showed that D-OccNet could generate better
3D objects, in terms of visual detail compared to the current state of the art
method.

## Acknowledgment

This research primarily took place in Pakistan Institute of Engineering and
Applied Sciences, Pakistan. Minhaj Uddin Ansari and Talha Bilal were both
undergraduate students at that time. Minhaj is now pursuing masters in Queen’s
University, Canada. Talha Bilal is now pursuing masters under the Erasmus
program. We (Minhaj and Talha) are grateful to Pakistan Institute of
Engineering and Applied Sciences for providing us the resources, and for the
continuous support and motivation given to us by our supervisor Dr. Naeem
Akhter.

## References

  * [1] C. B. Choy, D. Xu, J. Gwak, K. Chen, and S. Savarese, “3d-r2n2: A unified approach for single and multi-view 3d object reconstruction,” _CoRR_ , vol. abs/1604.00449, 2016.
  * [2] H. Fan, H. Su, and L. Guibas, “A point set generation network for 3d object reconstruction from a single image,” _2017 IEEE Conference on Computer Vision and Pattern Recognition (CVPR)_ , pp. 2463–2471, 2017.
  * [3] T. Groueix, M. Fisher, V. G. Kim, B. C. Russell, and M. Aubry, “Atlasnet: A papier-mâché approach to learning 3d surface generation,” _CoRR_ , vol. abs/1802.05384, 2018.
  * [4] K. He, X. Zhang, S. Ren, and J. Sun, “Deep residual learning for image recognition,” _CoRR_ , vol. abs/1512.03385, 2015.
  * [5] D. Jack, J. K. Pontes, S. Sridharan, C. Fookes, S. Shirazi, F. Maire, and A. P. Eriksson, “Learning free-form deformations for 3d object reconstruction,” _CoRR_ , vol. abs/1803.10932, 2018.
  * [6] C. Kong, C.-H. Lin, and S. Lucey, “Using locally corresponding cad models for dense 3d reconstructions from a single image,” in _Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition (CVPR)_ , July 2017\.
  * [7] Y. Liao, S. Donne, and A. Geiger, “Deep marching cubes: Learning explicit surface representations,” in _IEEE Conference on Computer Vision and Pattern Recognition (CVPR)_. IEEE Computer Society, 2018.
  * [8] C. Lin, C. Kong, and S. Lucey, “Learning efficient point cloud generation for dense 3d object reconstruction,” _CoRR_ , vol. abs/1706.07036, 2017.
  * [9] W. E. Lorensen and H. E. Cline, “Marching cubes: A high resolution 3d surface construction algorithm,” p. 163–169, 1987.
  * [10] P. Mandikal, N. K. L., M. Agarwal, and R. V. Babu, “3d-lmnet: Latent embedding matching for accurate and diverse 3d point cloud reconstruction from a single image,” _CoRR_ , vol. abs/1807.07796, 2018.
  * [11] L. M. Mescheder, M. Oechsle, M. Niemeyer, S. Nowozin, and A. Geiger, “Occupancy networks: Learning 3d reconstruction in function space,” _CoRR_ , vol. abs/1812.03828, 2018.
  * [12] C. R. Qi, H. Su, K. Mo, and L. J. Guibas, “Pointnet: Deep learning on point sets for 3d classification and segmentation,” _CoRR_ , vol. abs/1612.00593, 2016.
  * [13] N. Wang, Y. Zhang, Z. Li, Y. Fu, W. Liu, and Y. Jiang, “Pixel2mesh: Generating 3d mesh models from single RGB images,” _CoRR_ , vol. abs/1804.01654, 2018\.
  * [14] J. Wu, Y. Wang, T. Xue, X. Sun, B. Freeman, and J. Tenenbaum, “Marrnet: 3d shape reconstruction via 2.5d sketches,” in _Advances in Neural Information Processing Systems_ , 2017, pp. 540–550.
  * [15] J. Wu, C. Zhang, T. Xue, W. T. Freeman, and J. B. Tenenbaum, “Learning a probabilistic latent space of object shapes via 3d generative-adversarial modeling,” ser. NIPS’16, 2016, p. 82–90.
  * [16] Z. Wu, S. Song, A. Khosla, F. Yu, L. Zhang, X. Tang, and J. Xiao, “3d shapenets: A deep representation for volumetric shapes.” in _CVPR_. IEEE Computer Society, 2015, pp. 1912–1920.
  * [17] H. Xie, H. Yao, X. Sun, S. Zhou, S. Zhang, and X. Tong, “Pix2vox: Context-aware 3d reconstruction from single and multi-view images,” _CoRR_ , vol. abs/1901.11153, 2019.
  * [18] X. Yan, J. Yang, E. Yumer, Y. Guo, and H. Lee, “Perspective transformer nets: Learning single-view 3d object reconstruction without 3d supervision,” _CoRR_ , vol. abs/1612.00814, 2016.
  * [19] Y. Zhao, T. Birdal, H. Deng, and F. Tombari, “3d point-capsule networks,” _CoRR_ , vol. abs/1812.10775, 2018.


# Differentiable Hybrid Traffic Simulation

Sanghyun Son  University of Maryland, College ParkUSA<EMAIL_ADDRESS>,  Yi-
Ling Qiao  University of Maryland, College ParkUSA<EMAIL_ADDRESS>,  Jason
Sewall  NVIDIAUSA<EMAIL_ADDRESS>and  Ming C. Lin  University of
Maryland, College ParkUSA<EMAIL_ADDRESS>
###### Abstract.

We introduce a novel differentiable hybrid traffic simulator, which simulates
traffic using a hybrid model of both macroscopic and microscopic models and
can be directly integrated into a neural network for traffic control and flow
optimization. This is the first differentiable traffic simulator for
macroscopic and hybrid models that can compute gradients for traffic states
across time steps and inhomogeneous lanes. To compute the gradient flow
between two types of traffic models in a hybrid framework, we present a novel
intermediate conversion component that bridges the lanes in a differentiable
manner as well. We also show that we can use analytical gradients to
accelerate the overall process and enhance scalability. Thanks to these
gradients, our simulator can provide more efficient and scalable solutions for
complex learning and control problems posed in traffic engineering than other
existing algorithms. Refer to https://sites.google.com/umd.edu/diff-hybrid-
traffic-sim for our project.

Traffic Simulation, Differentiable Programming, Machine Learning

††submissionid: 372††copyright: acmcopyright††journal: TOG††journalyear:
2022††journalvolume: 41††journalnumber: 6††article: 258††publicationmonth:
12††doi: 10.1145/3550454.3555492††journal: TOG††ccs: Computing methodologies
Agent / discrete models††ccs: Computing methodologies Continuous
simulation††ccs: Computing methodologies Continuous models††ccs: Computing
methodologies Real-time simulation††ccs: Computing methodologies Multiscale
systems††ccs: Computing methodologies Simulation by animation††ccs: Computing
methodologies Simulation environments

##  1\. Introduction

Automobile traffic is a dynamic phenomenon that emerges from the interactions
of a multitude of distinct entities, and its complexity has led analysts to
rely on simulations to solve complex real-world traffic problems, such as
signal control, congestion, road design, and urban planning (Lieberman and
Rathi, 1997). Traffic simulation has become increasingly more important due to
population growth, fast-growing technologies, and demands for autonomous
driving. Autonomous driving agents need accurate and efficient models to
simulate all types of plausible traffic scenarios to capture corner cases and
accelerate learning-based training that requires data not easily available
from real-world capturing (Akhauri et al., 2020).

Techniques for traffic simulation can be broadly classified as either
_macroscopic_ or _microscopic_ based on their modeling assumptions.
Macroscopic models describe traffic evolution as a system of partial
differential equations (PDE), and the traffic state is represented as a
continuum of values across the road network. Since we treat vehicles as
particles translated by convection under this assumption, it is
computationally efficient but provides coarse simulation results. In contrast,
microscopic models represent traffic through individual vehicles, or agents,
that are evolved individually and which together characterize the traffic
state of the road network (Sewall et al., 2011). Therefore, this approach
gives us fine-grained details but often requires much more computation
resources than macroscopic models.

Figure 1.  Traffic simulation in an urban environment. Traffic simulation can
be used to analyze complex traffic dynamics. Among two popular traffic models,
by simulating areas of interest with the microscopic one and using macroscopic
simulation elsewhere, we can reduce the overall computational cost without
compromising significant details. For instance, the center of intersection can
be simulated with discrete, agent-based models for higher-fidelity inter-
vehicle dynamic interaction. (This scene is rendered in CARLA (Dosovitskiy et
al., 2017).)

In this work, we propose to adopt both microscopic and macroscopic models to
create a more general traffic simulator that is differentiable (Figure 1).
This kind of hybrid approach has already been proven to simulate large-scale
traffic environments much more efficiently than either of the approaches while
maintaining integrity (Sewall et al., 2011). Here we leverage the power of the
hybrid approach to support large-scale traffic scenes, which will be
demonstrated in our experiments. During simulation, our framework computes
gradient information that provides abundant insight into traffic dynamics.
This approach results in enhanced sample efficiency in our simulation, as each
sample of traffic simulation comes with gradient information explaining how
and why such an event has occurred. Moreover, this differentiability allows us
to integrate our framework with neural networks to support a wide range of
traffic control, planning, management, and flow optimization problems.

We first derive the analytical formulation of the differentiable traffic
models for both macroscopic and microscopic paradigms. The resulting
representations can simulate traffic on road networks for any metropolitan
region more efficiently. Furthermore, we present a probabilistic method to
convert traffic states between macroscopic and microscopic ones, which is also
differentiable. With this novel differentiable traffic simulation model, we
can find near-optimal solutions to some traffic problems that were otherwise
not previously solvable, and we can do so efficiently. To summarize, in this
paper we introduce the following main results:

  * •

Derivation of a differentiable macroscopic traffic model based on ARZ model
(Aw and Rascle, 2000; Zhang, 2002) (Sec. 3);

  * •

Analytical formulation of differentiable Intelligent Driver Model

(Treiber et al., 2000) for agent-based vehicle dynamics (Sec. 4);

  * •

Differentiable conversion between deterministic and probabilistic
instantiation of two traffic models (Sec. 5);

  * •

Application of the first differentiable hybrid traffic simulation framework to
traffic control (Sec. 6).

In our experiments, we observe up to an order of magnitude speedup in runtime
performance, in both forward simulation and back propagation process, compared
to a baseline differentiable simulator based on automatic differentiation.

##  2\. Previous Work

###  2.1. Traffic models

The preponderance of macroscopic models are based on the _density_ and _flux_
of the traffic flow; (Lighthill and Whitham, 1955; Richards, 1956) suggested
one of the earliest models based on this idea; the LWR model is a non-linear
PDE on the density of vehicles. To overcome the limitations of this model,
(Payne, 1971; Whitham, 2011) added a momentum term to the LWR model, which
allowed for more complex traffic dynamics than before. However, modeling
traffic flow as isotropic (Cassidy and Windover, 1995; Daganzo, 1995) led to
non-physical behavior, such as negative velocity. Based on the observation
that traffic flow is anisotropic, (Aw and Rascle, 2000; Zhang, 2002) modified
the momentum term in the previous model. The ARZ model has been used to
visually simulate the traffic flow (Sewall et al., 2010), and to solve traffic
congestion problems (Yu and Krstic, 2019). Because of the basic assumption,
the computational cost of these macroscopic models is proportional to the
length of the simulated road, not the number of vehicles therein. Therefore,
they are often more computationally efficient than the other models for large-
scale environments. However, they are not suitable for simulations where fine
details, such as the behavior of individual vehicles, are needed.

Microscopic models describe the motion of individual vehicles, typically
tracking velocity, the bumper-to-bumper distance to its leading vehicle, and
the relative velocity between them (Kesting et al., 2007). The research
history dates back to (Gazis et al., 1959, 1961) and (Newell, 1961). Modern
variations of the model often include behavioral traits of individual
vehicles, which make the simulation more realistic and descriptive. (Gipps,
1981; Bando et al., 1995; Treiber et al., 2000; Jiang et al., 2001) are such
models. These models are widely adopted in current agent-based traffic
simulators, such as (Lopez et al., 2018). Since every individual vehicle
observes its surrounding environments and decides its actions in these models,
we can get more fine-grained details of vehicle motions than in the
macroscopic models. However, their computational cost is proportional to the
number of vehicles in the scene, which makes them harder to be adopted for
large-scale simulations.

Hybrid models seek to combine these models in various fashions to take
advantage of their complementary properties; (Magne et al., 2000; Bourrel and
Lesort, 2003; Mammar et al., 2006; Sewall et al., 2011) fall into this
category and show the computational gains we can get from these models. This
hybridization is commonly a spatial one, with different parts of a network
running under distinct regimes. That is, these models often simulate only the
regions of interest with microscopic models for accuracy and use macroscopic
models elsewhere to maintain the overall flow correctness while attaining
computational efficiency. This typically necessitates a technique for
converting between different traffic flow representations (Bourrel and Lesort,
2003). In this work, we adopt this hybrid approach to support large-scale
scenarios and let our simulator be more general so that users can select the
simulation modes based on their needs.

###  2.2. Differentiable models

Differentiable models have been widely used in graphics and robotics
applications like visualization (Nimier-David et al., 2019; Li et al., 2018),
design (Du et al., 2020; Cascaval et al., 2021), and control (Heiden et al.,
2021). This paradigm enables gradient information to flow across complicated
functions, and facilitates the machine learning process by improving sample
efficiency (Mora et al., 2021; Shen et al., 2021). They often show better
results in finding optimal solutions for various high-dimensional problems (Ma
et al., 2021; Hu et al., 2020).

Related to traffic simulation, there especially have been a number of
differentiable simulations for a variety of systems including rigid bodies (de
Avila Belbute-Peres et al., 2018; Qiao et al., 2020), ariticulated bodies
(Qiao et al., 2021b), soft bodies (Geilinger et al., 2020; Du et al., 2021;
Qiao et al., 2021a), cloth (Li et al., 2022), and fluids (Holl et al., 2020;
Takahashi et al., 2021).

Because many traffic models are differentiable in nature, we can apply this
technique to traffic simulation and control. Recently, (Andelfinger, 2021) has
presented a differentiable agent-based traffic simulation, and showed how it
can be used to solve traffic signal control problems. However, the scope of
this work was limited to microscopic models for traffic, and the
implementation relies on automatic differentiation, which can be
computationally inefficient. In addition to traffic, this hybrid
differentiable simulation framework can potentially generalize to other large-
scale, multi-agent systems, such as crowds, insects, fluids, grains, etc. (He
et al., 2020; Colas et al., 2022; Hädrich et al., 2021; Ishiwaka et al., 2021)
as well.

##  3\. Macroscopic Model

In a macroscopic traffic model, traffic evolution is described by a system of
Partial Differential Equations (PDEs) on the road network.

###  3.1. Formulation

The ARZ (Aw and Rascle, 2000; Zhang, 2002) model describes vehicle flow in a
single lane of traffic through the following system of differential equations
in one dimension:

(1a)
$\mathbf{q}_{t}+\mathbf{f(q)}_{x}=0,\quad\mathbf{f(q)}=\mathbf{q}u=\begin{bmatrix}\rho
u\\\ yu\end{bmatrix}$ where
$\mathbf{q}(x,t)=\left[\rho(x,t),y(x,t)\right]^{\intercal}$; $\rho$ denotes
the density of traffic (cars per car length), $y$ is the relative flow of
traffic, and $u$ is the velocity of traffic. As per the notation, these terms
are observed at a certain position $x$ at a time $t$, and they are related to
each other as follows: (1b) $\displaystyle y(\rho,u)$
$\displaystyle=\rho\left(u-u_{\mathrm{eq}}(\rho)\right)$ (1c) $\displaystyle
u_{\mathrm{eq}}(\rho)$ $\displaystyle=u_{\max}\left(1-\rho^{\gamma}\right)$

where $u_{\max}$ denotes the maximum velocity, or speed limit, of the given
lane, and $0<\gamma<1$ is a constant tuning parameter. Conceptually,
$u_{\mathrm{eq}}$ represents a ‘comfortable’ velocity based on density. We
used $\gamma=0.5$ in all of the experiments.

Note that the system (1a) is a conservative system, where the sum of $\rho$
and $y$ over the entire lane is conserved, excluding the incoming and outgoing
fluxes at the borders of the lane.

###  3.2. Numerical Solution

Conservation laws like the ARZ system of PDEs are typically discretized with
the Finite Volume Method (FVM) (LeVeque et al., 2002) and integrated
explicitly in time. The discretized quantities in a lane are
$\mathbf{Q}_{i}^{n}$, where $n$ refers to the time step and $i$ refers to the
$i^{\mathrm{th}}$ cell. The solution procedure follows:

  1. (1)

Compute wave speeds $\lambda$ and fluxes $\mathbf{f}$ (from Eq. (1a)) by
solving the Riemann problem at the interface of each adjacent cell.

  2. (2)

Compute global time step $\Delta t\leq\max\left|\lambda\right|$; in this work
we choose a conservative, constant $\Delta t$.

  3. (3)

Integrate: $\mathbf{Q}^{n+1}_{i}=\mathbf{Q}^{n}_{i}-\frac{\Delta t}{\Delta
x}\left[\mathbf{f}(\mathbf{q}(b))-\mathbf{f}(\mathbf{q}(a))\right]$.
$\mathbf{q}(a)$ and $\mathbf{q}(b))$ are the intermediate states identified
between $\left(\mathbf{Q}^{n}_{i-1},\mathbf{Q}^{n}_{i}\right)$ and
$\left(\mathbf{Q}^{n}_{i},\mathbf{Q}^{n}_{i+1}\right)$ respectively. $\Delta
x$ is the length of cell $i$.

Details can be found in (Sewall, 2011); Appendix A.1 reproduces the six
relevant cases for the ARZ Riemann problems.

###  3.3. Differentiation

Our differential technique built on this requires analytical gradients of the
above; these gradients form the basis of still other gradients that are used
in solving more complex traffic control problems.

In our update scheme, we compute the gradients of $\mathbf{Q}^{n+1}_{i}$ using
$\mathbf{Q}^{n}_{i-1}$,$\mathbf{Q}^{n}_{i}$, and $\mathbf{Q}^{n}_{i+1}$:

(2a)
$\displaystyle\frac{\partial\mathbf{Q}^{n+1}_{i}}{\partial\mathbf{Q}^{n}_{i-1}}$
$\displaystyle=-\frac{\Delta t}{\Delta
x}\left[-\mathbf{f}^{\prime}(\mathbf{q}(a))\frac{\partial\mathbf{q}(a)}{\partial\mathbf{Q}^{n}_{i-1}}\right]$
(2b)
$\displaystyle\frac{\partial\mathbf{Q}^{n+1}_{i}}{\partial\mathbf{Q}^{n}_{i}}$
$\displaystyle=\mathbf{I}-\frac{\Delta t}{\Delta
x}\left[\mathbf{f}^{\prime}(\mathbf{q}(b))\frac{\partial\mathbf{q}(b)}{\partial\mathbf{Q}^{n}_{i}}-\mathbf{f}^{\prime}(\mathbf{q}(a))\frac{\partial\mathbf{q}(a)}{\partial\mathbf{Q}^{n}_{i}}\right]$
(2c)
$\displaystyle\frac{\partial\mathbf{Q}^{n+1}_{i}}{\partial\mathbf{Q}^{n}_{i+1}}$
$\displaystyle=-\frac{\Delta t}{\Delta
x}\left[\mathbf{f}^{\prime}(\mathbf{q}(b))\frac{\partial\mathbf{q}(b)}{\partial\mathbf{Q}^{n}_{i+1}}\right]$

We thus need the Jacobian $\mathbf{f}(\mathbf{q})$ and the partial derivatives
$\frac{\partial\mathbf{q}(a)}{\partial\mathbf{Q}^{n}_{i-1}}$,
$\frac{\partial\mathbf{q}(a)}{\partial\mathbf{Q}^{n}_{i}}$,
$\frac{\partial\mathbf{q}(b)}{\partial\mathbf{Q}^{n}_{i}}$,
$\frac{\partial\mathbf{q}(b)}{\partial\mathbf{Q}^{n}_{i+1}}$; this arises from
differentiating solutions to Riemann problems as shown in Appendix A.2. The
analytical gradients computed here play an important role in accelerating our
traffic simulator, as shown in Section 6.

####  3.3.1. Time step sizes

We have derived these gradients under the assumption that $\Delta t$ is
constant and satisfies the CFL and stability conditions (LeVeque et al., 2002)
for the entire simulation. It would be possible to have a dynamic $\Delta t$,
since this is also differentiable (the speed $\lambda$ is determined by the
eigenvalues of the flux function $\mathbf{f}$). However, for the sake of the
hybrid approach we ultimately use, we use a constant $\Delta t$ to ensure that
the microscopic simulation is stable.

####  3.3.2. Continuity issues

We have computed the gradients above assuming that the update scheme is
differentiable, and in fact, we can see that it is continuous and piece-wise
differentiable across most cases suggested in Appendix A.

We do not prove it thoroughly here, as it is generally trivial when we keep in
mind that $\lim_{u_{l}\rightarrow u_{r}}\mathbf{q}_{m}=\mathbf{q}_{l}$.
However, Case 1 is exceptional; there is a shock between the phase states
$q_{l}$ and $q_{m}$, which makes the two states discontinuous. Specifically,
we cannot guarantee that $q_{m}$ converges to $q_{l}$ when $\lambda_{s}<0$
converges to $0$. Therefore, the gradients we compute when $\lambda_{s}$ is
near zero could become unstable. However, in most cases in our simulation,
these cases rarely happened and did not affect the quality of our solution
even when we used the possibly unstable gradients as they are; it is possible
that numerical viscosity plays a role here.

##  4\. Microscopic Model

In a microscopic traffic model, traffic flow is described by interactions
between multiple discrete vehicles that follow certain rules. Here we use the
Intelligent Driver Model (IDM) (Treiber et al., 2000) to simulate such
vehicles.

###  4.1. Formulation

Under the microscopic viewpoint, we can describe the state of a discrete
vehicle as

(3) $\mathbf{q}(t)=\begin{bmatrix}p(t)\\\ v(t)\end{bmatrix},$

where $p$ and $v$ denote the position and velocity of the vehicle respectively
for a given time $t$.

Using these states to describe discrete vehicles occupying a lane, IDM
describes the motion of each individual vehicle based on its relationship with
the vehicle directly ahead of it; the _leading vehicle_. For the
$i^{\mathrm{th}}$ vehicle, let $h(i)^{\mathrm{th}}$ its leading vehicle.
According to the IDM, the acceleration of the $i^{\mathrm{th}}$ vehicle is
determined as follows:

(4a) $\displaystyle{\Delta p}$
$\displaystyle=p_{h(i)}(t)-p_{i}(t)-{length}_{h(i)}$ (4b)
$\displaystyle{\Delta v}$ $\displaystyle=v_{i}(t)-v_{h(i)}(t)$ (4c)
$\displaystyle s_{\mathrm{opt}}$
$\displaystyle=s_{\min}+v_{i}(t)T_{\mathrm{pref}}+\frac{v_{i}(t){\Delta
v}}{2\sqrt{a_{\max}a_{\mathrm{pref}}}}$ (4d) $\displaystyle a_{i}(t)$
$\displaystyle=a_{\max}\left[1-{\left(\frac{v_{i}(t)}{v_{\mathrm{targ}}}\right)}^{\delta}-{\left(\frac{s_{\mathrm{opt}}}{\Delta
p}\right)}^{2}\right]$

There are various hyperparameters included in the model. They characterize a
vehicle’s motion; see Appendix B.1 for more details. In our experiments, we
randomly initialized those hyperparameters for every single vehicle. Note that
we compute $s_{\mathrm{opt}}$ before we compute the acceleration term. It
represents the _optimal_ space that the vehicle should have to avoid collision
to its leading vehicle.

With the acceleration $a_{i}(t)$, we can update the vehicle’s state using
Euler’s method.

###  4.2. Differentiation

Based on this update scheme, we can compute the analytical gradients of the
IDM. Since the derivation is straightforward, we offer the analytical
gradients in Appendix B.2. Note that compared to the macroscopic model, there
is no significant theoretical challenge that we have to deal with in the IDM.

##  5\. Macro-Micro-Macro Conversion

Our traffic simulation can transition between the continuous, flow-based
macroscopic model and the discrete, agent-based microscopic model. We need to
consider the differentiation across the interfaces.

###  5.1. Macro to Micro

Given a lane in the macroscopic regime that flows into a microscopic lane, we
must convert the continuous representations therein to discrete vehicles at
their junction. This can be done in a deterministic or stochastic fashion.

####  5.1.1. Deterministic Instantiation

Assume that $\rho(t)$ is the density of the traffic flow at the end of a macro
lane on time $t$, and $v(t)$ is the velocity. We can set a _flux capacitor_
(c.f. (Sewall, 2011)) at the interface; this counts the total number of
vehicles $N_{d}(t)$ that have reached this point:

(5) $N_{d}(t)=\left\lfloor\int_{0}^{t}\rho(t)v(t)dt\right\rfloor$

where $\left\lfloor x\right\rfloor$ is the floor function that gives the
greatest integer less than or equal to the real number $x$. $N_{d}(t)$ will
then start from 0, and we will instantiate a vehicle when it increases by 1.
The velocity of the instantiated vehicle is set to be the macroscopic traffic
flow’s velocity.

####  5.1.2. Stochastic Instantiation

Another possible interpretation of the macro to micro conversion is as a
Poisson process. The flux of the traffic flow can be viewed as the intensity
of a Poisson process. In each time step of length $dt$, the number of
instantiated vehicles $X(t)$ follows a Poisson distribution $X(t)\sim
P(\rho(t)v(t)dt)$, $P(X=k)=e^{-\rho(t)v(t)dt}(\rho(t)v(t)dt)^{k}k!^{-1}$.
Similarly, the total number of vehicles $N_{s}(t)$ can be expressed as

(6) $N_{s}(t)=\int_{0}^{t}X(t)$

As above, the velocity of the emitted discrete vehicle is the velocity of the
macroscopic state at the junction of the two lanes.

###  5.2. Micro to Macro

Compared to the macro to micro conversion, agent-based information can be
converted to a continuum representation more simply. The density of a cell at
the interval $(l,r]$ in a macro lane can be defined as

(7)
$\displaystyle\rho(l,r)=\frac{1}{r-l}\sum_{i=1}^{n}\mathbf{1}_{(l,r]}(p_{i})$

where $\mathbf{1}(\cdot)$ is the indicator function that identifies vehicles
that lie in the interval, $n$ is the number of vehicles and $p_{i}$ is the
position of $i^{th}$ vehicle. The velocity of this cell can be the average
velocities of all the vehicles therein.

###  5.3. Differentiation

We have described two types of discrete processes: instantiation of vehicles,
and the indicator function. Automatic differentiation is able to compute the
gradients of velocities $v$ but cannot handle the gradients of densities
$\rho$. To backpropagate derivatives to the densities $\rho$, we first create
an ancillary variable $a_{i}=1$ for each discrete vehicle. We rewrite Equation
7 as:

(8)
$\displaystyle\rho(l,r)=\frac{1}{r-l}\sum_{i=1}^{n}\mathbf{1}_{(l,r]}(p_{i})\cdot
a_{i}$

Since we always set $a_{i}=1$, the equation constantly holds for forward
simulation. But $a_{i}$ can also receive the gradients for discrete vehicles
from the following macro lanes. Let $L$ be the loss function; then for each
vehicle in $(l,r]$:

(9)
$\frac{\partial{L}}{\partial{a_{i}}}=\frac{\partial{L}}{\partial{\rho(l,r)}}$

The next challenge for differentiation is the vehicle instantiation process in
the macro to micro conversion. We will start with the deterministic strategy:
assume that the total number of vehicles reach $n$ and $n+1$ at time $t_{1}$
and $t_{2}$, respectively. We therefore have the following equation

(10) $n+1=n+\int_{t_{1}}^{t_{2}}\rho(t)v(t)dt=n+a_{i}$

where the $i^{th}$ vehicle is the one instantiated at $t_{2}$. So for
$t\in(t_{1},t_{2}]$,

(11)
$\frac{\partial{L}}{\partial{\rho(t)}}=\frac{\partial{L}}{\partial{a_{i}}}v(t)dt$

For the stochastic case, we can first compute the expectation value of the
total number of vehicles according to the Poisson distribution

(12) $E[N_{s}(t)]=\int_{0}^{t}E[X(t)]=\int_{0}^{t}\rho(t)v(t)dt.$

We notice that the expectation value of $N_{s}(t)$ is the same as the flux
capacitor in the deterministic strategy. Similarly, we can find time $t_{1}$
and $t_{2}$ when the total number of vehicles reach $n$ and $n+1$. For all
time step $t\in(t_{1},t_{2}]$,
$\frac{\partial{L}}{\partial{\rho(t)}}=\frac{\partial{L}}{\partial{a_{i}}}v(t)dt$.
Intuitively, if $\frac{\partial{L}}{\partial{a_{i}}}\leq 0$, it means the
density is higher than the desired level and we need to decrease the intensity
of the Poisson process. We will also validate in our numerical experiments
that such estimation of gradients can effectively optimize our objective
function.

##  6\. Experimental Results

We show the effectiveness of our differentiable hybrid traffic simulator with
application to solve a variety of traffic problems. We first show that we can
accelerate the differentiable traffic simulator with the analytical gradients
that we have computed above by comparing it against a baseline simulator using
automatic differentiation. At the same time, we will justify the use of the
hybrid model by comparing its computation time with other approaches in large-
scale scenarios. Then, we prove the correctness and efficacy of our analytical
gradients by solving parameter estimation problems. Lastly, we conduct
experiments with traffic control problems to illustrate how we can integrate
our simulator with neural networks and how we can use our simulator to solve
real-world traffic problems.

We have implemented our traffic simulator with Python, and used PyTorch 1.9
(Paszke et al., 2019) for automatic differentiation. All experiments were run
on an Intel® Xeon® W-2255 CPU @ 3.70GHz, and traffic rendering was mainly
based on our own implementation on Unity Engine.

###  6.1. Acceleration

####  6.1.1. Analytical Gradients

Here we present how our analytical gradients can contribute to the
acceleration of both forward simulation and backward propagation in our
simulator. We have compared our simulator to the baseline simulator that
relies on automatic differentiation to compute the gradients in the system.
For all of the simulations, we have computed the gradients of the states at
the last time step with respect to the states at the initial time step. For
macroscopic and microscopic simulations, we used a single lane. For hybrid
simulation, we used 3 lanes, which is the minimum case to cover both macro-to-
micro and micro-to-macro conversion. The first and third lanes are set as
macroscopic lanes and the second one is a microscopic lane. To show how the
computation time changes with the scale of the simulation, we have
experimented with different scales of settings.

Table 1 shows that our analytical gradient allows us to run both forward
simulation (FW) and back propagation (BP) of the system much faster than the
baseline framework. In macroscopic and hybrid simulations, the speedup was up
to 4x times for forward simulation, and 18x for back propagation. In the
microscopic simulation, the speedup was up to 6x times for forward simulation,
and 7x times for back propagation. This shows the computational efficiency of
our framework, which is necessary to implement large-scale traffic simulation.

####  6.1.2. Hybrid Model

In this paper, we selected the hybrid model to run large-scale traffic
simulations in interactive time, while not losing their fidelity. To justify
this design choice, here we present the average frame rate for simulating a
large environment, which includes approximately 10K vehicles. Note that we can
specify the ratio $\epsilon$ of the vehicles which would be simulated using a
microscopic model; $\epsilon$ could be regarded as an interpolant between the
two models, where $\epsilon=0$ means macroscopic and $\epsilon=1$ means
microscopic approach. We measured the average frame rate for different
$\epsilon$ using single core; see Figure 2. Note that the average frame rates
decreases as $\epsilon$ increases and thus the simulation becomes more
”microscopic”.

Table 1. Comparison of our technique against automatic differentiation for
various types of simulations. Scale denotes number of time steps taken and the
number of cells (for macroscopic, and hybrid)/ number of vehicles
(microscopic). FW: forward simulation; BP: backpropagation. Ours outperforms
autodifferentiation in both forward passes and backpropagation, by up to a
factor of 18.77, 7.84, and 15.27 in macroscopic, microscopic, and hybrid
traffic simulation, respectively.
| Scale | 10/1K | 50/5K
---|---|---|---
Macroscopic | FW(Auto) | 5.28s$\pm$0.16s | 125.88s$\pm$5.66s
FW(Ours) | 1.42s$\pm$0.04s | 30.06s$\pm$0.88s
Speedup | 3.71x | 4.19x
BP(Auto) | 2.42s$\pm$0.20s | 71.08s$\pm$6.01s
BP(Ours) | 0.18s$\pm$0.03s | 3.79s$\pm$0.25s
Speedup | 13.73x | 18.77x
Microscopic | FW(Auto) | 1.54s$\pm$0.08s | 66.74s$\pm$1.67s
FW(Ours) | 0.45s$\pm$0.02s | 10.86s$\pm$0.24s
Speedup. | 3.42x | 6.15x
BP(Auto) | 1.58s$\pm$0.06s | 45.36s$\pm$0.95s
BP(Ours) | 0.20s$\pm$0.03s | 7.10s$\pm$0.42s
Speedup | 7.84x | 6.39x
Hybrid | FW(Auto) | 11.81s$\pm$0.42s | 310.67s$\pm$15.02s
FW(Ours) | 3.41s$\pm$0.06s | 66.52s$\pm$0.47s
Speedup | 3.47x | 4.67x
BP(Auto) | 2.70s$\pm$0.11s | 111.06s$\pm$3.09s
BP(Ours) | 0.38s$\pm$0.02s | 7.27s$\pm$0.13s
Speedup | 7.13x | 15.27x
0.010.020.050.10.20.5$1$$2$$3$$\epsilon$ = Ratio of microscopic vehiclesFrame
Per Second (FPS) Figure 2. Comparison on the average frame rate (FPS) under
the hybrid model, with different $\epsilon$. The total number of vehicles is
set to be 10K, and the $\epsilon$ value denotes the ratio of microscopic
vehicles among them.

|
---|---
|

(a) Macroscopic

(b) Microscopic

(c) Hybrid (Macro-Micro-Macro)

(d) Hybrid (Micro-Macro-Micro)

Figure 3. Parameter estimation problem. For given final states of the system,
one has to estimate the initial states that lead to the final states. (a)–(d)
Our results rendered in blue exhibit much faster convergence to the correct
solutions than the other state-of-the-art gradient-free optimization methods
by up to one order of magnitude.

###  6.2. Parameter Estimation Problem

One of the most fundamental problems that we can solve with our framework is
the parameter estimation problem. There are many parameters involved in our
system; the speed limit is the most typical example. In our formulation, we
focus on retrieving preceding states in the system based on given subsequent
ones. The simplest form of this problem would be, assuming $0\leq t\leq 1$,
estimating the initial state $\mathbf{q}(0)$ given the last state
$\mathbf{q}(1)$. This problem is considered to be the most basic formulation
to solve, as it only requires fundamental gradients computed analytically.

In our experiments, the simulation length was set as 10 seconds for every
setting. Theoretically, we can use any simulation lengths for this experiment,
but we found out that too long simulation lengths make this problem much
harder than it should be and 10 seconds is a proper length to prove our
approach’s efficacy. Under this setting, we have to estimate the initial
traffic states that would end up in the given final traffic states in 10
seconds. For a proposed estimate $\mathbf{q}_{\mathrm{est}}(0)$, we can
compute the estimation error with the following loss function, where $n$
denotes the size of the state vector:

$L=\frac{1}{n}\left|\mathbf{q}(1)-\mathbf{q}_{\mathrm{est}}(1)\right|^{2}$

We have compared our framework with other gradient-free optimization
algorithms: CMA-ES (Hansen, 2006), SLSQP (Kraft et al., 1988), and Nelder-Mead
(Nelder and Mead, 1965) algorithm. Experiments were run five times with
randomly initialized initial states for each of the algorithms and for each of
the simulation modes. In the hybrid setting, we used 3 sequentially connected
lanes. One of the experiments was done for macro-micro-macro lanes, and the
other was done for micro-macro-micro lanes. The loss graph for each setting is
shown in Figure 3. We can observe that our framework retrieves more precise
initial states faster than other algorithms for every setting.

###  6.3. Traffic Control Problems

Our differentiable traffic simulator can be integrated with neural networks to
enhance its ability to learn and control a given task. To test our framework’s
capability, we have experimented with the “intersection signal control
problem” (ITSCP) for macroscopic and hybrid simulation, and the pace car
problem for microscopic simulation.

For all the experiments, we used two fully connected layers as our neural
network, where each of the layers has 256 nodes. The neural network is trained
to emit proper control input, which would be fed into our simulator. Then we
can compute the gradient of the aggregated reward of each episode with respect
to each control input therein, and apply gradient descent to directly maximize
the reward. Note that this approach is same as analytical policy gradient
(APG) method. Then we compare our results to other baseline gradient-free RL
algorithms, as they are widely used to solve these kinds of control problems.

Figure 4. 4x4 Grid of Intersections. Hybrid traffic simulation at a 4x4 grid
of intersections: the density of each cell is rendered in red, which shows the
volume of traffic agents to cross the intersection. One has to optimize time
allocations for traffic lights at each intersection to maximize overall
traffic flow.

|
---|---

(a) Macroscopic

(b) Hybrid

Figure 5.  ITSCP. (a), (b) Learning graphs for macroscopic and hybrid
simulations. Our framework shows better convergence rates and results than
other RL algorithms.

---

(a) Time = 0 sec, Initial States

(b) Time = 10 sec, Target Speed Limit = 30 m/sec

(c) Time = 20 sec, Target Speed Limit = 10 m/sec

Figure 6.  Pace Car Problem. Sequence of simulation images in time order, from
top to bottom. (a) A pace car rendered in green guides following vehicles to
maintain a target speed limit, which varies over time. It fails when it runs
too fast, too slow, or goes out of the lane. (b) When a high target speed
limit (30 m/sec) is given, the pace car accelerates fast, and keeps a wide
distance from its following vehicles to let them run fast. (c) However, when
the target speed limit is reduced (10 m/sec), it swiftly decelerates and keeps
short distance from the following vehicles, to prevent them from running over
the target speed limit.

####  6.3.1. ITSCP

ITSCP is one of the most widely studied traffic control problems, because of
its significance in controlling traffic congestion in urban environment (Eom
and Kim, 2020). We have devised a scenario that falls into this category,
where an agent has to compute the optimal time allocations for traffic lights
(Figure 4).

In this scenario, there are two different traffic lights at a given
intersection: one is green across the West-East direction (WE-light), and the
other one is across the North-South direction (NS-light). When one of the
lights is turned on, the other one has turned off automatically. In a single
signal phase, the WE-light is turned on first, and then the NS-light is turned
on next. There are multiple signal phases during the entire simulation, and
the agent has to determine optimal time allocations for each light in each
signal phase.

The input to the neural network is given as the schedule of incoming traffic
flow for every lane. Our experiments were inspired by a synthetic traffic
scenario devised for another ITSCP experiment (Wei et al., 2018). We also
followed a conventional scheme to compute the reward (Eom and Kim, 2020). The
reward for an action is computed as the weighted sum of traffic flow that
crosses the intersection ($R_{f}$) and the length of waiting queue in each
lane ($R_{q}$), which is formulated as follows:

$R=c_{1}\cdot R_{f}+c_{2}\cdot R_{q}.$

For our experiments, $c_{1}$ is set as 1 and $c_{2}$ is set as -1, to maximize
the traffic flow and minimize the length of waiting queue. In macroscopic
environments, the length of the queue is computed as the number of cells that
have speed less than certain threshold. In hybrid environments, the traffic
flow is computed as the number of discrete vehicles crossing the intersection,
as the central part of an intersection is simulated with a microscopic model.

For comparison, we used three baseline model-free RL algorithms; DDPG
(Lillicrap et al., 2015), PPO (Schulman et al., 2017), and SAC (Haarnoja et
al., 2018). Note that these algorithms do not use gradient information that
our simulator provides. In contrast, our gradient-based optimization algorithm
uses this information to optimize the objective function. Experiments were run
five times for each algorithm, and Figure 5 shows the learning graphs of both
macroscopic and hybrid environment. In the macroscopic setting, our framework
converged to the near-optimal solution very quickly, while DDPG failed to
learn at all and PPO and SAC learned, but not to the extent of ours. In hybrid
setting, our framework and SAC both succeeded in converging to the best
solution, but ours converged faster than SAC. Also, Table 2 shows that our
best reward is better than any other algorithms for both settings.

---
Figure 7.  Pace Car Learning Graph. Our result rendered in blue shows far better performance (as indicated by higher reward values) than other algorithms, which make no significant gains in learning even after long periods. Table 2. Comparison of maximum rewards from each control problem. Our framework is able to find better solutions (indicated by higher reward values) than the other baseline, model-free RL algorithms (DDPG, PPO, or SAC) for every experiment setting (Macro, Micro, or Hybrid). | Ours | DDPG | PPO | SAC
---|---|---|---|---
Macro | 124.38 | 113.62 | 104.53 | 115.73
Micro | 4482.87 | 4.97 | 3218.19 | 3373.99
Hybrid | 208.48 | 174.18 | 158.18 | 191.51

####  6.3.2. Pace Car Problem

In this problem, we assume there is a discrete pace car, which runs in front
of other vehicles that follow the IDM. The pace car leads the vehicles to
enforce a specified target speed limit, which is often smaller than the lane’s
speed limit (Figure 6).

The optimal solution for this problem is clear to the human observer; all
vehicles should maintain the desired speed limit. However, the pace car has to
determine both acceleration and steering at every frame to remain on the road
and lead the vehicles at the same time. To introduce even more complexity, in
our experiments, we used a speed limit target that varies over time.
Therefore, the neural network receives other vehicle’s states and the speed
limit as input, and has to compute the optimal acceleration and steering for
next several frames, which is quite challenging. The reward is formulated as
follows:

(13) $R=\sum_{i=1}^{t}\sum_{j=1}^{n}C_{\max}-(v_{\mathrm{targ}}-v_{i,j})^{2}$

where $C_{\max}$ denotes an arbitrary constant value, $t$ and $n$ denote
number of frames and number of vehicles respectively. $v_{\mathrm{targ}}$ and
$v_{i,j}$ represent the target speed limit and $j^{\mathrm{th}}$ vehicle’s
speed at the given frame $i$.

We also used three baseline RL algorithms for comparison, each of which were
run 5 times. The learning graph in Figure 7 and the maximum rewards in Table 2
show that our framework achieved far better overall reward than the other
algorithms. PPO and SAC also succeeded in accomplishing initial high rewards
but their performance soon plunged. Our framework continuously learned and
achieved better and improved results.

##  7\. Conclusions

We have proposed a novel differentiable hybrid traffic simulator. In our
simulator, gradients of traffic states across the time steps are computed
analytically and propagated across lanes modeled under complimentary but
deeply different regimes. This was made possible with our novel discrete-
continuous differentiable conversion module. Further, our analytic formulation
is much more efficient than one based on automatic differentiation; therefore,
our technique can offer real-time traffic simulation, as shown in Figure 8,
that would not have been possible otherwise.

We have also shown that we can use the gradients to solve classic traffic
problems. For the parameter estimation problems, our framework is able to find
far better estimates than other gradient-free optimization algorithms. For
control problems, our framework succeeded in finding near-optimal policies,
which was not possible with other model-free RL algorithms.

###  7.1. Generalization

The concept of hybrid simulation is applicable to many different types of
dynamical systems, such as fluids (Mohamed and Mohamad, 2010; Golas et al.,
2012; Narain et al., 2010), crowd (Treuille et al., 2006; Narain et al.,
2009), and many other multi-agent systems (Zheng et al., 2017). And, our
method to differentiate the conversion between macro and micro traffic flow is
applicable to a broad range of dynamical systems that are simulated by hybrid
techniques, transitioning between continuous and discrete simulation domains.
In general, other deterministic or stochastic instantiation process can also
take advantage of our technique introduced in this paper. For example, the
micro and macro representations for crowd share many similarities with
traffic. The conversion between individual pedestrians and the crowd flow is
analogous to the process of instantiating and removing vehicles in the roads.
We also hypothesize that the algorithmic and computational framework on
differentiable hybrid simulation, as proposed here, can also be extended to
hybrid control of complex and autonomous systems (Branicky et al., 1998;
Fierro et al., 2001).

###  7.2. Limitations and Future Directions

Given these promising results, we plan to further augment the capability and
application of our simulator. First, while hybrid techniques are naturally
suited to metropolitan-scale scenarios, we have not yet applied this
differentiable hybrid technique to these; we expect to achieve even higher
performance improvements from these mega-scale experiments. We can also
integrate realistic vehicle motions, such as lane changing, to further capture
complex traffic dynamics in higher fidelity. Finally, we used simulated data
for solving traffic control problems in this paper, but we expect that our
simulator is applicable to real-world data as well.

Figure 8.  Large-scale Hybrid Traffic simulation. These traffic scenes consist
of many hundreds of vehicles across an urban scene simulated and controlled in
real time using a differentiable hybrid traffic model.

###### Acknowledgements.

This research is supported in part by Dr. Barry Mersky and Capital One
E-Nnovate Endowed Professorships, and ARL Cooperative Agreement W911NF2120076.

## References

  * (1)
  * Akhauri et al. (2020) Shivam Akhauri, Laura Y Zheng, and Ming C Lin. 2020.  Enhanced transfer learning for autonomous driving with systematic accident simulation. In _2020 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_. IEEE, 5986–5993.
  * Andelfinger (2021) Philipp Andelfinger. 2021\.  Differentiable Agent-Based Simulation for Gradient-Guided Simulation-Based Optimization. In _Proceedings of the 2021 ACM SIGSIM Conference on Principles of Advanced Discrete Simulation_. 27–38.
  * Aw and Rascle (2000) AATM Aw and Michel Rascle. 2000.  Resurrection of” second order” models of traffic flow.  _SIAM journal on applied mathematics_ 60, 3 (2000), 916–938.
  * Bando et al. (1995) Masako Bando, Katsuya Hasebe, Akihiro Nakayama, Akihiro Shibata, and Yuki Sugiyama. 1995\.  Dynamical model of traffic congestion and numerical simulation.  _Physical review E_ 51, 2 (1995), 1035.
  * Bourrel and Lesort (2003) Emmanuel Bourrel and Jean-Baptiste Lesort. 2003.  Mixing microscopic and macroscopic representations of traffic flow: Hybrid model based on Lighthill–Whitham–Richards theory.  _Transportation Research Record_ 1852, 1 (2003), 193–200.
  * Branicky et al. (1998) Michael S Branicky, Vivek S Borkar, and Sanjoy K Mitter. 1998\.  A unified framework for hybrid control: Model and optimal control theory.  _IEEE transactions on automatic control_ 43, 1 (1998), 31–45.
  * Cascaval et al. (2021) Dan Cascaval, Mira Shalah, Phillip Quinn, Rastislav Bodik, Maneesh Agrawala, and Adriana Schulz. 2021\.  Differentiable 3D CAD Programs for Bidirectional Editing.  _arXiv preprint arXiv:2110.01182_ (2021).
  * Cassidy and Windover (1995) Michael J Cassidy and John R Windover. 1995.  vÍethodology for Assessing Dynamics of Freeway Traffic Flow.  (1995).
  * Colas et al. (2022) Adèle Colas, Wouter van Toll, Katja Zibrek, Ludovic Hoyet, Anne-Hélène Olivier, and Julien Pettré. 2022.  Interaction Fields: Intuitive Sketch-based Steering Behaviors for Crowd Simulation. In _Computer Graphics Forum_.
  * Daganzo (1995) Carlos F Daganzo. 1995\.  Requiem for second-order fluid approximations of traffic flow.  _Transportation Research Part B: Methodological_ 29, 4 (1995), 277–286.
  * de Avila Belbute-Peres et al. (2018) Filipe de Avila Belbute-Peres, Kevin Smith, Kelsey Allen, Josh Tenenbaum, and J Zico Kolter. 2018\.  End-to-end differentiable physics for learning and control.  _Advances in neural information processing systems_ 31 (2018), 7178–7189.
  * Dosovitskiy et al. (2017) Alexey Dosovitskiy, German Ros, Felipe Codevilla, Antonio Lopez, and Vladlen Koltun. 2017.  CARLA: An Open Urban Driving Simulator. In _Proceedings of the 1st Annual Conference on Robot Learning_. 1–16.
  * Du et al. (2021) Tao Du, Kui Wu, Pingchuan Ma, Sebastien Wah, Andrew Spielberg, Daniela Rus, and Wojciech Matusik. 2021.  DiffPD: Differentiable Projective Dynamics with Contact.  _arXiv:2101.05917_ (2021).
  * Du et al. (2020) Tao Du, Kui Wu, Andrew Spielberg, Wojciech Matusik, Bo Zhu, and Eftychios Sifakis. 2020\.  Functional Optimization of Fluidic Devices with Differentiable Stokes Flow.  _ACM Trans. Graph._ (2020).
  * Eom and Kim (2020) Myungeun Eom and Byung-In Kim. 2020.  The traffic signal control problem for intersections: a review.  _European transport research review_ 12, 1 (2020), 1–20.
  * Fierro et al. (2001) Rafael Fierro, Aveek K Das, Vijay Kumar, and James P Ostrowski. 2001.  Hybrid control of formations of robots. In _Proceedings 2001 ICRA. IEEE International Conference on Robotics and Automation (Cat. No. 01CH37164)_ , Vol. 1. IEEE, 157–162.
  * Gazis et al. (1959) Denos C Gazis, Robert Herman, and Renfrey B Potts. 1959\.  Car-following theory of steady-state traffic flow.  _Operations research_ 7, 4 (1959), 499–505.
  * Gazis et al. (1961) Denos C Gazis, Robert Herman, and Richard W Rothery. 1961\.  Nonlinear follow-the-leader models of traffic flow.  _Operations research_ 9, 4 (1961), 545–567.
  * Geilinger et al. (2020) Moritz Geilinger, David Hahn, Jonas Zehnder, Moritz Bächer, Bernhard Thomaszewski, and Stelian Coros. 2020.  ADD: Analytically differentiable dynamics for multi-body systems with frictional contact.  _ACM Transactions on Graphics (TOG)_ 39, 6 (2020).
  * Gipps (1981) Peter G Gipps. 1981\.  A behavioural car-following model for computer simulation.  _Transportation Research Part B: Methodological_ 15, 2 (1981), 105–111.
  * Golas et al. (2012) Abhinav Golas, Rahul Narain, Jason Sewall, Pavel Krajcevski, Pradeep Dubey, and Ming Lin. 2012\.  Large-scale fluid simulation using velocity-vorticity domain decomposition.  _ACM Transactions on Graphics (TOG)_ 31, 6 (2012), 1–9.
  * Haarnoja et al. (2018) Tuomas Haarnoja, Aurick Zhou, Pieter Abbeel, and Sergey Levine. 2018.  Soft actor-critic: Off-policy maximum entropy deep reinforcement learning with a stochastic actor. In _International conference on machine learning_. PMLR, 1861–1870.
  * Hädrich et al. (2021) Torsten Hädrich, Daniel T Banuti, Wojtek Pałubicki, Sören Pirk, and Dominik L Michels. 2021\.  Fire in paradise: mesoscale simulation of wildfires.  _ACM Transactions on Graphics (TOG)_ 40, 4 (2021), 1–15.
  * Hansen (2006) Nikolaus Hansen. 2006\.  The CMA evolution strategy: a comparing review.  _Towards a new evolutionary computation_ (2006), 75–102.
  * He et al. (2020) Feixiang He, Yuanhang Xiang, Xi Zhao, and He Wang. 2020\.  Informative scene decomposition for crowd analysis, comparison and simulation guidance.  _ACM Transactions on Graphics (TOG)_ 39, 4 (2020), 50–1.
  * Heiden et al. (2021) Eric Heiden, Miles Macklin, Yashraj S Narang, Dieter Fox, Animesh Garg, and Fabio Ramos. 2021\.  DiSECt: A Differentiable Simulation Engine for Autonomous Robotic Cutting. In _Proceedings of Robotics: Science and Systems_. Virtual.  https://doi.org/10.15607/RSS.2021.XVII.067
  * Holl et al. (2020) Philipp Holl, Vladlen Koltun, Kiwon Um, and Nils Thuerey. 2020\.  phiflow: A differentiable pde solving framework for deep learning via physical simulations. In _NeurIPS Workshop_.
  * Hu et al. (2020) Yuanming Hu, Luke Anderson, Tzu-Mao Li, Qi Sun, Nathan Carr, Jonathan Ragan-Kelley, and Frédo Durand. 2020.  DiffTaichi: Differentiable Programming for Physical Simulation. In _ICLR_.
  * Ishiwaka et al. (2021) Yuko Ishiwaka, Xiao S Zeng, Michael Lee Eastman, Sho Kakazu, Sarah Gross, Ryosuke Mizutani, and Masaki Nakada. 2021.  Foids: bio-inspired fish simulation for generating synthetic datasets.  _ACM Transactions on Graphics (TOG)_ 40, 6 (2021), 1–15.
  * Jiang et al. (2001) Rui Jiang, Qingsong Wu, and Zuojin Zhu. 2001.  Full velocity difference model for a car-following theory.  _Physical Review E_ 64, 1 (2001), 017101.
  * Kesting et al. (2007) Arne Kesting, Martin Treiber, and Dirk Helbing. 2007\.  General lane-changing model MOBIL for car-following models.  _Transportation Research Record_ 1999, 1 (2007), 86–94.
  * Kraft et al. (1988) Dieter Kraft et al. 1988\.  A software package for sequential quadratic programming.  (1988).
  * LeVeque et al. (2002) Randall J LeVeque et al. 2002\.  _Finite volume methods for hyperbolic problems_. Vol. 31.  Cambridge university press.
  * Li et al. (2018) Tzu-Mao Li, Miika Aittala, Frédo Durand, and Jaakko Lehtinen. 2018.  Differentiable Monte Carlo ray tracing through edge sampling.  _ACM Trans. Graph._ 37, 6 (2018).
  * Li et al. (2022) Yifei Li, Tao Du, Kui Wu, Jie Xu, and Wojciech Matusik. 2022.  DiffCloth: Differentiable Cloth Simulation with Dry Frictional Contact.  _ACM Trans. Graph._ (2022).
  * Lieberman and Rathi (1997) Edward Lieberman and Ajay K Rathi. 1997.  Traffic simulation.  _Traffic flow theory_ (1997).
  * Lighthill and Whitham (1955) Michael James Lighthill and Gerald Beresford Whitham. 1955.  On kinematic waves II. A theory of traffic flow on long crowded roads.  _Proceedings of the Royal Society of London. Series A. Mathematical and Physical Sciences_ 229, 1178 (1955), 317–345.
  * Lillicrap et al. (2015) Timothy P Lillicrap, Jonathan J Hunt, Alexander Pritzel, Nicolas Heess, Tom Erez, Yuval Tassa, David Silver, and Daan Wierstra. 2015.  Continuous control with deep reinforcement learning.  _arXiv preprint arXiv:1509.02971_ (2015).
  * Lopez et al. (2018) Pablo Alvarez Lopez, Michael Behrisch, Laura Bieker-Walz, Jakob Erdmann, Yun-Pang Flötteröd, Robert Hilbrich, Leonhard Lücken, Johannes Rummel, Peter Wagner, and Evamarie Wießner. 2018.  Microscopic traffic simulation using sumo. In _2018 21st International Conference on Intelligent Transportation Systems (ITSC)_. IEEE, 2575–2582.
  * Ma et al. (2021) Pingchuan Ma, Tao Du, John Z Zhang, Kui Wu, Andrew Spielberg, Robert K Katzschmann, and Wojciech Matusik. 2021.  Diffaqua: A differentiable computational design pipeline for soft underwater swimmers with shape interpolation.  _ACM Transactions on Graphics (TOG)_ 40, 4 (2021), 1–14.
  * Magne et al. (2000) Laurent Magne, Sylvestre Rabut, and Jean-François Gabard. 2000\.  Towards an hybrid macro-micro traffic flow simulation model. In _INFORMS spring 2000 meeting_.
  * Mammar et al. (2006) Salim Mammar, Saïd Mammar, and Jean-Patrick Lebacque. 2006\.  Highway traffic hybrid macro-micro simulation model.  _IFAC Proceedings Volumes_ 39, 12 (2006), 627–632.
  * Mohamed and Mohamad (2010) Khaled M Mohamed and AA Mohamad. 2010.  A review of the development of hybrid atomistic–continuum methods for dense fluids.  _Microfluidics and Nanofluidics_ 8, 3 (2010), 283–302.
  * Mora et al. (2021) Miguel Angel Zamora Mora, Momchil P Peychev, Sehoon Ha, Martin Vechev, and Stelian Coros. 2021\.  PODS: Policy Optimization via Differentiable Simulation. In _International Conference on Machine Learning_. PMLR, 7805–7817.
  * Narain et al. (2009) Rahul Narain, Abhinav Golas, Sean Curtis, and Ming C Lin. 2009\.  Aggregate dynamics for dense crowd simulation.  In _ACM SIGGRAPH Asia 2009 papers_. 1–8.
  * Narain et al. (2010) Rahul Narain, Abhinav Golas, and Ming C Lin. 2010.  Free-flowing granular materials with two-way solid coupling.  In _ACM SIGGRAPH Asia 2010 papers_. 1–10.
  * Nelder and Mead (1965) John A Nelder and Roger Mead. 1965.  A simplex method for function minimization.  _The computer journal_ 7, 4 (1965), 308–313.
  * Newell (1961) Gordon Frank Newell. 1961\.  Nonlinear effects in the dynamics of car following.  _Operations research_ 9, 2 (1961), 209–229.
  * Nimier-David et al. (2019) Merlin Nimier-David, Delio Vicini, Tizian Zeltner, and Wenzel Jakob. 2019.  Mitsuba 2: A retargetable forward and inverse renderer.  _ACM Transactions on Graphics (TOG)_ 38, 6 (2019).
  * Paszke et al. (2019) Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, et al. 2019\.  Pytorch: An imperative style, high-performance deep learning library.  _Advances in neural information processing systems_ 32 (2019), 8026–8037.
  * Payne (1971) Harold J Payne. 1971\.  Model of freeway traffic and control.  _Mathematical Model of Public System_ (1971), 51–61.
  * Qiao et al. (2020) Yi-Ling Qiao, Junbang Liang, Vladlen Koltun, and Ming C. Lin. 2020.  Scalable Differentiable Physics for Learning and Control. In _ICML_.
  * Qiao et al. (2021a) Yi-Ling Qiao, Junbang Liang, Vladlen Koltun, and Ming C. Lin. 2021a.  Differentiable Simulation of Soft Multi-body Systems. In _Conference on Neural Information Processing Systems (NeurIPS)_.
  * Qiao et al. (2021b) Yi-Ling Qiao, Junbang Liang, Vladlen Koltun, and Ming C. Lin. 2021b.  Efficient Differentiable Simulation of Articulated Bodies. In _ICML_.
  * Richards (1956) Paul I Richards. 1956\.  Shock waves on the highway.  _Operations research_ 4, 1 (1956), 42–51.
  * Schulman et al. (2017) John Schulman, Filip Wolski, Prafulla Dhariwal, Alec Radford, and Oleg Klimov. 2017.  Proximal policy optimization algorithms.  _arXiv preprint arXiv:1707.06347_ (2017).
  * Sewall et al. (2011) Jason Sewall, David Wilkie, and Ming C Lin. 2011\.  Interactive hybrid simulation of large-scale traffic. In _Proceedings of the 2011 SIGGRAPH Asia Conference_. 1–12.
  * Sewall et al. (2010) Jason Sewall, David Wilkie, Paul Merrell, and Ming C Lin. 2010\.  Continuum traffic simulation. In _Computer Graphics Forum_ , Vol. 29. Wiley Online Library, 439–448.
  * Sewall (2011) Jason Douglas Sewall. 2011\.  _Efficient, scalable traffic and compressible fluid simulations using hyperbolic models_.  Ph.D. Dissertation. The University of North Carolina at Chapel Hill.
  * Shen et al. (2021) Siyuan Shen, Yang Yin, Tianjia Shao, He Wang, Chenfanfu Jiang, Lei Lan, and Kun Zhou. 2021.  High-order differentiable autoencoder for nonlinear model reduction.  _arXiv preprint arXiv:2102.11026_ (2021).
  * Takahashi et al. (2021) Tetsuya Takahashi, Junbang Liang, Yi-Ling Qiao, and Ming C Lin. 2021\.  Differentiable Fluids with Solid Coupling for Learning and Control. In _Proceedings of the AAAI Conference on Artificial Intelligence_ , Vol. 35. 6138–6146.
  * Treiber et al. (2000) Martin Treiber, Ansgar Hennecke, and Dirk Helbing. 2000\.  Congested traffic states in empirical observations and microscopic simulations.  _Physical review E_ 62, 2 (2000), 1805.
  * Treuille et al. (2006) Adrien Treuille, Seth Cooper, and Zoran Popović. 2006\.  Continuum crowds.  _ACM Transactions on Graphics (TOG)_ 25, 3 (2006), 1160–1168.
  * Wei et al. (2018) Hua Wei, Guanjie Zheng, Huaxiu Yao, and Zhenhui Li. 2018\.  Intellilight: A reinforcement learning approach for intelligent traffic light control. In _Proceedings of the 24th ACM SIGKDD International Conference on Knowledge Discovery & Data Mining_. 2496–2505.
  * Whitham (2011) Gerald Beresford Whitham. 2011\.  _Linear and nonlinear waves_. Vol. 42.  John Wiley & Sons.
  * Yu and Krstic (2019) Huan Yu and Miroslav Krstic. 2019.  Traffic congestion control for Aw–Rascle–Zhang model.  _Automatica_ 100 (2019), 38–51.
  * Zhang (2002) H Michael Zhang. 2002\.  A non-equilibrium traffic model devoid of gas-like behavior.  _Transportation Research Part B: Methodological_ 36, 3 (2002), 275–290.
  * Zheng et al. (2017) Yuanshi Zheng, Jingying Ma, and Long Wang. 2017.  Consensus of hybrid multi-agent systems.  _IEEE transactions on neural networks and learning systems_ 29, 4 (2017), 1359–1365.

##  Appendix A Solution to the Riemann problem in ARZ model

###  A.1. Formulation

For the given left state $\mathbf{q_{l}}$ and right state $\mathbf{q_{r}}$, we
can determine the intermediate state $\mathbf{q_{0}}$ as follows.

Case 0 : $u_{l}=u_{r}$

(14) $q_{0}=q_{l}$

Case 1 : $u_{l}>u_{r}$

(15) $q_{0}=\begin{cases}q_{l}&\mbox{if }\lambda_{s}\geq 0\\\ q_{m}&\mbox{if
}\lambda_{s}<0\\\
\end{cases}(\lambda_{s}=\frac{\rho_{m}u_{m}-\rho_{l}u_{l}}{\rho_{m}-\rho_{l}})$

Case 2 : $u_{r}-u_{max}\rho_{l}^{\gamma}<u_{l}<u_{r}$

(16) $q_{0}=\begin{cases}q_{l}&\mbox{if }\lambda_{0l}\geq 0\\\ q_{m}&\mbox{if
}\lambda_{0m}\leq 0\\\ \tilde{q}(0)&\mbox{if
}\lambda_{0l}<0,\lambda_{0m}>0\end{cases}$

where

$\displaystyle\lambda_{0l}$
$\displaystyle=u_{l}-u_{max}\gamma\rho_{l}^{\gamma}$
$\displaystyle\lambda_{0m}$
$\displaystyle=u_{m}-u_{max}\gamma\rho_{m}^{\gamma}$
$\displaystyle=u_{r}-u_{max}\gamma\rho_{l}^{\gamma}+\gamma(u_{r}-u_{l}).$

Case 3 : $u_{l}\leq u_{r}-u_{max}\rho_{l}^{\gamma}$

(17) $q_{0}=\begin{cases}q_{l}&\mbox{if }\lambda_{0l}\geq 0\\\
\tilde{q}(0)&\mbox{if }\lambda_{0l}<0\end{cases}$

Case 4 : $\rho_{l}=0,\rho_{r}>0$

(18) $q_{0}=q_{l}=\begin{bmatrix}0\\\ 0\end{bmatrix}$

Case 5 : $\rho_{l}>0,\rho_{r}=0$

Same as Case 3.

From above, $q_{m}$ is given as

$\displaystyle\rho_{m}$
$\displaystyle={(\rho_{l}^{\gamma}+\frac{u_{l}-u_{r}}{u_{max}})}^{\frac{1}{\gamma}}$
$\displaystyle u_{m}$ $\displaystyle=u_{r}$ $\displaystyle y_{m}$
$\displaystyle=\rho_{m}(u_{m}-u_{eq}(\rho_{m})),$

and $\tilde{q}(0)$ is given as

$\displaystyle\tilde{\rho}(0)$
$\displaystyle=(\frac{u_{l}+u_{max}\rho_{l}^{\gamma}}{(\gamma+1)u_{max}})^{\frac{1}{\gamma}}$
$\displaystyle\tilde{u}(0)$
$\displaystyle=\frac{\rho_{l}}{\rho_{l}+1}(u_{l}+u_{max}\rho_{l}^{\gamma})$
$\displaystyle\tilde{y}(0)$
$\displaystyle=\tilde{\rho}(0)(\tilde{u}(0)-u_{eq}(\tilde{\rho}(0))).$

###  A.2. Differentiation

We can see from the above formulation that we have to compute the partial
derivatives of $q_{l},q_{m},$ and $\tilde{q}(0)$ with respect to $q_{l}$ and
$q_{r}$, as those values comprise the solutions of the Riemann problem. We can
compute the partial derivatives as shown below.

The Jacobian of $\mathbf{f}$ is:

(19) $\mathbf{f}^{\prime}(\mathbf{q})=\begin{bmatrix}u_{eq}+\rho
u_{eq}^{\prime}&1\\\
yu_{eq}^{\prime}-\frac{y^{2}}{\rho^{2}}&\frac{2y}{\rho}+u_{eq}\end{bmatrix}$

Partial derivatives of $\mathbf{q_{l}}$

$\displaystyle\frac{\partial q_{l}}{\partial q_{l}}$ $\displaystyle=I_{2\times
2}$ $\displaystyle\frac{\partial q_{l}}{\partial q_{r}}$
$\displaystyle=\begin{bmatrix}0&0\\\ 0&0\end{bmatrix}$

Partial derivatives of $\mathbf{q_{m}}$

$\displaystyle\frac{\partial q_{m}}{\partial q_{l}}$
$\displaystyle=\begin{bmatrix}\frac{\partial\rho_{m}}{\partial\rho_{l}}&\frac{\partial\rho_{m}}{\partial
y_{l}}\\\ \frac{\partial y_{m}}{\partial\rho_{l}}&\frac{\partial
y_{m}}{\partial y_{l}}\\\ \end{bmatrix}$

where

$\displaystyle\frac{\partial\rho_{m}}{\partial\rho_{l}}$
$\displaystyle=\frac{1}{\gamma}(\rho_{l}^{\gamma}+\frac{u_{l}-u_{r}}{u_{max}})^{\frac{1-\gamma}{\gamma}}(\gamma\rho_{l}^{\gamma-1}+\frac{1}{u_{max}}\frac{\partial
u_{l}}{\partial\rho_{l}}),$ $\displaystyle\frac{\partial\rho_{m}}{\partial
y_{l}}$
$\displaystyle=\frac{1}{\gamma}(\rho_{l}^{\gamma}+\frac{u_{l}-u_{r}}{u_{max}})^{\frac{1-\gamma}{\gamma}}(\frac{1}{u_{max}}\frac{\partial
u_{l}}{\partial y_{l}}),$ $\displaystyle\frac{\partial
y_{m}}{\partial\rho_{l}}$
$\displaystyle=\frac{\partial\rho_{m}}{\rho_{l}}(u_{m}-u_{eq}(\rho_{m}))+\rho_{m}(-u^{\prime}_{eq}(\rho_{m})\frac{\partial\rho_{m}}{\partial\rho_{l}}),$
$\displaystyle\frac{\partial y_{m}}{\partial y_{l}}$
$\displaystyle=\frac{\partial\rho_{m}}{y_{l}}(u_{m}-u_{eq}(\rho_{m}))+\rho_{m}(-u^{\prime}_{eq}(\rho_{m})\frac{\partial\rho_{m}}{\partial
y_{l}}).$

$\displaystyle\frac{\partial q_{m}}{\partial q_{r}}$
$\displaystyle=\begin{bmatrix}\frac{\partial\rho_{m}}{\partial\rho_{r}}&\frac{\partial\rho_{m}}{\partial
y_{r}}\\\ \frac{\partial y_{m}}{\partial\rho_{r}}&\frac{\partial
y_{m}}{\partial y_{r}}\\\ \end{bmatrix}$

where

$\displaystyle\frac{\partial\rho_{m}}{\partial\rho_{r}}$
$\displaystyle=\frac{1}{\gamma}(\rho_{l}^{\gamma}+\frac{u_{l}-u_{r}}{u_{max}})^{\frac{1-\gamma}{\gamma}}(-\frac{1}{u_{max}}\frac{\partial
u_{r}}{\partial\rho_{r}}),$ $\displaystyle\frac{\partial\rho_{m}}{\partial
y_{r}}$
$\displaystyle=\frac{1}{\gamma}(\rho_{l}^{\gamma}+\frac{u_{l}-u_{r}}{u_{max}})^{\frac{1-\gamma}{\gamma}}(-\frac{1}{u_{max}}\frac{\partial
u_{r}}{\partial y_{r}}),$ $\displaystyle\frac{\partial
y_{m}}{\partial\rho_{r}}$
$\displaystyle=\frac{\partial\rho_{m}}{\rho_{r}}(u_{m}-u_{eq}(\rho_{m}))+\rho_{m}(\frac{\partial
u_{r}}{\partial\rho_{r}}-u^{\prime}_{eq}(\rho_{m})\frac{\partial\rho_{m}}{\partial\rho_{r}}),$
$\displaystyle\frac{\partial y_{m}}{\partial y_{r}}$
$\displaystyle=\frac{\partial\rho_{m}}{y_{r}}(u_{m}-u_{eq}(\rho_{m}))+\rho_{m}(\frac{\partial
u_{r}}{\partial
y_{r}}-u^{\prime}_{eq}(\rho_{m})\frac{\partial\rho_{m}}{\partial y_{r}}).$

Partial derivatives of $\mathbf{\tilde{q}(0)}$

$\displaystyle\frac{\partial\tilde{q}(0)}{\partial q_{l}}$
$\displaystyle=\begin{bmatrix}\frac{\partial\tilde{\rho}(0)}{\partial\rho_{l}}&\frac{\partial\tilde{\rho}(0)}{\partial
y_{l}}\\\
\frac{\partial\tilde{y}(0)}{\partial\rho_{l}}&\frac{\partial\tilde{y}(0)}{\partial
y_{l}}\\\ \end{bmatrix}$

where

$\displaystyle\frac{\partial\tilde{\rho}(0)}{\partial\rho_{l}}$
$\displaystyle=\frac{1}{\gamma}(\frac{u_{l}+u_{max}\rho_{l}^{\gamma}}{(\gamma+1)u_{max}})^{\frac{1-\gamma}{\gamma}}(\frac{1}{(\gamma+1)u_{max}})(\frac{\partial
u_{l}}{\partial\rho_{l}}+u_{max}\gamma\rho_{l}^{\gamma-1}),$
$\displaystyle\frac{\partial\tilde{\rho}(0)}{\partial y_{l}}$
$\displaystyle=\frac{1}{\gamma}(\frac{u_{l}+u_{max}\rho_{l}^{\gamma}}{(\gamma+1)u_{max}})^{\frac{1-\gamma}{\gamma}}(\frac{1}{(\gamma+1)u_{max}})(\frac{\partial
u_{l}}{\partial y_{l}}),$
$\displaystyle\frac{\partial\tilde{y}(0)}{\partial\rho_{l}}$
$\displaystyle=\frac{\partial\tilde{\rho}(0)}{\partial\rho_{l}}(\tilde{u}(0)-u_{eq}(\tilde{\rho}(0)))+\tilde{\rho}(0)(\frac{\partial\tilde{u}(0)}{\partial\rho_{l}}-u^{\prime}_{eq}(\tilde{\rho}(0))\frac{\partial\tilde{\rho}(0)}{\partial\rho_{l}}),$
$\displaystyle\frac{\partial\tilde{y}(0)}{\partial y_{l}}$
$\displaystyle=\frac{\partial\tilde{\rho}(0)}{\partial
y_{l}}(\tilde{u}(0)-u_{eq}(\tilde{\rho}(0)))+\tilde{\rho}(0)(\frac{\partial\tilde{u}(0)}{\partial
y_{l}}-u^{\prime}_{eq}(\tilde{\rho}(0))\frac{\partial\tilde{\rho}(0)}{\partial
y_{l}}).$ $\displaystyle\frac{\partial\tilde{q}(0)}{\partial q_{r}}$
$\displaystyle=\begin{bmatrix}0&0\\\ 0&0\end{bmatrix}$

##  Appendix B Differentiating IDM

###  B.1. Hyperparameters

In the IDM, following hyperparameters are used to describe different
behavioral traits of discrete vehicles (Treiber et al., 2000).

  * •

$s_{min}$: Minimum desired distance to the leading vehicle.

  * •

$T_{pref}$: Desired time to move forward with current speed.

  * •

$a_{max}$: Upper bound of the computed acceleration.

  * •

$a_{pref}$: Comfortable braking deceleration, positive value.

  * •

$v_{targ}$: Target velocity it wants to maintain.

  * •

$length_{i}$: Length of the $i$-th vehicle.

###  B.2. Differentiation

First we compute the gradients of $\mathbf{a_{i}}(t)$ with respect to
$\mathbf{p_{i}}(t)$, $\mathbf{v_{i}}(t)$, $\mathbf{p_{h(i)}}(t)$, and
$\mathbf{v_{h(i)}}(t)$, where ${h(i)}$-th vehicle is the leading vehicle of
$i$-th vehicle.

$\displaystyle\frac{\partial\mathbf{a_{i}}(t)}{\partial\mathbf{p_{i}}(t)}$
$\displaystyle=-2a_{max}\frac{{s_{opt}}^{2}}{(p_{h(i)}(t)-p_{i}(t)-length_{h(i)})^{3}},$
$\displaystyle\frac{\partial\mathbf{a_{i}}(t)}{\partial\mathbf{p_{h(i)}}(t)}$
$\displaystyle=2a_{max}\frac{{s_{opt}}^{2}}{(p_{h(i)}(t)-p_{i}(t)-length_{h(i)})^{3}},$
$\displaystyle\frac{\partial\mathbf{a_{i}}(t)}{\partial\mathbf{v_{i}}(t)}$
$\displaystyle=a_{max}(-\delta\frac{v_{i}(t)^{\delta-1}}{{v_{targ}}^{\delta}}-2\frac{s_{opt}}{(p_{h(i)}(t)-p_{i}(t)-length_{h(i)})^{2}}\frac{\partial
s_{opt}}{\partial v_{i}(t)}),$
$\displaystyle\frac{\partial\mathbf{a_{i}}(t)}{\partial\mathbf{v_{h(i)}}(t)}$
$\displaystyle=a_{max}(-2\frac{s_{opt}}{(p_{h(i)}(t)-p_{i}(t)-length_{h(i)})^{2}}\frac{\partial
s_{opt}}{\partial v_{h(i)}(t)}),$

where

$\displaystyle\frac{\partial s_{opt}}{\partial v_{i}(t)}$
$\displaystyle=T_{pref}+\frac{2v_{i}(t)-v_{h(i)}(t)}{2\sqrt{a_{max}a_{pref}}},$
$\displaystyle\frac{\partial s_{opt}}{\partial v_{h(i)}(t)}$
$\displaystyle=\frac{-v_{i}(t)}{2\sqrt{a_{max}a_{pref}}}.$

Now for arbitrary $i$ and $j$, we can compute the gradient of the state of the
$i$-th vehicle with respect to that of the $j$-th vehicle by plugging in the
above values. Note that the gradient is non-zero only when $j$ equals to $i$
or $h(i)$.

(20)
$\frac{\partial\mathbf{p_{i}}(t+1)}{\partial\mathbf{p_{j}}(t)}=\begin{cases}1&\mbox{if
}j=i\\\ 0&\mbox{if }j=h(i)\\\ 0&\mbox{if }else\\\ \end{cases}$

(21)
$\frac{\partial\mathbf{p_{i}}(t+1)}{\partial\mathbf{v_{j}}(t)}=\begin{cases}\Delta{t}&\mbox{if
}j=i\\\ 0&\mbox{if }j=h(i)\\\ 0&\mbox{if }else\\\ \end{cases}$

(22)
$\frac{\partial\mathbf{v_{i}}(t+1)}{\partial\mathbf{p_{j}}(t)}=\begin{cases}\Delta{t}(-2a_{max}\frac{{s_{opt}}^{2}}{(p_{h(i)}(t)-p_{i}(t)-length_{h(i)})^{3}})&\mbox{if
}j=i\\\
\Delta{t}(2a_{max}\frac{{s_{opt}}^{2}}{(p_{h(i)}(t)-p_{i}(t)-length_{h(i)})^{3}})&\mbox{if
}j=h(i)\\\ 0&\mbox{if }else\\\ \end{cases}$

(23)
$\frac{\partial\mathbf{v_{i}}(t+1)}{\partial\mathbf{v_{j}}(t)}=\begin{cases}1+\Delta{t}a_{max}(-\delta\frac{{v_{i}(t)}^{\delta-1}}{{v_{targ}}^{\delta}}-2\frac{s_{opt}}{(p_{h(i)}(t)-p_{i}(t)-length_{h(i)})^{2}}(T_{pref}+\frac{2v_{i}(t)-v_{j}(t)}{2\sqrt{a_{max}a_{pref}}}))&\mbox{if
}j=i\\\
\Delta{t}a_{max}(-2\frac{s_{opt}}{(p_{h(i)}(t)-p_{i}(t)-length_{h(i)})^{2}}\frac{-v_{i}(t)}{2\sqrt{a_{max}a_{pref}}})&\mbox{if
}j=h(i)\\\ 0&\mbox{if }else\\\ \end{cases}$