Navier Stokes Turbulence / WCNSFVTurbulencePhysics

Define a turbulence model for a incompressible or weakly-compressible Navier Stokes flow with a finite volume discretization

Mixing length turbulence model

See the mixing length theory page for additional information.

If the turbulence model is selected to be the mixing-length model, a field variable representing the mixing length is computed using the WallDistanceMixingLengthAux.

The turbulent dynamic viscosity is then computed using a MixingLengthTurbulentViscosityFunctorMaterial. The following kernels are then added:

commentnote

These kernels are only added if each of these equations are being defined using their respective Physics.

K-Epsilon turbulence model

The WCNSFVTurbulencePhysics can be set to create the k-epsilon two-equation model.

The turbulent viscosity is then computed with:

The k equation is created with:

The epsilon equation is created with:

The boundary conditions are not set in this object for the TKE and TKED variables, as they are computed by the wall-functions in the relevant kernels. A boundary condition is set for the turbulent viscosity when using an auxiliary variable, with a INSFVTurbulentViscosityWallFunction.

Coupling with other Physics

A turbulence model can be added to a heat advection solve by using both a WCNSFVTurbulencePhysics and a Navier Stokes Fluid Heat Transfer / WCNSFVFluidHeatTransferPhysics. The following input performs this coupling for weakly compressible flow for the mixing length turbulence model in a 2D channel. No system parameters are passed, so the equations are solved in a fully coupled manner in the same nonlinear system.

[Physics]
  [NavierStokes]
    [Flow]
      [flow]
        compressibility = 'weakly-compressible'

        velocity_variable = 'vel_x vel_y'

        density = 'rho'
        dynamic_viscosity = 'mu'

        initial_velocity = '${inlet_v} 1e-15 0'
        initial_pressure = '${outlet_pressure}'

        inlet_boundaries = 'left'
        momentum_inlet_types = 'fixed-velocity'
        momentum_inlet_functors = '${inlet_v} 0'

        wall_boundaries = 'top bottom'
        momentum_wall_types = 'noslip noslip'

        outlet_boundaries = 'right'
        momentum_outlet_types = 'fixed-pressure'
        pressure_functors = '${outlet_pressure}'

        mass_advection_interpolation = 'average'
        momentum_advection_interpolation = 'average'
      []
    []
    [FluidHeatTransfer]
      [energy]
        coupled_flow_physics = flow

        thermal_conductivity = 'k'
        specific_heat = 'cp'

        initial_temperature = '${inlet_temp}'

        energy_inlet_types = 'fixed-temperature'
        energy_inlet_functors = '${inlet_temp}'
        energy_wall_types = 'heatflux heatflux'
        energy_wall_functors = '0 0'

        external_heat_source = 'power_density'
        energy_advection_interpolation = 'average'
      []
    []
    [Turbulence]
      [turbulence]
        coupled_flow_physics = flow
        fluid_heat_transfer_physics = energy
      []
    []
  []
[]
(contrib/moose/modules/navier_stokes/test/tests/finite_volume/wcns/channel-flow/2d-transient-physics.i)

A turbulence model can be added to a scalar advection solve by using both a WCNSFVTurbulencePhysics and a Navier Stokes Scalar Transport / WCNSFVScalarTransportPhysics. The following input performs this coupling for incompressible flow for the mixing length turbulence model in a 2D channel. No system parameters are passed, so the equations are solved in a fully coupled manner in the same nonlinear system.

[Physics]
  [NavierStokes]
    [Flow]
      [flow]
        compressibility = 'incompressible'

        density = ${rho}
        dynamic_viscosity = ${mu}

        initial_velocity = '1 1 0'
        initial_pressure = 0.0

        inlet_boundaries = 'left'
        momentum_inlet_types = 'fixed-velocity'
        momentum_inlet_function = '1 0'

        wall_boundaries = 'top bottom'
        momentum_wall_types = 'noslip symmetry'

        outlet_boundaries = 'right'
        momentum_outlet_types = 'fixed-pressure'
        pressure_function = '0'

        momentum_advection_interpolation = ${advected_interp_method}
        mass_advection_interpolation = ${advected_interp_method}
      []
    []
    [ScalarTransport]
      [scalars]
        add_scalar_equation = true
        passive_scalar_names = 'scalar'
        passive_scalar_source = 0.1
        passive_scalar_inlet_types = 'fixed-value'
        passive_scalar_inlet_function = '1'
        passive_scalar_advection_interpolation = ${advected_interp_method}
      []
    []

    [Turbulence]
      [mixing-length]
        turbulence_handling = 'mixing-length'
        coupled_flow_physics = flow
        scalar_transport_physics = scalars

        passive_scalar_schmidt_number = 1.0

        von_karman_const = ${von_karman_const}
        mixing_length_delta = 1e9
        mixing_length_walls = 'top bottom'
        mixing_length_aux_execute_on = 'initial'
      []
    []
  []
[]
(contrib/moose/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/2d-mixing-length-physics.i)

Input Parameters

  • C_pl10Production Limiter Constant Multiplier.

    Default:10

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:Production Limiter Constant Multiplier.

  • Pr_tPr_tTurbulent Prandtl number. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

    Default:Pr_t

    C++ Type:MooseFunctorName

    Unit:(no unit assumed)

    Controllable:No

    Description:Turbulent Prandtl number. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

  • blockBlocks (subdomains) that this Physics is active on.

    C++ Type:std::vector<SubdomainName>

    Unit:(no unit assumed)

    Controllable:No

    Description:Blocks (subdomains) that this Physics is active on.

  • initial_mu_tInitial value for the turbulence viscosity

    C++ Type:FunctionName

    Unit:(no unit assumed)

    Controllable:No

    Description:Initial value for the turbulence viscosity

  • neglect_advection_derivativesFalseWhether to remove the off-diagonal velocity term in the TKE and TKED advection term

    Default:False

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Whether to remove the off-diagonal velocity term in the TKE and TKED advection term

  • output_mu_tTrueWhether to add mu_t to the field outputs

    Default:True

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Whether to add mu_t to the field outputs

  • tke_advection_interpolationupwindThe numerical scheme to interpolate the TKE to the face when in the advection kernel.

    Default:upwind

    C++ Type:MooseEnum

    Unit:(no unit assumed)

    Options:average, upwind

    Controllable:No

    Description:The numerical scheme to interpolate the TKE to the face when in the advection kernel.

  • tke_nametkeName of the turbulent kinetic energy variable. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

    Default:tke

    C++ Type:MooseFunctorName

    Unit:(no unit assumed)

    Controllable:No

    Description:Name of the turbulent kinetic energy variable. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

  • tked_advection_interpolationupwindThe numerical scheme to interpolate the TKED to the face when in the advection kernel.

    Default:upwind

    C++ Type:MooseEnum

    Unit:(no unit assumed)

    Options:average, upwind

    Controllable:No

    Description:The numerical scheme to interpolate the TKED to the face when in the advection kernel.

  • tked_nameepsilonName of the turbulent kinetic energy dissipation variable. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

    Default:epsilon

    C++ Type:MooseFunctorName

    Unit:(no unit assumed)

    Controllable:No

    Description:Name of the turbulent kinetic energy dissipation variable. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

  • transientsame_as_problemWhether the physics is to be solved as a transient

    Default:same_as_problem

    C++ Type:MooseEnum

    Unit:(no unit assumed)

    Options:true, false, same_as_problem

    Controllable:No

    Description:Whether the physics is to be solved as a transient

  • turbulence_handlingnoneThe way turbulent diffusivities are determined in the turbulent regime.

    Default:none

    C++ Type:MooseEnum

    Unit:(no unit assumed)

    Options:mixing-length, k-epsilon, none

    Controllable:No

    Description:The way turbulent diffusivities are determined in the turbulent regime.

  • verboseFalseFlag to facilitate debugging a Physics

    Default:False

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Flag to facilitate debugging a Physics

Optional Parameters

  • C1_epsC1 coefficient for the turbulent kinetic energy dissipation equation

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:C1 coefficient for the turbulent kinetic energy dissipation equation

  • C2_epsC2 coefficient for the turbulent kinetic energy dissipation equation

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:C2 coefficient for the turbulent kinetic energy dissipation equation

  • initial_tke0Initial value for the turbulence kinetic energy

    Default:0

    C++ Type:FunctionName

    Unit:(no unit assumed)

    Controllable:No

    Description:Initial value for the turbulence kinetic energy

  • initial_tked0Initial value for the turbulence kinetic energy dissipation

    Default:0

    C++ Type:FunctionName

    Unit:(no unit assumed)

    Controllable:No

    Description:Initial value for the turbulence kinetic energy dissipation

  • sigma_epsScaling coefficient for the turbulent kinetic energy dissipation diffusion term. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

    C++ Type:MooseFunctorName

    Unit:(no unit assumed)

    Controllable:No

    Description:Scaling coefficient for the turbulent kinetic energy dissipation diffusion term. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

  • sigma_kScaling coefficient for the turbulent kinetic energy diffusion term. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

    C++ Type:MooseFunctorName

    Unit:(no unit assumed)

    Controllable:No

    Description:Scaling coefficient for the turbulent kinetic energy diffusion term. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

K-Epsilon Model Parameters

  • C_mu0.09Coupled turbulent kinetic energy closure.

    Default:0.09

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:Coupled turbulent kinetic energy closure.

  • bulk_wall_treatmentTrueWhether to treat the wall cell as bulk

    Default:True

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Whether to treat the wall cell as bulk

  • wall_treatment_TneqThe method used for computing the temperature wall functions

    Default:neq

    C++ Type:MooseEnum

    Unit:(no unit assumed)

    Options:eq_newton, eq_incremental, eq_linearized, neq

    Controllable:No

    Description:The method used for computing the temperature wall functions

  • wall_treatment_epsneqThe method used for computing the epsilon wall functions

    Default:neq

    C++ Type:MooseEnum

    Unit:(no unit assumed)

    Options:eq_newton, eq_incremental, eq_linearized, neq

    Controllable:No

    Description:The method used for computing the epsilon wall functions

K-Epsilon Wall Function Parameters

  • Sc_tTurbulent Schmidt numbers used for the passive scalar fields.

    C++ Type:std::vector<double>

    Unit:(no unit assumed)

    Controllable:No

    Description:Turbulent Schmidt numbers used for the passive scalar fields.

  • coupled_flow_physicsWCNSFVFlowPhysics generating the velocities

    C++ Type:PhysicsName

    Unit:(no unit assumed)

    Controllable:No

    Description:WCNSFVFlowPhysics generating the velocities

  • coupled_turbulence_physicsTurbulence Physics coupled with this Physics

    C++ Type:PhysicsName

    Unit:(no unit assumed)

    Controllable:No

    Description:Turbulence Physics coupled with this Physics

  • fluid_heat_transfer_physicsNavierStokesFVWCNSFVFluidHeatTransferPhysics generating the heat advection equations

    Default:NavierStokesFV

    C++ Type:PhysicsName

    Unit:(no unit assumed)

    Controllable:No

    Description:WCNSFVFluidHeatTransferPhysics generating the heat advection equations

  • scalar_transport_physicsNavierStokesFVWCNSFVScalarTransportPhysics generating the scalar advection equations

    Default:NavierStokesFV

    C++ Type:PhysicsName

    Unit:(no unit assumed)

    Controllable:No

    Description:WCNSFVScalarTransportPhysics generating the scalar advection equations

  • turbulent_prandtl1Turbulent Prandtl number for energy turbulent diffusion

    Default:1

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:Turbulent Prandtl number for energy turbulent diffusion

Coupled Physics Parameters

  • active__all__ If specified only the blocks named will be visited and made active

    Default:__all__

    C++ Type:std::vector<std::string>

    Unit:(no unit assumed)

    Controllable:No

    Description:If specified only the blocks named will be visited and made active

  • control_tagsAdds user-defined labels for accessing object parameters via control logic.

    C++ Type:std::vector<std::string>

    Unit:(no unit assumed)

    Controllable:No

    Description:Adds user-defined labels for accessing object parameters via control logic.

  • define_variablesTrueWhether to define variables if the variables with the specified names do not exist. Note that if the variables are defined externally from the Physics, the initial conditions will not be created in the Physics either.

    Default:True

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Whether to define variables if the variables with the specified names do not exist. Note that if the variables are defined externally from the Physics, the initial conditions will not be created in the Physics either.

  • ghost_layers2Number of layers of elements to ghost near process domain boundaries

    Default:2

    C++ Type:unsigned short

    Unit:(no unit assumed)

    Controllable:No

    Description:Number of layers of elements to ghost near process domain boundaries

  • inactiveIf specified blocks matching these identifiers will be skipped.

    C++ Type:std::vector<std::string>

    Unit:(no unit assumed)

    Controllable:No

    Description:If specified blocks matching these identifiers will be skipped.

Advanced Parameters

  • initial_from_file_timestepLATESTGives the time step number (or "LATEST") for which to read the Exodus solution

    Default:LATEST

    C++ Type:std::string

    Unit:(no unit assumed)

    Controllable:No

    Description:Gives the time step number (or "LATEST") for which to read the Exodus solution

  • initialize_variables_from_mesh_fileFalseDetermines if the variables that are added by the action are initializedfrom the mesh file (only for Exodus format)

    Default:False

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Determines if the variables that are added by the action are initializedfrom the mesh file (only for Exodus format)

Restart From Exodus Parameters

  • k_t_as_aux_variableFalseWhether to use an auxiliary variable for the turbulent conductivity

    Default:False

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Whether to use an auxiliary variable for the turbulent conductivity

  • linearize_sink_sourcesFalseWhether to linearize the source term

    Default:False

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Whether to linearize the source term

  • mu_t_as_aux_variableFalseWhether to use an auxiliary variable instead of a functor material property for the turbulent viscosity

    Default:False

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Whether to use an auxiliary variable instead of a functor material property for the turbulent viscosity

  • tke_face_interpolationaverageThe numerical scheme to interpolate the TKE to the face (separate from the advected quantity interpolation).

    Default:average

    C++ Type:MooseEnum

    Unit:(no unit assumed)

    Options:average, skewness-corrected

    Controllable:No

    Description:The numerical scheme to interpolate the TKE to the face (separate from the advected quantity interpolation).

  • tke_scaling1The scaling factor for the turbulent kinetic energy equation.

    Default:1

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:The scaling factor for the turbulent kinetic energy equation.

  • tke_two_term_bc_expansionTrueIf a two-term Taylor expansion is needed for the determination of the boundary valuesof the turbulent kinetic energy.

    Default:True

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:If a two-term Taylor expansion is needed for the determination of the boundary valuesof the turbulent kinetic energy.

  • tked_face_interpolationaverageThe numerical scheme to interpolate the TKED to the face (separate from the advected quantity interpolation).

    Default:average

    C++ Type:MooseEnum

    Unit:(no unit assumed)

    Options:average, skewness-corrected

    Controllable:No

    Description:The numerical scheme to interpolate the TKED to the face (separate from the advected quantity interpolation).

  • tked_scaling1The scaling factor for the turbulent kinetic energy dissipation equation.

    Default:1

    C++ Type:double

    Unit:(no unit assumed)

    Controllable:No

    Description:The scaling factor for the turbulent kinetic energy dissipation equation.

  • tked_two_term_bc_expansionTrueIf a two-term Taylor expansion is needed for the determination of the boundary valuesof the turbulent kinetic energy dissipation.

    Default:True

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:If a two-term Taylor expansion is needed for the determination of the boundary valuesof the turbulent kinetic energy dissipation.

  • turbulent_viscosity_interp_methodharmonicFace interpolation method for the turbulent viscosity

    Default:harmonic

    C++ Type:MooseEnum

    Unit:(no unit assumed)

    Options:average, harmonic

    Controllable:No

    Description:Face interpolation method for the turbulent viscosity

  • turbulent_viscosity_two_term_bc_expansionTrueIf a two-term Taylor expansion is needed for the determination of the boundary valuesof the turbulent viscosity.

    Default:True

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:If a two-term Taylor expansion is needed for the determination of the boundary valuesof the turbulent viscosity.

K-Epsilon Model Numerical Parameters

  • mixing_length_aux_execute_onWhen the mixing length aux kernels should be executed.

    C++ Type:ExecFlagEnum

    Unit:(no unit assumed)

    Options:NONE, INITIAL, LINEAR, NONLINEAR_CONVERGENCE, NONLINEAR, POSTCHECK, TIMESTEP_END, TIMESTEP_BEGIN, MULTIAPP_FIXED_POINT_END, MULTIAPP_FIXED_POINT_BEGIN, FINAL, CUSTOM

    Controllable:No

    Description:When the mixing length aux kernels should be executed.

  • mixing_length_delta1Tunable parameter related to the thickness of the boundary layer.When it is not specified, Prandtl's original unbounded wall distance mixing length model isretrieved. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

    Default:1

    C++ Type:MooseFunctorName

    Unit:(no unit assumed)

    Controllable:No

    Description:Tunable parameter related to the thickness of the boundary layer.When it is not specified, Prandtl's original unbounded wall distance mixing length model isretrieved. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

  • mixing_length_namemixing_lengthName of the mixing length auxiliary variable

    Default:mixing_length

    C++ Type:AuxVariableName

    Unit:(no unit assumed)

    Controllable:No

    Description:Name of the mixing length auxiliary variable

  • mixing_length_two_term_bc_expansionTrueIf a two-term Taylor expansion is needed for the determination of the boundary valuesof the mixing length field.

    Default:True

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:If a two-term Taylor expansion is needed for the determination of the boundary valuesof the mixing length field.

  • turbulence_wallsWalls where the mixing length model should be utilized.

    C++ Type:std::vector<BoundaryName>

    Unit:(no unit assumed)

    Controllable:No

    Description:Walls where the mixing length model should be utilized.

  • von_karman_const0.41Von Karman parameter for the mixing length model. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

    Default:0.41

    C++ Type:MooseFunctorName

    Unit:(no unit assumed)

    Controllable:No

    Description:Von Karman parameter for the mixing length model. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

  • von_karman_const_00.09'Escudier' model parameter. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

    Default:0.09

    C++ Type:MooseFunctorName

    Unit:(no unit assumed)

    Controllable:No

    Description:'Escudier' model parameter. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.

Mixing Length Model Parameters