$\def\pd{\partial}$ $\def\l{\left}\def\r{\right}$ $\def\mdot{{\dot{m}}}$ $\def\eps{\varepsilon}$ Consider a tube with longitudinal coordinate $x$ from $0$ to $l$ and varying cross-section $A(x)$. Derek S. Bale (a former student of Prof. Randall J. LeVeque) gives in his PhD (http://faculty.washington.edu/rjl/students/dbale/thesis.ps.gz) the Euler equations for density $\rho$, impulse density $(\rho u)$ and total energy density $e$ in such a tube in the form $$ \begin{array}{ccccl} \pd_t \rho &+& \frac1{A}\pd_x(A(\rho u)) &=& 0\\ \pd_t (\rho u) &+& \frac1{A}\pd_x\l(A\l(\frac{(\rho u)^2}\rho+p\r)\r) &=& \frac{A'}{A}p\\ \pd_t e &+& \frac1{A}\pd_x\l(A(e+p)\frac{(\rho u)}{\rho}\r)&=& 0 \end{array} $$ as an example for a spatially varying flux function.

Thereby, the pressure is determined by the gamma-law $p =(\gamma-1)\l(e-\frac12 \frac{(\rho u)^2}{\rho}\r)$. Note, that I have slightly reformulated the system such that beside the pressure $p$ only conserved quantities are used.

But, is it not better to change the conserved variables to the mass per length $m':=\rho A$, the mass flow $\mdot:=A\rho u$ and the energy per length $E':=Ae$?

As an intermediate variable one could introduce a pressure force $P:=A p=(\gamma-1)\l(E'-\frac12\frac{\mdot^2}{m'}\r)$.

Euler's equations reformulated with these new conserved quantities are: $$ \begin{array}{ccccl} \pd_t m' &+& \pd_x\mdot &=& 0\\ \pd_t \mdot &+& \pd_x\l(\frac{\mdot^2}{m'}+P\r) &=& \frac{A'}{A}P\\ \pd_t E' &+& \pd_x\l((E'+P)\frac{\mdot}{m'}\r)&=& 0 \end{array} $$ In this way one obtains a system with the structure of the original Euler equations (including a source term) where the flux function does not depend on the spatial variable. Therefore, the numerics should work fine. I could use Roe's method on this problem.

Do I miss here something?


Bale uses the first form to calculate supersonic waves through a tube with a narrowing halfway. The speed-up in the narrowing causes a supersonic wave and the corresponding shock wave.

The variation of $A(x)$ is smooth: $$ A(x) = \begin{cases} 1& \text{ for }x<1\\ 1+\frac\eps2\bigl(\cos(\pi(x-1))-1\bigr)&\text{for }1\leq x \leq 3\\ 1&\text{for }x>3 \end{cases} $$

I am aware of the fact that smooth transformations of the conserved quantities, such as $(m',\mdot,E')=(A\rho,A(\rho u),A e)$ have the potential to change the structure of discontinusous solutions.

But, is this really the case here?

  • $\begingroup$ I do use the second form that you have written. In fact, the second form is the more fundamental one since it is a consequence of applying conservation principle on a control volume. The first form is obtained after some differential calculus manipulations of the second form, assuming enough smoothness in the solution. $\endgroup$ – cfdlab Mar 14 '14 at 16:37
  • $\begingroup$ Yes, that is exactly what I am nervous about. I have edited the question in this regard. $\endgroup$ – Tobias Mar 14 '14 at 16:58
  • $\begingroup$ You want to use the second form, right ? That is indeed the "conservation form". The first form is obtained by transforming the second form. So in my opinion you have to worry about the correctness of the first form. But if LeVeque is using this form, maybe there is some other issue involved !!! $\endgroup$ – cfdlab Mar 14 '14 at 17:18
  • $\begingroup$ @PraveenChandrashekar: I think both forms are correct. See my answer below. $\endgroup$ – Tobias Mar 14 '14 at 19:56
  • $\begingroup$ @PraveenChandrashekar Note, that your reassuring comments were important to me. Thanks. Do you have some literature reference with the second form? It would be nice if I could cite it in my work. $\endgroup$ – Tobias Mar 15 '14 at 7:19

$\def\rmC{{\mathrm C}}$ $\def\vr{{\vec r}}$ $\def\vA{{\vec A}}$ $\def\ve{{\vec e}}$ $\def\l{\left}\def\r{\right}$ In the following I try to show that Bale's form is correct. If you complain, leave a comment.

The equations we can trust are the integral equations.

Mass Balance:

$$ \l[\int_0^x \rho(x',t') A(x') dx'\r]_{t'=0}^t + \int_{t'=0}^t \l[(\rho u)(x',t')A(x')\r]_{x'=0}^x d t' =0 $$ with the density $\rho(x,t)$ and the momentum flow $(\rho u)(x,t)$ as averages over $A(x)$.

Momentum balance:

$$ \l[\int_0^x (\rho u)(x',t') A(x')dx'\r]_{t'=0}^t + \ve_x\cdot\int_{t'=0}^t\l( \int_{\partial V}(\rho \vec{v})\vec{v}\cdot d\vec{A} + \int_{\partial V} p d\vec{A}\r)dt =0 $$ with the volume $V=\bigcup_{x'\in[0,x]}A(x')$.

There is no flow out of complementary part $A_\rmC:=\pd V\setminus \bigl(\{0\}\times A(0)\cup \{x\}\times A(x)\bigr)$ of the tube cross-sections $\{0\}\times A(0)$ and $\{x\}\times A(x)$. Furthermore, we approximate $\int_{A(x)} (\rho u)(\vr,t) \vec{v}(\vr,t)\cdot d\vA \approx (\rho u)(x,t) u(x,t) A(x)$ and remark that for integration in x-direction we have $\ve_x\cdot d\vA = d A_x(x')$.

That gives: $$ \begin{array}{l} \l[\int_0^x (\rho u)(x',t') A(x')dx'\r]_{t'=0}^t+ \int_{t'=0}^t \l[(\rho u)(x',t') u(x',t') A(x') + pA(x')\r]_{x'=0}^{x} dt =\\ =\quad-\int_{t'=0}^t \int_{x'=0}^x p(x',t') d A(x') dt' \end{array} $$

Thereby $\int_{x'=0}^x p(x',t') d A(x')$ can be interpreted as Riemann-Stieltjes integral which can be re-written as $\int_{x'=0}^x p(x',t') A'(x') dx'$ in the case that $A(x)$ depends smoothly on $x$.

The Riemann-Stieltjes integral suggests to use $$ p_k\frac{A_{k+\frac12}-A_{k-\frac12}}{A_k} $$ as a 1st-order update for the source term to capture discontinuities of $A(x)$. The 1st-order update of the flux function in the FVM-method is already designed such that discontinuities of $A(x)$ are captured.

Energy balance:

With a similar approximation for the averages over the cross-section as for the momentum balance and with the zero flux out of $A_\rmC$ we get: $$ \l[\int_0^x e(x',t')A(x')dx'\r]_{t'=0}^t + \int_{t'=0}^t \biggl[(e(x',t') + p(x',t'))u(x')A(x')\biggr]_{x'=0}^x dt = 0 $$


LeVeque describes in chapter 6.16 "Capacity-Form Differencing" of his FVM-book that conservation laws of the form $$ \frac{d}{dt}\int_{x_1}^{x_2} \kappa(x) q(x,t) dx = f(q(x_1,t))-f(q(x_2,t)) $$ with a capacity function $\kappa(x)$ can be treated with the help of capacity-form differencing.

The more general integral form is $$ \l[\int_{x_1}^{x_2} \kappa(x) q(x,t) dx\r]_{t_1}^{t_2} = \int_{t_1}^{t_2}(f(q(x_1,t))-f(q(x_2,t))) dt. $$ In the above conservation laws $A(x')$ can be seen as capacity function.

We see that the equations in Bale's form are correct even for shock waves with capacity-form differencing.

The other form I have given in my question is also correct. For that case the conserved quantities stand right under the spatial integral on the left hand side. Everything else is just combination of terms without differentiation.

Now I believe that Bale used his form of the equations to have a simple example for a system with spatially varying flux function.

| cite | improve this answer | |
  • $\begingroup$ You are right. There are identical formulations. $\endgroup$ – cfdlab Mar 15 '14 at 4:35
  • $\begingroup$ @PraveenChandrashekar: I have added a note about the interpretation of the source term as Riemann-Stieltjes integral and the corresponding numerical approximation. This discretization is also given in eprints.maths.ox.ac.uk/1343. Thanks for the reference. Do you know wether it works well for discontinuities in $A$ with explicite FVM-methods? (I will try it anyway. I am already coding...) $\endgroup$ – Tobias Mar 17 '14 at 15:04

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.