On Existence and Admissibility of Singular Solutions for Systems of Conservation Laws

A weak notion of solution for systems of conservation laws in one dimension is put forward. In the framework introduced here, it can be shown that the Cauchy problem for any \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$n\times n$$\end{document}n×n system of conservation laws has a solution. The solution concept is an extension of the notion of singular \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\delta $$\end{document}δ-shocks which have been used to provide solutions for Riemann problems in various systems, for example in cases where strict hyperbolicity or the genuine-nonlinearity condition are not satisfied, or in cases where initial conditions have large variation. We also introduce admissibility conditions which eliminate a wide range of unreasonable solutions. Finally, we provide an example from the shallow water system which justifies introduction of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\delta $$\end{document}δ-distributions as a part of solutions to systems of conservation laws.


Introduction
We are interested in the general n × n system of one-dimensional conservation laws . . .
The standard theory of hyperbolic conservation laws is concerned with solutions which are at worst locally integrable. More precisely, if the system is strictly hyperbolic and genuinely nonlinear 1 , and if the total variation of u 0 is small enough, then the Cauchy problem has a solution [5,6,21,39,41,42]. If the initial functions are step functions, then the Cauchy problem is called the Riemann problem. If the differences in function values at the jumps are small then the Riemann problem can be solved uniquely using Lax-admissible solutions consisting of rarefaction waves, compressive shock waves and contact discontinuities.
On the other hand, if the above conditions are not fulfilled, then the Cauchy, or even the Riemann problem may not admit a Lax admissible weak solution or even any weak solution. Indeed, the standard theory is far from complete, as even some fairly simple problems cannot be resolved in a satisfactory manner (see e.g. [13,16,35,37,45]). Nevertheless, one might expect that in the case of physically relevant systems (such as e.g. the p-system arising in gas dynamics and nonlinear elasticity) one should be able to prove existence of weak solutions of bounded variation (BV) of the corresponding Cauchy problem. On the other hand, in fairly recent contributions [8,64] where the p-system was considered, it was shown that the Glimm scheme and wave front tracking, respectively, for a Cauchy problem with large initial data blow up (in the sense that uniform BV bounds on the approximate solutions are not available). However, note that for certain functions p, Cauchy problems for the p system can be resolved [4].
In order to deal with situations where the standard theory fails, expanding the space of possible solutions to the space of Radon measures has been useful in some cases. In particular, if solutions are allowed to contain Dirac δ-distributions, some systems to which the standard theory does not apply can be resolved. The first result in this direction can be found in [37]. However, the interest in such an approach (as well as its necessity) was particularly raised after it was shown in [35] that strictly hyperbolic and genuinely nonlinear system of conservation laws with compact Hugoniot locus feature these non-standard solutions.
The question which naturally appears is how to incorporate the δ-distribution as a part of solution to (1.1). If we substitute the distribution directly into (1.1), the problem of multiplication of singular distributions arises. It turns out that when one deals with systems which are linear with respect to one of the unknowns (and since we can approximate any locally integrable initial data by a step function), it is sufficient to understand the multiplication of a δ distribution and a Heaviside function (see [61] and references therein).
Let us focus therefore on the issue of defining the multiplication of a δ-distribution and a Heaviside function. There are several reasonable ways to define such a product, but we shall mention only two main approaches. The first one is to define H (x)δ(x) dx = H (0) which means that one needs to define the value of the Heaviside function at the jump point. Such an approach has been used in [13,28] to introduce the notion of measure-type solution for (1.1). The framework introduced in [28] yields uniqueness of solutions if an additional condition of Oleinik-type is required. However, that approach can be applied only to equations which are linear with respect to one of the unknowns, as already alluded to above.
The other main approach to making sense of singular solutions is the weak asymptotic method which entails defining smooth approximations (δ ε ) and (H ε ), and considering the product (δ ε H ε ). The limit as ε → 0 of the family (δ ε H ε ) can then be defined to be the product of the δ-distribution in the context of the weak asymptotic method, and one obtains approximate solutions to (1.1) for which one concludes that they converge weakly to a measure containing a δ-distribution. This approach has been used in [14,16,65] among many others, and is similar to the vanishing viscosity and the vanishing pressure methods used for example in [29,35,45]. One may also note that the family (δ ε H ε ) can be considered as an element of a Colombeau algebra, such as defined in [11]. This approach has proved fruitful in some cases (see [46] for example). In this direction, we also mention the α-product by Sarricco (see e.g. [51,52,57,58] and references therein) which is applied in the case of various systems of conservation laws and similar equations and which works in the case of fully nonlinear equation (i.e. nonlinear with respect to all unknowns).
Close examination of the weak asymptotic method reveals that is does not give precise meaning to whether and how the limiting distributions satisfy the original differential equations. This was first achieved in the work by Danilov and Shelkovich [17] where a variational formulation for the δ-type framework was proposed for systems which are linear with respect to one of the unknown functions and for δ-solutions where components containing the δ-distributions are supported on isolated curves. Indeed, the same idea can be used to define δ-solutions for general 2 × 2 systems of conservation laws, such as explained in [30]. In the present contribution, we extend the variational formulation to the case of n × n systems such as (1.1) and we prove that the Riemann problem for this system always has a solution in the framework of the definition given. We then extend the existence theorem to the general Cauchy problem by using a front-tracking approach. The key step is to estimate the BV norm of the solution after each interaction, and this is facilitated by introducing the singular parts of the solution.
The variational formulation of the Cauchy problem is very weak, and uniqueness of solutions cannot be expected. In fact, it can be shown that even the Riemann problem for a 2 × 2 system cannot be solved uniquely in the context of the variational formulation of singular solutions (see [30]). We partially deal with this issue by formulating an admissibility criterion which while not providing uniqueness at least eliminates a wide range of possible trivial solutions. The admissibility criterion is physically motivated and requires that at least one of the equations of the system satisfies the Rankine-Hugoniot condition at every shock.
Finally, it is shown that at least in the special case of the Riemann problem and a polynomial flux function, the solution can be obtained from the weak asymptotic method satisfying the admissibility conditions. This result establishes that the solution stems from an approximation procedure within the framework of smooth approximate solution converging towards δ-type solutions which is similarly in nature the above mentioned vanishing viscosity or vanishing pressure approaches. This tool enables us to make the connection between Rankine-Hugoniot deficits as appearing in the variational formulation, and the δ-type solutions to the original system.
As the results proved here do not rely on the assumption of hyperbolicity, the impression may arise that the solution concept provided here is too broad to be useful. We partially agree with such an assessment, but we also note that so far there are no general results for solvability of systems of the form (1.1), (1.2) in the literaure (except for some very special situations). In addition, our contribution is not an attempt to solve the problem in the standard space of locally integrable functions, but to investigate an alternative approach which would eventually lead to a satisfactory solution to (1.1), (1.2) from both the physical and the mathematical point of view. Finally, as shown in Sect. 3, we are able to find approximate smooth solutions to Riemann problems for a wide range of conservation laws which converge towards the δ-type solutions introduced here.
The paper is organized as follows. The definition of the variational formulation, admissibility conditions, and the existence theorem are given in Sect. 2. Sect. 3 contains a detailed account of the application of the weak asymptotic method to the Riemann problem for a 2 × 2 system. Finally, in Sect. 4, we review the physical validity of singular solutions of conservation laws and show how the present approach arises naturally in the context of the the shallow-water system.

Variational Formulation of Delta-Shock Solutions
In this section, we shall introduce a solution framework in which any Cauchy problem corresponding to (1.1) with BV -initial data has a solution. It is based on the variational formulation of δ-shock solutions introduced in [17] in the case of 2 × 2 systems which are linear with respect to one of the unknowns. Note that in such systems the δ distribution arises as a part of the solution defined along characteristics [16] as well as in the limit of the vanishing viscosity approximations [17,35]. As shown in [31], the variational formulation may also be applied to to arbitrary 2×2 systems of conservation laws, and this extension may be justified using the complex weak asymptotic method. However, in general this variational notion of solution is so weak that uniqueness can only be proved in special cases. In the following, we will generalize the variational formulation to n × n systems of conservation laws and Cauchy data, such as given by (1.1), (1.2).
As a first step, we shall define δ-type solutions for Riemann problems corresponding to systems of conservation laws. Suppose is a graph in the closed upper half plane containing a single Lipschitz continuous arc γ . Let 0 = {x 0 } be the initial point of the arc γ . Let u k , k = 1, . . . , n, be distributions of the form where U k ∈ L ∞ (R + × R) and α k (t, x) are real-valued functions. Let ∂ϕ(t,x) ∂l denote the tangential derivative of a function ϕ on the graph γ , and let γ denote the line integral over the arc γ .
First note that if we assume γ = {(t, γ (t)) : t ∈ [0, T ]}, we can also assume that ϕ(T , x) = 0 in which case one should add the terms − γ α k (T , γ (T ))ϕ(T , γ (T )) on the left-hand sides of (2.2). Also notice that (2.2) represents a generalization of the classical weak solution concept for (1.1). Indeed, if we omit the integral over the curve γ , we reach the standard variational formulation of (1.1). However, Definition 2.1 is very weak in the sense that a wide range of distributions satisfy its requirements. This is not surprising since even in the case of the standard weak solutions for scalar conservation laws, non-uniqueness issues arise [6]. Let us first look at the system (1.1) with singular Riemann data Using Definition 2.1, it is not difficult to see that for any c ∈ R, and any given initial states (2.4) and the amplitude α(t) = (α 1 (t), . . . , α d (t)) of the singular part of the shock is given by In other words, the following theorem holds.

4), and α(t) is given by (2.5). Then u(t, x) is a solution of the Riemann problem
Proof The proof of the theorem follows by substituting u into (2.2) (for simplicity, we shall assume here and in the sequel T = ∞ unless stated otherwise). After standard transformations, we reach the identities From here, the statement of the theorem follows immediately.
Before introducing admissibility conditions which will eliminate trivial situations given by the previous theorem, we will define a solution framework which can be applied to any Cauchy problem associated to (1.1). Notice that if we assume that the curves from the Definition 2.1 are straight lines of the form simplifies in the sense that instead of the tangential derivative appearing there, we simply have the expression (assume that supp' (2.7) Moreover, if we assume that we are solving the Riemann problem, then we can also assume that α are constants (at least in the intervals between possible waves interactions). With this example in mind, we introduce a subspace of the dual to the Bochner space L 1 ([0, T ]; C 0 (R)) for a given T > 0 as follows.

Definition 2.2
Let T > 0 be given, and let C 0 (R) be the space of continuous functions on R which decay to 0 at infinity. We denote by G e (T ) the subspace of the dual of the Bochner It can be proved that the closure of G e (T ) with respect to the weak topology actually coincides with the entire dual space is the space of real-valued finite Radon measures. Indeed, we have the following proposition.
where the closure Cl is understood with respect to the weak topology.
Proof Let μ ∈ L 1 ([0, T ]; C 0 (R)) be arbitrary. We construct an approximation as follows. For a given nN, we introduce a uniform partition of the interval [−n, n] by setting x n k = k n , k = −n 2 , · · · , n 2 . We then define Since μ(t, ·) has finite variation, and indeed, sup ). Thus, to prove weak convergence, it is sufficient to consider test functions from a dense subspace. In particular, we can show that and since ϕ has compact support, the last term in the expression above is zero for large enough n. Moreover, such ϕ(t, x) is uniformly continuous in (t, x), and can be approximated by a function taking a finite number of values (a simple function). Given a challenge number ε, we can take n large enough such that Since μ has bounded variation, we have for large enough n, and this concludes the proof.
We next want to define an operation similar to (2.7) for elements of We use a partition {x n k } as in the previous proposition. We then take the following approximation of the measure μ by elements of G e : We can now introduce the following definition.

Definition 2.3 We say that the functional
where α n k are given by (2.9).

Remark 1
We note that in the case when the measure has the form then the corresponding tangential derivative has the form Indeed, by the definition of the generalized tangential derivative, we need to approximate the measure μ by the atomic measures of the form μ n (t, from where we see that and (2.11) follows directly.
Clearly, m l will always exist since a bounded sequence in L ∞ w * ([0, T ]; M(R)) is weakly precompact. On the other hand, such m l might not be unique since different subsequences could converge to different limits. So while the tangential derivative in general is not well defined, fixing the approximation in the from (2.9) reduces the number of possible generalized tangential derivatives significantly.
A simple example is the δ-distribution of the form α(t)δ(x − ct). In this case, it is straightforward to see that the generalized tangential derivative should be α (t)δ(x − ct) and that it is the unique tangential derivative of α(t)δ(x − ct). As a more interesting example, let us compute the tangential derivative for a measure defined by a regular function. Since the δ-distribution is the distributional derivative of a shock wave, an interesting case will be the measure defined by the derivative of a rarefaction wave. We stress that we are not solving the conservation law given in the example below, but merely compute the tangential derivative of a measure. A solution to the scalar conservation law does not contain any special structures as part of the solution (see Definition 2.4 below), they are, as well known, ordinary locally integrable functions.
Example 1 Rarefaction waves are given by the functions of the form u( x t ) which, in the scalar case, solve the equation Denote by m the measure defined by the function ∂ x u( m(t, x) = u x t . Let us compute ∂ l m. According to the definition of the generalized tangential derivative, we consider the approximation The tangential derivative is a weak limit (along a subsequence) of ∂ t m n (t, x). We have which converges weakly toward ∂ x f (u( x t )). From the above, we conclude that the generalized tangential derivative of m is With these preliminaries in place, we can introduce L ∞ w * ([0, T ]; M(R))-solutions to an arbitrary Cauchy problem corresponding to (1.1).

Definition 2.4
We say that the pair of distributions u(t, x) = (U(t, x), m(t, x)) := U(t, x) + m(t, x), (2.12) where U(t, x) = (U 1 (t, x), is a measure-type solution to (1.1) with the initial data u| t=0 = u 0 if the following relations hold for any ϕ ∈ C 1 c (R + × R) for some generalized tangential derivatives ∂ l m = (∂ l m 1 , . . . , ∂ l m d ) of the measure m.
As is obvious from the formulation of the functionals m and ∂m ∂l , this definition is a generalized version of Definition 2.1. Indeed, Definition 2.4 reduces to Definition 2.1 if the δ-distributions involved in the solution to (1.1) are supported on isolated smooth curves (see (2.7)).
However, Definition 2.4 is quite general as for any initial data u 0 ∈ BV (R) and any fixed speed vector c, the distribution represents a measure-type solution to (1.1). To reduce the number of possibilities, we introduce basic admissibility conditions which will not provide uniqueness of this type of solution, but which eliminates pathological solutions such as (2.14).

Definition 2.5
We say that the measure-type solution to (1.1) satisfies basic admissibility conditions if it represents a weak limit to the following sequence of measures where 0 ≤ t n j ≤ T n j for every n, j ∈ N, α n j (t) = (α n 1 j , . . . , α n d j )(t) are affine functions, and U n (t, x) = (U n  1 (t, x), . . . , U n d (t, x)) are step functions such that at each step at least one of the equations in system (1.1) is satisfied by U n (t, x) in the standard weak sense (i.e. the corresponding Rankine-Hugoniot conditions are satisfied).
The latter definition is motivated by the wave front tracking procedure. Actually, using such an approach, we are able to prove that there exists the measure-type solution to any Cauchy problem corresponding to (1.1) with BV-initial data satisfying the basic admissibility conditions 2 . For simplicity, we shall provide the proof in the case of 2 × 2 system. The proof is analogous in the case of an n × n system, and it is based on BV-estimates which do not increase during the interactions. Instead, we have to increment the strength of δ-functions appearing during the interaction, but the increment is also finite. We note also that we use Remark 1 in the course of the proof of the next theorem. Theorem 2.2 Let n = 2, and assume that f 1 , f 2 ∈ C 1 (R 2 ) and u 10 , u 20 ∈ L ∞ (R)∩ BV (R). Then, the system (1.1) with the initial data u 1 | t=0 = U 10 , u 2 | t=0 = U 20

admits a solution in the sense of Definition 2.4 satisfying the basic admissibility conditions.
Proof We shall apply the method of wave front tracking [6,27]. In the framework of that method, one approximates the initial data by piecewise constant functions and then solves Riemann problems at each step. We thus obtain new waves which mutually interact. In the case of the standard wave front tracking approach which implies solving a Riemann problem using Lax admissible shock waves, the BV-character of the solution can be lost as the result of the interactions of different waves (the solution can blow up [7,64]).
Using our definition, it is not difficult to construct approximate solution for an arbitrary Cauchy problem using the wave front tracking algorithm. Instead of a possible increase of the local BV-bounds, we will obtain δ-type solutions, and the coefficient accompanying the δ-function will fix the Rankine-Hugoniot deficit.
As we have seen in the proof of Theorem 2.1, when solving a Riemann problem, we can choose practically any speed of the shock which solves the problem. We correct the mistake originating from the Rankine-Hugoniot conditions using the δ-shocks. What is reasonable to do is to minimize the strength of the δ-shocks (since in this way we get a smaller error with respect to the standard situation when we use the weak solution concept). The correction will be smaller if the speed of the shock is smaller. In other words, we shall use the following procedure: (A) We approximate the functions U 01 and U 02 by the piecewise constant functions where χ {(x j ,x j+1 )} is the characteristic function on the interval (x j , x j+1 ) defined in the proof of Proposition 2.1, and we require U 01 − U N 01 L 1 (K ) = o(1), K ⊂⊂ R, for i = 1, 2. For j ∈ Z, denote by The corresponding speed c j , appearing in the place of c from (2.6), satisfies Items A) and B) are standard in wave front tracking. As for the item C), we use Theorem 2.1. More precisely, we solve each of the Riemann problems so that we adjoin a δ-distribution either to u 1 or to u 2 (or to both u 1 and u 2 ) so that the speed of the corresponding δ-shock is minimal and that, at the same time, it is given by the Rankine-Hugoniot conditions of one of the equations.
We thus obtain the solution (u N 1 , u N 2 ) of the approximate problem until the set of first interactions. The interaction means that the shock waves (together with the accompanying δdistributions) are at the same point in R + ×R. At that moment, the shocks and δ-distributions will merge (the state between two shocks will disappear and δ will have the strength equal to the strength of the interacting δ-s at the moment of interaction; see Fig. 1).
At that moment, we stop the time and solve a new set of Riemann problems as done in Theorem 2.1 (this time involving the δ-distributions as part of initial data). Thus, for any t ∈ [0, T ], we conclude that the distributions u ε 1 and u ε 2 have the form (2.16) for some piecewise constant functions U N 1 and U N 2 , and where t j denotes the time of creation of the singularity, and T j denotes the time of the next interaction. The solution satisfies the relations (2.17) and . Notice that, according to the construction (see Fig. 1), the interactions do not raise the BV-bound of the regular part of the solution, while the total sum of the strength-increment of the δ-distributions can be bounded by the following expression (see (2.5)): where BV (w) = sup i |w(x i ) − w(x i−1 )| and M is the L ∞ norm of the initial data. From here, we also conclude that ∂ t U N is bounded in L 1 ([0, T ]; W −1,1 loc (R)) since for any compactly supported ϕ ∈ C 1 c ([0, T ] × R) From here and since the BV-bound of the L ∞ -parts U N 1 and U N 2 of the families u N 1 and u N 2 , respectively, are finite, according to the Aubin-Lions lemma, we conclude that U N 1 and U N 2 are strongly L 1 loc -precompact. Thus, they admit an L 1 loc limit along a subsequence as N → ∞, and we denote the limit by U 1 and U 2 .
Moreover, since the sum of the strength of all δ-functions from (2.16) remains uniformly bounded with respect to ε, for every T > 0 there exists the weak limit to the functions from (2.16) in Li p([0, T ] × M loc (R)). In view of the previous paragraph, we denote these limits by We can assume that the limit exists along the same subsequences of (U N 1 ) and (U N 2 ) as well as for the terms of the form (2.9).
Passing to the limit in (2.17) and (2.18) along the subsequence defining simultaneously the generalized tangential derivative of m 1 and m 2 , and the converging subsequences of (u N 1 ) and (u N 2 ), we conclude that they satisfy (2.13).

Weak Asymptotics
Close examination of the variational formulation (2.2) of weak singular solutions reveals that it is given in terms of the Rankine-Hugoniot deficit, and it is not entirely clear how this deficit is connected to a singular solution. Therefore, it is important to be able to construct the singular solutions as limits of approximate solutions of a regularized form of the equation.
Following the method laid out in [17] and [30], we use the complex-valued weak asymptotic method in order to make the connection between the Rankine-Hugoniot deficit and the singular solutions containing δ-distributions.
In order to outline the weak asymptotic method, let us first define a vanishing family of distributions. a family of distributions depending on ε ∈ (0, 1). We The estimate on the right-hand side is understood in the usual Landau sense. Thus, we may say that a family of distributions approach zero in the sense defined above if for a given test function φ, the pairing f ε , φ converges to zero as ε approaches zero.

Definition 3.2
We say that the family of complex-valued distributions (u ε ) = (u 1ε , . . . , u dε ) represents a weak asymptotic solution to (1.1) if there exist real-valued distributions u = (u 1 , . . . , u d ) ∈ C(R + ; D (R)), such that for every fixed t ∈ R + = (0, ∞) u ε u as ε → 0, in D (R), and In the following, we will show that any 2×2 system of conservation laws with polynomial flux function admits a sequence of approximate solutions converging towards a δ-type solution of the Riemann problem. To this end, let us consider the system where a 1 , . . . , a m and b 1 , . . . , b n are constants while p j , q j , j = 1, . . . , m, and r j , s j , j = 1, . . . , n, are non-negative integers. First, we shall formally extend the system by an equation which is trivially satisfied so that we can control our approximations. Without loss of generality, we assume that s n = max{s 1 , . . . , s n } > 0 and b n = 1, and we write with the condition w(0, x) = 1 and, for simplicity, r n = 0. The case r n = 0 can be handled similarly. Clearly, w ≡ 1 and the systems (3.1), (3.2), and (3.3), (3.4), (3.5) are equivalent. But since we need to control approximations on infinitesimal intervals, the function w (or more precisely its approximation) will play a substantial role. We have the following theorem: Then, there exist (complex valued) families (u ε ), (v ε ) and (w ε ) such that (3.8) such that w ε t=0 1 in D (R), and w ε 1 in D (R + × R), Moreover, where c is given by the Rankine-Hugoniot conditions from Eq.
while α is the Rankine-Hugoniot deficit.
Proof Fix large M > 0 (say M > 100) and consider the following approximation w ε : Then, let From the above assumptions, we see that (3.7) reduces to (3.9) and this is clearly satisfied according to the choice of the speed c.
As for the family (v ε ), we take a standard approximation of the Dirac distribution: (δ ε ) is non-negative family of functions compactly supported at (−Mε, 0) (3.10) Observe that the functions (v ε ) may not be real-valued. If we notice that Eq. (3.5) becomes in the limit Thus, choosing α (t) to be the Rankine-Hugoniot deficit, we conclude the theorem.
We note that by introducing more auxiliary equations of the type ∂ t w = 0, w(0, x) = 1 and positioning them in appropriate places in the original equation, we can get a large number of different approximate solutions consistent with the conclusions of Theorem 2.1.

Energy loss
Momentum loss Fig. 2 Left panel: schematic picture of a hydraulic jump. Momentum is conserved, but energy is lost due to surface undulations and turbulence. Right panel: flow and under a sluice gate. Momentum loss occurs, but energy is conserved

Discussion
One natural question to ask is whether the appearance of δ-type solutions is physically reasonable. Many of the systems which were shown to have these singular solutions were found using mathematical rather than physical considerations. Nevertheless, it may be argued that δ-type solutions appear rather naturally especially from a modeling point of view. First of all, δ-distributions are experimentally confirmed as solutions of systems of conservation laws. This has recently been observed in [44] in nonlinear chromatography. Namely, during an experiment involving chemical interactions of different substances, the authors of [44] noticed abrupt increment of concentrations on specific isolated sets. The results are confirmed mathematically in several papers where the system of chromatography equations is considered (cf. [36,62,65]). Similar concentration effects were also observed in a system modeling the polymer flooding of a porous medium [26], in pressureles gas dynamics [67] and in thin-film equations [40,59].
The necessity to consider Rankine-Hugoniot deficits and corresponding singular solutions can also be understood in terms of conservation of mass, momentum and energy. Consider the the shallow-water system If the solutions of the system are smooth, they also satisfy a corresponding conservation equation for horizontal momentum given by and an energy equation which takes the form Let us consider the situations sketched in Fig. 2. In both cases given there, the solutions develop discontinuities. Since mass should always be conserved, we need to decide whether to pair mass conservation with momentum or energy conservation in order to find appropriate solutions.
The phenomenon occurring in the case shown in the left panel is called a hydraulic jump which is frequently observed in open channel flow such as rivers and spillways. In this case, it is well known that the conservation of mass and momentum is required [25]. If the hydraulic jump is traveling at a speed c, the corresponding Rankine-Hugoniot conditions take the form It is also well known that there is a loss of energy across the transition region [2,66], even if the flow is not turbulent.
On the other hand, in the flow under a sluice gate such as shown in the right panel of Fig. 2, mass and energy are conserved. Momentum is lost because the sluice gate exerts a force F on the fluid. However energy is conserved since when the fluid touches the gate, its velocity is zero. Therefore, no work is done (dW /dt = Fu where W is the energy, and u = 0). Thus the Rankine-Hugoniot conditions take the form where most appropriately, c = 0.
Let us now assume that we have both situations simultaneously, i.e. a flow containing both a traveling hydraulic jump and a passage under a sluice gate in an adjacent section of the channel. If both phenomena are to be described by a single system of two equations (since they happen simultaneously), a natural choice is the system (4.1) and (4.2). Locally, proper discontinuous solutions are obtained using the Rankine-Hugoniot conditions corresponding to either equations (4.1) and (4.3) for the hydraulic jump, or to (4.1) and (4.4) for the flow under the sluice gate. These local solutions will not in general satisfy the Rankine-Hugoniot conditions corresponding to (4.1) and (4.2), and it can be shown in particular in the case of a traveling hydraulic jump that a nonzero Rankine-Hugoniot deficit must occur in (4.2). However, this deficit can be handled by incorporating a δ-distribution into the solution u.
The necessity to go beyond the standard BV-solutions can also be ascribed to the use of mathematically ineligible operations during the model derivation; for instance the tacit assumption on smoothness of solutions. Indeed, it is well understood that the introduction of the notion of weak solution causes non-uniqueness. In the case of scalar conservation laws, this imperfection is recovered using the Kruzhkov entropy concept [38] but only if the flux is Lipschitz continuous. Clearly, the number of possibly ineligible operations is greater in the case of systems of conservation laws, and may cause a variety of problems such as non-existence and non-uniqueness.
Indeed, uniqueness is an issue which is quite rarely considered in the framework of δshock solutions. Uniqueness was obtained in [28] using entropy inequalities of Oleinik-type [49], and in [47] with the help of entropy inequalities of Kruzhkov-type [38]. These methods may yield results in special cases, but they do not offer a clear reason why the uniqueness is gained, and what the entropy inequality means physically or even mathematically.
The weak asymptotic method appears to lead to the correct choice of the shock speed for singular solutions. In fact, in cases where the equations have singular solutions given in exact form, the weak asymptotic method can be shown to give the right solution (cf. [34,56] for example). However, there is not a firm proof that the weak asymptotic method works in every case.
It is our belief that one of the reasons for the lack of a rigorous uniqueness concept, beside obvious mathematical difficulties, essentially lies in the oversimplification of the physical system, leading to ambiguities in the mathematical treatment. One possible remedy would be to consider the shock structure as part of the system, such as for example explained in [66]. However, then the relative simplicity of the conservation laws would also be lost.
Funding Open access funding provided by Austrian Science Fund (FWF). The authors have not disclosed any funding.
Data Availability Data sharing is not applicable to this article as no new data were created or analyzed in this study.

Declarations
Competing Interests The authors have not disclosed any competing interests.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.