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 65 66### Advection 67 68This problem solves the convection (advection) equation for the total (scalar) energy density, 69transported by the (vector) velocity field. 70 71This is 3D advection given in two formulations based upon the weak form. 72 73State Variables: 74 75 *q = ( rho, U<sub>1</sub>, U<sub>2</sub>, U<sub>3</sub>, E )* 76 77 *rho* - Mass Density 78 79 *U<sub>i</sub>* - Momentum Density , *U<sub>i</sub> = rho ui* 80 81 *E* - Total Energy Density, *E = rho Cv T + rho (u u) / 2 + rho g z* 82 83Advection Equation: 84 85 *dE/dt + div( E _u_ ) = 0* 86 87#### Initial Conditions 88 89Mass Density: 90 Constant mass density of 1.0 91 92Momentum Density: 93 Rotational field in x,y with no momentum in z 94 95Energy Density: 96 Maximum of 1. x0 decreasing linearly to 0. as radial distance increases 97 to 1/8, then 0. everywhere else 98 99#### Boundary Conditions 100 101Mass Density: 102 0.0 flux 103 104Momentum Density: 105 0.0 106 107Energy Density: 108 0.0 flux 109 110### Density Current 111 112This problem solves the full compressible Navier-Stokes equations, using 113operator composition and design of coupled solvers in the context of atmospheric 114modeling. This problem uses the formulation given in Semi-Implicit Formulations 115of the Navier-Stokes Equations: Application to Nonhydrostatic Atmospheric Modeling, 116Giraldo, Restelli, and Lauter (2010). 117 118The 3D compressible Navier-Stokes equations are formulated in conservation form with state 119variables of density, momentum density, and total energy density. 120 121State Variables: 122 123 *q = ( rho, U<sub>1</sub>, U<sub>2</sub>, U<sub>3</sub>, E )* 124 125 *rho* - Mass Density 126 127 *U<sub>i</sub>* - Momentum Density , *U<sub>i</sub> = rho u<sub>i</sub>* 128 129 *E* - Total Energy Density, *E = rho c<sub>v</sub> T + rho (u u) / 2 + rho g z* 130 131Navier-Stokes Equations: 132 133 *drho/dt + div( U ) = 0* 134 135 *dU/dt + div( rho (u x u) + P I<sub>3</sub> ) + rho g khat = div( F<sub>u</sub> )* 136 137 *dE/dt + div( (E + P) u ) = div( F<sub>e</sub> )* 138 139Viscous Stress: 140 141 *F<sub>u</sub> = mu (grad( u ) + grad( u )^T + lambda div ( u ) I<sub>3</sub>)* 142 143Thermal Stress: 144 145 *F<sub>e</sub> = u F<sub>u</sub> + k grad( T )* 146 147Equation of State: 148 149 *P = (gamma - 1) (E - rho (u u) / 2 - rho g z)* 150 151Temperature: 152 153 *T = (E / rho - (u u) / 2 - g z) / c<sub>v</sub>* 154 155Constants: 156 157 *lambda = - 2 / 3*, From Stokes hypothesis 158 159 *mu* , Dynamic viscosity 160 161 *k* , Thermal conductivity 162 163 *c<sub>v</sub>* , Specific heat, constant volume 164 165 *c<sub>p</sub>* , Specific heat, constant pressure 166 167 *g* , Gravity 168 169 *gamma = c<sub>p</sub> / c<sub>v</sub>*, Specific heat ratio 170 171#### Initial Conditions 172 173Potential Temperature: 174 175 *theta = thetabar + deltatheta* 176 177 *thetabar = theta0 exp( N * * 2 z / g )* 178 179 *deltatheta = 180 r <= rc : theta0(1 + cos(pi r)) / 2 181 r > rc : 0* 182 183 *r = sqrt( (x - xc) * * 2 + (y - yc) * * 2 + (z - zc) * * 2 )* 184 with *(xc,yc,zc)* center of domain 185 186Exner Pressure: 187 188 *Pi = Pibar + deltaPi* 189 190 *Pibar = g * * 2 (exp( - N * * 2 z / g ) - 1) / (cp theta0 N * * 2)* 191 192 *deltaPi = 0* (hydrostatic balance) 193 194Velocity/Momentum Density: 195 196 *U<sub>i</sub> = u<sub>i</sub> = 0* 197 198Conversion to Conserved Variables: 199 200 *rho = P0 Pi**(c<sub>v</sub>/R<sub>d</sub>) / (R<sub>d</sub> theta)* 201 202 *E = rho (c<sub>v</sub> theta Pi + (u u)/2 + g z)* 203 204Constants: 205 206 *theta0* , Potential temperature constant 207 208 *thetaC* , Potential temperature perturbation 209 210 *P0* , Pressure at the surface 211 212 *N* , Brunt-Vaisala frequency 213 214 *c<sub>v</sub>* , Specific heat, constant volume 215 216 *c<sub>p</sub>* , Specific heat, constant pressure 217 218 *R<sub>d</sub>* = c<sub>p</sub> - c<sub>v</sub>, Specific heat difference 219 220 *g* , Gravity 221 222 *r<sub>c</sub>* , Characteristic radius of thermal bubble 223 224 *l<sub>x</sub>* , Characteristic length scale of domain in x 225 226 *l<sub>y</sub>* , Characteristic length scale of domain in y 227 228 *l<sub>z</sub>* , Characteristic length scale of domain in z 229 230 231#### Boundary Conditions 232 233Mass Density: 234 0.0 flux 235 236Momentum Density: 237 0.0 238 239Energy Density: 240 0.0 flux 241 242### Time Discretization 243 244For all different problems, the time integration is performed with an explicit formulation, therefore 245it can be subject to numerical instability, if run for large times or with large time steps. 246 247### Space Discretization 248 249The geometric factors and coordinate transformations required for the integration of the weak form 250are described in the file [`common.h`](common.h) 251