Runtime options

Common Options

The Navier-Stokes HONEE app is controlled via command-line options. The following options are common among all problem types:

Table 1 Common Runtime Options

Option

Description

Default value

-ceed

CEED resource specifier

/cpu/self/opt/blocked

-problem

Problem to solve (advection, density_current, euler_vortex, shocktube, blasius, channel, gaussian_wave, and taylor_green)

density_current

-implicit

Use implicit time integrator formulation

-degree

Polynomial degree of tensor product basis (must be >= 1)

1

-q_extra

Number of extra quadrature points

0

-ts_monitor_solution

PETSc output format, such as cgns:output-%d.cgns (requires PETSc --download-cgns)

-ts_monitor_solution_interval

Number of time steps between visualization output frames.

1

-viewer_cgns_batch_size

Number of frames written per CGNS file if the CGNS file name includes a format specifier (%d).

20

-checkpoint_interval

Number of steps between writing binary checkpoints. 0 has no output, -1 outputs final state only

10

-checkpoint_vtk

Checkpoints include VTK (*.vtu) files for visualization. Consider -ts_monitor_solutioninstead.

false

-viz_refine

Use regular refinement for VTK visualization

0

-output_dir

Output directory for binary checkpoints and VTK files (if enabled).

.

-output_add_stepnum2bin

Whether to add step numbers to output binary files

false

-continue

Continue from previous solution (input is step number of previous solution)

0

-continue_filename

Path to solution binary file from which to continue from

[output_dir]/ns-solution.bin

-continue_time_filename

Path to time stamp binary file (only for legacy checkpoints)

[output_dir]/ns-time.bin

-bc_wall

Use wall boundary conditions on this list of faces

-wall_comps

An array of constrained component numbers for wall BCs

-bc_slip

Use weak slip boundary condition on this list of faces

-bc_symmetry_x

Use symmetry boundary conditions, for the x component, on this list of faces

-bc_symmetry_y

Use symmetry boundary conditions, for the y component, on this list of faces

-bc_symmetry_z

Use symmetry boundary conditions, for the z component, on this list of faces

-bc_inflow

Use inflow boundary conditions on this list of faces

-bc_outflow

Use outflow boundary conditions on this list of faces

-bc_freestream

Use freestream boundary conditions on this list of faces

-ts_monitor_turbulence_spanstats_collect_interval

Number of timesteps between statistics collection

1

-ts_monitor_turbulence_spanstats_viewer

Sets the PetscViewer for the statistics file writing, such as cgns:output-%d.cgns (requires PETSc --download-cgns). Also turns the statistics collection on.

-ts_monitor_turbulence_spanstats_viewer_interval

Number of timesteps between statistics file writing (-1 means only at end of run)

-1

-ts_monitor_turbulence_spanstats_viewer_cgns_batch_size

Number of frames written per CGNS file if the CGNS file name includes a format specifier (%d).

20

-ts_monitor_wall_force

Viewer for the force on each no-slip wall, e.g., ascii:force.csv:ascii_csv to write a CSV file.

-mesh_transform

Transform the mesh, usually for an initial box mesh.

none

-snes_view

View PETSc SNES nonlinear solver configuration

-log_view

View PETSc performance log

-help

View comprehensive information about run-time options

-test_type

Run in test mode and specify whether solution (solver) or turbulent statistics (turb_spanstats) output should be verified

none

-compare_final_state_atol

Test absolute tolerance

1E-11

-compare_final_state_filename

Test filename

For the case of a square/cubic mesh, the list of face indices to be used with -bc_wall, bc_inflow, bc_outflow, bc_freestream and/or -bc_symmetry_x, -bc_symmetry_y, and -bc_symmetry_z are:

Table 2 2D Face ID Labels

PETSc Face Name

Cartesian direction

Face ID

faceMarkerBottom

-z

1

faceMarkerRight

+x

2

faceMarkerTop

+z

3

faceMarkerLeft

-x

4

Table 3 3D Face ID Labels

PETSc Face Name

Cartesian direction

Face ID

faceMarkerBottom

-z

1

faceMarkerTop

+z

2

faceMarkerFront

-y

3

faceMarkerBack

+y

4

faceMarkerRight

+x

5

faceMarkerLeft

-x

6

Boundary conditions

Boundary conditions for compressible viscous flows are notoriously tricky. Here we offer some recommendations.

Inflow

If in a region where the flow velocity is known (e.g., away from viscous walls), use bc_freestream, which solves a Riemann problem and can handle inflow and outflow (simultaneously and dynamically). It is stable and the least reflective boundary condition for acoustics.

If near a viscous wall, you may want a specified inflow profile. Use bc_inflow and see Blasius boundary layer and discussion of synthetic turbulence generation for ways to analytically generate developed inflow profiles. These conditions may be either weak or strong, with the latter specifying velocity and temperature as essential boundary conditions and evaluating a boundary integral for the mass flux. The strong approach gives sharper resolution of velocity structures. We have described the primitive variable formulation here; the conservative variants are similar, but not equivalent.

Outflow

If you know the complete exterior state, bc_freestream is the least reflective boundary condition, but is disruptive to viscous flow structures. If thermal anomalies must exit the domain, the Riemann solver must resolve the contact wave to avoid reflections. The default Riemann solver, HLLC, is sufficient in this regard while the simpler HLL converts thermal structures exiting the domain into grid-scale reflecting acoustics.

If acoustic reflections are not a concern and/or the flow is impacted by walls or interior structures that you wish to resolve to near the boundary, choose bc_outflow. This condition (with default outflow_type: riemann) is stable for both inflow and outflow, so can be used in areas that have recirculation and lateral boundaries in which the flow fluctuates.

The simpler bc_outflow variant, outflow_type: pressure, requires that the flow be a strict outflow (or the problem becomes ill-posed and the solver will diverge). In our experience, riemann is slightly less reflective but produces similar flows in cases of strict outflow. The pressure variant is retained to facilitate comparison with other codes, such as PHASTA-C, but we recommend riemann for general use.

Periodicity

PETSc provides two ways to specify periodicity:

  1. Topological periodicity, in which the donor and receiver dofs are the same, obtained using:

dm_plex:
  shape: box
  box_faces: 10,12,4
  box_bd: none,none,periodic

The coordinates for such cases are stored as a new field with special cell-based indexing to enable wrapping through the boundary. This choice of coordinates prevents evaluating boundary integrals that cross the periodicity, such as for the outflow Riemann problem in the presence of spanwise periodicity.

  1. Isoperiodicity, in which the donor and receiver dofs are distinct in local vectors. This is obtained using zbox, as in:

dm_plex:
  shape: zbox
  box_faces: 10,12,4
  box_bd: none,none,periodic

Isoperiodicity enables standard boundary integrals, and is recommended for general use. At the time of this writing, it only supports one direction of periodicity. The zbox method uses Z-ordering to construct the mesh in parallel and provide an adequate initial partition, which makes it higher performance and avoids needing a partitioning package.

Advection-Diffusion

There is a reduced mode for pure advection, which holds density \(\rho\) and momentum density \(\rho \bm u\) constant while advecting “total energy density” \(E\). The advection problems can be run in both 2D and 3D, based on the DM defined for the problem. The following additional command-line options are available:

Table 4 Advection Runtime Options

Option

Description

Default value

Unit

-units_meter

1 meter in scaled length units

1E-2

-units_second

1 second in scaled time units

1E-2

-units_kilogram

1 kilogram in scaled mass units

1E-6

-strong_form

Strong (1) or weak/integrated by parts (0) residual

0

-stab

Stabilization method (none, su, or supg)

none

-stab_tau

Formulation for \(\tau\) in stabilization (ctau, advdiff_shakib)

ctau

-Ctau_t

Scaling factor on the temporal portion of the \(\tau\) formulation

-Ctau_a

Scaling factor on the advection portion of the \(\tau\) formulation

\(P^2\)

-Ctau_d

Scaling factor on the diffusion portion of the \(\tau\) formulation

\(P^4\)

-CtauS

Scale coefficient for stabilization tau (nondimensional)

0

-wind_type

Wind type in Advection (rotation, translation, boundary_layer)

rotation

-wind_translation

Constant wind vector when -wind_type translation

1,0,0

-diffusion_coeff

Diffusion coefficient

0

-E_wind

Total energy of inflow wind when -wind_type translation

1E6

J

-advection_ic_type

Initial condition type, (sphere, cylinder, cosine_hill, skew, wave, boundary_layer)

sphere

-advection_ic_bubble_rc

For sphere or cylinder IC, characteristic radius of thermal bubble

1000

m

-advection_ic_bubble_continuity

For sphere or cylinder IC, different shapes of bubble, (smooth, back_sharp, thick, cosine)

smooth

-advection_ic_wave_type

For wave IC, the wave form used for -advection_ic_type wave (sine, square)

sine

-advection_ic_wave_frequency

For wave IC, frequency of the wave

\(2\pi\)

1/s

-advection_ic_wave_phase

For wave IC, phase angle of the wave

\(2\pi\)

-advection_ic_bl_height_factor

For boundary_layer IC, sets the height of the linear boundary layer initial condition in proportion to the domain height

\(1\)

For 3D advection, an example of the rotation mode can be run with:

./navierstokes -problem advection -dm_plex_box_faces 10,10,10 -dm_plex_dim 3 -dm_plex_box_lower 0,0,0 -dm_plex_box_upper 8000,8000,8000 -bc_wall 1,2,3,4,5,6 -wall_comps 4 -wind_type rotation -implicit -stab su

and the translation mode with:

./navierstokes -problem advection -dm_plex_box_faces 10,10,10 -dm_plex_dim 3 -dm_plex_box_lower 0,0,0 -dm_plex_box_upper 8000,8000,8000 -wind_type translation -wind_translation .5,-1,0 -bc_inflow 1,2,3,4,5,6

For 2D advection, an example of the rotation mode can be run with:

./navierstokes -problem advection -dm_plex_box_faces 20,20 -dm_plex_box_lower 0,0 -dm_plex_box_upper 1000,1000 -bc_wall 1,2,3,4 -wall_comps 4 -wind_type rotation -implicit -stab supg

and the translation mode with:

./navierstokes -problem advection -dm_plex_box_faces 20,20 -dm_plex_box_lower 0,0 -dm_plex_box_upper 1000,1000 -units_meter 1e-4 -wind_type translation -wind_translation 1,-.5 -bc_inflow 1,2,3,4

Note the lengths in -dm_plex_box_upper are given in meters, and will be nondimensionalized according to -units_meter.

Inviscid Ideal Gas

Isentropic Euler vortex

For the Isentropic Vortex problem, the following additional command-line options are available:

Table 5 Isentropic Vortex Runtime Options

Option

Description

Default value

Unit

-center

Location of vortex center

(lx,ly,lz)/2

(m,m,m)

-units_meter

1 meter in scaled length units

1E-2

-units_second

1 second in scaled time units

1E-2

-mean_velocity

Background velocity vector

(1,1,0)

-vortex_strength

Strength of vortex < 10

5

-c_tau

Stabilization constant

0.5

This problem can be run with:

./navierstokes -problem euler_vortex -dm_plex_box_faces 20,20,1 -dm_plex_box_lower 0,0,0 -dm_plex_box_upper 1000,1000,50 -dm_plex_dim 3 -bc_inflow 4,6 -bc_outflow 3,5 -bc_symmetry_z 1,2 -mean_velocity .5,-.8,0.

Sod shock tube

For the Shock Tube problem, the following additional command-line options are available:

Table 6 Shock Tube Runtime Options

Option

Description

Default value

Unit

-units_meter

1 meter in scaled length units

1E-2

-units_second

1 second in scaled time units

1E-2

-yzb

Use YZB discontinuity capturing

none

-stab

Stabilization method (none, su, or supg)

none

This problem can be run with:

./navierstokes -problem shocktube -yzb -stab su -bc_symmetry_z 3,4 -bc_symmetry_y 1,2 -bc_wall 5,6 -dm_plex_dim 3 -dm_plex_box_lower 0,0,0 -dm_plex_box_upper 1000,100,100 -dm_plex_box_faces 200,1,1 -units_second 0.1

Newtonian viscosity, Ideal Gas

For the Density Current, Channel, and Blasius problems, the following common command-line options are available:

Table 7 Newtonian Ideal Gas problems Runtime Options

Option

Description

Default value

Unit

-units_meter

1 meter in scaled length units

1

-units_second

1 second in scaled time units

1

-units_kilogram

1 kilogram in scaled mass units

1

-units_Kelvin

1 Kelvin in scaled temperature units

1

-stab

Stabilization method (none, su, or supg)

none

-c_tau

Stabilization constant, \(c_\tau\)

0.5

-Ctau_t

Stabilization time constant, \(C_t\)

1.0

-Ctau_v

Stabilization viscous constant, \(C_v\)

36, 60, 128 for degree = 1, 2, 3

-Ctau_C

Stabilization continuity constant, \(C_c\)

1.0

-Ctau_M

Stabilization momentum constant, \(C_m\)

1.0

-Ctau_E

Stabilization energy constant, \(C_E\)

1.0

-div_diff_flux_projection_method

Method used to calculate divergence of diffusive flux projection (none, direct, or indirect)

none

-div_diff_flux_projection_ksp*

Control the KSP object for the projection of the divergence of diffusive flux

N/A

-cv

Heat capacity at constant volume

717

J/(kg K)

-cp

Heat capacity at constant pressure

1004

J/(kg K)

-gravity

Gravitational acceleration vector

0,0,0

m/s^2

-lambda

Stokes hypothesis second viscosity coefficient

-2/3

-mu

Shear dynamic viscosity coefficient

1.8e-5

Pa s

-k

Thermal conductivity

0.02638

W/(m K)

-newtonian_unit_tests

Developer option to test properties

false

boolean

-state_var

State variables to solve solution with. conservative (\(\rho, \rho \bm{u}, \rho e\)), primitive (\(P, \bm{u}, T\)), or entropy (\(\frac{\gamma - s}{\gamma - 1} - \frac{\rho}{P} (e - c_v T),\ \frac{\rho}{P} \bm{u},\ -\frac{\rho}{P}\)) where \(s = \ln(P\rho^{-\gamma})\)

conservative

string

-idl_decay_time

Characteristic timescale of the pressure deviance decay. The timestep is good starting point

-1 (disabled)

s

-idl_start

Start of IDL in the x direction

0

m

-idl_length

Length of IDL in the positive x direction

0

m

-idl_pressure

Pressure used for IDL reference pressure

-reference_pressure

Pa

-sgs_model_type

Type of subgrid stress model to use. Currently only data_driven is available

none

string

-sgs_model_dd_leakyrelu_alpha

Slope parameter for Leaky ReLU activation function. 0 corresponds to normal ReLU

0

-sgs_model_dd_parameter_dir

Path to directory with data-driven model parameters (weights, biases, etc.)

./dd_sgs_parameters

string

-sgs_model_dd_model_implementation

Which computational implementation to use for SGS DD model (fused, sequential_ceed, sequential_torch)

fused

string

-sgs_model_dd_torch_model_path

Path to the PyTorch *.pt file containing the DD inference model

string

-sgs_model_dd_torch_model_device

What hardware to perform the model inference on (cpu, cuda, hip, xpu)

Default matches the libCEED backend

string

-diff_filter_monitor

Enable differential filter TSMonitor

false

boolean

-diff_filter_grid_based_width

Use filter width based on the grid size

false

boolean

-diff_filter_width_scaling

Anisotropic scaling for filter width in wall-aligned coordinates (snz)

1,1,1

m

-diff_filter_kernel_scaling

Scaling to make differential kernel size equivalent to other filter kernels

0.1

m^2

-diff_filter_wall_damping_function

Damping function to use at the wall for anisotropic filtering (none, van_driest)

none

string

-diff_filter_wall_damping_constant

Constant for the wall-damping function. \(A^+\) for van_driest damping function.

25

-diff_filter_friction_length

Friction length associated with the flow, \(\delta_\nu\). Used in wall-damping functions

0

m

-sgs_train_enable

Whether to enable in situ training of data-driven SGS model. Require building with SmartRedis.

false

boolean

-sgs_train_write_data_interval

Number of timesteps between writing training data into SmartRedis database

1

-sgs_train_overwrite_data

Whether new training data should overwrite old data on database

true

boolean

-sgs_train_filter_widths

List of scalar values for different filter widths to calculate for training data

m

-smartsim_collocated_num_ranks

Number of MPI ranks associated with each collocated database (i.e. ranks per node)

1

Gaussian Wave

The Gaussian wave problem has the following command-line options in addition to the Newtonian Ideal Gas options:

Table 8 Gaussian Wave Runtime Options

Option

Description

Default value

Unit

-freestream_riemann

Riemann solver for boundaries (HLL or HLLC)

hllc

-freestream_velocity

Freestream velocity vector

0,0,0

m/s

-freestream_temperature

Freestream temperature

288

K

-freestream_pressure

Freestream pressure

1.01e5

Pa

-epicenter

Coordinates of center of perturbation

0,0,0

m

-amplitude

Amplitude of the perturbation

0.1

-width

Width parameter of the perturbation

0.002

m

This problem can be run with the examples/gaussianwave.yaml file via:

./build/navierstokes -options_file examples/gaussianwave.yaml
problem: gaussian_wave
mu: 0 # Effectively solving Euler momentum equations

dm_plex_box_faces: 40,40,1
dm_plex_box_upper: 1,1,0.025
dm_plex_box_lower: 0,0,0
dm_plex_dim: 3
bc_freestream: 4,6,3,5
bc_symmetry_z: 1,2

reference:
  temperature: 0.25
  pressure: 71.75
freestream:
  # riemann: hll # causes thermal bubble to reflect acoustic waves from boundary
  velocity: 2,2,0

epicenter: 0.33,0.75,0
amplitude: 2
width: 0.05

ts:
  adapt_type: none
  max_steps: 100
  dt: 2e-3
  type: alpha
  #monitor_solution: cgns:nwave.cgns
  #monitor_solution_interval: 10

implicit: true
stab: supg
state_var: primitive

snes_rtol: 1e-4
ksp_rtol: 1e-2
snes_lag_jacobian: 20
snes_lag_jacobian_persists:

## Demonstrate acoustic wave dissipation using an internal damping layer
# idl:
#   decay_time: 2e-3
#   start: 0
#   length: .25

Vortex Shedding - Flow past Cylinder

The vortex shedding, flow past cylinder problem has the following command-line options in addition to the Newtonian Ideal Gas options:

Table 9 Vortex Shedding Runtime Options

Option

Description

Default value

Unit

-freestream_velocity

Freestream velocity vector

0,0,0

m/s

-freestream_temperature

Freestream temperature

288

K

-freestream_pressure

Freestream pressure

1.01e5

Pa

The initial condition is taken from -reference_temperature and -reference_pressure. To run this problem, first generate a mesh:

$ make -C examples/meshes

Then run by building the executable and running:

$ make -j
$ mpiexec -n 6 build/navierstokes -options_file examples/vortexshedding.yaml -{ts,snes}_monitor_

The vortex shedding period is roughly 5.6 and this problem runs until time 100 (2000 time steps). The above run writes a file named force.csv (see ts_monitor_wall_force in examples/vortexshedding.yaml), which can be postprocessed by running to create a figure showing lift and drag coefficients over time.

$ python postprocess/vortexshedding.py
problem: newtonian

# Time Stepping Settings
implicit: true
stab: supg

checkpoint_interval: 10

ts:
  adapt_type: 'none'
  type: alpha
  dt: .05
  max_time: 100
  alpha_radius: 0.5
  monitor_solution: cgns:vortexshedding-q3-g1-n08.cgns
  monitor_solution_interval: 5
  monitor_wall_force: ascii:force.csv:ascii_csv

# Reference state is used for the initial condition, zero velocity by default.

# This choice of pressure and temperature have a density of 1 and acoustic speed
# of 100. With velocity 1, this flow is Mach 0.01.
reference:
  pressure: 7143
  temperature: 24.92

# If the the outflow is placed close to the cylinder, this will recirculate cold
# fluid, demonstrating how the outflow BC is stable despite recirculation.
outflow:
  temperature: 20

# Freestream inherits reference state as default
freestream:
  velocity: 1,0,0
# Small gravity vector to break symmetry so shedding can start
gravity: 0,-.01,0

# viscosity corresponds to Reynolds number 100
mu: 0.01
k: 14.34 # thermal conductivity, Pr = 0.71 typical of air

## DM Settings:
degree: 3
dm_plex_filename: examples/meshes/cylinder-q1-n08.msh

# Boundary Settings
bc_symmetry_z: 6
bc_wall: 5
bc_freestream: 1
bc_outflow: 2
bc_symmetry_y: 3,4
wall_comps: 1,2,3

# Primitive variables are preferred at low Mach number
state_var: primitive

dm_view:
ts_monitor:
snes_lag_jacobian: 20
snes_lag_jacobian_persists:

#pmat_pbdiagonal:
#ksp_type: bcgsl
#pc_type: vpbjacobi
amat_type: shell

Density current

The Density Current problem has the following command-line options in addition to the Newtonian Ideal Gas options:

Table 10 Density Current Runtime Options

Option

Description

Default value

Unit

-center

Location of bubble center

(lx,ly,lz)/2

(m,m,m)

-dc_axis

Axis of density current cylindrical anomaly, or (0,0,0) for spherically symmetric

(0,0,0)

-rc

Characteristic radius of thermal bubble

1000

m

-theta0

Reference potential temperature

300

K

-thetaC

Perturbation of potential temperature

-15

K

-P0

Atmospheric pressure

1E5

Pa

-N

Brunt-Vaisala frequency

0.01

1/s

This problem can be run with:

./navierstokes -problem density_current -dm_plex_box_faces 16,1,8 -degree 1 -dm_plex_box_lower 0,0,0 -dm_plex_box_upper 2000,125,1000 -dm_plex_dim 3 -rc 400. -bc_wall 1,2,5,6 -wall_comps 1,2,3 -bc_symmetry_y 3,4 -mu 75

Channel flow

The Channel problem has the following command-line options in addition to the Newtonian Ideal Gas options:

Table 11 Channel Runtime Options

Option

Description

Default value

Unit

-umax

Maximum/centerline velocity of the flow

10

m/s

-theta0

Reference potential temperature

300

K

-P0

Atmospheric pressure

1E5

Pa

-body_force_scale

Multiplier for body force (-1 for flow reversal)

1

This problem can be run with the examples/channel.yaml file via:

./build/navierstokes -options_file examples/channel.yaml
problem: 'channel'
mu: .01

umax: 40
implicit: true
ts:
  type: 'beuler'
  adapt_type: 'none'
  dt: 5e-6

q_extra: 2

dm_plex_box_lower: 0,0,0
dm_plex_box_upper: .01,.01,.001
dm_plex_dim: 3
degree: 1
dm_plex_box_faces: 10,10,1
bc_symmetry_z: 1,2
bc_wall: 3,4
wall_comps: 1,2,3
dm_plex_box_bd: 'periodic,none,none'

Blasius boundary layer

The Blasius problem has the following command-line options in addition to the Newtonian Ideal Gas options:

Table 12 Blasius Runtime Options

Option

Description

Default value

Unit

-velocity_infinity

Freestream velocity

40

m/s

-temperature_infinity

Freestream temperature

288

K

-pressure_infinity

Atmospheric pressure, also sets IDL reference pressure

1.01E5

Pa

-temperature_wall

Wall temperature

288

K

-delta0

Boundary layer height at the inflow

4.2e-3

m

-platemesh_modify_mesh

Whether to modify the mesh using the given options below.

false

-platemesh_refine_height

Height at which -platemesh_Ndelta number of elements should refined into

5.9E-4

m

-platemesh_Ndelta

Number of elements to keep below -platemesh_refine_height

45

-platemesh_growth

Growth rate of the elements in the refinement region

1.08

-platemesh_top_angle

Downward angle of the top face of the domain. This face serves as an outlet.

5

degrees

-platemesh_y_node_locs_path

Path to file with y node locations. If empty, will use mesh warping instead.

""

-stg_use

Whether to use STG for the inflow conditions

false

-n_chebyshev

Number of Chebyshev terms

20

-chebyshev_

Prefix for Chebyshev snes solve

This problem can be run with the examples/blasius.yaml file via:

./build/navierstokes -options_file examples/blasius.yaml
problem: 'blasius'

implicit: true
ts:
  adapt_type: 'none'
  type: 'beuler'
  dt: 2e-6
  max_time: 1.0e-3
  #monitor_solution: cgns:blasius-%d.cgns
  #monitor_solution_interval: 10
checkpoint_interval: 10

## Linear Settings:
degree: 1
dm_plex_box_faces: 40,60,1
mesh_transform: platemesh
platemesh_nDelta: 45

# # Quadratic Settings:
# degree: 2
# dm_plex_box_faces: 20,30,1
# platemesh:
#   modify_mesh: true
#   nDelta: 22
#   growth: 1.1664 # 1.08^2

stab: 'supg'

dm_plex_box_lower: 0,0,0
dm_plex_box_upper: 4.2e-3,4.2e-3,5.e-5
dm_plex_dim: 3
# Faces labeled 1=z- 2=z+ 3=y- 4=y+ 5=x+ 6=x-
bc_symmetry_z: 1,2
bc_wall: 3
wall_comps: 1,2,3
bc_inflow: 6
bc_outflow: 5,4
gravity: 0,0,0

# stg:
#   use: false
#   inflow_path: "./STGInflow_blasius.dat"
#   mean_only: true

# ts_monitor_turbulence_spanstats:
#   collect_interval: 1
#   viewer_interval: 5
#   viewer: cgns:stats-%d.cgns
#   viewer_cgns_batch_size: 1

STG Inflow for Flat Plate

Using the STG Inflow for the blasius problem adds the following command-line options:

Table 13 Blasius Runtime Options

Option

Description

Default value

Unit

-stg_inflow_path

Path to the STGInflow file

./STGInflow.dat

-stg_rand_path

Path to the STGRand file

./STGRand.dat

-stg_alpha

Growth rate of the wavemodes

1.01

-stg_u0

Convective velocity, \(U_0\)

0.0

m/s

-stg_mean_only

Only impose the mean velocity (no fluctutations)

false

-stg_strong

Strongly enforce the STG inflow boundary condition

false

-stg_fluctuating_IC

“Extrude” the fluctuations through the domain as an initial condition

false

-stg_dx

Set the element size in the x direction. Default is calculated for box meshes, assuming equispaced elements.

m

-stg_h_scale_factor

Scale element size for cutoff frequency calculation

\(1/p\)

This problem can be run with the examples/blasius.yaml file via:

./build/navierstokes -options_file examples/blasius.yaml -stg_use true

Note the added -stg_use true flag This overrides the stg: use: false setting in the examples/blasius.yaml file, enabling the use of the STG inflow.