hub.darcs.net :: basilisk -> basilisk -> annotate -> src/compressible.h.page

Basilisk source code (http://basilisk.fr/src/)

Explosion test case for compressible gas dynamics Jose-Maria Fullana <fullana@ida.upmc.fr> Mon Mar 3 18:23:04 UTC 2014
1
2
3
/**
# Compressible gas dynamics
Documentation fixes for compressible solver Stephane Popinet <popinet@basilisk.fr> Mon Mar 3 19:41:04 UTC 2014
4
5
The Euler system of conservation laws for a compressible gas can be
written
Explosion test case for compressible gas dynamics Jose-Maria Fullana <fullana@ida.upmc.fr> Mon Mar 3 18:23:04 UTC 2014
6
7
8

$$
\partial_t\left(\begin{array}{c}
Documentation fixes for compressible solver Stephane Popinet <popinet@basilisk.fr> Mon Mar 3 19:41:04 UTC 2014
9
10
11
12
\rho \\
E    \\
w_x  \\
w_y  \\
Explosion test case for compressible gas dynamics Jose-Maria Fullana <fullana@ida.upmc.fr> Mon Mar 3 18:23:04 UTC 2014
13
14
15
\end{array}\right) 
+ 
\nabla_x \cdot\left(\begin{array}{c}
Documentation fixes for compressible solver Stephane Popinet <popinet@basilisk.fr> Mon Mar 3 19:41:04 UTC 2014
16
17
18
19
w_x \\
\frac{w_x}{\rho} ( E + p ) \\
\frac{w_x^2}{\rho} + p \\
\frac{w_y w_x}{\rho} \\
Explosion test case for compressible gas dynamics Jose-Maria Fullana <fullana@ida.upmc.fr> Mon Mar 3 18:23:04 UTC 2014
20
21
22
\end{array}\right)
+ 
\nabla_y \cdot\left(\begin{array}{c}
Documentation fixes for compressible solver Stephane Popinet <popinet@basilisk.fr> Mon Mar 3 19:41:04 UTC 2014
23
24
25
26
w_y \\
\frac{w_y}{\rho} ( E + p ) \\
\frac{w_y w_x}{\rho} \\
\frac{w_y^2}{\rho} + p \\
Explosion test case for compressible gas dynamics Jose-Maria Fullana <fullana@ida.upmc.fr> Mon Mar 3 18:23:04 UTC 2014
27
28
29
30
31
\end{array}\right)
= 0
$$

with $\rho$ the gas density, $E$ the total energy, $\mathbf{w}$ the
Documentation fixes for compressible solver Stephane Popinet <popinet@basilisk.fr> Mon Mar 3 19:41:04 UTC 2014
32
gas momentum and $p$ the pressure given by the equation of state
Explosion test case for compressible gas dynamics Jose-Maria Fullana <fullana@ida.upmc.fr> Mon Mar 3 18:23:04 UTC 2014
33
34
35
36
37

$$
p = (\gamma - 1)(E - \rho\mathbf{u}^2/2)
$$
Documentation fixes for compressible solver Stephane Popinet <popinet@basilisk.fr> Mon Mar 3 19:41:04 UTC 2014
38
39
with $\gamma$ the polytropic exponent. This system can be solved using
the generic solver for systems of conservation laws. */
Explosion test case for compressible gas dynamics Jose-Maria Fullana <fullana@ida.upmc.fr> Mon Mar 3 18:23:04 UTC 2014
40
41
42
43
44

#include "conservation.h"

/**
The conserved scalars are the gas density $\rho$ and the total energy
Documentation fixes for compressible solver Stephane Popinet <popinet@basilisk.fr> Mon Mar 3 19:41:04 UTC 2014
45
46
$E$. The only conserved vector is the momentum $\mathbf{w}$.  The constant
$\gamma$ is represented by *gammao* here, with a default value of 1.4. */
Explosion test case for compressible gas dynamics Jose-Maria Fullana <fullana@ida.upmc.fr> Mon Mar 3 18:23:04 UTC 2014
47
48
49
50
51
52
53
54
55
56
57
58

scalar rho[], E[];
vector w[];
scalar * scalars = {rho, E};
vector * vectors = {w};
double gammao = 1.4 ;

/**
The system is entirely defined by the `flux()` function called by the
generic solver for conservation laws. The parameter passed to the
function is the array `s` which contains the state variables for each
conserved field, in the order of their definition above (i.e. scalars
Documentation fixes for compressible solver Stephane Popinet <popinet@basilisk.fr> Mon Mar 3 19:41:04 UTC 2014
59
then vectors). */
Explosion test case for compressible gas dynamics Jose-Maria Fullana <fullana@ida.upmc.fr> Mon Mar 3 18:23:04 UTC 2014
60
61
62

void flux (const double * s, double * f, double e[2])
{
Documentation fixes for compressible solver Stephane Popinet <popinet@basilisk.fr> Mon Mar 3 19:41:04 UTC 2014
Proper indentation of documentation comments Stephane Popinet <popinet@basilisk.fr> Mon Mar 31 10:09:39 UTC 2014
64
65
66
67
  /**
  We first recover each value ($\rho$, $E$, $w_x$ and $w_y$) and then
  compute the corresponding fluxes (`f[0]`, `f[1]`, `f[2]` and
  `f[3]`). */
Documentation fixes for compressible solver Stephane Popinet <popinet@basilisk.fr> Mon Mar 3 19:41:04 UTC 2014
Compressible scheme in 3D and explosion test case Stephane Popinet <popinet@basilisk.fr> Sun Apr 12 16:06:42 UTC 2015
69
70
71
72
  double rho = s[0], E = s[1], wn = s[2], w2 = 0.;
  for (int i = 2; i < 2 + dimension; i++)
    w2 += sq(s[i]);
  double un = wn/rho, p = (gammao - 1.)*(E - 0.5*w2/rho);
Explosion test case for compressible gas dynamics Jose-Maria Fullana <fullana@ida.upmc.fr> Mon Mar 3 18:23:04 UTC 2014
73
74
75
76

  f[0] = wn;
  f[1] = un*(E + p);
  f[2] = un*wn + p;
Compressible scheme in 3D and explosion test case Stephane Popinet <popinet@basilisk.fr> Sun Apr 12 16:06:42 UTC 2015
77
78
  for (int i = 3; i < 2 + dimension; i++)
    f[i] = un*s[i];
Explosion test case for compressible gas dynamics Jose-Maria Fullana <fullana@ida.upmc.fr> Mon Mar 3 18:23:04 UTC 2014
Proper indentation of documentation comments Stephane Popinet <popinet@basilisk.fr> Mon Mar 31 10:09:39 UTC 2014
80
81
82
  /**
  The minimum and maximum eigenvalues for the Euler system are the
  characteristic speeds $u \pm \sqrt(\gamma p / \rho)$. */
Documentation fixes for compressible solver Stephane Popinet <popinet@basilisk.fr> Mon Mar 3 19:41:04 UTC 2014
Explosion test case for compressible gas dynamics Jose-Maria Fullana <fullana@ida.upmc.fr> Mon Mar 3 18:23:04 UTC 2014
84
85
86
87
  double c = sqrt(gammao*p/rho);
  e[0] = un - c; // min
  e[1] = un + c; // max
}