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