I wrote a program in Fortran which calculates Sod shock tube numerically. Now I want to compare it with exact (analytical) solution, but I don't know how to get it.

Can I find it somewhere or do I have to write a new program which calculates it? If it's the latter, where can I find the equations I need?


For the Sod shock tube, there are 5 states, two of which come from initial conditions (states I & V in the linked Wiki entry). Two of the remaining 3 can be calculated with relative ease given the linked Wikipedia entry; the most difficult state (III) requires a root solver to find the value of the pressure.

What's not covered in the Wiki entry is state II, which is the rarefaction wave region. In this region, we have that: \begin{align} u_{2}(x)&=\frac{2}{\gamma+1}\cdot\left(c_1+\frac{\left(x-x_\text{mid}\right)}{t}\right)\\ \rho_{2}(x)&=\rho_1\cdot\left(1-\frac{\gamma-1}{2}\cdot\frac{u_{2}(x)}{c_1}\right)^{2/(\gamma-1)}\\ P_{2}(x)&=P_1\cdot\left(1-\frac{\gamma-1}{2}\cdot\frac{u_{2}(x)}{c_1}\right)^{2\gamma/(\gamma-1)} \end{align} where $x_\text{mid}$ is the point that divides the two states (often taken to be 0.5 with boundaries at 0 & 1), $t$ the time since the simulation started and $c_1$ the speed of sound in state I (the left state)--note that this assumes that the left state is the over-pressure region--all other terms take their normal meaning.

Building a function or two to find the states should be pretty straight-forward, once the above is known.

  • 1
    $\begingroup$ I know Wikipedia mentions it, but these equations are for a calorically perfect gas. Thermally perfect or other EOS's (stiffened gas, Mie-Gruneisen, etc.) are also possible, but require more complex iterative solvers in the various states. We've still called them "Sod shock tubes." I think shallow-water equation folks call it a dam breaking problem. $\endgroup$ – tpg2114 Aug 20 '18 at 18:39
  • $\begingroup$ Thank you! I found pretty straight-forward manual also here: www3.nd.edu/~gtryggva/CFD-Course2017/Lecture-10-2017.pdf but I am not sure about what 's' is (on page 3, right upper slide). Also, by the iteration they talk about on page 2, left lower slide, do they mean fixed point iteration? I am pretty new to this area, so I am not sure. $\endgroup$ – Andrej Aug 21 '18 at 10:32
  • $\begingroup$ I meant 'page 2, right lower slide', sorry. $\endgroup$ – Andrej Aug 21 '18 at 10:39
  • $\begingroup$ @Andrej I'll have to read it (on mobile currently), but generally any root finding algorithm is considered iterative. $\endgroup$ – Kyle Kanos Aug 21 '18 at 11:08
  • $\begingroup$ @Andrej The value of the pressure in region 3 depends on itself (i.e., $P_3=f(P_3,...)$), so any iterative method (fixed point, bisection, Brent, etc) will work, that's what they mean with that statement. $\endgroup$ – Kyle Kanos Aug 21 '18 at 11:40

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.