Micro- and macroscopic modeling of crowding and pushing in corridors

Michael Fischer; Gaspard Jankowiak; Marie-Therese Wolfram; Michael Fischer; Gaspard Jankowiak; Marie-Therese Wolfram

doi:10.3934/nhm.2020025

Networks and Heterogeneous Media

2020, Volume 15, Issue 3: 405-426. doi: 10.3934/nhm.2020025

Previous Article Next Article

Micro- and macroscopic modeling of crowding and pushing in corridors

1.
Radon Institute for Applied and Computational Mathematics, Altenbergerstr. 69, 4040, Linz, Austria
2.
University of Warwick, Mathematics Institute, Gibbet Hill Road, CV47AL Coventry, UK

Received: 01 November 2019 Revised: 01 July 2020 Published: 09 September 2020
91C20, 68Q80, 35Q84

Experiments with pedestrians revealed that the geometry of the domain, as well as the incentive of pedestrians to reach a target as fast as possible have a strong influence on the overall dynamics. In this paper, we propose and validate different mathematical models at the micro- and macroscopic levels to study the influence of both effects. We calibrate the models with experimental data and compare the results at the micro- as well as macroscopic levels. Our numerical simulations reproduce qualitative experimental features on both levels, and indicate how geometry and motivation level influence the observed pedestrian density. Furthermore, we discuss the dynamics of solutions for different modeling approaches and comment on the analysis of the respective equations.

Keywords:

Microscopic and macroscopic pedestrian models,
Fokker-Planck equation,
parameter estimation,
scalar conservation laws,
eikonal equation

Citation: Michael Fischer, Gaspard Jankowiak, Marie-Therese Wolfram. Micro- and macroscopic modeling of crowding and pushing in corridors[J]. Networks and Heterogeneous Media, 2020, 15(3): 405-426. doi: 10.3934/nhm.2020025

Related Papers:

[1]	Yannick Holle, Michael Herty, Michael Westdickenberg . New coupling conditions for isentropic flow on networks. Networks and Heterogeneous Media, 2020, 15(4): 605-631. doi: 10.3934/nhm.2020016
[2]	Alexandre M. Bayen, Alexander Keimer, Nils Müller . A proof of Kirchhoff's first law for hyperbolic conservation laws on networks. Networks and Heterogeneous Media, 2023, 18(4): 1799-1819. doi: 10.3934/nhm.2023078
[3]	Gabriella Bretti, Roberto Natalini, Benedetto Piccoli . Numerical approximations of a traffic flow model on networks. Networks and Heterogeneous Media, 2006, 1(1): 57-84. doi: 10.3934/nhm.2006.1.57
[4]	Yogiraj Mantri, Michael Herty, Sebastian Noelle . Well-balanced scheme for gas-flow in pipeline networks. Networks and Heterogeneous Media, 2019, 14(4): 659-676. doi: 10.3934/nhm.2019026
[5]	Steinar Evje, Kenneth H. Karlsen . Hyperbolic-elliptic models for well-reservoir flow. Networks and Heterogeneous Media, 2006, 1(4): 639-673. doi: 10.3934/nhm.2006.1.639
[6]	Jens Lang, Pascal Mindt . Entropy-preserving coupling conditions for one-dimensional Euler systems at junctions. Networks and Heterogeneous Media, 2018, 13(1): 177-190. doi: 10.3934/nhm.2018008
[7]	Michael Herty . Modeling, simulation and optimization of gas networks with compressors. Networks and Heterogeneous Media, 2007, 2(1): 81-97. doi: 10.3934/nhm.2007.2.81
[8]	Adriano Festa, Simone Göttlich, Marion Pfirsching . A model for a network of conveyor belts with discontinuous speed and capacity. Networks and Heterogeneous Media, 2019, 14(2): 389-410. doi: 10.3934/nhm.2019016
[9]	Raul Borsche, Axel Klar, T. N. Ha Pham . Nonlinear flux-limited models for chemotaxis on networks. Networks and Heterogeneous Media, 2017, 12(3): 381-401. doi: 10.3934/nhm.2017017
[10]	Raimund Bürger, Harold Deivi Contreras, Luis Miguel Villada . A Hilliges-Weidlich-type scheme for a one-dimensional scalar conservation law with nonlocal flux. Networks and Heterogeneous Media, 2023, 18(2): 664-693. doi: 10.3934/nhm.2023029

Abstract

1. Introduction

We consider networks modeled by a directed graph where the dynamics on each edge are described by one-dimensional conservation laws. The dynamics are coupled at the vertices of the graph, called junctions. We are especially interested in (isentropic) gas dynamics, but there are many other applications for example in traffic, supply chains, data networks or blood circulation. This field became of interest to many researchers in the last two decades, see for example the overview by Bressan et al. [8]. A main challenge is posed by prescribing suitable coupling conditions at the junction. We consider novel conditions for the system of isentropic gas, also referred to as p-system.

The isentropic gas equations at a junction with $k = 1, \dots, d\in\mathbb N$ adjacent pipelines are given by

$\begin{equation} \begin{cases} \partial_t \rho_k+\partial_x(\rho_k u_k)& = 0\\ \partial_t (\rho_k u_k)+\partial_x(\rho_k u_k^2+\kappa \rho_k^\gamma)& = 0 \end{cases}\quad\text{for a.e. } t > 0, x > 0, \end{equation}$

(1.1)

where $\rho_k\ge0$ denotes the gas density, $u_k\in\mathbb R$ the mean velocity, and $p = \kappa \rho^\gamma$ the pressure given by the $\gamma$ -pressure law with $\kappa>0$ and $1<\gamma<3$ . The equation is supplemented by an entropy condition, an initial condition $(\hat\rho_k, \hat u_k)(x) = (\rho_k, u_k)(0+, x)$ , $x>0$ and a suitable coupling condition on the traces $(\bar\rho_k, \bar u_k)(t) = (\rho_k, u_k)(t, 0+)$ , $t>0$ . Furthermore, we are interested in weak solutions to (1.1). For a general introduction to the theory of conservation laws see the book by Dafermos [17].

1.1. Previous results

The most challenging problem in modeling (gas) networks is to find physically correct coupling conditions. A first condition is usually conservation of mass at the junction

$\begin{equation} \sum\limits_{k = 1}^d A_k \bar\rho_k\bar u_k = 0, \quad\text{for a.e. }t > 0, \end{equation}$

(1.2)

where $A_k>0$ denotes the cross-sectional area of the $k$ -th pipeline. To ensure uniqueness of solutions, we impose more conditions at the junction. The number of additional conditions depends on the sign of the characteristic speeds at the junction. Several coupling conditions have been proposed, for example equality of pressure [1,2]

$\begin{equation} \kappa \bar\rho_k^\gamma = \mathcal H_p(t), \quad \text{for }k = 1, \dots, d, \quad\text{a.e. }t > 0, \end{equation}$

(1.3)

equality of momentum flux [10,11]

$\begin{equation} \bar\rho_k\bar u_k^2+\kappa \bar\rho_k^\gamma = \mathcal H_{\mathrm{MF}}(t), \quad \text{for }k = 1, \dots, d, \quad\text{a.e. }t > 0, \end{equation}$

(1.4)

and equality of stagnation enthalpy/Bernoulli invariant [28]

$\begin{equation} \frac{\bar u_k^2}{2}+\frac{\kappa}{\gamma-1} \bar\rho_k^{\gamma-1} = \mathcal H_{\mathrm{SE}}(t), \quad \text{for }k = 1, \dots, d, \quad\text{a.e. }t > 0. \end{equation}$

(1.5)

Notice that the first two conditions are non-physically in the sense that energy may be produced at the junction [28]. Equality of stagnation enthalpy implies conservation of energy at the junction. We derive a coupling condition with dissipated energy at the junction. This is consistent with the isentropic gas equations where energy is also dissipated. The special situation with only two pipelines were studied in [13,21]

Usually existence and uniqueness of solutions to a generalization of the Riemann problem at the junction are studied locally in state space first. Reigstad [27] introduced a method to study existence and uniqueness almost globally in the subsonic region under a technical assumption. The results for the Riemann problem are used, to construct approximate solutions to the generalized Cauchy problem, usually by wave front tracking. See the book by Bressan [7] for a general introduction to the wave front tracking method. Colombo, Herty and Sachers [12] proved a general existence and uniqueness theorem for the generalized Cauchy problem by using this method. This theorem requires a transversality condition, subsonic data and sufficiently small total variation of the initial data.

Another approach to prove existence of solutions is the method of compensated compactness. This has been applied to a scalar traffic model [9] and the isentropic gas equations [22]. The method requires less assumptions on the regularity of the initial data but is restricted to systems with a large class of entropies. Moreover, less regularity of the solutions is obtained and the traces at the junction have to be considered carefully, see e.g. [5,22].

1.2. A new coupling condition

To supplement conservation of mass (1.2), we use an approach based on the kinetic model for isentropic gas and a maximum energy/entropy dissipation principle at the junction.

A kinetic model for the isentropic gas equations was introduced by Lions, Perthame and Tadmor [25]. The corresponding vector-valued BGK model were introduced by Bouchut [6]. For $f = f(t, x, \xi)$ we impose

$\begin{equation} \partial_t f^k_\epsilon+\xi\partial_x f^k_\epsilon = \frac{M[f^k_\epsilon]-f^k_\epsilon}{\epsilon}, \quad \text{for a.e. }t > 0, x > 0, \xi\in\mathbb R, k = 1, \dots, d, \end{equation}$

(1.6)

with Maxwellian $M[f]$ (will be defined later). The half-space solutions are coupled at the junction by a kinetic coupling condition

$\begin{equation} \Psi^k[f_\epsilon^1(t, 0, \cdot), \dots, f_\epsilon^d(t, 0, \cdot)](\xi) = f_\epsilon^k(t, 0, \xi), \quad \text{for a.e. }t > 0, \xi > 0. \end{equation}$

(1.7)

To select the function $\Psi$ , we follow an idea of Dafermos [16] and use maximum entropy/energy dissipation as a selection criteria for the physically correct kinetic coupling condition. More precisely, we determine $\Psi^k$ such that as much energy is dissipated as possible under the condition of conserved mass. We obtain that for every pipeline the outgoing data is given by a Maxwellian with an artificial density $\rho_*^\epsilon$ and zero speed, i.e.

$\begin{equation} \Psi^k[f_\epsilon^1(t, 0, \cdot), \dots, f_\epsilon^d(t, 0, \cdot)](\xi) = M(\rho_*^\epsilon(t), 0, \xi), \quad \text{for a.e. }t > 0, \xi > 0. \end{equation}$

(1.8)

A formal limit argument leads to the definition of a generalized Riemann problem for $(\rho_k, u_k)$ . As in [20], each half-space solution is given by the restriction of the solution to a standard Riemann problem. The left Riemann initial state is given again by an artificial state with a suitable density and zero speed.

Definition 1. Let $(\hat\rho_k, \hat u_k)\in D$ , $k = 1, \dots, d$ . Then, we call $(\rho_k, u_k)\colon (0, \infty)_t\times(0, \infty)_x\to D$ a weak solution to the generalized Riemann problem if the following assertions hold true:

$\rm RP0:$ The solution satisfies the constant initial condition

$\begin{equation*} (\rho_k, u_k)(0+, x) = (\hat{\rho}_k, \hat{u}_k)\in D, \quad \text{for all }x > 0, \, k = 1, \dots, d. \end{equation*}$

$\rm RP1:$ There exists $\rho_*\ge 0$ such that $(\rho_k, u_k)$ is equal to the restriction to $x>0$ of the weak entropy solution in the sense of Lax with initial condition

$\begin{equation*} (\rho_k, u_k)(0+, x) = \begin{cases} (\hat{\rho}_k, \hat{u}_k), \quad& x > 0, \\ (\rho_*, 0), \quad& x < 0, \end{cases} \end{equation*}$

for all $k = 1, \dots, d$ .

$\rm RP2:$ Mass is conserved at the junction

$\begin{equation*} \sum\limits_{k = 1}^d A_k\bar{\rho}_k\bar{u}_k = 0, \end{equation*}$

where $(\bar{\rho}_k, \bar{u}_k) = (\rho_k, u_k)(t, 0+)\in D,$ for a.e. $t>0$ , $k = 1, \dots, d$ .

The set $D$ is the state space and will be defined in Section 2. Notice that the condition RP1 can be reformulated using a Riemann problem formulation for boundary conditions $\mathcal V(\rho_*, 0)$ . This formulation was used by Dubois and LeFloch [18] and will be defined later. It is illustrated in Figure 3. The set $\mathcal V(\rho_*, 0)$ will be used to define solutions to the generalized Cauchy problem (see Definition 6).

Figure 3. The boundary Riemann set

$\mathcal V(\rho_*, 0)$ .

DownLoad: Full-Size Img PowerPoint

Since the new coupling condition is based on restrictions of standard Riemann problems, we get a simple wave structure for the solutions in the sense of Definition 1. This structure allows us to prove existence and uniqueness of solutions globally in state space. We can use techniques by Reigstad [27] in the subsonic regime and extend them to the full state space. A general local existence and uniqueness result [12] for the Cauchy problem applies to the new condition. As a by-product we obtain Lipschitz continuous dependence on the initial data.

The coupling condition satisfies several properties regarding the energy/entropy dissipation. First, we obtain that entropy is non-increasing at the junction for a large class of symmetric entropies and in particular for the physical energy. A corollary of this property is a maximum principle on the Riemann invariants. More precisely, if the Riemann invariants of the initial data are bounded, then the Riemann invariants of the solution are bounded for all times. Furthermore, we prove a relation between the traces of the stagnation enthalpy at the junction.

We present an example in which the new coupling condition is the only condition leading to the physically correct wave types. The solutions to the generalized Riemann problem are computed numerically by using Newton's method. Furthermore, we study level sets associated to different coupling conditions and consider their geometry.

Our approach can be easily generalized to other hyperbolic systems. We extend it to full gas dynamics and obtain a similar coupling condition with an artificial density, zero speed and an artificial temperature at the junction. For more details and a brief literature overview see Section 8.

1.3. Organization of the paper

In Section 2, we recall several properties of the isentropic gas equations and the initial boundary value problem. In Section 3, we give a detailed motivation and a formal derivation of the coupling condition. Existence and uniqueness of solutions to the generalized Riemann problem will be proven in Section 4 and the corresponding results for the generalized Cauchy problem will be proven in Section 5. In Section 6, we derive several physical properties of the coupling condition, e.g. non-increasing energy, a maximum principle on the Riemann invariants and a relation for the traces of the stagnation enthalpy. Section 7 is devoted to numerical considerations. In Section 8, the extension of our approach to full gas dynamics is given. In Section 9, we finish with a conclusion.

2. The isentropic gas equations and basic definitions

2.1. Entropy solutions, Riemann invariants and Lax curves

The isentropic gas equations in one space dimension are given by

$\begin{equation} \begin{cases} \partial_t \rho+\partial_x(\rho u)& = 0\\ \partial_t (\rho u)+\partial_x(\rho u^2+\kappa \rho^\gamma)& = 0 \end{cases}\quad \text{a.e. }t\ge 0, x\in \mathbb R. \end{equation}$

(2.1)

Furthermore, we impose the entropy condition

$\begin{equation} \partial_t\eta_S(\rho, u)+\partial_x G_S(\rho, u)\le0\quad\text{a.e. }t, x, \end{equation}$

(2.2)

for all (weak) entropy pairs $(\eta_S, G_S)$ , where $\eta_S$ is a convex function with

$\begin{equation} G_S' = \eta_S' F', \quad \eta_S(\rho = 0, u) = 0, \quad \eta_{S, \rho}'(\rho = 0, u) = S(u), \quad\text{ for all }u\in\mathbb R, \end{equation}$

(2.3)

and a suitable convex function $S\colon\mathbb R\to\mathbb R$ . The involved derivatives are taken with respect to the conserved quantities $(\rho, \rho u)$ . We recall some basic definitions and notation:

$\begin{gather} D = \{(x_0, x_1)\in \mathbb R^2|x_0 > 0 \text{ or }x_0 = x_1 = 0\}, \end{gather}$

(2.4)

$\begin{gather} \chi(\rho, \xi) = c_{\gamma, \kappa}(a_\gamma^2 \rho^{\gamma-1}-\xi^2)_+^\lambda, \end{gather}$

(2.5)

$\begin{gather} \theta = \frac{\gamma-1}{2}, \quad\lambda = \frac{1}{\gamma-1}-\frac{1}{2}, \quad c_{\gamma, \kappa} = \frac{a_\gamma^{2/(\gamma-1)}}{J_\lambda}, \end{gather}$

(2.6)

$\begin{gather} J_\lambda = \int_{-1}^1 (1-z^2)^\lambda\, \mathrm{d}z, \quad a_\gamma = \frac{2\sqrt{\gamma \kappa}}{\gamma-1}. \end{gather}$

(2.7)

The isentropic gas equations admit the Riemann invariants

$\begin{align} \omega_1& = u-a_\gamma\rho^\theta, \quad\omega_2 = u+a_\gamma\rho^\theta, \end{align}$

(2.8)

for $(\rho, u)\in D$ . The eigenvalues are given by

$\begin{align} \lambda_1(\rho, u) = u-\sqrt{\kappa\gamma}\rho^\theta, \quad\lambda_2(\rho, u) = u+\sqrt{\kappa\gamma}\rho^\theta, \end{align}$

(2.9)

and the eigenvectors by

$\begin{align} r_1(\rho, u) = \begin{pmatrix}1 \\u-\sqrt{\kappa\gamma}\rho^\theta \end{pmatrix}, \quad r_2(\rho, u) = \begin{pmatrix}1 \\u+\sqrt{\kappa\gamma}\rho^\theta \end{pmatrix}, \end{align}$

(2.10)

for $(\rho, u)\in D$ . We call a state $(\rho, u)\in D$

● subsonic if $\lambda_1(\rho, u)<0<\lambda_2(\rho, u)$ ;

● sonic if $\lambda_1(\rho, u) = 0$ or $\lambda_2(\rho, u) = 0$ ;

● supersonic if $0<\lambda_1(\rho, u)<\lambda_2(\rho, u)$ or $\lambda_1(\rho, u)<\lambda_2(\rho, u)<0$ .

Next, we define several quantities corresponding to the kinetic (BGK) model for isentropic gas dynamics (1.6). The vector-valued Maxwellian $M[f]$ for $f\colon \mathbb R\to D$ is defined by

$\begin{equation} M[f](\xi) = M(\rho_f, u_f, \xi) \end{equation}$

(2.11)

where

$\begin{equation} \rho_f = \int_{\mathbb R}f_0(\xi)\, \mathrm{d}\xi, \quad\rho_f u_f = \int_{\mathbb R}f_1(\xi)\, \mathrm{d}\xi \end{equation}$

(2.12)

and

$\begin{gather} M(\rho, u, \xi) = (\chi(\rho, \xi-u), ((1-\theta)u+\theta \xi)\chi(\rho, \xi-u)). \end{gather}$

(2.13)

The kinetic entropies are defined by

$\begin{align} H_S(f, \xi)& = \int_\mathbb R \Phi(\rho(f, \xi), u(f, \xi), \xi, v)S(v)\;\text{d}v\quad \text{for }f\in D\backslash \{ 0\}, \end{align}$

(2.14)

$\begin{align} H_S(0, \xi)& = 0, \end{align}$

(2.15)

where

$\begin{align} u(f, \xi)& = \frac{f_1/f_0-\theta \xi}{1-\theta}, \end{align}$

(2.16)

$\begin{align} \rho(f, \xi)& = a_\gamma^{-\tfrac{2}{\gamma-1}}\left(\left(\frac{f_1/f_0-\xi}{1-\theta}\right)^2+\left(\frac{f_0}{c_{\gamma, \kappa}}\right)^{1/\lambda}\right)^{\tfrac{1}{\gamma-1}}. \end{align}$

(2.17)

The kernel $\Phi$ is defined by

$\begin{gather} \Phi(\rho, u, \xi, v) = \frac{(1-\theta)^2}{\theta}\frac{c_{\gamma, \kappa}}{J_\lambda}\mathbb{1}_{\omega_1 < \xi < \omega_2}\mathbb{1}_{\omega_1 < v < \omega_2}|\xi-v|^{2\lambda-1}\Upsilon_{\lambda-1}(z), \end{gather}$

(2.18)

$\begin{gather} z = \frac{(\xi+v)(\omega_1+\omega_2)-2(\omega_1\omega_2+\xi v)}{(\omega_2-\omega_1)|\xi-v|}, \end{gather}$

(2.19)

$\begin{gather} \Upsilon_{\lambda-1}(z) = \int_1^z (y^2-1)^{\lambda -1}\;\text{d}y, \quad z\ge 1. \end{gather}$

(2.20)

The kinetic Riemann invariants are given by

$\begin{align} \omega_1& = u(f, \xi)-a_\gamma\rho(f, \xi)^\theta, \quad\omega_2 = u(f, \xi)+a_\gamma\rho(f, \xi)^\theta, \end{align}$

(2.21)

for $f\in D\backslash \{ 0\}$ . The macroscopic entropy and entropy flux are given by

$\begin{align} \eta_S(\rho, u)& = \int_{\mathbb R}\chi(\rho, v-u) S(v)\;\text{d}v = \int_{\mathbb R} H_S(M(\rho, u, \xi), \xi)\;\text{d}\xi, \end{align}$

(2.22)

$\begin{align} G_S(\rho, u)& = \int_\mathbb R [(1-\theta)u+\theta v]\ \chi(\rho, v-u)S(v)\;\text{d}v \end{align}$

(2.23)

$\begin{align} & = \int_\mathbb R\xi H_S(M(\rho, u, \xi), \xi)\;\text{d}\xi, \end{align}$

(2.24)

for $(\rho, u)\in D$ . If additionally $S\in C^1(\mathbb R, \mathbb R)$ , then the gradient of $\eta$ with respect to the conserved variables is given by

$\begin{equation} \eta_S'(\rho, u) = \frac{1}{J_\lambda}\int_{-1}^1 (1-z^2)^\lambda \begin{pmatrix} S(u+a_\gamma\rho^\theta z)+(\theta a_\gamma\rho^\theta z- u)S'(u+a_\gamma\rho^\theta z)\\ S'(u+a_\gamma\rho^\theta z)\end{pmatrix}\mathrm{d}z, \end{equation}$

(2.25)

for $(\rho, u)\in D$ . The kinetic entropy parametrized by $S(v) = v^2/2$ is given by

$\begin{equation} H(f, \xi) = \frac{\theta}{1-\theta}\frac{\xi^2}{2} f_0+\frac{\theta}{2c_{\gamma, \kappa}^{1/\lambda}}\frac{f_0^{1+1/\lambda}}{1+1/\lambda} +\frac{1}{1-\theta}\frac{1}{2}\frac{f_1^2}{f_0}-\frac{\theta}{1-\theta}\xi f_1, \end{equation}$

(2.26)

and the corresponding macroscopic entropy pair is given by the physical energy and energy flux

$\begin{align} \eta(\rho, u) = \frac{\rho u^2}{2}+\frac{\kappa}{\gamma-1}\rho^\gamma, \quad G(\rho, u) = \frac{\rho u^3}{2}+\frac{\gamma\kappa}{\gamma-1}\rho^\gamma u. \end{align}$

(2.27)

To construct solutions to the generalized Riemann problem, we need the (forward) Lax wave curves $\mathcal W_1(\rho_0, u_0)$ and $\mathcal W_2(\rho_0, u_0)$ which are the composition of the corresponding rarefaction and shock curves. The rarefaction curves are given by

$\begin{align} \mathcal R_1(\rho_0, u_0):\quad u = u_0+a_\gamma\rho_0^\theta-a_\gamma\rho^\theta, \quad \text{for }\rho < \rho_0, \end{align}$

(2.28)

$\begin{align} \mathcal R_2(\rho_0, u_0):\quad u = u_0-a_\gamma\rho_0^\theta+a_\gamma\rho^\theta, \quad \text{for }\rho > \rho_0, \end{align}$

(2.29)

and the shock curves are given by

$\begin{align} \mathcal S_1(\rho_0, u_0):\quad u = u_0-\sqrt{\frac{\kappa(\rho^\gamma-\rho_0^\gamma)(\rho-\rho_0)}{\rho\rho_0}}, \quad \text{for }\rho > \rho_0, \end{align}$

(2.30)

$\begin{align} \mathcal S_2(\rho_0, u_0):\quad u = u_0-\sqrt{\frac{\kappa(\rho^\gamma-\rho_0^\gamma)(\rho-\rho_0)}{\rho\rho_0}}, \quad \text{for }\rho < \rho_0. \end{align}$

(2.31)

We will use the notation $\mathcal S_2^-(\rho, u), \mathcal R_2^-(\rho, u), \dots$ for the reversed wave curves. They satisfy the same condition as the forward wave curves but the fixed variables are $(\rho, u)$ instead of $(\rho_0, u_0)$ . We always consider the self-similar solutions to Riemann problems in the sense of Lax and denote them by $\mathcal{RP}(\rho_l, u_l, \rho_r, u_r)(t/x)$ for initial data

$\begin{equation} \begin{cases} (\rho_l, u_l), \quad &x < 0, \\ (\rho_r, u_r), \quad &x > 0, \end{cases} \end{equation}$

(2.32)

where $(\rho_l, u_l), (\rho_r, u_r)\in D$ .

2.2. The initial boundary value problem

In this subsection, we recall some basic properties of the initial boundary value problem. The sets $\mathcal E(\rho_b, u_b)$ and $\mathcal V(\rho_b, u_b)$ of admissible boundary values were introduced in [18]. We recall their definitions.

Definition 2. Let $(\rho_b, u_b)\in D$ . $\mathcal V(\rho_b, u_b)$ is the set of states $(\rho, u)\in D$ with

$\begin{equation*} (\rho, u) = \mathcal{RP}(\rho_b, u_b, \rho_r, u_r)(0+), \quad\text{for a state }(\rho_r, u_r)\in D. \end{equation*}$

Definition 3. Let $(\rho_b, u_b)\in D$ . $\mathcal E(\rho_b, u_b)$ is the set of states $(\rho, u)\in D$ with

$\begin{equation*} G_S(\rho, u)-G_S(\rho_b, u_b)-\eta_S'(\rho_b, u_b)\left(F(\rho, u)-F(\rho_b, u_b)\right)\le 0 \end{equation*}$

for all entropy pairs $(\eta_S, G_S)$ of class $C^1$ (i.e. $S\in C^1(\mathbb R, \mathbb R$ )).

We recall the following result.

Proposition 1. ([23,Theorem 3.4]). Let $(\rho_b, u_b)\in D$ , then $\mathcal V(\rho_b, u_b)\subset \mathcal E(\rho_b, u_b)$ . The reversed set inclusion does not hold true in general.

Next, we define subsets of $D$ to consider different situations in the initial boundary value problem. We are especially interested in the case $(\rho_b, u_b) = (\rho_*, 0)$ .

Definition 4. Let $(\rho_b, u_b)\in D$ . Then,

● $\mathcal A$ is the set of states which are connected to $\mathcal W_1(\rho_b, u_b)\cap\{\lambda_1\ge 0\}$ by its reversed 2-wave curve;

● $\mathcal B$ is the set of states which are connected to $\mathcal W_1(\rho_b, u_b)\cap\{\lambda_1<0<\lambda_2\}$ by its reversed 2-wave curve with positive wave speed;

● $\mathcal C$ is the set of states which are connected to $\mathcal W_1(\rho_b, u_b)\cap\{\lambda_2\le 0\}$ by its reversed 2-wave curve or are connected to $\mathcal W_1(\rho_b, u_b)\cap\{\lambda_1<0<\lambda_2\}$ by its reversed 2-wave curve with non-positive wave speed;

● $\mathcal J$ is the set of states which are connected to $\mathcal W_1(\rho_b, u_b)\cap\{\lambda_1<0<\lambda_2\}$ by its reversed 2-wave curve with zero wave speed.

We write $\mathcal A_*, \, \mathcal B_*, \, \mathcal C_*, \, \mathcal J_*$ if $(\rho_b, u_b) = (\tilde\rho_*, 0)$ with $\tilde\rho_*\ge 0$ .

The sets $\mathcal A_*, \, \mathcal B_*, \, \mathcal C_*, \, \mathcal J_*$ are shown in Figure 2. The sets $\mathcal A_*$ and $\mathcal B_*$ are separated by the 2-wave curve $\mathcal W_2(\rho_\alpha, u_\alpha)$ , where $(\rho_\alpha, u_\alpha)$ is the unique state in $\{\lambda_1 = 0\}\cap \mathcal R_1(\tilde\rho_*, 0)$ . The sets $\mathcal B_*$ and $\mathcal C_*$ are separated by $\mathcal J_*$ and $\mathcal R_2(\rho_\beta, u_\beta)$ , where $(\rho_\beta, u_\beta)$ is the unique state in $\{\lambda_2 = 0\}\cap \mathcal S_1(\tilde\rho_*, 0)$ . For the construction of $(\rho_\alpha, u_\alpha)$ and $(\rho_\beta, u_\beta)$ , we refer to Figure 1.

Figure 1. Construction of

$(\rho_\alpha, u_\alpha)$ and

$(\rho_\beta, u_\beta)$ .

DownLoad: Full-Size Img PowerPoint

Figure 2. The sets

$\mathcal A_*$ ,

$\mathcal B_*$ ,

$\mathcal C_*$ and

$\mathcal J_*$ .

DownLoad: Full-Size Img PowerPoint

Next, we construct a solution which satisfies all properties in Definition 1 except of conservation of mass RP2. They will be used to construct the desired solution to Definition 1 later. To clarify that RP2 does not necessarily hold true, we denote the artificial density by $\tilde\rho_*$ .

Lemma 1. Let $(\hat \rho_k, \hat u_k)\in D$ and $\tilde\rho_*\ge 0$ . Then, there exists a unique function $(\rho_k, u_k)\colon(0, \infty)_t\times(0, \infty)_x\to D$ which coincides with the self-similar Lax solution to the standard Riemann problem with initial condition

$\begin{equation*} (\rho_k, u_k)(0, x) = \begin{cases} (\hat\rho_k, \hat u_k), \quad &x > 0, \\ (\tilde\rho_*, 0), \quad &x < 0, \end{cases} \end{equation*}$

for a.e. $t>0$ , $x>0$ . Furthermore, we have the following properties for the trace $(\bar\rho_k, \bar u_k) = (\rho_k, u_k)(t, 0+)$ illustrated in Figure 3:

(i) $(\hat\rho_k, \hat u_k) = (\bar\rho_k, \bar u_k)$ if and only if $(\hat\rho_k, \hat u_k)\in \mathcal{V}(\tilde\rho_*, 0)$ .

(ii) $(\bar\rho_k, \bar u_k)$ cannot be supersonic with $\lambda_1(\bar\rho_k, \bar u_k)> 0$ .

(iii) $(\bar\rho_k, \bar u_k)$ is sonic with $\lambda_1(\bar\rho_k, \bar u_k) = 0$ if and only if $(\hat\rho_k, \hat u_k)\in\mathcal A$ . Furthermore, $(\bar\rho_k, \bar u_k)$ is the unique element in $\{\lambda_1 = 0\}\cap \mathcal R_1(\tilde\rho_*, 0)$ .

(iv) $(\bar\rho_k, \bar u_k)$ is subsonic if and only if $(\hat\rho_k, \hat u_k)\in\mathcal B$ .

(v) $(\bar\rho_k, \bar u_k)$ is sonic with $\lambda_2(\bar\rho_k, \bar u_k) = 0$ if and only if $(\bar\rho_k, \bar u_k)$ is connected to $(\hat\rho_k, \hat u_k)$ by a 2-rarefaction curve and $(\hat\rho_k, \hat u_k)\in\mathcal C\backslash \{\lambda_2\ge 0\}$ .

(vi) $(\bar\rho_k, \bar u_k)$ is supersonic with $\lambda_2(\bar\rho_k, \bar u_k)< 0$ if and only if $(\hat\rho_k, \hat u_k)\in\mathcal C\cap\{\lambda_2<0\}$ .

Proof. The existence and uniqueness of self-similar Lax solutions to Riemann problems is well-known. For Riemann problems with vacuum initial data see [26]. It remains to prove the properties for $(\bar\rho_k, \bar u_k)$ . They follow from the considerations in [18].

3. Motivation and formal derivation

In this section we give a physical motivation and formal derivation for the new coupling condition. Both are based on the kinetic model for isentropic gas and a maximum energy/entropy dissipation principle. First, we specify the kinetic coupling condition which conserves mass and dissipates as much energy as possible. In the second step, we consider the macroscopic limit of the kinetic coupling condition. This relaxation works only on a formal level since the currently available results for passing to the limit at the junction are not strong enough. Nevertheless, we are able to take the formal limit towards the macroscopic coupling condition. Finally, we also give an interpretation of the resulting conditions.

3.1. A maximum energy dissipation principle applied to kinetic coupling conditions

Since Dafermos [16] introduced the entropy rate admissibility criterion it is a natural approach to maximize the entropy dissipation in the field of hyperbolic conservation laws. This technique can be used to single out the physically correct solutions. We adapt this approach and aim to find the most dissipative kinetic coupling condition (with the constrained of conservation of mass). Since the physical energy is an entropy for the system of isentropic gas and is the physically relevant entropy, we maximize the energy dissipation.

We consider the kinetic BGK model of the isentropic gas equations which is given by

$\begin{equation} \partial_t f_\epsilon+ \xi\partial_x f_\epsilon = \frac{M[f_\epsilon]-f_\epsilon}{\epsilon}, \quad \text{for a.e. } t > 0, x\in\mathbb R, \xi\in\mathbb R, \end{equation}$

(3.1)

where $f_\epsilon = f_\epsilon(t, x, \xi)\in D$ .

Remark 1. The BGK model, its relaxation limit and boundary conditions were studied by Berthelin and Bouchut [3,4,5]. These results were extended to networks in [22]. Notice that our considerations are independent of the right hand side of the kinetic equation as long as the kinetic solution converges to an entropy solution of the macroscopic equation.

To couple the half-space solutions, we have to define a kinetic coupling condition

$\begin{equation} \Psi\colon L^1_\mu((-\infty, 0)_\xi, D)^d\to L^1_\mu((0, \infty)_\xi, D)^d;\quad g\mapsto\Psi[g]. \end{equation}$

(3.2)

The half-space solutions $f_\epsilon^k\colon(0, \infty)_t\times(0, \infty)_x\times\mathbb R_\xi\to D$ to the kinetic model are coupled by

$\begin{equation} \Psi^k[f_\epsilon^1(t, 0, \cdot), \dots, f_\epsilon^d(t, 0, \cdot)](\xi) = f_\epsilon^k(t, 0, \xi), \quad\text{for a.e. }t > 0, \xi > 0. \end{equation}$

(3.3)

We are interested in kinetic coupling conditions which conserve mass. More precisely, we require that

$\begin{align} \sum\limits_{k = 1}^d A_k \left(\int_0^\infty \xi \Psi^k_0[g](\xi)\, \mathrm{d}\xi+\int_{-\infty}^0 \xi g^k_0(\xi)\, \mathrm{d}\xi\right) = 0, \end{align}$

(3.4)

holds for all $g\in L^1_\mu((-\infty, 0)_\xi, D)^d$ . Our aim is to find the coupling condition which dissipates as much energy as possible and conserves mass. The energy dissipation at the junction is given by

$\begin{equation} \sum\limits_{k = 1}^d A_k \left(\int_0^\infty \xi H(\Psi^k[g](\xi), \xi)\, \mathrm{d}\xi+\int_{-\infty}^0 \xi H(g^k(\xi), \xi)\, \mathrm{d}\xi\right). \end{equation}$

(3.5)

To find the unique minimizer of this functional, we use the convexity of the kinetic energy. More precisely, we use the sub-differential inequality (see e.g. [4])

$\begin{equation} H(g, \xi)\ge H(M(\rho, u, \xi), \xi)+\eta'(\rho, u)(g-M(\rho, u, \xi)), \end{equation}$

(3.6)

for every $g\in D$ , $(\rho, u)\in D$ , $\xi\in \mathbb R$ . Applying the sub-differential inequality leads to

$\begin{align} &\sum\limits_{k = 1}^d A_k \int_0^\infty \xi H(\Psi^k[g](\xi), \xi)\, \mathrm{d}\xi\\ \ge& \sum\limits_{k = 1}^d A_k \int_0^\infty \xi H(M(\rho_*, 0, \xi), \xi)\, \mathrm{d}\xi\\ &+ \eta'(\rho_*, 0)\left[ \sum\limits_{k = 1}^d A_k \left(\int_0^\infty \xi \Psi^k_0[g](\xi)\, \mathrm{d}\xi -\int_0^\infty \xi M(\rho_*, 0, \xi)\, \mathrm{d}\xi\right)\right]\\ = &\sum\limits_{k = 1}^d A_k \int_0^\infty \xi H(M(\rho_*, 0, \xi), \xi)\, \mathrm{d}\xi, \end{align}$

(3.7)

where $\rho_*\ge 0$ is uniquely defined by

$\begin{equation} \sum\limits_{k = 1}^d A_k \left(\int_0^\infty \xi M_0(\rho_*, 0, \xi)\, \mathrm{d}\xi +\int_{-\infty}^0 \xi g^k_0(\xi)\, \mathrm{d}\xi\right) = 0. \end{equation}$

(3.8)

Such an $\rho_*$ always exists since

$\begin{equation} \int_0^\infty \xi M_0(\rho_*, 0, \xi)\, \mathrm{d}\xi = c_{\gamma, \kappa}\left(a_\gamma \rho_*^\theta\right)^{\frac{\gamma+1}{\gamma-1}}\int_0^1 z\, (1-z^2)^\lambda\, \mathrm{d}z. \end{equation}$

(3.9)

The uniqueness follows from the strict convexity of $H$ . We obtain the following result.

Theorem 1. The unique kinetic coupling condition

$\begin{equation*} \Psi\colon L^1_\mu((-\infty, 0)_\xi, D)^d\to L^1_\mu((0, \infty)_\xi, D)^d;\quad g\mapsto\Psi[g] \end{equation*}$

which conserves mass (3.4) and minimizes (3.5) is given by

$\begin{equation} \Psi^k[g](\xi): = M(\rho_*, 0, \xi), \end{equation}$

(3.10)

where $\rho_*\ge 0$ is defined by equation (3.8).

According to [22], weak solutions to the kinetic BGK model on networks exist under suitable conditions on the initial data.

Instead of minimizing the kinetic energy, other kinetic entropies could be considered. Notice that every kinetic entropy $H_S$ parametrized by a strictly convex, symmetric function $S\in C^1(\mathbb R, \mathbb R)$ leads to the kinetic coupling condition obtained in (3.10). This can be proven by applying the sub-differential inequality for $H_S$ and using $G_S(\rho_*, 0) = \partial_{\rho u}\eta_{S}(\rho_*, 0) = 0$ .

3.2. Formal relaxation limit

We take the formal limit at the junction with $\Psi^k[g](\xi) = M(\rho_*, 0, \xi)$ . Since (3.6) holds for every convex kinetic entropy $H_S$ , we get

$\begin{align} &\int_{\mathbb R}\xi H_S(f^k_\epsilon(t, 0, \xi), \xi)\, \mathrm{d}\xi-G_S(\rho_*^\epsilon(t), 0)\\ &\quad-\eta_S'(\rho_*^\epsilon(t), 0)\left(\int_{\mathbb R}\xi f^k_\epsilon(t, 0, \xi)\, \mathrm{d}\xi-F(\rho_*^\epsilon(t), 0)\right)\le 0, \end{align}$

(3.11)

for $t>0$ , $k = 1, \dots, d$ , where $\rho_*^\epsilon(t)\ge 0$ is defined by (3.8) with $g^k(\xi) = f^k_\epsilon(t, 0, \xi)$ and fixed $\epsilon>0$ . Assuming that $\rho^\epsilon_*\to\rho_*$ converges strongly in $L^1_{\mathrm loc}$ as $\epsilon\to 0$ and using the arguments in [22], we get

$\begin{equation} \overline{G_S(\rho_k, u_k)}(t)-G_S(\rho_*(t), 0)-\eta_S'(\rho_*(t), 0)(\overline{F(\rho_k, u_k)}(t)-F(\rho_*(t), 0))\le 0, \end{equation}$

(3.12)

for a.e. $t>0$ , $k = 1, \dots, d$ and every convex $S\in C^1(\mathbb R, \mathbb R)$ . Notice that it is open if the strong limit $\rho_*^\epsilon(t)\to\rho_*(t)$ can be justified. Furthermore, (3.12) is the entropy formulation of boundary conditions induced by $\mathcal E(\rho_*(t), 0)$ . It was proven in [22], that mass remains conserved at the junction after taking the limit. More precisely, we have

$\begin{equation} \sum\limits_{k = 1}^d A_k \overline{\rho_k u_k}(t) = 0, \quad\text{a.e. }t > 0, \end{equation}$

(3.13)

for the weak traces $\overline{\rho_k u_k}$ at $x = 0$ .

We summarize that after formally taking the (strong) limit, the traces at $x = 0$ satisfy the entropy formulation of the boundary condition $\mathcal E(\rho_*(t), 0)$ with boundary data $(\rho_*(t), 0)$ and mass is conserved at the junction. Assuming that the stronger formulation $\mathcal V(\rho_*(t), 0)$ of the boundary condition holds true, we obtain immediately Definition 1 and Definition 6.

3.3. Interpretation of the macroscopic coupling condition and necessary physical properties

In this subsection, we restrict ourselves to the generalized Riemann problem since it is a building block for the Cauchy problem.

The new macroscopic coupling condition is an implicit condition compared to the known coupling conditions in the literature. The idea of the new coupling condition is to assume the existence of left hand states of zero speed, independent of $k$ and such that mass is conserved. This is different to the known coupling conditions which are based on a coupling of traces of physical quantities.

We made a particular choice by choosing $u_* = 0$ for all left states. This choice can be interpreted in the following way. On the kinetic level the particles are stopped immediately after arriving at the junction. Then, the particles are instantaneously redistributed equally and into all pipelines. This artificial process leads to a coupling condition which does not prefer any pipeline and ignores the momentum of the incoming particles.

We interpret the macroscopic coupling condition by gas being stopped at the junction. Therefore, it is reasonable that we state a relation between the traces of the stagnation enthalpy at the junction. The stagnation enthalpy determines the enthalpy at a stagnation point after the gas is brought to a stop. This relation is given by inequalities depending on the signs of $\bar{u}_k$ at the junction (see Corollary 1). Furthermore, the macroscopic coupling satisfies several properties which seem to be necessary for a physically correct coupling condition. These properties are non-increasing energy at the junction and a maximum principle on the Riemann invariants (see Section 6 for more details). Furthermore, the same derivation technique applied to the full Euler equations leads to very similar results (see Section 8).

We emphasis that the coupling condition does not coincide with the Rankine-Hugoniot conditions in the case $d = 2$ . This can be easily checked since momentum is not conserved at the junction. Notice that this fact is not a disadvantage of the coupling condition since we want to model the coupling condition with maximum energy dissipation and conservation of mass but we neglect conservation of momentum. An interpretation of our coupling condition for $d = 2$ can be given by an infinitesimal small point were turbulence occurs due to a geometric effect at the junction.

Summarizing, the derivation of the coupling condition by the kinetic model and the maximum energy dissipation principle at the junction lead to a choice of an artificial state of zero speed at the junction. Furthermore, the interpretation of the macroscopic coupling condition by particles stopped at the junction and redistributed is only possible with a state of zero speed. From a formal mathematical point of view, the proofs of the physical properties in Section 6 work only if the artificial state has zero speed. This observation is due to the structure of (3.12) and the fact that $\eta'_{S, \rho u}(\rho_*, 0) = 0$ for symmetric $S$ .

4. The generalized Riemann problem

In this section, we prove existence and uniqueness of solutions to the generalized Riemann problem. Our strategy is as follows. Due to Lemma 1, there exists a solution for a fixed artificial density $\tilde\rho_*$ but possibly without conservation of mass at the junction. We define the mass production at the junction as a function of the artificial density $\tilde\rho_*$ and prove its continuity and monotonicity. We conclude with the intermediate value theorem. The proof is similar to the proof in [27]. Notice that the artificial density is a monotone momentum related coupling constant in the sense of [27]. The structure of the generalized Riemann problem allows us to extend the result to the supersonic region.

Proposition 2. Assume that initial data $(\hat\rho_k, \hat u_k)\in D$ are given. Let $(\rho_k, u_k)$ be the function obtained in Lemma 1 with artificial density $\tilde\rho_*\ge 0$ . Then, the trace $\bar\rho_k\bar u_k = (\rho_k u_k)(t, 0+)$ is continuous with respect to $\tilde\rho_*$ .

Proof. First, we prove the continuity with respect to the artificial density at fixed $\tilde\rho_*>0$ . Let $\mathcal A_*, \mathcal B_*, \mathcal C_*, \mathcal J_*$ be as in Definition 4. If $(\hat\rho_k, \hat u_k)$ lies in the interior of $\mathcal A_*$ , $\mathcal B_*$ and $\mathcal C_*$ , the continuity follows from the fact that the wave curves and the curves defined by $\lambda_i(\rho, u) = 0$ , $i = 1, 2$ are continuous. Therefore, it remains to prove continuity at the boundaries.

Step 1: First, we consider the boundary between $\mathcal A_*$ and $\mathcal B_*$ . More precisely, the 2-wave curve $\mathcal W_2(\rho_\alpha, u_\alpha)$ with $\{(\rho_\alpha, u_\alpha)\} = \{\lambda_1 = 0\}\cap \mathcal R_1(\tilde\rho_*, 0)$ . We have

$\begin{equation} \lim\limits_{\rho\searrow\tilde\rho_*}(\bar\rho_k, \bar u_k) = (\rho_\alpha, u_\alpha) = \lim\limits_{\rho\nearrow\tilde\rho_*}(\bar\rho_k, \bar u_k), \end{equation}$

(4.1)

since the 2-wave curve, the 1-rarefaction curve and the curve defined by $\{\lambda_1 = 0\}$ are continuous.

Step 2: Next, we consider the boundary between $\mathcal B_*$ and $\mathcal C_*$ . The continuity along $\mathcal R_2(\rho_\beta, u_\beta)$ is trivial since all involved curves are continuous. It remains to prove the continuity on $\mathcal J_*$ . For $\rho>\tilde\rho_*$ , we have

$\begin{equation} (\bar\rho_k, \bar u_k)(\rho) = (\hat\rho_{k}, \hat u_{k}), \end{equation}$

(4.2)

where $(\bar\rho_k, \bar u_k)(\rho)$ is the trace at $x = 0$ of the function obtained in Lemma 1 with artificial density $\rho$ . For $\rho<\tilde\rho_*$ and $|\rho-\tilde\rho_*|$ sufficiently small, the state $(\hat\rho_{k}, \hat u_{k})$ is connected to the boundary state $(\bar\rho_k, \bar u_k)(\rho)\in\mathcal S_1(\tilde\rho_*, 0)$ by a 2-shock with (small) positive speed. We get

$\begin{equation} \lim\limits_{\rho\nearrow\tilde\rho_*}\bar\rho_k \bar u_k = \hat\rho_{k} \hat u_{k}, \end{equation}$

(4.3)

since the speed of the 2-shock tends to zero as $\rho\nearrow\tilde\rho_*$ . The continuity on $\mathcal J_*$ follows from (4.2 – 4.3). Notice, that $\bar\rho_k$ and $\bar u_k$ itself are not continuous on $\mathcal J_*$ .

The continuity at $\tilde\rho_* = 0$ follows from similar considerations, see also [26].

Lemma 2. ([27,Remark 1]). Along the reversed 2-wave curves monotonicity in $\rho_0$ is equivalent to monotonicity in $u_0$ . More precisely,

$\begin{align} \frac{\mathrm{d}{u_0}}{\mathrm{d}{\rho_0}}\Big|_{\mathcal W_2^-} > 0, \quad\mathit{\text{for all}}(\rho, \rho u)\in D\backslash\{0\}. \end{align}$

(4.4)

The subscript denotes differentiation along the reversed $2$ -wave curve $\mathcal W^-_2(\rho, u)$ .

Proof. By the formula for the reversed 2-rarefaction wave, we have

$\begin{equation*} \frac{\mathrm{d}{u_0}}{\mathrm{d}{\rho_0}}\Big|_{\mathcal R_2^-} = \sqrt{\kappa\gamma}\rho_0^{\frac{\gamma-3}{2}} > 0. \end{equation*}$

Along the reversed 2-shock curve, we get

$\begin{align*} \frac{\mathrm{d}{u_0}}{\mathrm{d}{\rho_0}}\Big|_{\mathcal S^-_2}& = \frac{\kappa}{2\rho_0\rho} \frac{(1-\gamma)\rho_0^{\gamma-1}\rho+\gamma\rho_0^\gamma-\frac{\rho^{\gamma+1}}{\rho_0}}{\sqrt{\frac{\kappa(\rho_0^\gamma-\rho^\gamma)(\rho_0-\rho)}{\rho_0\rho}}}\\ & = \frac{\kappa}{2\rho_0\rho (u_0-u)}\left(\frac{\rho}{\rho_0}(\rho_0^\gamma-\rho^\gamma)+\gamma\rho_0^{\gamma-1}(\rho_0-\rho)\right) > 0, \end{align*}$

since $\rho<\rho_0$ and $u<u_0$ .

In the subsonic regime, we can determine the artificial density by a function $\tilde\rho_* = R_*(\bar\rho_k, \bar u_k)$ .

Definition 5. Let $R_*\colon D\to (0, \infty)$ be defined by

$\begin{align} R_*(\rho, u) = \begin{cases} \left(\rho^\theta+\frac{u}{a_\gamma}\right)^{1/\theta}, \quad &\text{if }u\ge 0, \\ R_*, \quad \text{with }u = -\sqrt{\frac{\kappa(\rho^\gamma-R_*^\gamma)(\rho-R_*)}{\rho R_*}}, R_* < \rho, \quad &\text{if }u < 0. \end{cases} \end{align}$

(4.5)

Notice that $R_*$ is well-defined, since for fixed $\rho>0$ , the function

$\begin{equation*} (0, \rho]\to(-\infty, 0], \quad R_*\mapsto -\sqrt{\frac{\kappa(\rho^\gamma-R_*^\gamma)(\rho-R_*)}{\rho R_*}}, \end{equation*}$

is bijective. By this definition we can reformulate RP1 in the subsonic regime by

$\begin{equation} R_*(\bar\rho_k, \bar u_k)(t) = \mathcal H_{R_*}(t), \quad \text{for } k = 1, \dots, d, \quad\text{for a.e. }t > 0. \end{equation}$

(4.6)

Compare the new condition (4.6) with the coupling conditions in (1.3 – 1.5) and note that they are different even in the subsonic regime.

Lemma 3. We have

$\begin{align} \frac{\mathrm{d}{R_*}}{\mathrm{d}{\rho_0}}\Big|_{\mathcal W_2^-} > 0, \quad\mathit{\text{for all}}(\rho, u)\in D\backslash\{0\}. \end{align}$

(4.7)

The subscript represents the differentiation along the reversed 2-wave curve.

Proof. For $u>0$ , differentiation along the 2-wave curve gives

$\begin{align} \frac{\mathrm{d}{R_*}}{\mathrm{d}{\rho_0}}& = \left(\rho_0^{\frac{\gamma-3}{2}}+\frac{1}{\sqrt{\kappa\gamma}}\frac{\mathrm{d}{u_0}}{\mathrm{d}{\rho_0}}\right) \left(\rho_0^\theta+\frac{u_0}{a_\gamma}\right)^{\frac{3-\gamma}{\gamma-1}} > 0, \end{align}$

(4.8)

since $\frac{\mathrm{d}{u_0}}{\mathrm{d}{\rho_0}}>0$ . For $u<0$ , differentiation of

$\begin{equation*} u_0 = -\sqrt{\frac{\kappa(\rho_0^\gamma-R_*^\gamma)(\rho_0-R_*)}{\rho_0 R_*}} \end{equation*}$

along the 2-wave curve gives

$\begin{equation*} \frac{\mathrm{d}{u_0}}{\mathrm{d}{\rho_0}} = \frac{\kappa}{u_0}\left[\gamma\left(\rho_0^{\gamma-1}-R_*^{\gamma-1}\frac{\mathrm{d}{R_*}}{\mathrm{d}{\rho_0}}\right)\left(\frac{1}{R_*}-\frac{1}{\rho_0}\right)+(\rho_0^\gamma-R^\gamma_*)\left(\frac{1}{\rho_0^2}-\frac{1}{R_*^2}\frac{\mathrm{d}{R_*}}{\mathrm{d}{\rho_0}}\right)\right], \end{equation*}$

or equivalently

$\begin{align*} \frac{\mathrm{d}{R_*}}{\mathrm{d}{\rho_0}} = \frac{\frac{2 u_0}{\kappa}\frac{\mathrm{d}{u_0}}{\mathrm{d}{\rho_0}}+\gamma\rho_0^{\gamma-1}\left(\frac{1}{\rho_0}-\frac{1}{R_*}\right)+\frac{1}{\rho_0^2}(R_*^\gamma-\rho_0^\gamma)}{\gamma R_*^{\gamma-1}\left(\frac{1}{\rho_0}-\frac{1}{R_*}\right)+\frac{1}{R^2_*}(R_*^\gamma-\rho_0^\gamma)}. \end{align*}$

The right hand side is strictly positive, since $\rho_0>R_*$ , $u_0<0$ and $\frac{\mathrm{d}{u_0}}{\mathrm{d}{\rho_0}}>0$ .

Lemma 4. We have

$\begin{equation*} \frac{\mathrm{d}{R_*}}{\mathrm{d}{\rho}}\Big|_{\lambda_1 = 0} > 0, \end{equation*}$

where the subscript denotes the differentiation along the curve $\{\lambda_1 = 0\}$ .

Proof. We differentiate $u$ along the curve defined by $\{\lambda_1 = 0\}$ and get

$\begin{equation*} \frac{\mathrm{d}{u}}{\mathrm{d}{\rho}} = \frac{\mathrm{d}{}}{\mathrm{d}{\rho}}\left(\sqrt{\kappa\gamma}\rho^\theta\right) = \sqrt{\kappa\gamma}\theta\rho^{\frac{\gamma-3}{2}} > 0. \end{equation*}$

The result follows with (4.8).

Proposition 3. Fix initial data $(\hat\rho_k, \hat u_k)\in D$ . Let $(\rho_k, u_k)$ be the solution obtained in Lemma 1 with artificial density $\tilde\rho_*\ge 0$ . Then, the trace of the momentum $\bar\rho_k\bar u_k$ is increasing in $\tilde\rho_*$ . It is strictly increasing if $(\hat\rho_k, \hat u_k)\in \mathcal A_* \cup\mathcal B_*$ , $\tilde\rho_*>0$ and constant if $(\hat\rho_k, \hat u_k)\in \operatorname{int}\mathcal C_*$ .

Proof. Step 1: We consider the case $(\hat\rho_k, \hat u_k)\in \operatorname{int}\mathcal A_*$ , $\tilde\rho_*>0$ . We already know that $(\bar\rho_k, \bar u_k)\in\mathcal R_1(\tilde\rho_*, 0)\cap\{\lambda_1 = 0\}$ . The formulas for $\mathcal R_1(\tilde\rho_*, 0)$ and $\lambda_1$ lead to

$\begin{equation*} \frac{\mathrm{d}(\bar\rho_k \bar u_k)}{\mathrm{d} \tilde\rho_*} = \frac{\mathrm d}{\mathrm{d}\tilde\rho_*}\left(\sqrt{\gamma\kappa}\left(\frac{2}{\gamma+1} \right)^{\frac{\gamma+1}{\gamma-1}}\tilde\rho_*^{\frac{\gamma+1}{2}}\right) = \sqrt{\kappa\gamma}\left(\frac{2}{\gamma+1}\right)^{\frac{2}{\gamma-1}}\tilde\rho_*^{\frac{\gamma-1}{2}} > 0. \end{equation*}$

Step 2: Next, we consider the case $(\hat\rho_k, \hat u_k)\in \mathcal B_*$ , $\tilde\rho_*>0$ . $(\bar\rho_k, \bar u_k)$ is the unique element in $\mathcal W_2(\hat\rho_k, \hat u_k)\cap \mathcal W_1(\tilde\rho_*, 0)$ . Since Lemma 3, we have

$\begin{equation*} \frac{\mathrm{d} R_*}{\mathrm{d}\rho_0}\Big|_{\mathcal W_2^-} > 0, \end{equation*}$

but this implies

$\begin{equation*} \frac{\mathrm{d}{\bar\rho_k}}{\mathrm{d}{\tilde\rho_*}} = \frac{\mathrm{d}{R_*^{-1}(\tilde\rho_*)}}{\mathrm{d}{\tilde\rho_*}} > 0, \end{equation*}$

where $R_*^{-1}(\tilde\rho_*)$ is the unique element in $\mathcal W_1(\tilde\rho_*, 0)\cap \mathcal W_2(\hat\rho_k, \hat u_k)$ . By chain rule, we get

$\begin{equation*} \frac{\mathrm{d}{(\bar\rho_k \bar u_k)}}{\mathrm{d}{\tilde\rho_*}} = \frac{\mathrm{d}{(\bar\rho_k \bar u_k)}}{\mathrm{d}{\bar\rho_k}}\cdot \frac{\mathrm{d}{\bar\rho_k}}{\mathrm{d}{\tilde\rho_*}}. \end{equation*}$

Therefore, it remains to prove

$\begin{equation*} \frac{\mathrm{d}{(\bar\rho_k \bar u_k)}}{\mathrm{d}{\bar\rho_k}}\Big|_{\mathcal W_2^-} > 0. \end{equation*}$

On $\mathcal R_2^-(\hat \rho_k, \hat u_k)$ , we have

$\begin{align*} \frac{\mathrm{d}{(\bar\rho_k \bar u_k)}}{\mathrm{d}{\bar\rho_k}} = &\bar u_k+\bar \rho_k\frac{\mathrm{d}{\bar u_k}}{\mathrm{d}{\bar\rho_k}} = \bar u_k+\sqrt{\kappa\gamma}\bar \rho_k^{\theta} = \lambda_2(\bar \rho_k, \bar u_k) > 0, \end{align*}$

since $(\bar \rho_k, \bar u_k)\in\mathcal B_*$ and Lemma 1. On $\mathcal S_2^-(\hat \rho_k, \hat u_k)$ , we use the fact that $\lambda_2(\bar\rho_k, \bar u_k)>0$ and get

$\begin{align*} \frac{\mathrm{d}{(\bar\rho_k \bar u_k)}}{\mathrm{d}{\bar\rho_k}}& = \bar u_k+\bar \rho_k\frac{\mathrm{d}{\bar u_k}}{\mathrm{d}{\bar\rho_k}}\\ & = \bar u_k+\frac{\sqrt{\kappa}}{2\bar\rho_k \hat\rho_k(\bar u_k-\hat u_k)}\left[\gamma\bar\rho_k^\gamma(\bar\rho_k-\hat\rho_k)+\hat\rho_k(\bar\rho_k^\gamma-\hat\rho_k^\gamma) \right]\\ & = \lambda_2(\bar\rho_k, \bar u_k)-\sqrt{\kappa\gamma}\bar\rho_k^\theta+\frac{\sqrt{\kappa}}{2\bar\rho_k \hat\rho_k(\bar u_k-\hat u_k)}\big[\gamma\bar\rho_k^\gamma(\bar\rho_k-\hat\rho_k)+\hat\rho_k(\bar\rho_k^\gamma-\hat\rho_k^\gamma) \big]\\ & = \lambda_2(\bar\rho_k, \bar u_k)+\frac{\sqrt{\kappa}}{2\bar\rho_k \hat\rho_k(\bar u_k-\hat u_k)}\\ &\cdot\left[ \gamma\bar\rho_k^\gamma(\bar\rho_k-\hat\rho_k)+\hat\rho_k(\bar\rho_k^\gamma-\hat\rho_k^\gamma)-2\sqrt{\gamma\bar\rho_k^\gamma(\bar\rho_k-\hat\rho_k)}\sqrt{\hat\rho_k (\bar\rho_k^\gamma-\hat\rho_k^\gamma)}\right]\\ & = \lambda_2(\bar\rho_k, \bar u_k)+\frac{\sqrt{\kappa}}{2\bar\rho_k \hat\rho_k(\bar u_k-\hat u_k)} \left[\sqrt{\gamma\bar\rho_k^\gamma(\bar\rho_k-\hat\rho_k)}-\sqrt{\hat\rho_k(\bar\rho_k^\gamma-\hat\rho_k^\gamma)}\right]^2\\ & > \lambda_2(\bar\rho_k, \bar u_k) > 0, \end{align*}$

for $\bar\rho_k\neq\hat\rho_k$ . For $(\bar \rho_k, \bar u_k)\in \mathcal W_1(\tilde\rho_*, 0)$ the left and right limit of the derivatives along $\mathcal W^-_2(\bar \rho_k, \bar u_k)$ are strictly positive and strict monotonicity in $\mathcal B_*$ follows.

Step 3: The strict monotonicity at the boundary between $\mathcal A_*$ and $\mathcal B_*$ can be shown by taking the left and right limit of the derivatives which are strictly positive.

Step 4: Finally, we consider the case $(\hat\rho_k, \hat u_k)\in \mathcal C_*$ . In the case $\lambda_2(\hat\rho_k, \hat u_k)\le 0$ , we have $(\bar \rho_k, \bar u_k) = (\hat \rho_k, \hat u_k)$ and observe that $(\bar \rho_k, \bar u_k)$ is locally constant as a function of $\tilde\rho_*$ . For $\lambda_2(\hat\rho_k, \hat u_k)\ge 0$ , we have

$\begin{equation*} \{(\bar\rho_k, \bar u_k)\} = \mathcal R_2(\hat\rho_k, \hat u_k)\cap \{\lambda_2 = 0\}. \end{equation*}$

Therefore, $(\bar\rho_k, \bar u_k)$ is locally constant with respect to $\tilde\rho_*$ .

Theorem 2. Assume that the initial states $(\hat\rho_k, \hat u_k)\in D$ are given. Then, there exists a unique solution $(\rho_k, u_k)$ to the generalized Riemann problem according to Definition 1 with a unique artificial density $\rho_*\ge 0$ .

Proof. For the solution $(\rho_k, u_k)$ obtained in Lemma 1 with artificial density $\tilde\rho_*\ge 0$ , the mass production at the junction is given by

$\begin{equation*} m(\tilde\rho_*) = \sum\limits_{k = 1}^d A_k \bar\rho_k \bar u_k. \end{equation*}$

Since Proposition 2 and 3, the function $m$ is continuous and increasing in $\tilde\rho_*$ . We will use these properties to prove existence and uniqueness by the intermediate value theorem. We divide the rest of the proof in two steps.

Step 1: We prove that there exist $0\le \rho_-<\rho_+$ such that

$\begin{equation*} m(\rho_-)\le 0\le m(\rho_+). \end{equation*}$

● We set

$\begin{equation*} \rho_- = \operatorname{argmin}\;\{\rho\ge 0\, |\, (\rho, 0)\in \mathcal W_2^-(\hat\rho_k, \hat u_k)\text{ for }k = 1, \dots, d\}. \end{equation*}$

Then, we have $(\bar\rho_k, \bar u_k)\in\mathcal V(\rho_-, 0)\cap \{\rho u\le 0\}$ for all $k = 1, \dots, d$ , but this implies $m(\rho_-)\le 0$ .

● We set

$\begin{equation*} \rho_+ = \operatorname{argmax}\;\{\rho\ge 0\, |\, (\rho, 0)\in \mathcal W_2^-(\hat\rho_k, \hat u_k)\text{ for }k = 1, \dots, d\}. \end{equation*}$

Then, we have $(\bar\rho_k, \bar u_k)\in \mathcal V(\rho_+, 0)\cap \{\rho u\ge 0\}$ for all $k = 1, \dots, d$ , but this implies $m(\rho_+)\ge 0$ .

By the intermediate value theorem, we conclude that there exists $\rho_*\ge 0$ such that $m(\rho_*) = 0$ .

Step 2: We prove that $m$ is strictly increasing at $\rho_*$ . Since $m(\rho_*) = 0$ , there exists $1\le k_0\le d$ such that $\bar\rho_{k_0}\bar u_{k_0}\ge 0$ . Due to Proposition 3 and $(\bar\rho_{k_0}, \bar u_{k_0})\in\mathcal A_*\cup\mathcal B_*$ , $\bar\rho_{k_0}\bar u_{k_0}$ is strictly increasing with respect to $\tilde\rho_*$ . Thus, $m$ is strictly increasing at $\rho_*$ . This implies the uniqueness of the artificial density $\rho_*\ge 0$ with $m(\rho_*) = 0$ .

Since Riemann problems admit unique self-similar Lax solutions and $\rho_*$ is uniquely defined by $(\bar\rho_k, \bar u_k)$ with $\bar\rho_k \bar u_k\ge 0$ , the solution to the generalized Riemann problem is unique.

5. The generalized Cauchy problem

In this section, we prove existence and uniqueness of solutions to the Cauchy problem. This result is based on a general existence theorem by Colombo, Herty and Sachers [12] and holds true in a neighborhood of a subsonic solution. We also obtain Lipschitz continuous dependence on the initial data.

Definition 6. Fix $(\hat\rho_1, \hat\rho_1\hat u_1, \dots, \hat\rho_d, \hat\rho_d\hat u_d)\in U^0+L^1(0, \infty)_x$ with $(\hat\rho_k, \hat\rho_k\hat u_k)\in D$ a.e. and $T\in(0, \infty]$ . Then, we call $(\rho_1, \rho_1 u_1, \dots, \rho_d, \rho_d u_d)\in C([0, T]_t, U^0+L^1(0, \infty)_x)$ a weak solution to the generalized Cauchy problem if $(\rho_k, \rho_k u_k)$ , $k = 1, \dots, d$ are weak entropy solutions to the isentropic gas equations and the following assertions hold true:

$\rm CP0:$ The solution satisfies the initial condition

$\begin{equation*} (\rho_k, u_k)(0+, x) = (\hat{\rho}_k, \hat{u}_k)(x)\in D, \end{equation*}$

for a.e. $x>0, \, k = 1, \dots, d$ .

$\rm CP1:$ For a.e. $t>0$ , there exists $\rho_*(t)\ge 0$ such that

$\begin{equation*} (\bar\rho_k, \bar u_k)(t)\in\mathcal V (\rho_*(t), 0) \end{equation*}$

for all $k = 1, \dots, d$ .

$\rm CP2:$ Mass is conserved at the junction

$\begin{equation*} \sum\limits_{k = 1}^d A_k\bar{\rho}_k\bar{u}_k = 0, \quad\text{for a.e. }t > 0. \end{equation*}$

Theorem 3. Fix a vector of subsonic states $U^0 = ( \rho^0_1, \rho_1^0 u^0_1, \dots, \rho^0_d, \rho_d^0 u_d^0)\in D^d$ such that the corresponding generalized Riemann problem admits a stationary solution. Then, there exist $\delta, L>0$ and a map $S\colon [0, \infty)\times \mathcal D\to\mathcal D$ such that

● $\mathcal D \supset \{ U\in U^0+L^1((0, \infty)_x, D^d)\, ;\, TV(U)\le \delta \}$ ;

● for $U\in \mathcal D$ , $S_0 U = U$ and for $s, t\ge 0$ , $S_s S_t U = S_{s+t}U$ ;

● for $U, V\in \mathcal D$ and $s, t\ge 0$ , $\| S_t U-S_s V\|_{L^1}\le L (\|U-V\|_{L^1}+|t-s|)$ ;

● if $U\in\mathcal D$ piecewise constant, then for $t>0$ sufficiently small, $S_t U$ coincides with the juxtaposition of the solution to Riemann problems centered at the points of jumps or at the junction.

Moreover, for every $U\in\mathcal D$ , the map $t\mapsto S_t U$ is a solution to the generalized Cauchy problem.

Proof. Since $(\hat\rho_k, \hat u_k)$ is subsonic, we can choose $\delta>0$ sufficiently small such that $\mathcal D$ is contained in the subsonic region. Therefore, CP1 is equivalent to

$\begin{equation*} R_*(\bar\rho_k, \bar u_k)(t) = \mathcal H_{R_*}(t), \quad \text{for } k = 1, \dots , d, \end{equation*}$

for a.e. $t>0$ . Note that $R_*$ is defined in Definition 5. Next, we apply Theorem 3.2 in [12]. Therefore, we define the function

$\begin{equation*} \Psi(U) = \begin{pmatrix} \sum _{k = 1}^d A_k \rho_k u_k\\ R_*(\rho_1, u_1)-R_*(\rho_2, u_2)\\ \vdots\\ R_*(\rho_{d-1}, u_{d-1})-R_*(\rho_d, u_d) \end{pmatrix}. \end{equation*}$

It remains to prove the transversality condition

$\begin{equation} \operatorname{det}\begin{bmatrix}D_1\Psi( U^0)\cdot r_2(\rho_1^0, u_1^0) &\cdots& D_d\Psi( U^0)\cdot r_d(\rho_d^0, u_d^0)\end{bmatrix}\neq 0, \end{equation}$

(5.1)

where $D_k = D_{(\rho_k, \rho_k u_k)}$ . By Lemma 3 and the proof of Proposition 3, we get

$\begin{align*} \frac{\partial \rho u}{\partial \rho}\Big|_{\mathcal W^-_2} > 0\quad \text{and} \quad\frac{\partial R_*}{\partial \rho}\Big|_{\mathcal W^-_2} > 0 \end{align*}$

in the subsonic region. We deduce that

$\begin{align*} D_k(\rho_k u_k)\cdot r_2(\rho_k , u_k) > 0\quad\text{and}\quad D_k R_*(\rho_k, u_k)\cdot r_2(\rho_k, u_k) > 0. \end{align*}$

This implies that the matrix involved in (5.1) has components with fixed sign which are given by

$\begin{equation*} \label{eq:SignsCellsMatrixTransversality} \begin{pmatrix} + & + & + &\cdots & \cdots & + \\ + & - & 0 & \cdots & \cdots & 0 \\ 0 & + & - & \ddots & &\vdots \\ \vdots & \ddots & \ddots & \ddots & \ddots&\vdots \\ \vdots & & \ddots & \ddots & \ddots&0 \\ 0 &\dots& \dots & 0 & + & - \end{pmatrix}. \end{equation*}$

A Laplace expansion implies that the determinant of this matrix has a fixed sign and is non-zero.

Remark 2. The existence and uniqueness result is restricted to subsonic initial data with sufficiently small total variation. The global result for the generalized Riemann problem and the large amount of inequalities for entropy fluxes at the junction (Propsition 4) motivate to prove a more general result. Notice that the method in [22] based on compensated compactness can be applied to the kinetic coupling condition (3.10). This result justifies the relaxation in the interior of the pipelines. Nevertheless, it is open how the traces relax at the junction and if the obtained macroscopic solution satisfies the coupling condition in Definition 6.

6. Energy/Entropy dissipation at the junction

In this section, we prove some physical properties of the coupling condition. In particular, we prove that energy is non-increasing at the junction, a relation for the stagnation enthalpy and a maximum principle on the Riemann invariants.

Proposition 4. Assume that initial states $(\hat\rho_k, \hat u_k)\in D$ are given. Let $(\rho_k, u_k)$ be the solution to the generalized Riemann or Cauchy problem. Then,

$\begin{equation*} \sum\limits_{k = 1}^d A_k G_S(\bar\rho_k, \bar u_k)\le 0, \end{equation*}$

for every convex $S\colon\mathbb R\to\mathbb R$ with $S(v) = S(-v)$ .

Proof. The case $\rho_* = 0$ is trivial. For $\rho_*\neq 0$ fix $S\in C^1(\mathbb R, \mathbb R)$ with $S(v) = S(-v)$ . Since $(\bar\rho_k, \bar u_k)\in\mathcal V(\rho_*, 0)\subset\mathcal E(\rho_*, 0)$ , we have

$\begin{equation*} G_S(\bar\rho_k, \bar u_k)-G_S(\rho_*, 0)-\eta_S'(\rho_*, 0)\left(F(\bar\rho_k, \bar u_k)-F(\rho_*, 0)\right)\le 0. \end{equation*}$

Notice that

$\begin{align*} G_S(\rho_*, 0)& = \int_{\mathbb R}\theta v \chi(\rho_*, v)S(v)\, \mathrm{d}v = 0, \quad \text{and}\\ \partial_{\rho u}\eta_{S}(\rho_*, 0)& = \frac{1}{J_\lambda}\int_{-1}^1 (1-z^2)^\lambda S'(a_\gamma\rho^\theta_* z)\mathrm{d}z = 0, \end{align*}$

since the integrands are anti-symmetric. These observations together with conservation of mass at the junction give

$\begin{align*} 0&\ge \sum\limits_{k = 1}^d A_k \left[ G_S(\bar\rho_k, \bar u_k)-G_S(\rho_*, 0)-\eta_S'(\rho_*, 0)(F(\bar \rho_k, \bar u_k)-F(\rho_*, 0))\right]\\ & = \sum\limits_{k = 1}^d A_k G_S(\bar\rho_k, \bar u_k)-\partial_\rho\eta_{S}(\rho_*, 0)\left(\sum\limits_{k = 1}^d A_k\bar \rho_k\bar u_k\right)\\ & = \sum\limits_{k = 1}^d A_k G_S(\bar\rho_k, \bar u_k). \end{align*}$

An approximation argument leads to the result for general convex functions $S$ .

Corollary 1 (Non-increasing energy). Assume that sufficiently regular initial data $(\hat\rho_k, \hat u_k)$ are given. Let $(\rho_k, u_k)$ be the solution to the generalized Riemann or Cauchy problem. Then, the following properties hold true:

(i) At the junction energy is non-increasing, i.e.

$\begin{equation*} \sum\limits_{k = 1}^d A_k G(\bar\rho_k, \bar u_k)\le 0, \quad \mathit{\text{for a.e.}}t > 0. \end{equation*}$

(ii) At the junction the traces of the stagnation enthalpy

$\begin{equation*} h(\rho, u) = \frac{u^2}{2}+\frac{\kappa\gamma}{\gamma-1}\rho^{\gamma-1} \end{equation*}$

are related by

$\begin{equation*} h(\bar\rho_k, \bar u_k)\le h(\rho_*, 0)\le h(\bar\rho_l, \bar u_l), \end{equation*}$

for $\bar u_l\le 0\le \bar u_k, 1\le k, l\le d$ , for a.e. $t>0$ .

Proof. Applying Proposition 4 to $S(v) = v^2/2$ gives

$\begin{equation*} \label{eq:proofEnergyDissipationJunction} \sum\limits_{k = 1}^d A_k G(\bar\rho_k, \bar u_k)\le 0. \end{equation*}$

Since $(\rho_k, u_k)\in\mathcal V(\rho_*, 0)\subset\mathcal E(\rho_*, 0)$ , we have

$\begin{equation*} G(\bar\rho_k, \bar u_k)-\partial_\rho \eta(\rho_*, 0)\bar\rho_k\bar u_k\le 0. \end{equation*}$

The result follows from dividing the inequality by $\bar\rho_k\bar u_k\neq 0$ and the fact $\partial_\rho \eta(\rho_*, 0) = h(\rho_*, 0)$ . The cases $\bar\rho_k\bar u_k = 0$ and $\rho_* = 0$ are trivial.

Corollary 2 (Maximum principle). Let $(\rho_k, u_k)\colon(0, \infty)_t\times(0, \infty)_x\to D$ be the solution to the generalized Riemann or Cauchy problem with initial condition $(\hat\rho_k, \hat u_k)$ . Assume that

$\begin{equation*} -\omega_M\le\omega_1(\hat\rho_k, \hat u_k)(x) < \omega_2(\hat\rho_k, \hat u_k)(x)\le\omega_M, \end{equation*}$

for a.e. $x>0$ , $k = 1, \dots, d$ . Then, we have

$\begin{equation*} -\omega_M\le\omega_1(\rho_k, u_k)(t, x) < \omega_2(\rho_k, u_k)(t, x)\le\omega_M, \end{equation*}$

for a.e. $t, x>0, \, k = 1, \dots, d$ .

Proof. We define the symmetric, positive function

$\begin{equation*} S_M(v) = (-\omega_M-v)^2_+ +(v-\omega_M)^2_+, \quad v\in\mathbb R. \end{equation*}$

As proven in [22], the definition of $\eta_{S_M}$ implies

$\begin{equation} \eta_{S_M}(\rho, u) = 0, \text{ if and only if }-\omega_M\le\omega_1(\rho, u) < \omega_2(\rho, u)\le\omega_M. \end{equation}$

(6.1)

Furthermore, the divergence theorem and the entropy condition give

$\begin{gather*} \int_0^\infty \eta_{S_M}(\rho_k, u_k)(T, x)\, \mathrm{d}x-\int_0^\infty \eta_{S_M}(\rho_k, u_k)(0, x)\, \mathrm{d}x\\ -\int_0^T G_{S_M}(\rho_k, u_k)(t, 0)\, \mathrm{d}t\le 0, \end{gather*}$

for $T>0$ . Taking the sum over $k$ , Proposition 4 and (6.1) lead to

$\begin{align*} 0\le \sum\limits_{k = 1}^d A_k \int_0^\infty \eta_{S_M}(\rho_k, u_k)(T, x)\, \mathrm{d}x\le \sum\limits_{k = 1}^d A_k \int_0^\infty \eta_{S_M}(\rho_k, u_k)(0, x)\, \mathrm{d}x = 0. \end{align*}$

The result follows from (6.1).

7. Numerical results

We give an example in which the new coupling condition produces the physically correct wave types. This observation is based on the assumption of the appearance of turbulence at the junction. The other coupling conditions produce different wave types. Furthermore, we study level sets associated to unphysical coupling conditions. We will consider the coupling conditions with equal pressure, equal momentum flux, equal stagnation enthalpy and the artificial density coupling condition.

7.1. A numerical example

We consider the shallow water equations and set $\gamma = 2$ , $\kappa = 5$ . We aim to compute solutions to the gerenalized Riemann problem in the sense of Lax.

The implementation is based on Newton's method applied to the coupling condition. The unknowns are the parameter of the reversed $2$ -wave curve of the initial data. We will consider examples which can be solved in the subsonic region. Therefore, we can use the function $R_*$ (see Definition 5) to implement the artificial density coupling condition.

We consider a three junction situation with one ingoing and two outgoing pipelines. The idea of this example is to assume that the sum of the momentum of the initial states is zero and the initial densities coincide. For this example we clearly get a stationary solution if we take the coupling condition with equal pressure. To make the example more concrete, we take the initial data in Table 1. The traces and energy dissipation of the numerical solutions are given in Table 2. The generalized Riemann problems are solved by the $2$ -waves in Table 3.

Table 1. Initial data.

pipeline	$\hat\rho_{k}$	$\hat\rho_{k} \hat u_{k}$
1	$+1.0000$	$-1.0000$
2	$+1.0000$	$+0.5000$
3	$+1.0000$	$+0.5000$

| Show Table

DownLoad: CSV

Table 2. Numerical results.

	Equal density		Equal momentum flux		Equal stagnation enthalpy		Equal artificial density
pipeline	$\bar\rho_{k}$	$\bar\rho_{k} \bar u_{k}$	$\bar\rho_{k}$	$\bar\rho_{k} \bar u_{k}$	$\bar\rho_{k}$	$\bar\rho_{k} \bar u_{k}$	$\bar\rho_{k}$	$\bar\rho_{k} \bar u_{k}$
1	$+1.0000$	$-1.0000$	$+0.8964$	$-1.1981$	$+0.8518$	$-1.2670$	$+1.1776$	$-0.5417$
2	$+1.0000$	$+0.5000$	$+1.0266$	$+0.5991$	$+1.0356$	$+0.6335$	$+0.9346$	$+0.2708$
3	$+1.0000$	$+0.5000$	$+1.0266$	$+0.5991$	$+1.0356$	$+0.6335$	$+0.9346$	$+0.2708$
Energy dissipation	$-7.5000\times 10^{-2}$		$-1.725\times 10^{-2}$		$\approx 0$		$-1.3852\times 10^{-1}$

| Show Table

DownLoad: CSV

Table 3. Wave types.

pipeline	Equal density	Equal momentum flux	Equal stagnation enthalpy	Equal artificial density
1	no waves	rarefaction wave	rarefaction wave	shock
2	no waves	shock	shock	rarefaction wave
3	no waves	shock	shock	rarefaction wave

| Show Table

DownLoad: CSV

We make the following observations:

● The $2$ -waves types obtained by solving with the artificial density coupling condition are different to the other coupling conditions. The wave types of the artificial density coupling conditions seem to be the physically correct ones. More precisely, we expect to have a shock in the incoming pipeline and rarefaction waves in the outgoing pipelines due to turbulence at the junction.

● The most energy is dissipated at the junction when the artificial density coupling condition is imposed.

● If the artificial density coupling condition is imposed, the momentum traces at the junction are smaller compared to the other coupling conditions. The momentum traces for the coupling conditions with equal pressure, equal momentum flux or equal stagnation enthalpy differ less strongly in comparison with each other.

7.2. Level sets associated the coupling conditions

In this section we consider the geometry of level sets corresponding to different coupling conditions. More precisely, we drop the condition on conservation of mass and compute the sets in which the attained boundary values may lie. For a solution to the generalized Riemann problem, the traces $\{(\bar\rho_k, \bar\rho_k\bar u_k)\, |\, k = 1, \dots, d \}$ are contained in one of these level sets. More precisely, we consider

$\begin{equation} \{(\rho, \rho u)\in D\, |\, \mathcal H(\rho, u) = \mathcal H(\rho_0, u_0)\}, \end{equation}$

(7.1)

where $(\rho_0, \rho_0 u_0)\in D$ is a fixed (subsonic) state and $\mathcal H$ denotes the pressure, momentum flux or stagnation enthalpy. For the artificial density coupling condition, we consider the set

$\begin{equation*} \mathcal V(\rho_*, 0), \quad\text{for suitable }\rho_*\ge 0, \end{equation*}$

which coincides with the definition in (7.1) in the subsonic case by taking $\mathcal H = R_*$ (see Definition 5). We set $\gamma = 1.4$ , $\kappa = 1$ and $(\rho_0, \rho_0 u_0) = (1, 0)$ . This choice leads to the level sets displayed in Figure 4.

Figure 4. Level sets for different coupling conditions.

DownLoad: Full-Size Img PowerPoint

We make the following observations:

● The coupling conditions with equal momentum flux and stagnation enthalpy induce bounded level sets.

● The level sets induced by equal pressure, momentum flux and stagnation enthalpy are symmetric w.r.t. the $\rho$ -axis and have the same tangent at $(1, 0)$ . This geometric property may lead to asymptotically similar behaviour of the coupling conditions near a stationary solution.

8. Full gas dynamics

In this section, we apply the derivation technique to full gas dynamics. We compute the maximum entropy dissipating kinetic coupling condition for kinetic models with the standard Maxwellian, e.g. the Boltzmann equation, the linear Boltzmann equation or the Boltzmann BGK model. We maximize the entropy dissipation and not the energy dissipation since energy is conserved. Again, we define a macroscopic coupling condition which can be formally obtained by a macroscopic limit. Notice that the formal macroscopic limit of the kinetic Boltzmann (type) equations is given by the full compressible Euler equations. The aim of this section is to underline that the presented approach is quite general and can be adopted easily to other hyperbolic systems equipped by a kinetic model and an entropy. Furthermore, the obtained results are very similar to the results for the isentropic gas equations.

8.1. Derivation of the kinetic coupling condition

First, we recall some basic definitions and explain the setting. As before, we consider a network of one-dimensional pipelines. We assume that one of the kinetic equations mentioned above is satisfied in the interior of the pipelines. The kinetic Boltzmann (type) equations admit the standard Maxwellian

$\begin{equation} M_{\rho, u, \theta}(\xi) = \frac{\rho}{\sqrt{2 \pi\theta}} \exp\left(\frac{-(u-\xi)^2}{2\theta}\right), \end{equation}$

(8.1)

for $\rho\ge 0$ , $u\in \mathbb R$ , $\theta > 0$ . Again, we define a kinetic coupling condition at the junction. Since energy is conserved in full gas dynamics, we obtain a second natural condition in addition to the conservation of mass. Nevertheless, two conditions are not enough to select a unique kinetic coupling condition. Therefore, we aim to select the kinetic coupling condition given by

$\begin{equation} \Psi\colon L^1_\mu((-\infty, 0)_\xi, [0, \infty))^d\to L^1_\mu((0, \infty)_\xi, [0, \infty))^d;\quad g\mapsto\Psi[g], \end{equation}$

(8.2)

which conserves mass and energy and dissipates as much entropy as possible. More precisely, we minimize

$\begin{equation} \sum\limits_{k = 1}^d A_k \int_0^\infty \xi\, \Psi^k[g](\xi) \log \Psi^k[g](\xi)\, \mathrm{d}\xi, \end{equation}$

(8.3)

with respect to

$\begin{align} \sum\limits_{k = 1}^d A_k\left(\int_0^\infty \xi\, \Psi^k[g](\xi) \, \mathrm{d}\xi+\int_{-\infty}^0 \xi\, g^k(\xi) \, \mathrm{d}\xi\right)& = 0, \end{align}$

(8.4)

$\begin{align} \sum\limits_{k = 1}^d A_k\left(\int_0^\infty \xi^3\, \Psi^k[g](\xi) \, \mathrm{d}\xi+\int_{-\infty}^0 \xi^3\, g^k(\xi) \, \mathrm{d}\xi\right)& = 0. \end{align}$

(8.5)

The unique minimizer of this problem is given by

$\begin{equation} \Psi^k[g](\xi) = M_{\rho_*, 0, \theta_*}(\xi), \quad\text{for }\xi > 0, \end{equation}$

(8.6)

where $\rho_*\ge 0, \, \theta_*>0$ are chosen such that (8.4 – 8.5) hold. The proof works similar to (3.7):

Since $v\mapsto v \log v$ is convex on $[0, \infty)$ and admits the derivative $v\mapsto \log v+1$ , we get

$\begin{align} &\sum\limits_{k = 1}^d A_k \int_0^\infty \xi\, \Psi^k[g](\xi) \log \Psi^k[g](\xi)\, \mathrm{d}\xi\\ \ge&\sum\limits_{k = 1}^d A_k\left( \int_0^\infty \xi\, M_{\rho_*, 0, \theta_*}(\xi) \log M_{\rho_*, 0, \theta_*}(\xi)\, \mathrm{d}\xi\right)\\ &+\sum\limits_{k = 1}^d A_k\left(\int_0^\infty \xi\, [\log M_{\rho_*, 0, \theta_*}(\xi)+1]\left(\Psi^k[g](\xi)-M_{\rho_*, 0, \theta_*}(\xi)\right)\, \mathrm{d}\xi\right)\\ = & \sum\limits_{k = 1}^d A_k\int_0^\infty \xi\, M_{\rho_*, 0, \theta_*}(\xi) \log M_{\rho_*, 0, \theta_*}(\xi)\, \mathrm{d}\xi. \end{align}$

(8.7)

The last step follows by

$\begin{equation} \log M_{\rho_*, 0, \theta_*}(\xi)+1 = \log \frac{\rho}{\sqrt{2 \pi\theta}} +1 -\frac{\xi^2}{2\theta} \end{equation}$

(8.8)

and (8.4 – 8.5). It can be easily checked that for every $g$ there exists $(\rho_*, 0, \theta_*)$ such that (8.4 – 8.5) hold for $\Psi^k[g](\xi) = M_{\rho_*, 0, \theta_*}(\xi)$ . By convexity of $v\mapsto v \log v$ , it follows that entropy is non-increasing at the junction, i.e.

$\begin{equation} \sum\limits_{k = 1}^d A_k\left(\int_0^\infty \xi\, M_{\rho_*, 0, \theta_*}(\xi) \log M_{\rho_*, 0, \theta_*}(\xi)\, \mathrm{d}\xi+\int_{-\infty}^0 \xi\, g(\xi) \log g(\xi)\, \mathrm{d}\xi\right)\le 0. \end{equation}$

(8.9)

8.2. The full Euler equations on networks and a new coupling condition

Next, we consider the macroscopic limit. As mentioned above, the kinetic equation converges formally towards the full Euler equations for ideal polytropic gas given by

$\begin{align} \begin{cases} \partial_t \rho +\partial_x(\rho u)& = 0, \\ \partial_t(\rho u)+\partial_x(\rho u^2+\rho\theta)& = 0, \\ \partial_t(\frac{\rho u^2}{2}+\frac{\rho\theta}{2})+\partial_x(\frac{\rho u^3}{2}+\frac{3}{2}\rho u \theta)& = 0, \end{cases}\quad \text{for a.e. }t > 0, x\in\mathbb R, \end{align}$

(8.10)

with density $\rho\ge 0$ , mean velocity $u\in\mathbb R$ , temperature $\theta>0$ and adiabatic index $\gamma = 3$ . The equations of full gas dynamics model conservation of mass, momentum and energy. As usually, we impose the additional entropy condition

$\begin{equation} \partial_t\left(\rho \log \left( \frac{\rho}{\theta^{1/2}}\right)\right)+\partial_x\left(\rho u \log \left( \frac{\rho}{\theta^{1/2}}\right) \right)\le 0, \quad \text{for a.e. }t > 0, x\in\mathbb R. \end{equation}$

(8.11)

The full Euler equations on networks were studied before by several authors [14,15,19,24]. We summarize the main ideas of the constructed coupling conditions. Analogous to isentropic gas dynamics, conservation of mass at the junction is imposed

$\begin{equation} \sum\limits_{k = 1}^d A_k\bar{\rho}_k\bar{u}_k = 0, \quad\text{a.e. }t > 0. \end{equation}$

(8.12)

Since energy is conserved in full gas dynamics, we additionally assume that energy is conserved at the junction, i.e.

$\begin{equation} \sum\limits_{k = 1}^d A_k\overline{\left(\frac{\rho u^3}{2}+\frac{3}{2}\rho u \theta \right)_k} = 0. \end{equation}$

(8.13)

There are more conditions needed to single out a unique solution. Most of them are a straight forward extension of a coupling condition for isentropic gas. We give a short overview of the coupling conditions in the literature:

Colombo and Mauri [15] introduced equality of momentum flux at the junction

$\begin{equation} \overline{\left(\rho u^2+\rho\theta \right)_k} = \mathcal H_{\mathrm{MF}}(t), \quad\text{for a.e. } t > 0, k = 1, \dots, d. \end{equation}$

(8.14)

Herty [19] used equality of pressure

$\begin{equation} \bar\rho_k\bar\theta_k = \mathcal H_p(t), \quad\text{for a.e. } t > 0, k = 1, \dots, d. \end{equation}$

(8.15)

Networks consisting of $d = 2$ pipelines with different cross-sectional area were studied by Colombo and Marcellini [14] with different coupling conditions. One of them is based on a smooth approximation of the discontinuity in the cross-section. Lang and Mindt [24] impose equality of stagnation enthalpy

$\begin{equation} \overline{\left(\frac{u^2}{2}+\frac{3}{2}\theta \right)_k} = \mathcal H_{\mathrm{SE}}(t), \quad\text{for a.e. } t > 0, k = 1, \dots, d, \end{equation}$

(8.16)

and equality of entropy for traces with outgoing flow

$\begin{gather} \log\left(\frac{\rho}{\theta^{1/2}}\right) = \mathcal H_{\mathrm{S}}(t), \quad\text{for a.e. }t > 0, \text{ for }\bar u_k > 0, \end{gather}$

(8.17)

$\begin{gather} \text{with}\quad \mathcal H_{\mathrm{S}}(t) = \frac{\sum_{\bar u_k < 0}A_k \overline{\left(\rho u \log\left(\frac{\rho}{\theta^{1/2}}\right) \right)_k}}{\sum_{\bar u_k < 0}A_k \bar \rho_k\bar u_k}. \end{gather}$

(8.18)

These two conditions imply conservation of energy and entropy at the junction. Notice that conservation of entropy at the junction is not consistent with the fact that entropy can be dissipated in full gas dynamics.

In full gas dynamics an additional phenomena appears since the number of ingoing/outgoing characteristics at the junction can change in the subsonic region. This fact makes it more complicated to prove existence and uniqueness results. Nevertheless, we can use the formal arguments in Section 3 and the derivation in the previous subsection to define the following new coupling condition for full gas dynamics.

Definition 7. Fix initial data $(\hat\rho_k, \hat u_k, \hat\theta_k)\in D_{3 \times 3} = \{(\rho, u, \theta)\, |\, \rho> 0, u\in\mathbb R, \theta>0 \text{ or }\rho = u = \theta = 0\}, k = 1, \dots, d$ . Then, we call $(\rho_k, u_k, \theta_k)\colon (0, \infty)_t\times(0, \infty)_x\to D$ a weak solution to the generalized Riemann problem if the following assertions hold true:

$\rm RP0:$ The solution satisfies the initial condition

$\begin{equation*} (\rho_k, u_k, \theta_k)(0+, x) = (\hat{\rho}_k, \hat{u}_k, \hat{\theta}_k)\in D_{3 \times 3}, \quad \text{for }x > 0, \, k = 1, \dots, d; \end{equation*}$

$\rm RP1:$ There exists $(\rho_*, 0, \theta_*)\in D_{3\times 3}$ such that $(\rho_k, u_k, \theta_k)$ is equal to the restriction to $x>0$ of the Lax solution to the standard Riemann problem with initial condition

$\begin{equation*} (\rho_k, u_k, \theta_k)(0+, x) = \begin{cases} (\hat{\rho}_k, \hat{u}_k, \hat{\theta}_k), \quad& x > 0, \\ (\rho_*, 0, \theta_*), \quad& x < 0, \end{cases} \end{equation*}$

for all $k = 1, \dots, d$ ;

$\rm RP2:$ Mass is conserved at the junction

$\begin{equation*} \sum\limits_{k = 1}^d A_k\bar{\rho}_k\bar{u}_k = 0, \quad\text{for all }t > 0; \end{equation*}$

$\rm RP3:$ Energy is conserved at the junction

$\begin{equation*} \sum\limits_{k = 1}^d A_k\overline{\left(\frac{\rho u^3}{2}+\frac{3}{2}\rho u \theta \right)_k} = 0, \quad\text{for all }t > 0. \end{equation*}$

Notice that this condition leads to conservation of mass and energy at the junction by definition. Furthermore, entropy is non-increasing at the junction by the entropy formulation of boundary conditions $\mathcal E(\rho_*, 0, \theta_*)$ , i.e.

$\begin{equation} \sum\limits_{k = 1}^d A_k\overline{\left(\rho u \log\left(\frac{\rho}{\theta^{1/2}}\right) \right)_k}\le 0. \end{equation}$

(8.19)

Therefore, the new coupling condition satisfies some necessary physical properties.

9. Conclusion

We introduced a new coupling condition for isentropic gas and proved existence and uniqueness of solutions to the generalized Riemann and Cauchy problem.

The derivation of the coupling condition is based on the kinetic model and the selection of the unique kinetic coupling condition which conserves mass and dissipates as much energy as possible. The obtained kinetic coupling condition distributes the incoming kinetic data into all pipelines by the same Maxwellian with suitable artificial density and zero speed. Formal arguments lead to a corresponding macroscopic definition to the generalized Riemann problem. In this definition the artificial state with zero speed appears as the (left) initial state for a standard Riemann problem.

In addition to the derivation, we proved physical properties of the coupling condition. The coupling condition ensures that energy is non-increasing at the junction and leads to a maximum principle on the Riemann invariants. Furthermore, a relation of the traces of the stagnation enthalpy at the junction was given. Notice that these properties hold true due to the choice of an artificial state with zero speed.

We gave an example in which the new coupling condition is the only condition producing the physically correct wave types. The solutions to the generalized Riemann problems were computed numerically. Furthermore, we studied level sets related to different coupling conditions and their geometry.

Finally, we considered the coupling condition in view of the model hierarchy of gas dynamics by applying the same approach to full gas dynamics and Boltzmann (type) equations. We took the kinetic coupling conditions with conservation of mass and energy at the junction and maximize the entropy dissipation. This consideration leads to very similar results. In particular, we obtained an artificial state with suitable density and temperature and again with zero speed.

In summary, we defined a new coupling condition, derived several physical and mathematical properties and gave a motivation. Future research may consider more detailed numerical aspects and the rigorous justification of the considerations in Section 3.

References

[1]	Crowding and queuing at entrances, (2018), http://ped.fz-juelich.de/da/doku.php?id=wuptmp.
[2]	Jupedsim, (2019), https://www.jupedsim.org.
[3]	J. Adrian, M. Boltes, S. Holl, A. Sieben and A. Seyfried, Crowding and queuing in entrance scenarios: Influence of corridor width in front of bottlenecks, arXiv e-prints, arXiv: 1810.07424.
[4]	When push comes to shove: Exclusion processes with nonlocal consequences. Physica A (2015) 437: 119-129.
[5]	M. Bardi and I. Capuzzo-Dolcetta, Optimal Control and Viscosity Solutions of Hamilton-Jacobi-Bellman Equations, Systems & Control: Foundations & Applications, Birkhäuser Boston, Inc., Boston, MA, 1997.
[6]	Lane formation by side-stepping. SIAM Journal on Mathematical Analysis (2016) 48: 981-1005.
[7]	Flow characteristics in a crowded transport model. Nonlinearity (2016) 29: 3528-3550.
[8]	Simulation of pedestrian dynamics using a two-dimensional cellular automaton. Physica A: Statistical Mechanics and its Applications (2001) 295: 507-525.
[9]	Distance functions and almost global solutions of eikonal equations. Communications in Partial Differential Equations (2010) 35: 391-414.
[10]	A. Corbetta, J. A. Meeusen, C.-m. Lee, R. Benzi and F. Toschi, Physics-based modeling and data representation of pairwise interactions among pedestrians, Physical Review E, 98 (2018), 062310.
[11]	E. Cristiani, B. Piccoli and A. Tosin, Multiscale Modeling of Pedestrian Dynamics, MS & A, Modeling, Simulation and Applications, 12. Springer, Cham, 2014.
[12]	On the Hughes' model for pedestrian flow: The one-dimensional case. Journal of Differential Equations (2011) 250: 1334-1362.
[13]	D. C. Duives, W. Daamen and S. Hoogendoorn, Trajectory analysis of pedestrian crowd movements at a Dutch music festival, Pedestrian and Evacuation Dynamics 2012, Springer, (2014), 151–166.
[14]	Parameter estimation for macroscopic pedestrian dynamics models from microscopic data. SIAM Journal on Applied Mathematics (2019) 79: 1475-1500.
[15]	B. Hein, Agent-Based Modelling for Crowding and Queuing in Front of Bottlenecks, Bachelor's Thesis, University of Wuppertal, 2019.
[16]	D. Helbing and P. Molnár, Social force model for pedestrian dynamics, Physical Review E, 51 (1995), 4282.
[17]	K. Hirai and K. Tarui, A simulation of the behavior of a crowd in panic, Systems and Control.
[18]	A continuum theory for the flow of pedestrians. Transportation Research Part B: Methodological (2002) 36: 507-535.
[19]	A. Johansson and D. Helbing, Analysis of empirical trajectory data of pedestrians, Pedestrian and Evacuation Dynamics 2008, (2009), 203–214.
[20]	A. Jüngel, The boundedness-by-entropy method for cross-diffusion systems, Nonlinearity, 28 (2015), 1963.
[21]	Simulation of evacuation processes using a bionics-inspired cellular automaton model for pedestrian dynamics. Physica A: Statistical Mechanics and its Applications (2002) 312: 260-276.
[22]	C. Koutschan, H. Ranetbauer, G. Regensburger and M.-T. Wolfram, Symbolic derivation of mean-field PDEs from lattice-based models, 2015 17th International Symposium on Symbolic and Numeric Algorithms for Scientific Computing (SYNASC), 27–33.
[23]	O. Ladyzhenskaya, Linear and Quasilinear Elliptic Equations, ISSN, Elsevier Science, 1968.
[24]	Explicit formula for scalar nonlinear conservation laws with boundary condition. Mathematical Methods in the Applied Sciences (1988) 10: 265-287.
[25]	Explicit formula for weighted scalar nonlinear hyperbolic conservation laws. Transactions of the American Mathematical Society (1988) 308: 667-683.
[26]	C. Lehrenfeld, On a Space-Time Extended Finite Element Method for the Solution of a Class of Two-Phase Mass Transport Problems, PhD thesis, RWTH Aachen, 2015, http://publications.rwth-aachen.de/record/462743.
[27]	A. Y. Leroux, Approximation de Quelques Problèmes Hyperboliques Non-Linéaires, Thèse d'état, Rennes, 1979.
[28]	B. Maury and S. Faure, Crowds in Equations. An Introduction to the Microscopic Modeling of Crowds, Advanced Textbooks in Mathematics. World Scientific Publishing Co. Pte. Ltd., Hackensack, NJ, 2019.
[29]	Experimental study of the behavioural mechanisms underlying self-organization in human crowds. Proceedings of The Royal Society B: Biological Sciences (2009) 276: 2755-2762.
[30]	S. Nowak and A. Schadschneider, Quantitative analysis of pedestrian counterflow in a cellular automaton model, Physical Review E, 85 (2012), 066128.
[31]	S. Okazaki, A study of pedestrian movement in architectural space, Trans. of A.I.J., 283.
[32]	Time-evolving measures and macroscopic modeling of pedestrian flow. Archive for Rational Mechanics and Analysis (2011) 199: 707-738.
[33]	Fast sweeping methods for eikonal equations on triangular meshes. SIAM J. Numerical Analysis (2007) 45: 83-107.
[34]	C. Rudloff, T. Matyus and S. Seer, Comparison of different calibration techniques on simulated data, Pedestrian and Evacuation Dynamics 2012, Springer, (2013), 657–672.
[35]	A. Schadschneider, C. Eilhardt, S. Nowak and R. Will, Towards a calibration of the floor field cellular automaton, Pedestrian and Evacuation Dynamics, (2011), 557–566.
[36]	A. Schadschneider, H. Klpfel, T. Kretz, C. Rogsch and A. Seyfried, Fundamentals of pedestrian and evacuation dynamics, IGI Global, (2009), 124–154.
[37]	Comparative study of macroscopic pedestrian models. Transportation Research Procedia (2014) 2: 477-485.
[38]	U. Weidmann, Transporttechnik der Fussgänger: Transporttechnische Eigenschaften des Fussgängerverkehrs, Schriftenreihe des IVT, IVT, 1993.
[39]	W. G. Weng, T. Chen, H. Y. Yuan and W. C. Fan, Cellular automaton simulation of pedestrian counter flow with different walk velocities, Physical Review. E, 74 (2006), 036102.
[40]	C. A. Yates, A. Parker and R. E. Baker, Incorporating pushing in exclusion-process models of cell migration, Physical Review E, 91 (2015), 052711.

This article has been cited by:

1.	Martin Gugat, Michael Herty, 2022, 23, 9780323850599, 59, 10.1016/bs.hna.2021.12.002
2.	Evgenii S. Baranovskii, Vyacheslav V. Provotorov, Mikhail A. Artemov, Alexey P. Zhabko, Non-Isothermal Creeping Flows in a Pipeline Network: Existence Results, 2021, 13, 2073-8994, 1300, 10.3390/sym13071300
3.	Sara Grundel, Michael Herty, Hyperbolic discretization of simplified Euler equation via Riemann invariants, 2022, 106, 0307904X, 60, 10.1016/j.apm.2022.01.006
4.	EIKE FOKKEN, SIMONE GÖTTLICH, MICHAEL HERTY, Efficient simulation of coupled gas and power networks under uncertain demands, 2022, 0956-7925, 1, 10.1017/S0956792522000079
5.	Michael Herty, Niklas Kolbe, Siegfried Müller, Central schemes for networked scalar conservation laws, 2022, 18, 1556-1801, 310, 10.3934/nhm.2023012
6.	Yannick Holle, Entropy Dissipation at the Junction for Macroscopic Traffic Flow Models, 2022, 54, 0036-1410, 954, 10.1137/21M1423920
7.	Michael Herty, Niklas Kolbe, Siegfried Müller, A Central Scheme for Two Coupled Hyperbolic Systems, 2024, 6, 2096-6385, 2093, 10.1007/s42967-023-00306-5
8.	Andrea Corli, Massimiliano D. Rosini, Ulrich Razafison, 2024, Mathematical Modeling of Chattering and the Optimal Design of a Valve*, 979-8-3503-1633-9, 76, 10.1109/CDC56724.2024.10886245
9.	Andrea Corli, Ulrich Razafison, Massimiliano D. Rosini, Coherence of Coupling Conditions for the Isothermal Euler System, 2025, 0170-4214, 10.1002/mma.10847

Reader Comments

Your name:*

Email:*
© 2020 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)