xref: /libCEED/examples/fluids/README.md (revision 565b71c80ca6063971908314385fa86604050bb3)
1## libCEED: Navier-Stokes Example
2
3This page provides a description of the Navier-Stokes example for the libCEED library, based on PETSc.
4
5The Navier-Stokes problem solves the compressible Navier-Stokes equations in three dimensions using an
6explicit time integration. The state variables are mass density, momentum density, and energy density.
7
8The main Navier-Stokes solver for libCEED is defined in [`navierstokes.c`](navierstokes.c)
9with different problem definitions according to the application of interest.
10
11Build by using
12
13`make`
14
15and run with
16
17`./navierstokes`
18
19Available runtime options are:
20
21|  Option                  | Meaning                                                                                         |
22| :----------------------- | :-----------------------------------------------------------------------------------------------|
23| `-ceed`                  | CEED resource specifier                                                                         |
24| `-test`                  | Run in test mode                                                                                |
25| `-problem`               | Problem to solve (`advection`, `advection2d`, or `density_current`)                             |
26| `-stab`                  | Stabilization method                                                                            |
27| `-implicit`              | Use implicit time integartor formulation                                                        |
28| `-bc_wall`               | Use wall boundary conditions on this list of faces                                              |
29| `-bc_slip_x`             | Use slip boundary conditions, for the x component, on this list of faces                        |
30| `-bc_slip_y`             | Use slip boundary conditions, for the y component, on this list of faces                        |
31| `-bc_slip_z`             | Use slip boundary conditions, for the z component, on this list of faces                        |
32| `-viz_refine`            | Use regular refinement for visualization                                                        |
33| `-degree`                | Polynomial degree of tensor product basis (must be >= 1)                                        |
34| `-units_meter`           | 1 meter in scaled length units                                                                  |
35| `-units_second`          | 1 second in scaled time units                                                                   |
36| `-units_kilogram`        | 1 kilogram in scaled mass units                                                                 |
37| `-units_Kelvin`          | 1 Kelvin in scaled temperature units                                                            |
38| `-theta0`                | Reference potential temperature                                                                 |
39| `-thetaC`                | Perturbation of potential temperature                                                           |
40| `-P0`                    | Atmospheric pressure                                                                            |
41| `-N`                     | Brunt-Vaisala frequency                                                                         |
42| `-cv`                    | Heat capacity at constant volume                                                                |
43| `-cp`                    | Heat capacity at constant pressure                                                              |
44| `-g`                     | Gravitational acceleration                                                                      |
45| `-lambda`                | Stokes hypothesis second viscosity coefficient                                                  |
46| `-mu`                    | Shear dynamic viscosity coefficient                                                             |
47| `-k`                     | Thermal conductivity                                                                            |
48| `-CtauS`                 | Scale coefficient for stabilization tau (nondimensional)                                        |
49| `-strong_form`           | Strong (1) or weak/integrated by parts (0) advection residual                                   |
50| `-lx`                    | Length scale in x direction                                                                     |
51| `-ly`                    | Length scale in y direction                                                                     |
52| `-lz`                    | Length scale in z direction                                                                     |
53| `-rc`                    | Characteristic radius of thermal bubble                                                         |
54| `-resx`                  | Resolution in x                                                                                 |
55| `-resy`                  | Resolution in y                                                                                 |
56| `-resz`                  | Resolution in z                                                                                 |
57| `-center`                | Location of bubble center                                                                       |
58| `-dc_axis`               | Axis of density current cylindrical anomaly, or {0,0,0} for spherically symmetric               |
59| `-output_freq`           | Frequency of output, in number of steps                                                         |
60| `-continue`              | Continue from previous solution                                                                 |
61| `-degree`                | Polynomial degree of tensor product basis                                                       |
62| `-qextra`                | Number of extra quadrature points                                                               |
63| `-of`                    | Output folder                                                                                   |
64
65For the case of a square/cubic mesh, the list of face indices to be used with `-bc_wall` and/or `-bc_slip_x`,
66`-bc_slip_y`, and `-bc_slip_z` are:
67
68* 2D:
69  - faceMarkerBottom = 1;
70  - faceMarkerRight  = 2;
71  - faceMarkerTop    = 3;
72  - faceMarkerLeft   = 4;
73*  3D:
74  - faceMarkerBottom = 1;
75  - faceMarkerTop    = 2;
76  - faceMarkerFront  = 3;
77  - faceMarkerBack   = 4;
78  - faceMarkerRight  = 5;
79  - faceMarkerLeft   = 6;
80
81### Advection
82
83This problem solves the convection (advection) equation for the total (scalar) energy density,
84transported by the (vector) velocity field.
85
86This is 3D advection given in two formulations based upon the weak form.
87
88State Variables:
89
90   *q = ( rho, U<sub>1</sub>, U<sub>2</sub>, U<sub>3</sub>, E )*
91
92   *rho* - Mass Density
93
94   *U<sub>i</sub>*  - Momentum Density    ,   *U<sub>i</sub> = rho ui*
95
96   *E*   - Total Energy Density,   *E  = rho Cv T + rho (u u) / 2 + rho g z*
97
98Advection Equation:
99
100   *dE/dt + div( E _u_ ) = 0*
101
102#### Initial Conditions
103
104Mass Density:
105    Constant mass density of 1.0
106
107Momentum Density:
108    Rotational field in x,y with no momentum in z
109
110Energy Density:
111    Maximum of 1. x0 decreasing linearly to 0. as radial distance increases
112    to 1/8, then 0. everywhere else
113
114#### Boundary Conditions
115
116Mass Density:
117    0.0 flux
118
119Momentum Density:
120    0.0
121
122Energy Density:
123    0.0 flux
124
125### Density Current
126
127This problem solves the full compressible Navier-Stokes equations, using
128operator composition and design of coupled solvers in the context of atmospheric
129modeling. This problem uses the formulation given in Semi-Implicit Formulations
130of the Navier-Stokes Equations: Application to Nonhydrostatic Atmospheric Modeling,
131Giraldo, Restelli, and Lauter (2010).
132
133The 3D compressible Navier-Stokes equations are formulated in conservation form with state
134variables of density, momentum density, and total energy density.
135
136State Variables:
137
138   *q = ( rho, U<sub>1</sub>, U<sub>2</sub>, U<sub>3</sub>, E )*
139
140   *rho* - Mass Density
141
142   *U<sub>i</sub>*  - Momentum Density   ,  *U<sub>i</sub> = rho u<sub>i</sub>*
143
144   *E*   - Total Energy Density,  *E  = rho c<sub>v</sub> T + rho (u u) / 2 + rho g z*
145
146Navier-Stokes Equations:
147
148   *drho/dt + div( U )                               = 0*
149
150   *dU/dt   + div( rho (u x u) + P I<sub>3</sub> ) + rho g khat = div( F<sub>u</sub> )*
151
152   *dE/dt   + div( (E + P) u )                       = div( F<sub>e</sub> )*
153
154Viscous Stress:
155
156   *F<sub>u</sub> = mu (grad( u ) + grad( u )^T + lambda div ( u ) I<sub>3</sub>)*
157
158Thermal Stress:
159
160   *F<sub>e</sub> = u F<sub>u</sub> + k grad( T )*
161
162Equation of State:
163
164   *P = (gamma - 1) (E - rho (u u) / 2 - rho g z)*
165
166Temperature:
167
168   *T = (E / rho - (u u) / 2 - g z) / c<sub>v</sub>*
169
170Constants:
171
172   *lambda = - 2 / 3*,  From Stokes hypothesis
173
174   *mu*              ,  Dynamic viscosity
175
176   *k*               ,  Thermal conductivity
177
178   *c<sub>v</sub>*              ,  Specific heat, constant volume
179
180   *c<sub>p</sub>*              ,  Specific heat, constant pressure
181
182   *g*               ,  Gravity
183
184   *gamma  = c<sub>p</sub> / c<sub>v</sub>*,  Specific heat ratio
185
186#### Initial Conditions
187
188Potential Temperature:
189
190   *theta = thetabar + deltatheta*
191
192   *thetabar   = theta0 exp( N * * 2 z / g )*
193
194   *deltatheta =
195        r <= rc : theta0(1 + cos(pi r)) / 2
196        r > rc : 0*
197
198   *r        = sqrt( (x - xc) * * 2 + (y - yc) * * 2 + (z - zc) * * 2 )*
199    with *(xc,yc,zc)* center of domain
200
201Exner Pressure:
202
203   *Pi = Pibar + deltaPi*
204
205   *Pibar      = g * * 2 (exp( - N * * 2 z / g ) - 1) / (cp theta0 N * * 2)*
206
207   *deltaPi    = 0* (hydrostatic balance)
208
209Velocity/Momentum Density:
210
211   *U<sub>i</sub> = u<sub>i</sub> = 0*
212
213Conversion to Conserved Variables:
214
215   *rho = P0 Pi**(c<sub>v</sub>/R<sub>d</sub>) / (R<sub>d</sub> theta)*
216
217   *E   = rho (c<sub>v</sub> theta Pi + (u u)/2 + g z)*
218
219Constants:
220
221   *theta0*          ,  Potential temperature constant
222
223   *thetaC*          ,  Potential temperature perturbation
224
225   *P0*              ,  Pressure at the surface
226
227   *N*               ,  Brunt-Vaisala frequency
228
229   *c<sub>v</sub>*              ,  Specific heat, constant volume
230
231   *c<sub>p</sub>*              ,  Specific heat, constant pressure
232
233   *R<sub>d</sub>*     = c<sub>p</sub> - c<sub>v</sub>,  Specific heat difference
234
235   *g*               ,  Gravity
236
237   *r<sub>c</sub>*              ,  Characteristic radius of thermal bubble
238
239   *l<sub>x</sub>*              ,  Characteristic length scale of domain in x
240
241   *l<sub>y</sub>*              ,  Characteristic length scale of domain in y
242
243   *l<sub>z</sub>*              ,  Characteristic length scale of domain in z
244
245
246#### Boundary Conditions
247
248Mass Density:
249    0.0 flux
250
251Momentum Density:
252    0.0
253
254Energy Density:
255    0.0 flux
256
257### Time Discretization
258
259For all different problems, the time integration is performed with an explicit
260or implicit formulation.
261
262### Space Discretization
263
264The geometric factors and coordinate transformations required for the integration of the weak form
265are described in the file [`common.h`](common.h)
266