11institutetext: German Research Center for Artificial Intelligence (DFKI),
Germany 22institutetext: RPTU Kaiserslautern-Landau, Germany
22email<EMAIL_ADDRESS>
# DWA: Differential Wavelet Amplifier for Image Super-Resolution

Brian B. Moser 11 2 2 Stanislav Frolov 1122 Federico Raue 11 Sebastian Palacio
11 Andreas Dengel 11 2 2

###### Abstract

This work introduces Differential Wavelet Amplifier (DWA), a drop-in module
for wavelet-based image Super-Resolution (SR). DWA invigorates an approach
recently receiving less attention, namely Discrete Wavelet Transformation
(DWT). DWT enables an efficient image representation for SR and reduces the
spatial area of its input by a factor of 4, the overall model size, and
computation cost, framing it as an attractive approach for sustainable ML. Our
proposed DWA model improves wavelet-based SR models by leveraging the
difference between two convolutional filters to refine relevant feature
extraction in the wavelet domain, emphasizing local contrasts and suppressing
common noise in the input signals. We show its effectiveness by integrating it
into existing SR models, e.g., DWSR and MWCNN, and demonstrate a clear
improvement in classical SR tasks. Moreover, DWA enables a direct application
of DWSR and MWCNN to input image space, reducing the DWT representation
channel-wise since it omits traditional DWT.

###### Keywords:

Differential Wavelet Amplifier Image Super-Resolution.

##  1 Introduction

Image Super-Resolution (SR) has an impressive legacy in Computer Vision (CV)
yet still presents an exhilarating challenge [31, 21]. SR is a task of
enhancing Low-Resolution (LR) images to High Resolution (HR). It is
challenging because many High Resolution (HR) images can correspond to a given
Low-Resolution (LR) image, rendering the task mathematically ill-posed.

In recent years, deep learning has fueled rapid development in SR, leading to
tremendous progress [6, 7]. While many techniques have improved the overall
quality of image reconstructions, there remains a pressing need for methods
capable of producing high-frequency details, particularly when dealing with
high magnification ratios [24]. Addressing this issue is crucial for the
continued advancement of SR. Influenced by achievements on other CV tasks,
recent research focused on trending approaches like Transformer-based networks
[28, 16, 27], Denoising Diffusion Probabilistic Models [15, 24, 23] or
Generative Adversarial Networks [30, 29]. Despite astonishing reconstruction
capabilities, they often lack an explicit focus on generating high-frequency
details, i.e., local variations.

This work aims to advance the field of SR by exploring wavelet-based networks.
Unfortunately, this technique has received less attention despite its
significant potential [21]. We seek to provide a fresh perspective and revive
research by re-evaluating these approaches. Discrete Wavelet Transformation
(DWT) enables an efficient image representation without losing information
compared to its naïve spatial representation, i.e., traditional RGB format. It
does so by separating high-frequency details in distinct channels and reducing
the spatial area of input image representation by a factor of 4. Therefore, a
smaller receptive field is required to capture the input during feature
extraction. Using DWT like in DWSR [9] and MWCNN [17] reduces the overall
model size and computational costs while performing similarly to state-of-the-
art image SR architectures.

This work introduces a new Differential Wavelet Amplifier (DWA) module
inspired by differential amplifiers from electrical engineering [2].
Differential amplifiers increase the difference between two input signals and
suppress the common voltage shared by the two inputs, called Common Mode
Rejection (CMR) [11]. In other words, it mitigates the impact of noise (e.g.,
electromagnetic interference, vibrations, or thermal noise) affecting both
source inputs while retaining valuable information and improving the integrity
of the measured input signal. Our proposed DWA layer adapts this idea to deep
learning and can be used as a drop-in module to existing SR models. This work
shows its effectiveness as exemplary for wavelet-based SR approaches. DWA
leverages the difference between two convolutional filters with a stride
difference to enhance relevant feature extraction in the wavelet domain,
emphasizing local contrasts and suppressing common noise in the input signals.
We demonstrate the effectiveness of DWA through extensive experiments and
evaluations, showing improved performance compared to existing wavelet-based
SR models without DWA: DWSR with DWA shows overall better performance w.r.t.
PSNR and SSIM, and MWCNN with DWA achieves better SSIM scores with comparable
PSNR values on the testing datasets Set5 [4], Set14 [32], and BSDS100 [20].

Taken together, our work makes the following key contributions:

  * •

Introduction of Differential Wavelet Amplifier (DWA): a novel module that
leverages the difference between two convolutional filters horizontally and
vertically in a wavelet-based image representation, which is applicable as
drop-in addition in existing network architectures.

  * •

Comprehensive evaluation demonstrating the improved performance by using DWA
on popular SR datasets such as Set5 [4], Set14 [32], and BSDS100 [20] by
adding DWA to existing wavelet-based SR models, namely, DWSR [9] and MWCNN
[17].

  * •

Experimental analysis showing that DWA enables a direct application of DWSR
and MWCNN to the input space by avoiding the DWT on the input image. This
application reduces the input channel-wise to 3 instead of 12 channels for RGB
images while keeping the spatial reduction benefit of DWT.

  * •

Visual examination of reconstructions showcasing that the DWSR with the DWA
module captures better distinct edges and finer details, which are also closer
to the ground truth residuals.

##  2 Background

This chapter provides comprehensive background information on 2D Discrete
Wavelet Transform (2D-DWT), how SR models (DWSR [9] and MWCNN [17]) use it,
and related work to Differential Wavelet Amplifiers (DWA). Additionally, we
introduce differential amplifiers from electrical engineering, which inspired
our proposed method DWA.

###  2.1 Discrete Wavelet Transform in SR

The 2D Discrete Wavelet Transform (2D-DWT) decomposes an image into four
unique sub-bands with distinct frequency components: a low-frequency
approximation sub-band and three high-frequency detail sub-bands representing
horizontal, vertical, and diagonal details. Let
$x\left[n\right]\in\mathbb{R}^{N}$ be a signal. The 1D Discrete Wavelet
Transformation (1D-DWT) with Haar wavelet passes the input signal first
through a half-band high-filter $h\left[n\right]$ and a low-pass filter
$l\left[n\right]$. Next, half of the sample is eliminated according to the
Nyquist rule [9]. The wavelet coefficients are calculated by repeating the
decomposition to each output coefficient iteratively [26]. In the case of
images, it applies $h\left[n\right]$ and $l\left[n\right]$ in different
combinations, resulting in four function applications.

The DWSR [9] SR model exploits the wavelet domain and gets the DWT
representation of the interpolated LR image as input. DWSR is composed of 10
convolution layers that are applied sequentially. It adds the interpolated LR
input as residual for the final reconstruction step, which results in learning
only the sparse residual information between the LR and HR domains.

MWCNN [17] exploits multi-level DWT (multiple applications of DWT) and
utilizes a U-Net architecture [22]. DWT replaces all downsizing steps, and the
inverse operation of DWT replaces all upsampling steps. Ultimately, it uses
the interpolated LR image as a residual connection for the final prediction.
The standard MWCNN setup consists of 24 convolution layers.

One caveat of DWSR and MWCNN in learning the residual is that they must
translate its rich information input to sparse representation, e.g., the
average band. To ease the burden, we present a Differential Wavelet Amplifier,
which directly transforms the input into sparse representations inspired by
differential amplifiers introduced next.

###  2.2 Differential Amplifier

An electronic amplifier is a standard electrical engineering device to
increase a signal’s power [2]. One type of electronic amplifier is the
differential amplifier that increases the difference between two input signals
and suppresses the common voltage shared by the two inputs [14]. Given two
inputs $V^{-}_{in},V^{+}_{in}\in\mathbb{R}^{N}$ and the differential gain of
the amplifier $A_{d}\in\mathbb{R}$, the output $V_{out}$ is calculated as

$V_{out}=A_{d}\left(V^{+}_{in}-V^{-}_{in}\right)$ (1)

The purpose of differential amplifiers is to suppress common signals or noise
sources that are present in multiple input channels while retaining valuable
information. In the literature, this is called Common Mode Rejection (CMR) and
is a critical property in many electrical engineering applications,
particularly in systems that measure small signals in the presence of noise or
interference, e.g., electromagnetic interference or thermal noise [11]. Hence,
using CMR improves the signal-to-noise ratio, overall system performance, and
signal integrity since the system can focus on the relevant differential
signals.

###  2.3 Differential Convolutions

Closest to our work is Sarıgül et al. [25], which applies differential
convolutions, i.e., the difference of two convolution layers, to emphasize
contrasts for image classification, which is inherently different to image
generation tasks such as image SR. Despite this, they do not consider a stride
difference vital for capturing variations. Knutsson et al. [13] theoretically
examine a normalized version of differential convolutions also with no stride
difference. Due to the time of publication, they did not try it in the case of
deep learning-based image SR. Newer applications like Canh et al. [5] consider
learnable parameters to turn the Difference of Gaussians (DoG) [18] into a
learnable framework, but has the same caveat: As Knutsson concluded, their
approaches can be interpreted as a standard convolution weighted with the
local energy minus the “mean” operator acting on the “mean” data, i.e., a more
elaborate convolution operation.

A similarity could also be seen in the approach of residual connections of
ResNets [10] when the kernel parameters have a negative sign. However,
residual connections are different since they force a convolution layer to
learn to extract the sparse details that are not apparent in the input. In
contrast, our proposed method with Differential Wavelet Amplifier (DWA)
explicitly produces sparse details by design due to the subtraction operator.
Therefore, DWA does not have to learn what input information should be removed
for the residual information. It can focus on relevant features that persist
when the stride convolution does not detect the same feature, thereby
emphasizing local contrast.

Figure 1:  Visualization of DWA. It takes the difference of two convolutional
filters with a stride difference of at least 1, vertically and horizontally.
Next, it concatenates the input with horizontal and vertical feature maps. In
the end, it applies a final convolution.

##  3 Differential Wavelet Amplifier (DWA)

This section presents our proposed Differential Wavelet Amplifier (DWA)
module. Inspired by differential amplifiers in electrical engineering, DWA is
designed to operate in the wavelet domain and exploits the difference between
two input signals to improve the performance of image SR methods based on
wavelet predictions.DWA is applied separately in the horizontal and vertical
axis of the input image. In each direction, we perform two convolutions with a
stride distance in one direction for both axis (from left to right, from top
to bottom, as in MDLSTMs [8]), allowing a fine-grained feature extraction and
emphasizing local contrasts while suppressing the common mode in the input,
similar to CMR in electrical engineering. Figure 1 visualizes all processes
involved in DWA.

Let $\mathbf{x}\in\mathbb{R}^{w\times h\times c_{in}}$ be an input image or
feature map with $c_{in}$ channels. We define
$\psi\left(\mathbf{x},(i,j)\right):\mathbb{R}^{w\times h\times
c_{in}}\times\mathbb{N}^{2}\rightarrow\mathbb{R}^{k\cdot k\times c_{in}}$ as a
function that extracts $k\cdot k$ points around a spatial position $(i,j)$. We
can then express the resulting feature maps for the horizontal
$\mathbf{H}\left(\mathbf{x}\right)$ and vertical
$\mathbf{V}\left(\mathbf{x}\right)$ axis as

$\begin{split}\mathbf{H}\left(\mathbf{x}\right)_{i,j}&=f\left(\psi\left(\mathbf{x},(i,j)\right);\theta_{1}\right)-f\left(\psi\left(\mathbf{x},(i+s,j)\right);\theta_{2}\right),\\\
\mathbf{V}\left(\mathbf{x}\right)_{i,j}&=f\left(\psi\left(\mathbf{x},(i,j)\right);\theta_{3}\right)-f\left(\psi\left(\mathbf{x},(i,j+s)\right);\theta_{4}\right),\end{split}$
(2)

where $f:\mathbb{R}^{k\cdot k\times c_{in}}\rightarrow\mathbb{R}^{c_{f}}$ is a
convolution operation with parameters $\theta_{n}$ for $0<n<4$ , $k\times k$
the kernel size and $s\in\mathbb{N}$ a pre-defined stride difference. As a
result, the local variance is captured in one direction for both axes, similar
to MDLSTMs [8]: from left to right with parameters $\theta_{1}$ and
$\theta_{2}$ and from top to bottom with parameters $\theta_{3}$ and
$\theta_{4}$. We obtain two distinct feature maps that capture complementary
input image information and provide richer feature representations for the
wavelet-based SR task. The input is directly translated to sparse
representations, which reduces the distance to residual target objectives in
networks that use residual connections for final prediction.

We concatenate the resulting feature maps alongside the input to ensure no
information is lost during the DWA processing. This combination creates a
comprehensive set of feature maps that retains the original input information
while incorporating the directional features obtained from both axes. More
formally:

$g\left(\mathbf{x}\right)=\mathbf{x}\odot\sigma\left(H\left(\mathbf{x}\right)\odot
V\left(\mathbf{x}\right)\right),$ (3)

where $\odot$ is a channel-wise concatenation operator and $\sigma$ is a non-
linear function like sigmoid, tanh or ReLU [1]. The concatenated feature map
is fed into an additional convolution layer $f_{final}:\mathbb{R}^{k\cdot
k\times(c_{in}+2\cdot c_{f})}\rightarrow\mathbb{R}^{c_{final}}$ and parameters
$\theta_{final}$, which maps the channel size after concatenation to a desired
target channel size $c_{final}$ such that our module can easily be
incorporated into existing models:

$\text{DWA}\left(\mathbf{x}\right)_{i,j}=f_{final}\left(\psi\left(g\left(\mathbf{x}\right),(i,j)\right);\theta_{final}\right)$
(4)

A SR model utilizing this DWA module exploits the comprehensive feature map to
learn the complex relationships between LR and HR images, ultimately
reconstructing the HR image with reduced noise. By employing the DWA, we aim
to harness the benefits of wavelet domain processing and the difference
between two convolutional filters. We demonstrate the effectiveness of our
approach through extensive experiments and evaluations in the following
sections.

###  3.1 Direct Application of DWA (DWA Direct)

One way to circumvent additional computation steps is to apply DWA directly on
the image space, omitting DWT and learning the transition between image and
frequency space implicitly via DWA. Thus, the interpolation of the input,
which effectively adds no additional information since it generates only
approximated values, can be reduced by half for networks like DWSR or MWCNN.
Consequently, the network is better adapted to the given values of the LR
input. In the experiments, we evaluate this alternative approach called DWA
Direct and show that it further enhances the performances of DWSR and MWCNN.

##  4 Experiments

We evaluate our proposed DWA module by integrating it into the wavelet-based
SR models DWSR and MWCNN. We begin this section by describing the experiments.
Next, we discuss the results quantitatively and qualitatively. We show the
effectiveness of DWA and that a direct application of wavelet-based SR models
with DWA to image space is feasible without forfeiting reconstruction quality.

###  4.1 Experimental Setup

We applied widely-used SR datasets to evaluate our method. In addition, we
utilized standard augmentation techniques such as rotation, horizontal and
vertical flipping. For testing, we employed the datasets Set5 [4], Set14 [32],
BSDS100 [20]. For training, we used different settings for DWSR and MWCNN to
match the original works for a fair comparison, as dissected in the following.
In all experiments, we train using the Adam optimizer [12] with a learning
rate of $10^{-4}$ with $L2$ regularization of $10^{-8}$ on a single A100 GPU.
Moreover, we use a learning rate decay schedule, which reduces the learning
rate by $20\%$ every 20 epochs.

Ablation Study: We use DIV2K [3] and follow the standard procedure by
extracting sub-images of $192\times 192$ for training. We iterate for 40
epochs over the training dataset. Since we compare with DWSR, we use $L1$-loss
as the learning objective, as reported by the authors of DWSR.

DWSR-Scenario: We use DIV2K [3] like in the ablation study, but we train for
100 epochs as reported in DWSR.

MWCNN-Scenario: We collect 800 images from DIV2K [3], 200 images from BSD [20]
and 4,744 images from WED [19] and train for 100 epochs. Contrary to DWSR, we
adapt the $L2$-loss like the authors of MWCNN. For sub-image extraction, we
use a size of $240\times 240$ to match the training settings of MWCNN.

Table 1:  Comparison of different striding settings on BSDS100 (2x and 4x scaling). Scale | Metric | no stride | s=1 | s=2 | s=3
---|---|---|---|---|---
2x |  PSNR $\uparrow$ | 31.8314 | 31.8660 | 31.8598 | 31.8588
SSIM $\uparrow$ | 0.9058 | 0.9061 | 0.9060 | 0.9059
4x |  PSNR $\uparrow$ | 27.2870 | 27.3048 | 27.2927 | 27.2872
SSIM $\uparrow$ | 0.7457 | 0.7471 | 0.7464 | 0.7466

##  5 Results

This section presents the quantitative and qualitative analysis of this work.
It shows that incorporating the DWA module into DWSR improves the performance
in every dataset and for all scaling factors. Moreover, we consistently
improve the SSIM scores by implementing DWA into MWCNN and achieve similar
PSNR results. This section starts with an ablation study to investigate
different striding settings and the effect of combining DWA with DWSR for the
direct application and the regular DWT case (see Section 3.1). Next, we
examine the performance scores of our DWA module on classical SR datasets with
DWSR and MWCNN. Finally, we visually compare the quality of the
reconstructions.

####  5.0.1 Ablation Study

Table 1 shows the impact of different striding settings for DWSR with DWA with
2x and 4x scaling. We observe an improvement for striding settings greater
than 0, significantly for PSNR and slightly for SSIM. The differences between
striding settings greater than 0 are minimal, with a slight decrease for
larger striding sizes. Nonetheless, they outperform DWA with no stride
difference consistently. Thus, having a stride difference to capture local
variations more effectively benefits the overall performance of DWSR.

We further investigate the impact of various model configurations, DWSR with
or without the DWA module, in a direct application or without (see Section
3.1). Figure 2 presents the results, where two graphs display the PSNR and
SSIM values [21], respectively, for each method. We apply the ablation study
with different model depths, ranging from 6 to 18, instead of using a standard
depth of 10 for DWSR. As a result, DWSR with DWA or DWA Direct consistently
outperforms the DWSR baseline across all model depths. This demonstrates the
effectiveness of incorporating the DWA module as the first layer in the DWSR
framework. Moreover, DWA Direct outperforms DWA applied to the DWT on the
input with greater model depths. Furthermore, we observe a considerable
performance drop in DWSR Direct without using the DWA module compared to all
other evaluated methods. This indicates that the DWA module is crucial in
enabling the Direct approach, as its absence significantly degrades
performance.

Figure 2:  Results of ablation study on BSDS100 with scaling factor 2x. We
tested different configurations: Baseline, Direct (application on the image
space), DWA, and DWA Direct (application on the image space).

####  5.0.2 Performance

Table 2 summarizes PSNR and SSIM scores when applying the DWA module to DWSR
and MWCNN for classical SR datasets on different scaling factors for a longer
training span. We observe that incorporating the DWA module into DWSR improves
the performance in every dataset and for all scaling factors. For MWCNN with
DWA, a similar observation can be made, especially for the SSIM scores, which
show overall the best performances. However, it has slightly decreased PSNR
values for some cases, e.g., for scaling factor 3. Both applications, DWSR
with DWA and MWCNN with DWA, are applied directly on the input image space,
omitting a DWT of the input.

Table 2: Comparison of DWSR, MWCNN, and DWA Direct with DWSR and MWCNN architecture on Set5, Set14, and BSDS100. Note that PSNR [dB] is a logarithmic scale, and SSIM reflects correlations (with values ranging from -1 to 1) [21].  Dataset | Scale | DWSR | MWCNN | DWA Direct | DWA Direct
---|---|---|---|---|---
[9] | [17] | [DWSR] | [MWCNN]
PSNR/SSIM | PSNR/SSIM | PSNR/SSIM | PSNR/SSIM
Set5 [4] | 2x | 37.43 / 0.9568 | 37.91 / 0.9600 | 37.79 / 0.9645 |  37.99 / 0.9652
3x | 33.82 / 0.9215 |  34.18 / 0.9272 | 33.85 / 0.9310 | 34.09 / 0.9329
4x | 31.39 / 0.8833 | 32.12 / 0.8941 | 31.76 / 0.8898 |  32.16 / 0.9054
Set14 [32] | 2x | 33.07 / 0.9106 |  33.70 / 0.9182 | 33.38 / 0.9237 |  33.70 / 0.9265
3x | 29.83 / 0.8308 |  30.16 / 0.8414 | 29.90 / 0.8504 | 30.12 / 0.8545
4x | 28.04 / 0.7669 | 28.41 / 0.7816 | 28.31 / 0.7928 |  28.70 / 0.8012
BSDS100 [20] | 2x | 31.80 / 0.8940 |  32.23 / 0.8999 | 32.01 / 0.9080 | 32.21 / 0.9102
3x | n.a. |  29.12 / 0.8060 | 28.79 / 0.8174 | 28.93 / 0.8211
4x | 27.25 / 0.7240 | 27.62 / 0.7355 | 27.38 / 0.7503 |  27.63 / 0.7573

####  5.0.3 Visual Comparison

Figure 3:  Comparison of an HR ground truth image (BSDS100, 2x scaling), DWSR,
and DWA. First row: the entire image space of the HR image and the
corresponding reconstructions obtained by the SR models. Second row: zoomed-in
regions within the images from the first row. Third row: residual image
representing the difference between the LR and HR images. As a result, the DWA
model captures edges and details closer to the ground truth residuals, as
opposed to the DWSR model (also regarding color).

Figure 3 displays the ground truth HR image alongside the DWSR and DWA
reconstructions. DWSR and DWA perform reasonably well in reconstructing the
images. However, the DWA reconstructions exhibit more accurate and sharp
details, particularly in the zoomed-in regions. Since the added bicubic
interpolation of the LR image in the reconstruction process provides a robust
base prediction, we also present the residual images, which are the
differences between the bicubic interpolations and the ground truth images, to
highlight the performance difference between both approaches.

Figure 4:  Feature maps of DWSR and DWA Direct. DWA Direct extracts local
contrasts and variations more effectively, closer than DWSR to the residual
target.

These residual images are the learning targets of the models to improve the
reconstruction quality beyond interpolation. By comparing the residual images,
we can see more clearly that the DWA model captures better distinct edges and
finer details, which are also closer to the ground truth residuals, as opposed
to the DWSR model. It has more substantial edges and finer points in the
residual images, which are also closer in color (see red colored lines of DWSR
reconstruction in Figure 3 as a comparison). This observation aligns with our
quantitative results, where DWA outperforms DWSR regarding various performance
metrics.

To provide deeper insights into our proposed models, Figure 4 presents feature
maps generated by the DWSR and DWA Direct models after the first layer. To
ensure diversity, we selected the top five channels from each method based on
the highest sum of distances between pairwise differences of all channels. Our
analysis reveals that although DWSR operates on the frequency space, it still
remains similar to the LR input and fails to capture the desired target
residual. In contrast, DWA Direct extracts local contrasts and variations more
effectively from the image space and performs better in mapping the target
residual.

##  6 Conclusion and Future Work

In this work, we presented a novel Differential Wavelet Amplifier (DWA)
module, which can be used as a drop-in module to existing wavelet-based SR
models. We showed experimentally on Set5, Set14, and BSDS100 for scaling
factors 2, 3, and 4 that it improves the reconstruction quality of the SR
models DWSR and MWCNN while enabling an application of them to the input image
space directly without harm to performance. This module captures more distinct
edges and finer details, which are closer to the ground truth residuals, which
wavelet-based SR models usually learn. This work is an opportunity to seek
further advancements for SR based on frequency-based representations.

For future work, an exciting research avenue would be to explore ways to
incorporate DWA on different DWT levels in MWCNN instead of only applying it
initially.

## Acknowledgments

This work was supported by the BMBF projects SustainML (Grant 101070408) and
by Carl Zeiss Foundation through the Sustainable Embedded AI project
(P2021-02-009).

## References

  * [1] Agarap, A.F.: Deep learning using rectified linear units (relu). arXiv preprint arXiv:1803.08375 (2018)
  * [2] Agarwal, A., Lang, J.: Foundations of analog and digital electronic circuits. Elsevier (2005)
  * [3] Agustsson, E., Timofte, R.: Ntire 2017 challenge on single image super-resolution: Dataset and study. In: Proceedings of the IEEE conference on computer vision and pattern recognition workshops. pp. 126–135 (2017)
  * [4] Bevilacqua, M., Roumy, A., Guillemot, C., Alberi-Morel, M.L.: Low-complexity single-image super-resolution based on nonnegative neighbor embedding (2012)
  * [5] Canh, T.N., Jeon, B.: Difference of convolution for deep compressive sensing. In: 2019 IEEE International Conference on Image Processing (ICIP). pp. 2105–2109. IEEE (2019)
  * [6] Dong, C., Loy, C.C., He, K., Tang, X.: Image super-resolution using deep convolutional networks. IEEE transactions on pattern analysis and machine intelligence 38(2), 295–307 (2015)
  * [7] Dong, C., Loy, C.C., Tang, X.: Accelerating the super-resolution convolutional neural network. In: European conference on computer vision. pp. 391–407. Springer (2016)
  * [8] Graves, A., Schmidhuber, J.: Offline handwriting recognition with multidimensional recurrent neural networks. Advances in neural information processing systems 21 (2008)
  * [9] Guo, T., Seyed Mousavi, H., Huu Vu, T., Monga, V.: Deep wavelet prediction for image super-resolution. In: Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition Workshops. pp. 104–113 (2017)
  * [10] He, K., Zhang, X., Ren, S., Sun, J.: Deep residual learning for image recognition. In: Proceedings of the IEEE conference on computer vision and pattern recognition. pp. 770–778 (2016)
  * [11] Horowitz, P., Hill, W., Robinson, I.: The art of electronics, vol. 2. Cambridge university press Cambridge (1989)
  * [12] Kingma, D.P., Ba, J.: Adam: A method for stochastic optimization. arXiv preprint arXiv:1412.6980 (2014)
  * [13] Knutsson, H., Westin, C.F.: Normalized and differential convolution. In: Proceedings of IEEE Conference on Computer Vision and Pattern Recognition. pp. 515–523. IEEE (1993)
  * [14] Laplante, P.A.: Comprehensive dictionary of electrical engineering. CRC Press (2018)
  * [15] Li, H., Yang, Y., Chang, M., Chen, S., Feng, H., Xu, Z., Li, Q., Chen, Y.: Srdiff: Single image super-resolution with diffusion probabilistic models. Neurocomputing 479, 47–59 (2022)
  * [16] Liang, J., Cao, J., Sun, G., Zhang, K., Van Gool, L., Timofte, R.: Swinir: Image restoration using swin transformer. In: Proceedings of the IEEE/CVF international conference on computer vision. pp. 1833–1844 (2021)
  * [17] Liu, P., Zhang, H., Zhang, K., Lin, L., Zuo, W.: Multi-level wavelet-cnn for image restoration. In: Proceedings of the IEEE conference on computer vision and pattern recognition workshops. pp. 773–782 (2018)
  * [18] Lowe, D.G.: Distinctive image features from scale-invariant keypoints. International journal of computer vision 60, 91–110 (2004)
  * [19] Ma, K., Duanmu, Z., Wu, Q., Wang, Z., Yong, H., Li, H., Zhang, L.: Waterloo exploration database: New challenges for image quality assessment models. IEEE Transactions on Image Processing 26(2), 1004–1016 (2016)
  * [20] Martin, D., Fowlkes, C., Tal, D., Malik, J.: A database of human segmented natural images and its application to evaluating segmentation algorithms and measuring ecological statistics. In: Proceedings Eighth IEEE International Conference on Computer Vision. ICCV 2001. vol. 2, pp. 416–423. IEEE (2001)
  * [21] Moser, B.B., Raue, F., Frolov, S., Palacio, S., Hees, J., Dengel, A.: Hitchhiker’s guide to super-resolution: Introduction and recent advances. IEEE Transactions on Pattern Analysis and Machine Intelligence pp. 1–21 (2023). https://doi.org/10.1109/TPAMI.2023.3243794
  * [22] Ronneberger, O., Fischer, P., Brox, T.: U-net: Convolutional networks for biomedical image segmentation. In: Medical Image Computing and Computer-Assisted Intervention–MICCAI 2015: 18th International Conference, Munich, Germany, October 5-9, 2015, Proceedings, Part III 18. pp. 234–241. Springer (2015)
  * [23] Sahak, H., Watson, D., Saharia, C., Fleet, D.: Denoising diffusion probabilistic models for robust image super-resolution in the wild. arXiv preprint arXiv:2302.07864 (2023)
  * [24] Saharia, C., Ho, J., Chan, W., Salimans, T., Fleet, D.J., Norouzi, M.: Image super-resolution via iterative refinement. IEEE Transactions on Pattern Analysis and Machine Intelligence (2022)
  * [25] Sarıgül, M., Ozyildirim, B.M., Avci, M.: Differential convolutional neural network. Neural Networks 116, 279–287 (2019)
  * [26] Stephane, M.: A wavelet tour of signal processing (1999)
  * [27] Sun, L., Dong, J., Tang, J., Pan, J.: Spatially-adaptive feature modulation for efficient image super-resolution. arXiv preprint arXiv:2302.13800 (2023)
  * [28] Vaswani, A., Shazeer, N., Parmar, N., Uszkoreit, J., Jones, L., Gomez, A.N., Kaiser, Ł., Polosukhin, I.: Attention is all you need. Advances in neural information processing systems 30 (2017)
  * [29] Wang, X., Xie, L., Dong, C., Shan, Y.: Real-esrgan: Training real-world blind super-resolution with pure synthetic data. In: Proceedings of the IEEE/CVF International Conference on Computer Vision (ICCV) Workshops. pp. 1905–1914 (October 2021)
  * [30] Wang, X., Yu, K., Wu, S., Gu, J., Liu, Y., Dong, C., Qiao, Y., Change Loy, C.: Esrgan: Enhanced super-resolution generative adversarial networks. In: Proceedings of the European conference on computer vision (ECCV) workshops. pp. 0–0 (2018)
  * [31] Wang, Z., Chen, J., Hoi, S.C.H.: Deep learning for image super-resolution: A survey. IEEE Transactions on Pattern Analysis and Machine Intelligence 43(10), 3365–3387 (2021). https://doi.org/10.1109/TPAMI.2020.2982166
  * [32] Zeyde, R., Elad, M., Protter, M.: On single image scale-up using sparse-representations. In: International conference on curves and surfaces. pp. 711–730. Springer (2010)


# Analysis of Heterogeneous Vehicular Traffic: Using Proportional Densities

Nanyondo Josephine<EMAIL_ADDRESS>Henry Kasumba
<EMAIL_ADDRESS>Department of Mathematics, Faculty of Science and
Education, Busitema University, P.O Box 236, Tororo, Uganda  Department of
Mathematics, School of Physical Sciences, Makerere University, P.O Box 7062,
Kampala, Uganda

###### Abstract

An extended multi-class Aw-Rascle (AR) model with pressure term described as a
function of area occupancy defined in form of proportional densities is
presented. Two vehicle classes that is; cars and motorcycles are considered
based on an assumption that proportions of these form total traffic density.
Qualitative properties of the proposed equilibrium velocity is established.
Conditions under which the proposed model is stable are determine by linear
stability analysis. To compute numerical flux, the model is discretized by the
original Roe decomposition scheme, where Roe matrix, averaged data variables
and wave strengths are explicitly derived. The Roe matrix is shown to be
hyperbolic, consistent and conservative. From the numerical results, the
effect of motorcycles proportion on the flow of vehicle classes is determined.
Results obtained remain within limits therefore, the proposed model is
realistic.

###### keywords:

Multi-class, Roe decomposition, heterogeneous traffic, area occupancy,
macroscopic, Aw-Rascle model, proportional densities

††journal: Appropriate journal to be identifiedt1t1footnotetext: This article
is a result of the PhD program funded by Makerere-Sida Bilateral Program under
Project 316 "Capacity Building in Mathematics and Its Applications".

##  1 Introduction

In developing countries, various vehicle classes use same roads at the same
time. This is why traffic flow models that combine vehicle classes are
necessary in the study of heterogeneous vehicular traffic flow. The traffic
flow models are classified as macroscopic, microscopic or mesoscopic.
Macroscopic models consider traffic flow as a compressible fluid formed by
vehicles [1] and capture general relationship between flow, density and
velocity of vehicles [2, 3]. Microscopic models describe interactions between
vehicles within a traffic stream [4, 5]. The microscopic models consider
characteristics such as; vehicle lengths, speeds, accelerations, time, space
headways, vehicle and engine capabilities and human characteristics that
describe the driving behaviors. Mesoscopic models are described by a
combination of microscopic and macroscopic characteristics through probability
distributions [2] but at a high level of aggregation usually by speed-density
relations and queuing theory approaches [3].

Traffic flow is either homogeneous, heterogeneous, equilibrium or non-
equilibrium. In homogeneous traffic, all vehicles follow lane discipline
whereas in heterogeneous traffic, vehicles violate lane discipline.
Equilibrium flow is characterized by constant velocity and spatial
homogeneity. Whereas, non-equilibrium traffic flow is characterized by changes
in both velocity and spatial homogeneity [2].

The macroscopic models are adopted because of their low computation
complexity. The first macroscopic model was proposed by [6] and [7], given by

$\dfrac{\partial\rho}{\partial t}+\dfrac{\partial(\rho v)}{\partial x}=0.$

To address the drawbacks of the first order methods in characterizing
overtaking and creeping some researchers (see for example, [8] and the
references therein) used various forms of definitions of the velocity function
(in static form) together with the equation for the conservation of mass.
However, the exact formulation of the velocity function for each class of
vehicles is a challenge. Using the porous flow approach, Nair et al [9]
formulated a multi-class heterogeneous traffic flow model using first order
methods. The same theory is applied also in [10, 11] to formulate a model for
heterogeneous traffic flow that takes overtaking and creeping into account.
Since these models are based on the first order method they may fail to
characterize some important traffic phenomena that occurs in heterogeneous
traffic flow, like, stop and go wave, hysteresis [5], platoon dispersion [12].

Second and higher order models combine the dynamics of the velocity function
together with the conservation equation of vehicles. For instance Payne [13,
14] has shown that the microscopic car-following model can only be generalized
by a second and higher order continuum models. Gupta and Dhiman [15] have
shown that a heterogeneous continuum second order model generalizes a
kinematic car-following model for a non lane-based system of traffic flow that
takes lateral separation into account. However, in their model, Gupta and
Dhiman used the same equilibrium velocity for all vehicle classes. In order to
characterize the phenomena of creeping (or gap filling), Mohan and Ramadurai
[12, 16] used lateral area occupancy method to analyze the flow of
heterogeneous traffic. However, they ignored the concept of proportional
densities and their effect on the flow of vehicle classes.

In this paper, the effect of motorcycles proportion on the flow of vehicle
classes is determined. The subsequent subsections of the paper are organized
as follows. Derivation of area occupancy defined in terms of proportional
densities and conservative form of the proposed model are shown in Section
$2$. In Section $3$, qualitative properties of equilibrium velocity are
discussed. A discretization scheme, Quasi-linear form, original Roe numerical
scheme in which Roe matrix, averaged data variables, wave strengths and
entropy condition are shown in Section $4.$ Stability analysis of the proposed
model equations is done in Section 5. Numerical results are shown in Section
$6$. Finally, Section $7$ draws conclusion on the paper.

##  2 Model formulation

In this section, the Aw-Rascle model is modified by introducing an area
occupancy that is expressed in terms of vehicle class proportions and other
lateral parameters in the pressure term.

###  2.1 Applying area occupancy for heterogeneous traffic

We then introduce area occupancy ($AO$) that is defined in terms of
proportional densities and derive it as follows. Suppose that the total
density, $\rho$ of vehicles that are occupying the road, at any time is given
by

$\rho=\rho_{m}+\rho_{c},$ (1)

where $\rho_{c}$ refers to density of cars and $\rho_{m}$ density of
motorcycles. Let

$\rho_{m}=\delta(\rho_{m}+\rho_{c})=\delta\rho,$ (2)

where $\delta$ represents the proportion of motorcycles that contribute to the
total vehicular density on a given road section. Then

$\rho_{m}=\dfrac{\delta\rho_{c}}{1-\delta}.$ (3)

Assuming that the width of a motorcycle is one third the width of a car and
considering $l,w$ to be length and width of vehicles, the equation for area
occupancy becomes

$AO=\dfrac{\sum_{i}\rho_{i}a_{i}}{W}=\dfrac{\rho_{c}a_{c}+\rho_{m}a_{m}}{W}=\dfrac{(1-\delta)\rho
a_{c}+\delta\rho a_{m}}{W},\leavevmode\nobreak\ i=c,m.$ (4)

Substitution of (3) into (4) leads to

$AO=\dfrac{1}{W}\left(\rho_{c}a_{c}+\dfrac{\delta\rho_{c}a_{m}}{1-\delta}\right)=\dfrac{\rho_{c}\left[(1-\delta)a_{c}+\delta
a_{m}\right]}{W(1-\delta)},$ (5)

where $a_{i}=l_{i}\times w_{i}.$ Further simplification yields

$\displaystyle AO$ $\displaystyle=$
$\displaystyle\dfrac{(1-\delta)\rho_{c}l_{c}w_{c}+\delta\rho_{c}l_{m}w_{m}}{W(1-\delta)},\leavevmode\nobreak\
w_{m}=\dfrac{1}{3}w_{c},$ $\displaystyle=$
$\displaystyle\dfrac{\rho_{c}w_{c}}{W(1-\delta)}\left((1-\delta)l_{c}+\dfrac{1}{3}\delta
l_{m}\right),$ $\displaystyle=$ $\displaystyle\rho_{c}\psi_{c},$

where

$\psi_{c}=\dfrac{w_{c}}{W(1-\delta)}\left((1-\delta)l_{c}+\dfrac{1}{3}\delta
l_{m}\right),$ (7)

and $AO$ refers to area occupancy expressed in terms of proportional
densities. Equivalently, $AO$ can be rewritten as

$\displaystyle AO$ $\displaystyle=$
$\displaystyle\dfrac{\rho_{m}w_{c}(l_{c}(1-\delta)+\dfrac{1}{3}\delta
l_{m})}{W\delta},$ $\displaystyle=$ $\displaystyle\rho_{m}\psi_{m},$

where

$\psi_{m}=\dfrac{w_{c}(l_{c}(1-\delta)+\dfrac{1}{3}\delta l_{m})}{W\delta}.$
(9)

It can be inferred from (3) that

$\displaystyle\rho_{c}=\dfrac{\rho_{m}(1-\delta)}{\delta}=(1-\delta)\rho.$
(10)

Consequently the total density in (1) then becomes

$\rho=\rho_{m}+\rho_{c}=\delta\rho+(1-\delta)\rho.$

Incorporating the two vehicle classes in AR model, we obtain the following
proposed model

$\displaystyle\begin{cases}\dfrac{\partial\rho_{m}}{\partial
t}+\dfrac{\partial\left(\rho_{m}v_{m}\right)}{\partial x}=0,\\\\[6.45831pt]
\dfrac{\partial\left(v_{m}+p_{m}\right)}{\partial
t}+v_{m}\dfrac{\partial\left(v_{m}+p_{m}\right)}{\partial
x}=\dfrac{1}{\tau_{m}}\left(v_{em}-v_{m}\right),\\\\[6.45831pt]
\dfrac{\partial\rho_{c}}{\partial
t}+\dfrac{\partial\left(\rho_{c}v_{c}\right)}{\partial x}=0,\\\\[6.45831pt]
\dfrac{\partial\left(v_{c}+p_{c}\right)}{\partial
t}+v_{c}\dfrac{\partial\left(v_{c}+p_{c}\right)}{\partial
x}=\dfrac{1}{\tau_{c}}\left(v_{ec}-v_{c}\right),\end{cases}$ (11)

where the respective pressure terms are given by

$p_{m}=\left(\psi_{m}\rho_{m}\right)^{\gamma_{m}}\leavevmode\nobreak\
\text{and}\leavevmode\nobreak\
p_{c}=\left(\psi_{c}\rho_{c}\right)^{\gamma_{c}}.$

Here $\psi_{m},\leavevmode\nobreak\ \psi_{c}$ are parameters defined as in
(7), (9) and $\gamma_{m},\leavevmode\nobreak\ \gamma_{c}$ represent the
relaxation time of motorcycles and cars, respectively.

We can express (11) in conservative form given by

$\begin{cases}\dfrac{\partial\rho_{m}}{\partial
t}+\dfrac{\partial\left(\rho_{m}v_{m}\right)}{\partial x}=0,\\\\[6.45831pt]
\dfrac{\partial\left(\rho_{m}\left(v_{m}+p_{m}\right)\right)}{\partial
t}+\dfrac{\partial\left(\rho_{m}v_{m}\left(v_{m}+p_{m}\right)\right)}{\partial
x}=\dfrac{\rho_{m}}{\tau_{m}}\left(v_{em}-v_{m}\right),\\\\[6.45831pt]
\dfrac{\partial\rho_{c}}{\partial
t}+\dfrac{\partial\left(\rho_{c}v_{c}\right)}{\partial x}=0,\\\\[6.45831pt]
\dfrac{\partial\left(\rho_{c}\left(v_{c}+p_{c}\right)\right)}{\partial
t}+\dfrac{\partial\left(\rho_{c}v_{c}\left(v_{c}+p_{c}\right)\right)}{\partial
x}=\dfrac{\rho_{c}}{\tau_{c}}\left(v_{ec}-v_{c}\right).\end{cases}$ (12)

where the conserved variables are $\rho_{c},\leavevmode\nobreak\
\rho_{c}\left(v_{c}+p_{c}\right),\leavevmode\nobreak\ \rho_{m}$ and
$\rho_{m}\left(v_{m}+p_{m}\right),$
Considering terms in system (12) with subscript $c,$ let the generalized
momentum be $X_{c},$ such that

$X_{c}=\rho_{c}\left(v_{c}+p_{c}\right).$ (13)

Rearranging terms in (13) we obtain

$\displaystyle\dfrac{X_{c}}{\rho_{c}}=\left(v_{c}+p_{c}\right),\leavevmode\nobreak\
\text{and}\leavevmode\nobreak\ \rho_{c}v_{c}=X_{c}-\rho_{c}p_{c}.$

This implies

$\rho_{c}v_{c}\left(v_{c}+p_{c}\right)=\left(X_{c}-\rho_{c}p_{c}\right)\dfrac{X_{c}}{\rho_{c}},$

which simplifies into

$\dfrac{X_{c}^{2}}{\rho_{c}}-p_{c}X_{c}.$

The same is done on terms with subscript $m$ to obtain the following conserved
model equations.

$\begin{cases}\dfrac{\partial\rho_{m}}{\partial
t}+\dfrac{\partial\left(X_{m}-\rho_{m}p_{m}\right)}{\partial
x}=0,\\\\[6.45831pt] \dfrac{\partial X_{m}}{\partial
t}+\dfrac{\partial\left(\dfrac{X_{m}^{2}}{\rho_{m}}-p_{m}X_{m}\right)}{\partial
x}=\dfrac{\rho_{m}}{\tau_{m}}\left(v_{em}-v_{m}\right),\\\\[6.45831pt]
\dfrac{\partial\rho_{c}}{\partial
t}+\dfrac{\partial\left(X_{c}-\rho_{c}p_{c}\right)}{\partial
x}=0,\\\\[6.45831pt] \dfrac{\partial X_{c}}{\partial
t}+\dfrac{\partial\left(\dfrac{X_{c}^{2}}{\rho_{c}}-p_{c}X_{c}\right)}{\partial
x}=\dfrac{\rho_{c}}{\tau_{c}}\left(v_{ec}-v_{c}\right).\end{cases}$ (14)

In vector form, system (14) becomes

$\dfrac{\partial U}{\partial t}+\dfrac{\partial f(U)}{\partial x}=S(U),$ (15)

where

$U=\begin{bmatrix}\rho_{m}\\\\[6.45831pt] X_{m}\\\\[6.45831pt]
\rho_{c}\\\\[6.45831pt] X_{c}\end{bmatrix},\leavevmode\nobreak\
f(U)=\begin{bmatrix}X_{m}-\rho_{m}p_{m}\\\\[6.45831pt]
\dfrac{X_{m}^{2}}{\rho_{m}}-p_{m}X_{m}\\\\[6.45831pt]
X_{c}-\rho_{c}p_{c}\\\\[6.45831pt]
\dfrac{X_{c}^{2}}{\rho_{c}}-p_{c}X_{c}\end{bmatrix},\leavevmode\nobreak\ \\\
S(U)=\begin{bmatrix}0\\\\[6.45831pt]
\dfrac{\rho_{m}}{\tau_{m}}\left(v_{em}-v_{m}\right)\\\\[6.45831pt]
0\\\\[6.45831pt]
\dfrac{\rho_{c}}{\tau_{c}}\left(v_{ec}-v_{c}\right)\end{bmatrix}.$

##  3 Qualitative properties of Equilibrium velocity

We adopt a modified version of the Greenshield’s equilibrium velocity

$\displaystyle\begin{cases}v^{e}_{m}=v_{m\max}\left(1-\dfrac{AO_{m}}{AO_{m}^{\max}}\right),\leavevmode\nobreak\
\text{if}\leavevmode\nobreak\ AO_{m}\leq AO_{m}^{\max},\\\\[6.45831pt]
v^{e}_{c}=v_{c\max}\left(1-\dfrac{AO_{c}}{AO_{c}^{\max}}\right),\leavevmode\nobreak\
\text{if}\leavevmode\nobreak\ AO_{c}\leq AO_{c}^{\max},\\\\[6.45831pt]
v^{e}_{i}=0\leavevmode\nobreak\ \text{otherwise,}\leavevmode\nobreak\
i=m,c,\end{cases}$

where $v_{m\max},\leavevmode\nobreak\ v_{c\max},\leavevmode\nobreak\
AO_{m}^{\max},\leavevmode\nobreak\ AO_{c}^{\max}$ are the respective maximum
velocity and maximum area occupancy of the vehicle classes. Velocities were
defined in terms of $AO$ since it is a standard measure of heterogeneous
traffic concentration that lacks in lane discipline. Note that by writing $AO$
in the form $AO_{i}=\psi_{i}\rho_{i},\leavevmode\nobreak\
i=m,\leavevmode\nobreak\ c,$ and using definitions

$w_{m}=\dfrac{w_{c}}{3}\leavevmode\nobreak\ \text{and}\leavevmode\nobreak\
a_{i}=l_{i}*w_{i},\leavevmode\nobreak\ i=m,c,$

in (2.1) gives

$v^{e}_{m}=v_{m\max}\left(1-\dfrac{\rho_{m}\left[a_{c}\left(1-\delta\right)+\delta
a_{m}\right]}{0.85W\delta}\right).$ (16)

Similarly, using (2.1) yields

$v^{e}_{c}=v_{c\max}\left(1-\dfrac{\rho_{c}\left[a_{c}\left(1-\delta\right)+\delta
a_{m}\right]}{0.74W(1-\delta)}\right).$ (17)

The dependence of equilibrium speed of the flow, $v^{e}_{i},$ on the area
occupied by all vehicles is shown in (18) and (19) by taking derivatives of
$dv^{e}_{c},$ with respect to $AO$ to obtain:

$\displaystyle\dfrac{dv^{e}_{m}}{dAO_{m}}$ $\displaystyle=$
$\displaystyle-\dfrac{v_{m\max}}{AO_{m}^{\max}},$ (18)
$\displaystyle\dfrac{dv^{e}_{c}}{dAO_{c}}$ $\displaystyle=$
$\displaystyle-\dfrac{v_{c\max}}{AO_{c}^{\max}}.$ (19)

From the derivatives above observe that the speed of traffic flow decreases as
the area occupied by vehicles increases.

By taking derivatives with respect to the motorcycles proportion $\delta,$ we
show how the speed of each vehicle class is affected by the variation in the
proportion. We make use of 16 to derive

$\displaystyle\dfrac{dv^{e}_{m}}{d\delta}$ $\displaystyle=$
$\displaystyle\dfrac{\rho_{m}a_{c}v_{m\max}}{0.85\delta^{2}W},$ (20)

It is concluded from (20) that motorcycles’ equilibrium speed increases with
of their proportion. For the case of cars, 17 is used in deriving

$\displaystyle\dfrac{dv^{e}_{c}}{d\delta}$ $\displaystyle=$
$\displaystyle-\dfrac{a_{m}\rho_{c}v_{c\max}}{W(1-\delta)^{2}}.$ (21)

It can be deduced from (20) that an increase in motorcycles proportion leads
to a decrease in velocity of cars and vise versa. Next, we look at the
dynamics of equilibrium velocity when density is varied. That is achieved
through derivatives

$\dfrac{dv_{m}^{e}}{d\rho}=-v_{m\max}\left(\dfrac{\left[a_{c}\left(1-\delta\right)+\delta
a_{m}\right]}{0.85W}\right),$ (22)

and

$\dfrac{dv_{c}^{e}}{d\rho}=-v_{c\max}\left(\dfrac{\left[a_{c}\left(1-\delta\right)+\delta
a_{m}\right]}{0.74W}\right),$ (23)

where $\rho_{m}=\delta\rho,\leavevmode\nobreak\ \text{and}\leavevmode\nobreak\
\rho_{c}=(1-\delta)\rho.$ It is observed from (22) and (23) that equilibrium
velocity of vehicles decreases with increase of density. This behavior is
further evidenced in Figure 1. In Figures 1(a) and 1(b), cars are observed to
stop much faster than motorcycles irrespective of motorcycles proportion. When
the proportion is very large, motorcycles keep moving even when cars have
stopped. It is because they can ably maneuver or filter through spaces
unfilled by cars. Figure 1(c) predicts jam density reached at much earlier if
cars dominate the road than when motorcycles do. It can be concluded that
total flow of vehicles increases with increase of motorcycles proportion.

(a) Equilibrium velocity of motorcycles vs total density

(b) Equilibrium velocity of cars vs total density

(c) Total flow vs total density

Figure 1: Equilibrium velocities and total flow of the two vehicle classes
versus total density, with varying values of motorcycles’ proportion,
$\delta.$

##  4 Discretization scheme

The proposed model (15) in a space-time domain
$\displaystyle{[0,J]\times[0,N]}$ is considered. Solving it would mean
evolving the solution $U(x,t)$ in time starting from the initial condition
$U(x,0)$ and subject to boundary conditions. The space and time grids are
divided into $J$ and $N$ equally spaced grid points $x_{j}=j\Delta
x,\leavevmode\nobreak\ j=0,\ldots,J$ and $t^{n}=n\Delta t,\leavevmode\nobreak\
n=0,\ldots,N,$ of length $\Delta x$ and $\Delta t,$ respectively. Approximate
or discrete values of $U(x,t),$ for the data given by (15) at a time step $n$
and spatial position $j,$ are represented by $U_{j}^{n}\equiv U\left(j\Delta
x,n\Delta t\right)\equiv U\left(x_{j},t^{n}\right).$ Next, periodic conditions
are applied to the left and right boundary points $U_{0}^{n}$ and $U_{J}^{n},$
respectively. Since at any time step $n,$ the numerical solutions are obtained
as either initial data values or values computed in a previous time step, a
conservative finite-difference scheme is then applied to (15). Hence, the
solution $U_{j}^{n+1},$ at the next time step is obtained as

$\displaystyle U_{j}^{n+1}$ $\displaystyle=$ $\displaystyle
U_{j}^{n}-\dfrac{\Delta t}{\Delta
x}\left(F(U_{j}^{n},U_{j+1}^{n})-F(U_{j-1}^{n},U_{j}^{n})\right)+\Delta
tS(U_{i}^{n}),$ (24)

where $F(U_{j}^{n},U_{j+1}^{n})$ in (24) denotes Roe’s numerical flux [17] and
is given by

$\displaystyle
F(U_{j}^{n},U_{j+1}^{n})=\dfrac{1}{2}\left(F(U_{j}^{n})+F(U_{j+1}^{n})\right)-\dfrac{1}{2}\left(B(U^{n}_{j+1/2})\left(U^{n}_{j+1}-U^{n}_{j}\right)\right),$
(25)

where $B(U^{n}_{j+1/2})$ represents the averaged Jacobian matrix of $f(U)$.
For each time step, scheme (24) is implemented in two intermediate steps.
Firstly, one finds

$U_{j}^{*}=U_{j}^{n}-\dfrac{\Delta t}{\Delta
x}\left(F(U_{j}^{n},U_{j+1}^{n})-F(U_{j-1}^{n},U_{j}^{n})\right),$ (26)

obtained by ignoring the non-homogeneous part. Secondly, the solution
$U_{j}^{n+1}$ at time step $t^{n+1}$ is found by adding the non-homogeneous
part to (26). That is to say;

$U_{j}^{n+1}=U_{j}^{*}+\Delta tS(U_{j}^{n}).$

Below we describe the Roe’s numerical scheme that is utilized in this work.

###  4.1 Roe decomposition

The Roe decomposition scheme is used because it is capable of capturing abrupt
changes or discontinuities in the data variables and providing accurate
numerical solutions for traffic flow models [2]. Roe’s scheme is employed
purposely to transform the proposed model (11) into quasi-linear form and
linearize locally by approximating the Jacobian matrix, $B(U)$ with Roe
averages and in every time step, the procedure is repeated [18]. The resulting
system can then approximate speed found from the eigenvalues of the averaged
Jacobian matrix, as described below [18].

####  4.1.1 Quasi-linear form of model equations

The model equations are expressed in Quasi-linear form which is suitable for
carrying out linearization. The homogeneous part of system (14) can be
rewritten in quasi-linear form as

$\dfrac{\partial U}{\partial t}+B(U)\dfrac{\partial U}{\partial x}=0,$

where

$\displaystyle B(U)$ $\displaystyle=$ $\displaystyle\dfrac{\partial
f}{\partial U},$ (27) $\displaystyle=$
$\displaystyle\begin{bmatrix}-(\gamma_{m}+1)p_{m}&1&0&0\\\
B_{21}&\dfrac{2X_{m}}{\rho_{m}}-p_{m}&0&0\\\ 0&0&-(\gamma_{c}+1)p_{c}&1\\\
0&0&-\left(\dfrac{X_{c}^{2}}{\rho_{c}^{2}}+\dfrac{\gamma_{c}p_{c}X_{c}}{\rho_{c}}\right)&\dfrac{2X_{c}}{\rho_{c}}-p_{c}\end{bmatrix}$

is the Jacobian matrix with
$B_{21}=-\left(\dfrac{X_{m}^{2}}{\rho_{m}^{2}}+\dfrac{\gamma_{m}p_{m}X_{m}}{\rho_{m}}\right).$
Next, we require eigenvalues and eigenvectors of (14) to be applied in the Roe
numerical flux and entropy fix. Note that the eigenvalues are obtained from
the characteristic equations of $B(U)$ and are given as follows:

$\lambda_{1}(U)=v_{m},\leavevmode\nobreak\
\lambda_{2}(U)=v_{m}-\gamma_{m}p_{m},\leavevmode\nobreak\
\lambda_{3}(U)=v_{c},\leavevmode\nobreak\
\lambda_{4}(U)=v_{c}-\gamma_{c}p_{c}.$ (28)

The right eigenvectors corresponding to (28) are obtained as:

$\displaystyle r_{1}=\left(\begin{array}[]{c}1\\\ (1+\gamma_{m})p_{m}+v_{m}\\\
0\\\ 0\end{array}\right),\leavevmode\nobreak\ r_{2}$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{c}1\\\ v_{m}+p_{m}\\\ 0\\\
0\end{array}\right),\leavevmode\nobreak\ $ (37) $\displaystyle
r_{3}=\left(\begin{array}[]{c}0\\\ 0\\\ 1\\\
(1+\gamma_{c})p_{c}+v_{c}\end{array}\right),\leavevmode\nobreak\ r_{4}$
$\displaystyle=$ $\displaystyle\left(\begin{array}[]{c}0\\\ 0\\\ 1\\\
v_{c}+p_{c}\end{array}\right).$ (46)

We then proceed by applying the original Roe method to (14) in order to
compute Roe averages. The detailed procedure of the method is described in the
following steps [17]. Firstly, a parameter vector $Z,$ is introduced such that
$U$ and $f(U)$ of (15) can be rewritten as

$U=U(Z),\leavevmode\nobreak\ f(U)=f(U(Z)).$

Next, the jumps

$\Delta U=U_{r}-U_{l},\leavevmode\nobreak\ \Delta f=f(U_{r})-f(U_{l}),$ (47)

are expressed in terms of the change

$\Delta Z=Z_{r}-Z_{l}.$

Letting the parameter vector

$Z=\left[\begin{array}[]{c}\alpha_{1}\\\ \beta_{1}\\\ \alpha_{2}\\\
\beta_{2}\end{array}\right]=\left[\begin{array}[]{c}\sqrt{\rho_{m}}\\\
\dfrac{X_{m}}{\sqrt{\rho_{m}}}\\\ \sqrt{\rho_{c}}\\\
\dfrac{X_{c}}{\sqrt{\rho_{c}}}\end{array}\right],$

implies that $\alpha_{1}^{2}=\rho_{m},\leavevmode\nobreak\
X_{m}=\alpha_{1}\beta_{1},\leavevmode\nobreak\
\alpha_{2}^{2}=\rho_{c}\leavevmode\nobreak\ \text{and}\leavevmode\nobreak\
X_{c}=\alpha_{2}\beta_{2}.$ So, the vectors $U=U(Z)\leavevmode\nobreak\
\text{and}\leavevmode\nobreak\ f=f(Z)$ become

$U=\left[\begin{array}[]{c}\alpha_{1}^{2}\\\ \alpha_{1}\beta_{1}\\\
\alpha_{2}^{2}\\\ \alpha_{2}\beta_{2}\end{array}\right],\leavevmode\nobreak\
f=\left[\begin{array}[]{c}\alpha_{1}\beta_{1}-\psi^{\gamma}\alpha_{1}^{2(1+\gamma_{m})}\\\
\beta_{1}^{2}-\psi^{\gamma}\beta_{1}\alpha_{1}^{(1+2\gamma_{m})}\\\
\alpha_{2}\beta_{2}-\psi^{\gamma}\alpha_{2}^{2(1+\gamma_{c})}\\\
\beta_{2}^{2}-\psi^{\gamma}\beta_{2}\alpha_{2}^{(1+2\gamma_{c})}\end{array}\right].$

The averages of vector $Z$ then become
$\overline{Z}=[\overline{\alpha}_{1},\leavevmode\nobreak\
\overline{\beta}_{1},\leavevmode\nobreak\
\overline{\alpha}_{2},\leavevmode\nobreak\ \overline{\beta}_{2}]^{T}.$ These
are found by simple arithmetic averaging and obtained as:

$\displaystyle\overline{Z}=\dfrac{1}{2}\left[\begin{array}[]{c}Z_{l}+Z_{r}\end{array}\right]=\dfrac{1}{2}\left[\begin{array}[]{c}\alpha_{1l}+\alpha_{1r}\\\
\beta_{1l}+\beta_{1r}\\\ \alpha_{2l}+\alpha_{2r}\\\
\beta_{2l}+\beta_{2r}\end{array}\right].$

Following work by [17], two matrices $\overline{B}=\overline{B}(\overline{Z})$
and $\overline{C}=\overline{C}(\overline{Z})$ are computed such that the jumps
$\Delta U$ and $\Delta f$ in (47) are expressed in terms of the jump $\Delta
Z.$ That is;

$\Delta U=\overline{B}\Delta Z\leavevmode\nobreak\
\text{and}\leavevmode\nobreak\ \Delta f=\overline{C}\Delta Z.$ (49)

Using (49) leads to $\Delta f=\overline{C}\overline{B}^{-1}\Delta U,$ and
comparing it with property (C) of Theorem 1, equation (54), results into the
Roe matrix $\overline{A}=\overline{C}\overline{B}^{-1}.$ Matrices
$\overline{B}$ and $\overline{C}$ satisfying (49) are:

$\displaystyle\overline{B}=\dfrac{dU}{dZ}_{\mid{\overline{Z}}}=\left(\begin{array}[]{cccc}2\overline{\alpha}_{1}&0&0&0\\\
\overline{\beta}_{1}&\overline{\alpha}_{1}&0&0\\\
0&0&2\overline{\alpha}_{2}&0\\\
0&0&\overline{\beta}_{2}&\overline{\alpha}_{2}\end{array}\right),$

whose inverse is computed by applying the usual definition
$\overline{B}^{-1}=\dfrac{1}{Det(\overline{B})}\times Adj(\overline{B})$ and

$\displaystyle\overline{C}$ $\displaystyle=$ $\displaystyle
f^{\prime}(U(Z))\mid_{\overline{Z}},$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{cccc}\overline{\beta}_{1}-2(\gamma_{m}+1)\psi_{m}^{\gamma_{m}}\overline{\alpha}_{1}^{(1+2\gamma_{m})}&\overline{\alpha}_{1}&0&0\\\
-(1+2\gamma_{m})\psi_{m}^{\gamma_{m}}\overline{\beta}_{1}\overline{\alpha}_{1}^{2\gamma_{m}}&\overline{C}_{22}&0&0\\\
0&0&\overline{\beta}_{2}-2(\gamma_{c}+1)\psi_{c}^{\gamma_{c}}\overline{\alpha}_{2}^{(1+2\gamma_{c})}&\overline{\alpha}_{2}\\\
0&0&-(1+2\gamma_{c})\psi_{c}^{\gamma_{c}}\overline{\beta}_{2}\overline{\alpha}_{2}^{2\gamma_{c}}&\overline{C}_{44}\end{array}\right),$

where
$\overline{C}_{22}=2\overline{\beta}_{1}-\psi^{\gamma}_{m}\overline{\alpha}_{1}^{(1+2\gamma_{m})},\leavevmode\nobreak\
\overline{C}_{44}=2\overline{\beta}_{2}-\psi^{\gamma}_{c}\overline{\alpha}_{2}^{(1+2\gamma_{c})}.$
Therefore, the required Roe averaged matrix is obtained as

$\overline{A}=\left(\begin{array}[]{cccc}\overline{A}_{11}&1&0&0\\\
\overline{A}_{21}&2\left(\dfrac{\overline{\beta}_{1}}{\overline{\alpha}_{1}}\right)-\psi_{m}^{\gamma_{m}}\overline{\alpha}_{1}^{2\gamma_{m}}&0&0\\\
0&0&-(1+\gamma_{c})\psi_{c}^{\gamma_{c}}\overline{\alpha}_{2}^{2\gamma_{c}}&1\\\
0&0&\overline{A}_{43}&2\left(\dfrac{\overline{\beta}_{2}}{\overline{\alpha}_{2}}\right)-\psi_{c}^{\gamma_{c}}\overline{\alpha}_{2}^{2\gamma_{c}}\end{array}\right),$
(52)

where
$\overline{A}_{11}=-(1+\gamma_{m})\psi_{m}^{\gamma_{m}}\overline{\alpha}_{1}^{2\gamma_{m}},\leavevmode\nobreak\
\overline{A}_{21}=-\left(\gamma_{m}\psi^{\gamma_{m}}\overline{\alpha}_{1}^{2\gamma_{m}}\dfrac{\overline{\beta}_{1}}{\overline{\alpha}_{1}}+\dfrac{\overline{\beta}_{1}^{2}}{\overline{\alpha}_{1}^{2}}\right),\\\
\overline{A}_{43}=-\left(\gamma_{c}\psi^{\gamma_{c}}\overline{\alpha}_{2}^{2\gamma_{c}}\dfrac{\overline{\beta}_{2}}{\overline{\alpha}_{2}}+\dfrac{\overline{\beta}_{2}^{2}}{\overline{\alpha}_{2}^{2}}\right)\leavevmode\nobreak\
\text{and}\leavevmode\nobreak\
\dfrac{\overline{\beta}_{1}}{\overline{\alpha}_{1}},$
$\dfrac{\overline{\beta}_{2}}{\overline{\alpha}_{2}}$ represent the sum of Roe
averaged velocity and pressure of motorcycles and cars, respectively. These
are given by

$\displaystyle\dfrac{\overline{\beta}_{1}}{\overline{\alpha}_{1}}=\dfrac{\dfrac{X_{mr}}{\sqrt{\rho_{mr}}}+\dfrac{X_{ml}}{\sqrt{\rho_{ml}}}}{\sqrt{\rho_{mr}}+\sqrt{\rho_{ml}}}=\overline{v}_{m}+\overline{p}_{m},\leavevmode\nobreak\
\dfrac{\overline{\beta}_{2}}{\overline{\alpha}_{2}}=\dfrac{\dfrac{X_{cr}}{\sqrt{\rho_{cr}}}+\dfrac{X_{cl}}{\sqrt{\rho_{cl}}}}{\sqrt{\rho_{cr}}+\sqrt{\rho_{cl}}}=\overline{v}_{c}+\overline{p}_{c}.$
(53)

Using the Roe averaged velocity, $\overline{v}$ and pressure, $\overline{p},$
the data variables of the proposed model can be approximated over a given road
section. Comparing (52) with the matrix in (27), no averaged $\overline{\rho}$
is required. Having derived the Roe matrix, $\overline{A}$ one ensures that it
satisfies properties in Theorem 1 below.

######  Theorem 1 (Toro [17])

The proposed model is hyperbolic if the Roe Jacobian matrix $\overline{A}$
satisfies the three properties namely:

  1. (A)

The proposed model is hyperbolic if $\overline{A}$ has real eigenvalues
$\overline{\lambda}_{i}=\overline{\lambda}_{i}(U_{l},\leavevmode\nobreak\
U_{r}),$ where
$\overline{\lambda}_{1}\leq\overline{\lambda}_{2}\leq\cdots\leq\overline{\lambda}_{m}$
and $\overline{r}_{1},\leavevmode\nobreak\
\overline{r}_{2},\ldots,\overline{r}_{m}$ linearly independent right
eigenvectors.

  2. (B)

Consistency with the exact Jacobian $\overline{A}(U,U)=B(U).$

  3. (C)

Conservation across discontinuities

$f(U_{r})-f(U_{l})=\overline{A}(U_{r}-U_{l}).$ (54)

Recalling the Jacobian matrix, eigenvalues and corresponding right
eigenvectors of the proposed model given in (27) - (46), the averaged
eigenvalues are deduced as

$\overline{\lambda}_{1}(\overline{U})=\overline{v}_{m}-\gamma_{m}\overline{p}_{m},\leavevmode\nobreak\
\overline{\lambda}_{2}(\overline{U})=\overline{v}_{m},\leavevmode\nobreak\
\overline{\lambda}_{3}(\overline{U})=\overline{v}_{c}-\gamma_{c}\overline{p}_{c},\leavevmode\nobreak\
\overline{\lambda}_{4}(\overline{U})=\overline{v}_{c}.$

These are shown to be real. Their corresponding averaged right eigenvectors

$\displaystyle\overline{r}_{1}$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{c}1\\\
\overline{v}_{m}+\overline{p}_{m}\\\ 0\\\
0\end{array}\right),\leavevmode\nobreak\
\overline{r}_{2}=\left(\begin{array}[]{c}1\\\
\overline{v}_{m}+(1+\gamma_{m})\overline{p}_{m}\\\ 0\\\
0\end{array}\right),\leavevmode\nobreak\ $ $\displaystyle\overline{r}_{3}$
$\displaystyle=$ $\displaystyle\left(\begin{array}[]{c}0\\\ 0\\\ 1\\\
\overline{v}_{c}+\overline{p}_{c}\end{array}\right),\leavevmode\nobreak\
\overline{r}_{4}=\left(\begin{array}[]{c}0\\\ 0\\\ 1\\\
\overline{v}_{c}+(1+\gamma_{c})\overline{p}_{c}\end{array}\right),$

are checked to be linearly independent. Hence, property (A) of Theorem 1 is
satisfied. Property (B), is proved by evaluating the derived Roe matrix (52)
at $U=\left[\rho_{m},\leavevmode\nobreak\ X_{m},\leavevmode\nobreak\
\rho_{c},\leavevmode\nobreak\ X_{c}\right]^{T}.$ Property $(C)$ is proved by
substituting the conserved variables into (52). To find the wave strengths,
$s_{i}$ we project the jump, $\Delta U=U_{r}-U_{l}$ onto the right
eigenvectors and obtain

$\displaystyle\Delta U$ $\displaystyle=$ $\displaystyle
U_{r}-U_{l}=\sum_{i=1}^{4}s_{i}\overline{r}_{i},$
$\displaystyle\Delta\rho_{m}$ $\displaystyle=$
$\displaystyle\rho_{mr}-\rho_{ml}=s_{1}+s_{2},$ $\displaystyle\Delta X_{m}$
$\displaystyle=$ $\displaystyle
X_{mr}-X_{ml}=\left(\overline{v}_{m}+\overline{p}_{m}\right)s_{1}+\left(\overline{v}_{m}+\left(1+\gamma_{m}\right)\overline{p}_{m}\right)s_{2},$
$\displaystyle\Delta\rho_{c}$ $\displaystyle=$
$\displaystyle\rho_{cr}-\rho_{cl}=s_{3}+s_{4},$ $\displaystyle\Delta X_{c}$
$\displaystyle=$ $\displaystyle
X_{cr}-X_{cl}=\left(\overline{v}_{c}+\overline{p}_{c}\right)s_{3}+\left(\overline{v}_{c}+\left(1+\gamma_{c}\right)\overline{p}_{c}\right)s_{4}.$

The solutions $s_{i}$ are verified to be

$\displaystyle\begin{cases}s_{1}=\Delta\rho_{m}-\dfrac{\Delta
X_{m}-\left(\overline{v}_{m}+\overline{p}_{m}\right)\Delta\rho_{m}}{\gamma_{m}\overline{p}_{m}},\\\
s_{2}=\dfrac{\Delta
X_{m}-\left(\overline{v}_{m}+\overline{p}_{m}\right)\Delta\rho_{m}}{\gamma_{m}\overline{p}_{m}},\\\
s_{3}=\Delta\rho_{c}-\dfrac{\Delta
X_{c}-\left(\overline{v}_{c}+\overline{p}_{c}\right)\Delta\rho_{c}}{\gamma_{c}\overline{p}_{c}},\\\
s_{4}=\dfrac{\Delta
X_{c}-\left(\overline{v}_{c}+\overline{p}_{c}\right)\Delta\rho_{c}}{\gamma_{c}\overline{p}_{c}},\end{cases}$
(57)

with the obvious definitions
$\Delta\rho_{c}=\rho_{cr}-\rho_{cl},\leavevmode\nobreak\
\Delta\rho_{m}=\rho_{mr}-\rho_{ml},\leavevmode\nobreak\ \Delta
X_{c}=X_{cr}-X_{cl},\leavevmode\nobreak\ \Delta X_{m}=X_{mr}-X_{ml}$. The
expressions of $\Delta X_{m},\leavevmode\nobreak\ \Delta X_{c}$ are given by

$\displaystyle\Delta X_{m}$ $\displaystyle=$
$\displaystyle\Delta\left(\rho_{m}\left(v_{m}+p_{m}\right)\right)=\Delta\left(\rho_{m}v_{m}\right)+\Delta\left(\rho_{m}p_{m}\right),$
$\displaystyle\Delta X_{c}$ $\displaystyle=$
$\displaystyle\Delta\left(\rho_{c}\left(v_{c}+p_{c}\right)\right)=\Delta\left(\rho_{c}v_{c}\right)+\Delta\left(\rho_{c}p_{c}\right),$

where

$\displaystyle\Delta\left(\rho_{m}v_{m}\right)$ $\displaystyle=$
$\displaystyle\overline{\rho}_{m}\Delta
v_{m}+\overline{v}_{m}\Delta\rho_{m}+O_{1}(\Delta^{2}),$
$\displaystyle\Delta\left(\rho_{m}p_{m}\right)$ $\displaystyle=$
$\displaystyle\overline{\rho}_{m}\Delta
p_{m}+\overline{p}_{m}\Delta\rho_{m}+O_{2}(\Delta^{2}),$
$\displaystyle\Delta\left(\rho_{c}v_{c}\right)$ $\displaystyle=$
$\displaystyle\overline{\rho}_{c}\Delta
v_{c}+\overline{v}_{c}\Delta\rho_{c}+O_{3}(\Delta^{2}),$
$\displaystyle\Delta\left(\rho_{c}p_{c}\right)$ $\displaystyle=$
$\displaystyle\overline{\rho}_{c}\Delta
p_{c}+\overline{p}_{c}\Delta\rho_{c}+O_{4}(\Delta^{2}),$ $\displaystyle
O_{1}(\Delta^{2})$ $\displaystyle=$
$\displaystyle\left(\rho_{mr}-\rho_{ml}\right)\left(v_{mr}-\overline{v}_{m}\right)-\left(\rho_{ml}-\overline{\rho}_{m}\right)\left(v_{ml}-\overline{v}_{m}\right),$
$\displaystyle O_{2}(\Delta^{2})$ $\displaystyle=$
$\displaystyle\left(\rho_{mr}-\rho_{ml}\right)\left(p_{mr}-\overline{p}_{m}\right)-\left(\rho_{ml}-\overline{\rho}_{m}\right)\left(p_{ml}-\overline{p}_{m}\right),$
$\displaystyle O_{3}(\Delta^{2})$ $\displaystyle=$
$\displaystyle\left(\rho_{cr}-\rho_{cl}\right)\left(v_{cr}-\overline{v}_{c}\right)-\left(\rho_{cl}-\overline{\rho}_{c}\right)\left(v_{cl}-\overline{v}_{c}\right),$
$\displaystyle O_{4}(\Delta^{2})$ $\displaystyle=$
$\displaystyle\left(\rho_{cr}-\rho_{cl}\right)\left(p_{cr}-\overline{p}_{c}\right)-\left(\rho_{cl}-\overline{\rho}_{c}\right)\left(p_{cl}-\overline{p}_{c}\right).$

Finally, using (53) – (57) into any of the expressions (25) or

$\overline{F}=\dfrac{1}{2}\left(F_{l}+F_{r}\right)-\dfrac{1}{2}\sum_{i=1}^{4}s_{i}\mid\overline{\lambda}_{i}\mid\overline{r}_{i},$

leads to the sought Roe numerical flux.

####  4.1.2 Entropy fix

In this subsection, discontinuities that develop particularly at road section
boundaries are dissolved by applying an entropy condition to the Roe
decomposition. By doing this, a continuous numerical solution is obtained. An
entropy solution $e\mid\Gamma\mid e^{-1},$ is then substituted for the
Jacobian matrix (27), where
$\mid\Gamma\mid=\left[\overline{\lambda}_{1},\leavevmode\nobreak\
\overline{\lambda}_{2},\leavevmode\nobreak\ \ldots,\leavevmode\nobreak\
\overline{\lambda}_{K}\right]$ refers to a diagonal matrix made up of averaged
eigenvalues $\overline{\lambda}_{k},\leavevmode\nobreak\ k=1,2,\ldots,K$ of
(27), $e$ is a matrix made up of averaged eigenvectors and $e^{-1}$ an inverse
of $e$ [2]. Next, the eigenvalues are modified by employing Harten and Hyman
entropy fix scheme in order to accurately characterize the numerical solutions
[2]. Hence, the modified eigenvalues become

$\displaystyle\overline{\lambda}_{k}=\begin{cases}\overline{\delta}_{k},\leavevmode\nobreak\
\mbox{for}\leavevmode\nobreak\
\mid\lambda_{k}\mid\leq\overline{\delta}_{k},\\\
\mid\overline{\delta}_{k}\mid,\leavevmode\nobreak\
\mbox{for}\leavevmode\nobreak\
\mid\lambda_{k}\mid\geq\overline{\delta}_{k},\end{cases}$

where $\overline{\delta}_{k}=\max\left(0,\leavevmode\nobreak\
\lambda_{j+1/2}-\lambda_{j},\leavevmode\nobreak\
\lambda_{j+1}-\lambda_{j+1/2}\right).$ This ensures that the eigenvalues
$\overline{\delta}_{k}$ are non negative and have different values at the road
section boundaries. However, for abrupt changes at the road section
boundaries, $\overline{\delta}_{k}$ is set to zero. Using similar procedure as
in [2], the transformed Jacobian matrix is given by $e\mid\Gamma\mid e^{-1},$
where

$\displaystyle e$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{cccc}r_{1}&r_{2}&r_{3}&r_{4}\end{array}\right)\mid_{\overline{\rho}_{m},\overline{v}_{m}\overline{\rho}_{c},\overline{v}_{c}m},$
$\displaystyle\mid\Gamma\mid$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{cccc}\mid\overline{\lambda_{1}}(U)\mid&0&0&0\\\
0&\mid\overline{\lambda_{2}}(U)\mid&0&0\\\
0&0&\mid\overline{\lambda_{3}}(U)\mid&0\\\
0&0&0&\mid\overline{\lambda_{4}}(U)\mid\end{array}\right),$ $\displaystyle
e^{-1}$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{cccc}\dfrac{-(\overline{v}_{m}+\overline{p}_{m})}{\gamma_{m}\overline{p}_{m}}&\dfrac{1}{\gamma_{m}\overline{p}_{m}}&0&0\\\
\dfrac{\gamma_{m}\overline{p}_{m}+(\overline{v}_{m}+\overline{p}_{m})}{\gamma_{m}\overline{p}_{m}}&\dfrac{-1}{\gamma_{m}\overline{p}_{m}}&0&0\\\
0&0&\dfrac{-(\overline{v}_{c}+\overline{p}_{c})}{\gamma_{c}\overline{p}_{c}}&\dfrac{1}{\gamma_{c}\overline{p}_{c}}\\\
0&0&\dfrac{\gamma_{c}\overline{p}_{c}+(\overline{v}_{c}+\overline{p}_{c})}{\gamma\overline{p}_{c}}&\dfrac{-1}{\gamma\overline{p}_{c}}\end{array}\right),$
$\displaystyle\overline{p}_{m}$ $\displaystyle=$
$\displaystyle\psi^{\gamma}\left(\dfrac{\rho_{mr}+\rho_{mr}}{2}\right)^{\gamma},\leavevmode\nobreak\
\overline{p}_{c}=\psi^{\gamma}\left(\dfrac{\rho_{cl}+\rho_{cr}}{2}\right)^{\gamma},$

The sought Roe numerical flux is finally derived by evaluating
$F(U_{j}),\leavevmode\nobreak\ F(U_{j+1})$ and substituting $e\mid\Gamma\mid
e^{-1}$ for $B(U_{j+1/2})$ into (25). Hence, the solutions at subsequent time
steps $n+1,$ are obtained using the derived flux (25) into (24). Stability of
the discretization schem is established by satisfying the Courant, Friedrich
and Lewy (CFL) [19].

##  5 Stability analysis of model equations

In this Section, behavior of the proposed traffic flow model when subjected to
small disturbances is determined following the procedure outlined in [15]. To
determine the model’s stability, initial densities
$\rho_{m0},\leavevmode\nobreak\ \rho_{c0}$ at time $t_{0}$ and corresponding
equilibrium velocities, $v_{m0},\leavevmode\nobreak\
v_{c0}=v_{em}(\psi\rho_{c0}),\leavevmode\nobreak\ v_{ec}(\psi\rho_{c0})$ are
considered [2]. The initial data is then subjected to changes (or
perturbations) denoted $\sigma_{\rho_{m}}(x,t),\leavevmode\nobreak\
\sigma_{\rho_{c}}(x,t)$ and $\sigma_{vm}(x,t),\leavevmode\nobreak\
\sigma_{vc}(x,t),$ and given by

$\begin{cases}\sigma_{\rho m}(x,t)=\rho_{m}(x,t)-\rho_{m0},\\\
\sigma_{vm}(x,t)=v_{m}(x,t)-v_{m0},\\\ \sigma_{\rho
c}(x,t)=\rho_{c}(x,t)-\rho_{c0},\\\
\sigma_{vc}(x,t)=v_{c}(x,t)-v_{c0}.\end{cases}$

Conditions under which changes in density and velocity grows, and eventually
causes traffic congestion are determined by assuming sufficiently small
changes, smooth functions and periodic functions [15]. Thus, stability of the
proposed model is implied by stability of a linear combination of such
functions. Hence, the changes are characterized as

$\begin{cases}\sigma_{\rho m}(x,t)=\rho_{m0}e^{ikx+rt},\\\
\sigma_{vm}(x,t)=v_{m0}e^{ikx+rt},\\\ \sigma_{\rho
c}(x,t)=\rho_{c0}e^{ikx+rt},\\\ \sigma_{vc}(x,t)=v_{c0}e^{ikx+rt},\end{cases}$
(61)

where $i=\sqrt{-1},$ $r$ refers to frequency of oscillations, $k$ number of
changes caused by spatial evolution and $kx$ the spatial evolution [2]. It can
be inferred from (61) that traffic flow is a periodic function of $kx$ since
$e^{ikx}=\cos(kx)+i\sin(kx).$ At time $t,$ changes in densities and velocities
of motorcycles and cars are represented by
$\rho_{m0}e^{rt},\leavevmode\nobreak\ \rho_{c0}e^{rt},$ and
$v_{m0}e^{rt},\leavevmode\nobreak\ v_{c0}e^{rt},$ respectively and $rt$ is the
growth rate of the changes. Therefore, the changes with respect to space and
time become

$\begin{array}[]{cc}\dfrac{\partial\sigma_{\rho c}(x,t)}{\partial
x}=ik\rho_{c0}e^{ikx+rt},&\dfrac{\partial\sigma_{\rho m}(x,t)}{\partial
x}=ik\rho_{m0}e^{ikx+rt},\\\ \dfrac{\partial\sigma_{\rho c}(x,t)}{\partial
t}=r\rho_{c0}e^{ikx+rt},&\dfrac{\partial\sigma_{\rho m}(x,t)}{\partial
t}=r\rho_{m0}e^{ikx+rt},\\\ \dfrac{\partial\sigma_{vc}(x,t)}{\partial
x}=ikv_{c0}e^{ikx+rt},&\dfrac{\partial\sigma_{vm}(x,t)}{\partial
x}=ikv_{m0}e^{ikx+rt},\\\ \dfrac{\partial\sigma_{vc}(x,t)}{\partial
t}=rv_{c0}e^{ikx+rt},&\dfrac{\partial\sigma_{vm}(x,t)}{\partial
t}=rv_{m0}e^{ikx+rt}.\end{array}$

Suppose that a sufficiently small and smooth perturbation $(\sigma_{\rho
c},\sigma_{vc},\sigma_{\rho m},\sigma_{vm})$ is given to a steady-state
solution $\rho_{c0},v_{ec}(\rho_{c0}),\rho_{m0},v_{em}(\rho_{m0}).$ We then
examine the effect of the parameter $\psi,$ on the propagation of small
disturbance. Substituting a small change (61) into the proposed model, using
Taylor’s series expansion at $\left(\rho_{c0},v_{c0},\rho_{m0},v_{m0}\right),$
collecting linear terms where necessary and dropping all non-linear terms [15]
the following linearized system is obtained:

$\displaystyle\begin{cases}\dfrac{\partial\sigma_{\rho m}}{\partial
t}+v_{m0}\dfrac{\partial\sigma_{\rho m}}{\partial
x}+\rho_{m0}\dfrac{\partial\sigma_{vm}}{\partial x}=0,\\\
\dfrac{\partial\sigma_{vm}}{\partial
t}+v_{m0}\dfrac{\partial\sigma_{vm}}{\partial x}+\dfrac{\partial
p_{m}(\rho_{m0})}{\partial\rho_{m}}\dfrac{\partial\sigma_{\rho_{m}}}{\partial
t}+v_{m0}\dfrac{\partial
p_{m}(\rho_{m0})}{\partial\rho_{m}}\dfrac{\partial\sigma_{\rho m}}{\partial
x}-A_{1}=0,\\\ \dfrac{\partial\sigma_{\rho c}}{\partial
t}+v_{c0}\dfrac{\partial\sigma_{\rho c}}{\partial
x}+\rho_{c0}\dfrac{\partial\sigma_{vc}}{\partial x}=0,\\\
\dfrac{\partial\sigma_{vc}}{\partial
t}+v_{c0}\dfrac{\partial\sigma_{vc}}{\partial x}+\dfrac{\partial
p_{c}(\rho_{c0})}{\partial\rho_{c}}\dfrac{\partial\sigma_{\rho_{c}}}{\partial
t}+v_{c0}\dfrac{\partial
p_{c}(\rho_{c0})}{\partial\rho_{c}}\dfrac{\partial\sigma_{\rho c}}{\partial
x}-A_{2}=0,\end{cases}$ (62)

where $A_{1}=\dfrac{\psi_{m}v_{em}^{\prime}(\psi_{m}\rho_{m0})\sigma_{\rho
m}-\sigma_{vm}}{\tau_{c}},\leavevmode\nobreak\
A_{2}=\dfrac{\psi_{c}v_{ec}^{\prime}(\psi_{c}\rho_{c0})\sigma_{\rho
c}-\sigma_{vc}}{\tau_{c}}.$ Following [15], we assume that

$\dfrac{\sigma_{\rho c}(x,t)}{\rho_{c0}}\ll 1,\leavevmode\nobreak\
\dfrac{\sigma_{\rho m}(x,t)}{\rho_{m0}}\ll 1,\leavevmode\nobreak\
\dfrac{\sigma_{vc}(x,t)}{v_{c0}}\ll 1,\leavevmode\nobreak\
\dfrac{\sigma_{vm}(x,t)}{v_{m0}}\ll 1.$

Hence, system (62) simplifies to

$\displaystyle\begin{cases}r\sigma_{\rho m}+v_{m0}ik\sigma_{\rho
m}+\rho_{m0}ik\sigma_{vm}=0,\\\
r\sigma_{vm}+v_{m0}ik\sigma_{vm}+r\dfrac{\partial
p_{m}(\rho_{c0})}{\partial\rho_{m}}\sigma_{\rho m}+v_{m0}ik\dfrac{\partial
p_{m}(\rho_{m0})}{\partial\rho_{m}}\sigma_{\rho m}-A_{2}=0,\\\ r\sigma_{\rho
c}+v_{c0}ik\sigma_{\rho c}+\rho_{c0}ik\sigma_{vc}=0,\\\
r\sigma_{vc}+v_{c0}ik\sigma_{vc}+r\dfrac{\partial
p_{c}(\rho_{c0})}{\partial\rho_{c}}\sigma_{\rho c}+v_{c0}ik\dfrac{\partial
p_{c}(\rho_{c0})}{\partial\rho_{c}}\sigma_{\rho c}-A_{1}=0.\end{cases}$

The resulting Jacobian matrix, $J$ is given by

$\begin{bmatrix}r+v_{c0}ik&\rho_{c0}ik&0&0\\\
\left(r+v_{c0}ik\right)\dfrac{\partial
p_{c}(\rho_{c0})}{\partial\rho_{c}}-\dfrac{\psi_{c}v_{ec}^{\prime}(\psi_{c}\rho_{c0})}{\tau_{c}}&r+v_{c0}ik+\dfrac{1}{\tau_{c}}&0&0\\\
0&0&r+v_{m0}ik&\rho_{m0}ik\\\ 0&0&J_{43}&J_{44}\end{bmatrix},$

where $J_{43}=\left(r+v_{m0}ik\right)\dfrac{\partial
p_{m}(\rho_{c0})}{\partial\rho_{m}}-\dfrac{\psi_{m}v_{em}^{\prime}(\psi_{m}\rho_{m0})}{\tau_{m}},\leavevmode\nobreak\
J_{44}=r+v_{m0}ik+\dfrac{1}{\tau_{m}}.$ For non-trivial solutions, $det(J)=0.$
Thus, the solutions of $J(\sigma_{\rho c},\sigma_{vc},\sigma_{\rho
m},\sigma_{vm})=0$ are denoted $r_{1\pm},\leavevmode\nobreak\ r_{2\pm},$
respectively and given by

$\dfrac{\dfrac{\partial
p_{m}}{\partial\rho_{m}}ik\rho_{m0}\tau_{m}-2ik\tau_{m}v_{m0}-1\pm\sqrt{G}}{2\tau_{m}},\leavevmode\nobreak\
\dfrac{\dfrac{\partial
p_{c}}{\partial\rho_{c}}ik\rho_{c0}\tau_{c}-2ik\tau_{c}v_{c0}-1\pm\sqrt{G}}{2\tau_{c}},$

where

$G={\left(\dfrac{\partial
p_{m}}{\partial\rho_{m}}\right)^{2}}i^{2}k^{2}\rho_{m0}^{2}\tau_{m}^{2}-4\psi
ik\rho_{m0}\tau_{m}v_{em}^{\prime}-2\dfrac{\partial
p_{m}}{\partial\rho_{m}}ik\rho_{m0}\tau_{m}+1.$

For stability of the proposed model, changes in densities and velocities need
to decrease as time progresses. Hence, the real parts
$Re(r_{1\pm,\leavevmode\nobreak\ 2\pm}),$ of $det(J)$ must be strictly
negative. That is to say; $Re(r_{1\pm}),Re(r_{2\pm})<0.$ However, the square
root contains a complex number. So, it becomes difficult to determine the
sign. Hence, we apply the formula [20]

$\sqrt{R\pm i|I|}=\sqrt{\dfrac{1}{2}\left(\sqrt{R^{2}+I^{2}}+R\right)}\pm
i\sqrt{\dfrac{1}{2}\left(\sqrt{R^{2}+I^{2}}-R\right)},$ (63)

where

$\displaystyle R_{c}$ $\displaystyle=$ $\displaystyle 1-{\left(\dfrac{\partial
p_{c}(\rho_{c0})}{\partial\rho_{c}}\right)^{2}}k^{2}\rho_{c0}^{2}\tau_{c}^{2},\leavevmode\nobreak\
I_{c}=4\psi_{c}k\rho_{c0}\tau_{c}v_{ec}^{\prime}+2\dfrac{\partial
p_{c}}{\partial\rho_{c}}k\rho_{c0}\tau_{c},$ $\displaystyle R_{m}$
$\displaystyle=$ $\displaystyle 1-{\left(\dfrac{\partial
p_{m}(\rho_{m0})}{\partial\rho_{m}}\right)^{2}}k^{2}\rho_{m0}^{2}\tau_{m}^{2},\leavevmode\nobreak\
I_{m}=4\psi_{m}k\rho_{m0}\tau_{m}v_{em}^{\prime}+2\dfrac{\partial
p_{m}}{\partial\rho_{m}}k\rho_{m0}\tau_{m}.$

Let us consider $Re(r_{1+}).$ To have $Re(r_{1+})<0,$ would mean

$-1+\sqrt{\dfrac{1}{2}\left(\sqrt{R_{m}^{2}+I_{m}^{2}}+R_{m}\right)}<0.$

After simplifying and substituting for $R_{m}$ and $I_{m}$ of (63) into the
previous the stability condition is derived as

$\psi_{m}v_{em}^{\prime}<-\dfrac{\partial p_{m}}{\partial\rho_{m}}.$ (64)

In the same way, the stability condition can be derived for $Re(r_{2+})<0$ and
obtained as

$\psi_{c}v_{em}^{\prime}<-\dfrac{\partial p_{c}}{\partial\rho_{c}}.$ (65)

Data for the proposed model is chosen such that it satisfies the stability
conditions in (64) and (65).

##  6 Numerical results

The following subsections describe vehicles’ density, velocity and flow
distributions when motorcycles’ proportion is set to $20\%$ and $90\%.$
Considering (2)and (10), we look at freeway and congested traffic on
roundabout. Roe’s numerical scheme is utilized in which periodic boundary
conditions are applied in order to represent a circular road. In an attempt to
dissolve discontinuities at boundaries, an entropy fix is applied to the Roe’s
numerical scheme. Stability of the scheme is established by applying CFL
stability condition. Table 1 shows parameter values that are used for
simulations. Note that when $\delta$ takes on values $0$ and $1,$ results
become unstable (undefined). Figures 2 \- 6 display simulation results of the
proposed model.

Table 1: Simulation parameters. Description | Value | Source
---|---|---
Length of the road | 200m |
Road step, h | 5 m | [2]
Time step | 0.05 s |
Relaxation time,$\tau_{m},\leavevmode\nobreak\ \tau_{c}$ | 2s, 2.5s | [2]
Equilibrium velocity $v_{e}(\rho_{m}),v_{e}(\rho_{c})$ | Greenshields | [2]
Maximum normalized density | $\rho=\rho_{m}+\rho_{c}=1$ | [2]
Maximum speed, $v_{m\max}$ | 11m/s | [21]
Maximum speed, $v_{c\max}$ | 13.8m/s | [22]
Maximum area occupancy, $AO_{m}^{\max}$ | 0.85 | [23]
Maximum area occupancy, $AO_{c}^{\max}$ | 0.74 | [23]
$\gamma_{i},\leavevmode\nobreak\ i=m,c$ | 2.23, 2.12 | [23]
Width of the road | 12m | [22]
Width of a car | 1.6m | [24]
Vehicle class proportion, $\delta$ | 20%, 90% |
Length of a car | 4m | [24]
Length of a motorcycle | 1.8m | [24]
Simulation time, T | 60 seconds | [2]

###  6.1 Freeway traffic on a roundabout

Under this subsection, we consider a free traffic flow with motorcycles
proportion $\delta,$ initially fixed at $20\%,$ for
$0s,1s,\leavevmode\nobreak\ 20s,\leavevmode\nobreak\ 40s$ and $60s.$ The
initial total density of traffic flow denoted $\rho_{0},$ is set to be

$\displaystyle\rho_{0}$ $\displaystyle=$
$\displaystyle\begin{cases}0.1,&\text{for}\leavevmode\nobreak\ x<100\\\
0.2,&\text{for}\leavevmode\nobreak\ x\geq 100.\end{cases}$ (66)

Since $\delta=20\%$ it follows that $\rho_{m}=0.2\rho_{0}$ and
$\rho_{c}=0.8\rho_{0}.$ The initial density distributions of the respective
vehicle classes are plotted in Figure 2(a), and corresponding velocities
depicted by Figure 2(b). In Figure 2, a shock wave develops at $1s$ for both
vehicle classes at $100m$ (when the low density stream meets high density
stream). After $20s,$ the shock wave is smoothed and both vehicle classes move
with velocities close to their respective maximum values. That is; $13.8m/s$
and $11m/s$ for cars and motorcycles, respectively (see Figures 2(f) \- 2(j)).
Next, $\delta$ is set to $90\%,$ such that $\rho_{m}=0.9\rho_{0}$ and
$\rho_{c}=0.1\rho_{0},$ for motorcycles and cars, respectively. The plot of
such initial density profile is shown in Figure 3(a) with corresponding
velocities in Figure 3(b). Same results are observed when $\delta=90\%$ (see
Figures 3(f) \- 3(j)) Therefore, it can be concluded that there is minimal
effect of motorcycles proportion during freeway traffic flow.

(a) Density vs x at $0s$

(b) Velocity vs x at $0s$

(c) Density vs x at $1s$

(d) Velocity vs x at $1s$

(e) Density vs x at $20s$

(f) Velocity vs x at $20s$

(g) Density vs x at $40s$

(h) Velocity vs x, at $40s$

(i) Density vs x at $60s$

(j) Velocity vs x, at $60s$

Figure 2: Densities and velocities of the proposed model on a $200m$ freeway
circular road when $\delta=20\%,$ at $T=0s,1s,20s,40s$ and $60s.$

(a) Density vs x at $0s$

(b) Velocity vs x at $0s$

(c) Density vs x at $1s$

(d) Velocity vs x at $1s$

(e) Density vs x, at $20s$

(f) Velocity vs x at $20s$

(g) Density vs x at $40s$

(h) Velocity vs x at $40s$

(i) Density vs x at $60s$

(j) Velocity vs x at $60s$

Figure 3: Densities and velocities of the proposed model on a $200m$ freeway
circular road when $\delta=90\%,$ at $T=0s,1s,20s,40s$ and $60s.$

###  6.2 Congested traffic on a roundabout

Now we consider a congested $200m$ circular road with total initial density,
$\rho_{0}$ given by

$\displaystyle\rho_{0}$ $\displaystyle=$
$\displaystyle\begin{cases}0.3,&\text{if}\leavevmode\nobreak\ x\leq 130\\\
0.6,&\text{if}\leavevmode\nobreak\ 130<x<180\\\
0.1,&\text{if}\leavevmode\nobreak\ x\geq 180.\end{cases}$ (67)

The proposed model densities and velocities at time,
$T=0s,1s,\leavevmode\nobreak\ 20s,\leavevmode\nobreak\ 40s$ and $60s,$ with
$\delta$ initially fixed at $20\%$ are displayed in Figure 4. Plots of initial
densities and velocities are displayed by Figures 4(a) and 4(b), respectively.
It is observed in Figures 4(c) \- 4(j) that at all selected times, both
vehicle classes develop shock waves as high and low densities of vehicles
meet. However, the shock waves are smoothed at $20s$ and beyond, since an
entropy fix is applied to Roe numerical scheme. We further observe that
velocities for both vehicle classes decrease far away from their maximum
values as time progresses. Different from the freeway traffic flow scenario,
velocities of both vehicle classes are non uniform at all times since they are
congested. Cars form a rarefaction wave at $T=1s$ and $180m,$ but propagates
to the left at $20s$ and $40s$ and later disappears at $60s.$ A rarefaction
wave develops for motorcycles at $T=1s$ and $180m$ but disappears at $20s.$ A
shock wave develops for both vehicle classes at $T=1s$ and $130m,$ but for
cars it propagates to the left at $20s$ and $40s$ and later moves back to the
right at $60s.$ At $20s$ and $60s,$ motorcycles hit jam before cars. At $40s,$
cars hit jam faster than motorcycles. At any time, jam is created in one
region for cars while motorcycles density is maintained at low values.

(a) Density vs x at $0s$

(b) Velocity vs x at $0s$

(c) Density vs x at $1s$

(d) Velocity vs x at $1s$

(e) Density vs x at $20s$

(f) Velocity vs x at $20s$

(g) Density vs x at $40s$

(h) Velocity vs x, at $40s$

(i) Density vs x at $60s$

(j) Velocity vs x at $60s$

Figure 4: Densities and velocities of the proposed model on a $200m$ congested
circular road when $\delta=20\%,$ at time $T=0s,1s,20s,40s$ and $60s.$

Next $\delta$ is set to be $90\%,$ that is; fewer cars and many motorcycles.
The initial density and velocity distributions are indicated in Figures 5(a)
and 5(b). Density for cars is observed to be almost uniform after $1s$ while
for motorcycles is non uniform at all times. A shock wave is shown to clear
faster for cars than for motorcycles. A stationary shock wave but with
decreasing amplitude develops for motorcycles in one region after $1s,$ where
they find difficulties squeezing through spaces unfilled by the many cars.
Figures 6(a) \- 6(h) display the space-time evolution of densities and
velocities for both vehicle classes, when $\delta=20\%$ and $90\%.$ The
densities and velocities become smooth as time progresses.

(a) Density vs x at $0s$

(b) Velocity vs x at $0s$

(c) Density vs x at $1s$

(d) Velocity vs x at $1s$

(e) Density vs x at $20s$

(f) Velocity vs x at $20s$

(g) Density vs x at $40s$

(h) Velocity vs x, at $40s$

(i) Density vs x at $60s$

(j) Velocity vs x at $60s$

Figure 5: Densities and velocities of the proposed model on a $200m$ congested
circular road when $\delta=90\%,$ at $0s,1s,20s,40s$ and $60s.$

(a) $\rho_{m}$ vs x vs t, $\delta=20\%$

(b) $v_{m}$ vs x vs t, $\delta=20\%$

(c) $\rho_{c}$ vs x vs t, at $\delta=20\%$

(d) $v_{c}$ vs x vs t, $\delta=20\%$

(e) $\rho_{m}$ vs x vs t, $\delta=90\%$

(f) $v_{m}$ vs x vs t, $\delta=90\%$

(g) $\rho_{c}$ vs x vs t, $\delta=90\%$

(h) $v_{c}$ vs x vs t, $\delta=90\%$

Figure 6: Densities and velocities of the proposed model on a $200m$ congested
circular road when $\delta=20\%$ and $90\%,$ at $60s.$

##  7 Conclusion

The aim of this paper is to determine the effect of proportional densities on
the flow of traffic. A new expression for area occupancy, defined in terms of
proportion densities is proposed. Roe’s scheme is explicitly applied, where
Roe matrix, averaged velocity, averaged pressure, wave strengths are derived.
An entropy condition is applied to Roe decomposition. Stability conditions for
each vehicle class are obtained by linear stability analysis whereas stability
of the numerical scheme is determined by the CFL condition.

Qualitative analysis shows that Roe decomposition scheme can be applied to AR
model and other variants. Numerical results show that when motorcycles are
more than cars, all vehicle classes move with higher velocities than when
motorcycles are fewer. A shock wave is observed to clear faster for cars than
for motorcycles since motorcycles are slower than cars. On the other hand if
motorcycles are fewer than cars, both vehicle classes develop shock waves but
this effect is more felt during congested traffic flow. Traffic jam is shown
to occur much earlier if cars dominate the road than when motorcycles do. At
any time, whenever cars develop shock wave motorcycles too develop the shock.
Densities and velocities for each vehicle class become smoother over time and
were shown to remain within limits. Thus, the results obtained from the
proposed model are realistic. Total traffic flow increases with increase of
motorcycles proportion. As a control for traffic jam, enforcement of lane
discipline is recommended.

## Acknowledgment

The first author acknowledges her other research supervisors Prof. J.Y.T
Mugisha of Department of Mathematics, Makerere University, P.O Box 7062,
Kampala, Uganda and Prof. Semu Mitiku Kassa of Department of Mathematics and
Statistical Sciences, Botswana International University of Science and
Technology (BIUST), P/Bag $16,$ Palapye, Botswana for providing technical
guidance and mentorship. The first author also acknowledges the Sida bilateral
program with Makerere University, 2015-2020, project 316 “Capacity building in
Mathematics and its applications” for providing financial support.

## References

  * [1] Y. Li, D. Sun, Microscopic car-following model for the traffic flow: the state of the art, Journal of Control Theory and Applications 10 (2) (2012) 133–143.
  * [2] Z. H. Khan, W. Imran, S. Azeem, K. S. Khattak, T. A. Gulliver, M. S. Aslam, A macroscopic traffic model based on driver reaction and traffic stimuli, Applied Sciences 9 (14) (2019) 2848.
  * [3] D. P. Moller, Introduction to Transportation Analysis, Modeling and Simulation: Computational Foundations and Multimodal Applications, Springer, London, 2014\.
  * [4] S. Maerivoet, B. De Moor, Transportation planning and traffic flow models, arXiv preprint physics/0507127.
  * [5] V. Knoop, Introduction to traffic flow theory: An introduction with exercises.
  * [6] M. J. Lighthill, G. B. Whitham, 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.
  * [7] P. I. Richards, Shock waves on the highway, Operations research 4 (1) (1956) 42–51.
  * [8] S. Fan, D. B. Work, A heterogeneous multiclass traffic flow model with creeping, SIAM J. Appl. Math. 75 (2015) 813–835.
  * [9] R. Nair, H. S. Mahmassani, E. Miller-Hooks, A porous flow approach to modeling heterogeneous traffic in disordered systems, Procedia-Social and Behavioral Sciences 17 (2011) 611–627.
  * [10] S. Gashaw, P. Goatin, J. Härri, Modeling and analysis of mixed flow of cars and powered two wheelers, Transportation Research Part C: Emerging Technologies 89 (2018) 148–167.  doi:10.1016/j.trc.2018.02.004.
  * [11] J. Hawkins, L. Kattan, Microsimulating heterogeneous traffic: An implementation of the porous flow approach, Transportation in Developing Economies 4 (2018) 1–12.
  * [12] R. Mohan, G. Ramadurai, Heterogeneous traffic flow modelling using macroscopic continuum model, Procedia-Social and Behavioral Sciences 104 (2013) 402–411.
  * [13] H. J. Payne, Models of freeway traffic and control, in: Mathematical Models of Public Systems, ed. by G.A. Bekey (Simulation Council, La Jolla), Vol. 1, Simulation Council, Inc., 1971, pp. 51 – 61.
  * [14] H. J. Payne, FREFLO: A macroscopic simulation model of freeway traffic, Transportation Research Record 722 (1979) 68 – 77.
  * [15] A. K. Gupta, I. Dhiman, Analyses of a continuum traffic flow model for a nonlane-based system, International Journal of Modern Physics C 25 (10) (2014) 1450045.
  * [16] R. Mohan, G. Ramadurai, Heterogeneous traffic flow modelling using second-order macroscopic continuum model, Physics Letters A 381 (3) (2017) 115–123.  doi:10.1016/j.physleta.2016.10.042.
  * [17] E. F. Toro, Riemann solvers and numerical methods for fluid dynamics: a practical introduction, Springer Science & Business Media, German, 2013.
  * [18] P. Kachroo, S. J. Al-Nasur, S. A. Wadoo, A. Shende, Pedestrian dynamics: Feedback control of crowd evacuation, Springer Science & Business Media, 2008\.
  * [19] R. J. LeVeque, R. J. Leveque, Numerical methods for conservation laws, Vol. 214, Springer, 1992.
  * [20] D. Helbing, A. Johansson, On the controversy around daganzo’s requiem for and aw–rascle’s resurrection of second-order traffic flow models, in: Modelling and Optimisation of Flows on Networks, Springer, 2013, pp. 271–302.
  * [21] R. Mohan, G. Ramadurai, Heterogeneous traffic flow modelling using second-order macroscopic continuum model, Physics Letters A 381 (3) (2017) 115–123.
  * [22] Ministry of works and transport, Uganda road design manual 1 (2010) 190–204.
  * [23] R. Mohan, Multi-class AR model: comparison with microsimulation model for traffic flow variables at network level of interest and the two-dimensional formulation, International Journal of Modelling and Simulation 41 (2) (2021) 81–91.
  * [24] V. T. Arasan, G. Dhivya, Measuring heterogeneous traffic density, in: Proceedings of International Conference on Sustainable Urban Transport and Environment, World Academy of Science, Engineering and technology, Bangkok, Vol. 36, 2008, p. 342.