Navier Stokes Scalar Transport / WCNSFVScalarTransportPhysics

Define the Navier Stokes weakly-compressible scalar field transport equation(s)

Equation

This Physics object creates the kernels and boundary conditions to solve the advection-diffusion-reaction equation for several scalar quantities advected by the flow.

ϕit+(ϕiv)(kiϕi)Qiλiϕi=0\dfrac{\partial \phi_i}{\partial t} + \nabla \cdot (\phi_i \mathbf{v}) - \nabla \cdot (k_i \nabla \phi_i) - Q_i - \lambda_i \phi_i = 0

where:

  • ϕi\phi_i is the i-th scalar quantity

  • \mathbf{v} is the advecting velocity

  • kik_i the i-th scalar diffusivity

  • QiQ_i is the i-th precursor source

  • λi\lambda_i is a reaction coefficient. It should be negative for a loss term

The kernels created are:

Coupling with other Physics

Scalar advection equations can be solved concurrently with the flow equations by combining the Navier Stokes Flow / WCNSFVFlowPhysics with the WCNSFVScalarTransportPhysics. The following input performs this coupling for incompressible flow 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}

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

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

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

        mass_advection_interpolation = 'average'
        momentum_advection_interpolation = 'average'
      []
    []
    [FluidHeatTransfer]
      [heat]
        thermal_conductivity = ${k}
        specific_heat = ${cp}

        energy_inlet_types = 'fixed-temperature'
        energy_inlet_function = '1'

        energy_wall_types = 'heatflux heatflux'
        energy_wall_function = '0 0'

        energy_advection_interpolation = 'average'
      []
    []

    [ScalarTransport]
      [heat]
        passive_scalar_names = 'scalar'

        passive_scalar_diffusivity = ${diff}
        passive_scalar_source = 0.1
        passive_scalar_coupled_source = U
        passive_scalar_coupled_source_coeff = 0.1

        passive_scalar_inlet_types = 'fixed-value'
        passive_scalar_inlet_function = '1'

        passive_scalar_advection_interpolation = 'average'
      []
    []
  []
[]
(contrib/moose/modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/2d-scalar-transport-physics.i)

Input Parameters

  • add_scalar_equationFalseWhether to add the scalar transport equation. This parameter is not necessary if using the Physics syntax

    Default:False

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Whether to add the scalar transport equation. This parameter is not necessary if using the Physics syntax

  • 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.

  • preconditioningnoneWhich preconditioning to use for this Physics

    Default:none

    C++ Type:MooseEnum

    Unit:(no unit assumed)

    Options:default, none

    Controllable:No

    Description:Which preconditioning to use for this Physics

  • 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

  • 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

  • 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

  • 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

Coupled Physics 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

  • initial_scalar_variablesInitial values of the passive scalar variables.

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

    Unit:(no unit assumed)

    Controllable:No

    Description:Initial values of the passive scalar variables.

  • passive_scalar_namesVector containing the names of the advected scalar variables.

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

    Unit:(no unit assumed)

    Controllable:No

    Description:Vector containing the names of the advected scalar variables.

Variable Parameters

  • passive_scalar_advection_interpolationupwindThe numerical scheme to use for interpolating passive scalar field, as an advected quantity, to the face.

    Default:upwind

    C++ Type:MooseEnum

    Unit:(no unit assumed)

    Options:average, upwind, sou, min_mod, vanLeer, quick, skewness-corrected

    Controllable:No

    Description:The numerical scheme to use for interpolating passive scalar field, as an advected quantity, to the face.

  • passive_scalar_face_interpolationaverageThe numerical scheme to interpolate the passive scalar field variables 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 passive scalar field variables to the face (separate from the advected quantity interpolation).

  • passive_scalar_scalingThe scaling factor for the passive scalar field variables.

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

    Unit:(no unit assumed)

    Controllable:No

    Description:The scaling factor for the passive scalar field variables.

  • passive_scalar_two_term_bc_expansionTrueIf a two-term Taylor expansion is needed for the determination of the boundary valuesof the advected passive scalar 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 advected passive scalar field.

Numerical Scheme Parameters

  • passive_scalar_coupled_sourceCoupled variable names for the sources used for the passive scalar fields. If multiple sources for each equation are specified, major (outer) ordering by equation.

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

    Unit:(no unit assumed)

    Controllable:No

    Description:Coupled variable names for the sources used for the passive scalar fields. If multiple sources for each equation are specified, major (outer) ordering by equation.

  • passive_scalar_coupled_source_coeffCoupled variable multipliers for the sources used for the passive scalar fields. If multiple sources for each equation are specified, major (outer) ordering by equation.

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

    Unit:(no unit assumed)

    Controllable:No

    Description:Coupled variable multipliers for the sources used for the passive scalar fields. If multiple sources for each equation are specified, major (outer) ordering by equation.

  • passive_scalar_diffusivityFunctor names for the diffusivities used for the passive scalar fields.

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

    Unit:(no unit assumed)

    Controllable:No

    Description:Functor names for the diffusivities used for the passive scalar fields.

  • passive_scalar_sourcePassive scalar sources

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

    Unit:(no unit assumed)

    Controllable:No

    Description:Passive scalar sources

Passive Scalar Control Parameters

  • passive_scalar_inlet_functorsFunctors for inlet boundaries in the passive scalar equations.

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

    Unit:(no unit assumed)

    Controllable:No

    Description:Functors for inlet boundaries in the passive scalar equations.

  • passive_scalar_inlet_typesTypes for the inlet boundaries for the passive scalar equation.

    C++ Type:MultiMooseEnum

    Unit:(no unit assumed)

    Options:fixed-value, flux-mass, flux-velocity

    Controllable:No

    Description:Types for the inlet boundaries for the passive scalar equation.

Inlet Boundary Parameters