NavierStokesFV System
Overview
The NavierStokesFV system is dedicated to decrease the effort required by the user to prepare simulations which need to solve the Navier Stokes equations. The action is capable of setting up:
Incompressible and weakly-compressible simulations for
Clean fluids and flows in porous media
The action handles boundary conditions, necessary materials for fluid properties, variables, variable initialization and various objects for turbulence modeling.
This action only supports Rhie-Chow interpolation for the determination of face velocities in the advection terms. The face interpolation of the advected quantities (e.g. upwind, average) can be controlled through the *_advection_interpolation
action parameters.
Automatically defined variables
The NavierStokesFV action automatically sets up the variables which are necessary for the solution of a given problem. These variables can then be used to couple fluid flow simulations with other physics. The list of variable names commonly used in the action syntax is presented below:
Velocities for non-porous-medium simulations:
(moose/modules/navier_stokes/include/base/NS.h)static const std::string velocity_x = "vel_x"; static const std::string velocity_y = "vel_y"; static const std::string velocity_z = "vel_z";
Velocities for porous medium simulations:
(moose/modules/navier_stokes/include/base/NS.h)static const std::string superficial_velocity_x = "superficial_vel_x"; static const std::string superficial_velocity_y = "superficial_vel_y"; static const std::string superficial_velocity_z = "superficial_vel_z";
Pressure and temperature:
(moose/modules/navier_stokes/include/base/NS.h)static const std::string pressure = "pressure";
(moose/modules/navier_stokes/include/base/NS.h)static const std::string T_fluid = "T_fluid";
For the default names of other variables used in this action, visit this site.
Bernoulli pressure jump treatment
Please see the Bernouilli pressure variable documentation for more information.
Examples
Incompressible fluid flow in a lid-driven cavity
In the following examples we present how the NavierStokesFV action can be utilized to simplify input files. We start with the simple lid-driven cavity problem with an Incompressible Navier Stokes formulation. The original description of the problem is available under the following link. First, we present the input file where the simulation is set up manually, by defining every kernel, boundary condition and material explicitly.
mu = 1
rho = 1
k = .01
cp = 1
velocity_interp_method = 'rc'
advected_interp_method = 'average'
[GlobalParams]
rhie_chow_user_object = 'rc'
[]
[UserObjects]
[rc]
type = INSFVRhieChowInterpolator
u = vel_x
v = vel_y
pressure = pressure
[]
[]
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = 1
ymin = 0
ymax = 1
nx = 32
ny = 32
[]
[]
[Variables]
[vel_x]
type = INSFVVelocityVariable
[]
[vel_y]
type = INSFVVelocityVariable
[]
[pressure]
type = INSFVPressureVariable
[]
[T_fluid]
type = INSFVEnergyVariable
[]
[lambda]
family = SCALAR
order = FIRST
[]
[]
[AuxVariables]
[U]
order = CONSTANT
family = MONOMIAL
fv = true
[]
[]
[AuxKernels]
[mag]
type = VectorMagnitudeAux
variable = U
x = vel_x
y = vel_y
[]
[]
[FVKernels]
[mass]
type = INSFVMassAdvection
variable = pressure
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = ${rho}
[]
[mean_zero_pressure]
type = FVIntegralValueConstraint
variable = pressure
lambda = lambda
[]
[u_advection]
type = INSFVMomentumAdvection
variable = vel_x
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = ${advected_interp_method}
rho = ${rho}
momentum_component = 'x'
[]
[u_viscosity]
type = INSFVMomentumDiffusion
variable = vel_x
mu = ${mu}
momentum_component = 'x'
[]
[u_pressure]
type = INSFVMomentumPressure
variable = vel_x
momentum_component = 'x'
pressure = pressure
[]
[v_advection]
type = INSFVMomentumAdvection
variable = vel_y
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = ${advected_interp_method}
rho = ${rho}
momentum_component = 'y'
[]
[v_viscosity]
type = INSFVMomentumDiffusion
variable = vel_y
mu = ${mu}
momentum_component = 'y'
[]
[v_pressure]
type = INSFVMomentumPressure
variable = vel_y
momentum_component = 'y'
pressure = pressure
[]
[temp_conduction]
type = FVDiffusion
coeff = 'k'
variable = T_fluid
[]
[temp_advection]
type = INSFVEnergyAdvection
variable = T_fluid
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = ${advected_interp_method}
[]
[]
[FVBCs]
[top_x]
type = INSFVNoSlipWallBC
variable = vel_x
boundary = 'top'
function = 'lid_function'
[]
[no_slip_x]
type = INSFVNoSlipWallBC
variable = vel_x
boundary = 'left right bottom'
function = 0
[]
[no_slip_y]
type = INSFVNoSlipWallBC
variable = vel_y
boundary = 'left right top bottom'
function = 0
[]
[T_hot]
type = FVDirichletBC
variable = T_fluid
boundary = 'bottom'
value = 1
[]
[T_cold]
type = FVDirichletBC
variable = T_fluid
boundary = 'top'
value = 0
[]
[]
[Materials]
[functor_constants]
type = ADGenericFunctorMaterial
prop_names = 'cp k'
prop_values = '${cp} ${k}'
[]
[ins_fv]
type = INSFVEnthalpyMaterial
temperature = 'T_fluid'
rho = ${rho}
[]
[]
[Functions]
[lid_function]
type = ParsedFunction
expression = '4*x*(1-x)'
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
nl_rel_tol = 1e-12
[]
[Outputs]
exodus = true
[]
(moose/modules/navier_stokes/test/tests/finite_volume/ins/lid-driven/lid-driven-with-energy.i)Careful! The utilization of central difference (average
) advected interpolation may lead to oscillatory behavior in certain scenarios. Even though it is not the case for this example, if this phenomenon arises, we recommend using first order upwind
or second order TVD chemes.
The same simulation can be set up using the action syntax which improves input file readability:
mu = 1
rho = 1
k = .01
cp = 1
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = 1
ymin = 0
ymax = 1
nx = 32
ny = 32
[]
[]
[Modules]
[NavierStokesFV]
compressibility = 'incompressible'
add_energy_equation = true
density = 'rho'
dynamic_viscosity = 'mu'
thermal_conductivity = 'k'
specific_heat = 'cp'
initial_pressure = 0.0
initial_temperature = 0.0
inlet_boundaries = 'top'
momentum_inlet_types = 'fixed-velocity'
momentum_inlet_function = 'lid_function 0'
energy_inlet_types = 'fixed-temperature'
energy_inlet_function = '0'
wall_boundaries = 'left right bottom'
momentum_wall_types = 'noslip noslip noslip'
energy_wall_types = 'heatflux heatflux fixed-temperature'
energy_wall_function = '0 0 1'
pin_pressure = true
pinned_pressure_type = average
pinned_pressure_value = 0
mass_advection_interpolation = 'average'
momentum_advection_interpolation = 'average'
energy_advection_interpolation = 'average'
[]
[]
[AuxVariables]
[U]
order = CONSTANT
family = MONOMIAL
fv = true
[]
[]
[AuxKernels]
[mag]
type = VectorMagnitudeAux
variable = U
x = vel_x
y = vel_y
[]
[]
[Materials]
[functor_constants]
type = ADGenericFunctorMaterial
prop_names = 'cp k rho mu'
prop_values = '${cp} ${k} ${rho} ${mu}'
[]
[]
[Functions]
[lid_function]
type = ParsedFunction
expression = '4*x*(1-x)'
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
nl_rel_tol = 1e-12
[]
[Outputs]
exodus = true
[]
(moose/modules/navier_stokes/test/tests/finite_volume/ins/lid-driven/lid-driven-with-energy-action.i)It is visible that in this case we defined the "compressibility" parameter to be incompressible
. Furthermore, the energy (enthalpy) equation is solved as well. The user can request this by setting "add_energy_equation" to true
. The boundary types are grouped into wall, inlet and outlet types. For more information on the available boundary types, see the list of parameters at the bottom of the page.
Incompressible fluid flow in porous medium
The following input file sets up a simulation of an incompressible fluid flow within a channel which contains a homogenized structure that is treated as a porous medium. The model accounts for the heat exchange between the fluid and the homogenized structure as well. For more description on the used model, visit the Porous medium Incompressible Navier Stokes page. First, the input file with the manually defined kernels and boundary conditions is presented:
mu = 1
rho = 1
k = 1e-3
cp = 1
u_inlet = 1
T_inlet = 200
advected_interp_method = 'average'
velocity_interp_method = 'rc'
[Mesh]
[mesh]
type = CartesianMeshGenerator
dim = 2
dx = '5 5'
dy = '1.0'
ix = '50 50'
iy = '20'
subdomain_id = '1 2'
[]
[]
[GlobalParams]
rhie_chow_user_object = 'rc'
[]
[UserObjects]
[rc]
type = PINSFVRhieChowInterpolator
u = superficial_vel_x
v = superficial_vel_y
pressure = pressure
porosity = porosity
[]
[]
[Variables]
inactive = 'T_solid'
[superficial_vel_x]
type = PINSFVSuperficialVelocityVariable
initial_condition = ${u_inlet}
[]
[superficial_vel_y]
type = PINSFVSuperficialVelocityVariable
initial_condition = 1e-6
[]
[pressure]
type = INSFVPressureVariable
[]
[T_fluid]
type = INSFVEnergyVariable
[]
[T_solid]
family = 'MONOMIAL'
order = 'CONSTANT'
fv = true
[]
[]
[AuxVariables]
[T_solid]
family = 'MONOMIAL'
order = 'CONSTANT'
fv = true
initial_condition = 100
[]
[porosity]
family = MONOMIAL
order = CONSTANT
fv = true
initial_condition = 0.5
[]
[]
[FVKernels]
inactive = 'solid_energy_diffusion solid_energy_convection'
[mass]
type = PINSFVMassAdvection
variable = pressure
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = ${rho}
[]
[u_advection]
type = PINSFVMomentumAdvection
variable = superficial_vel_x
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = ${rho}
porosity = porosity
momentum_component = 'x'
[]
[u_viscosity]
type = PINSFVMomentumDiffusion
variable = superficial_vel_x
mu = ${mu}
porosity = porosity
momentum_component = 'x'
[]
[u_pressure]
type = PINSFVMomentumPressure
variable = superficial_vel_x
momentum_component = 'x'
pressure = pressure
porosity = porosity
[]
[v_advection]
type = PINSFVMomentumAdvection
variable = superficial_vel_y
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = ${rho}
porosity = porosity
momentum_component = 'y'
[]
[v_viscosity]
type = PINSFVMomentumDiffusion
variable = superficial_vel_y
mu = ${mu}
porosity = porosity
momentum_component = 'y'
[]
[v_pressure]
type = PINSFVMomentumPressure
variable = superficial_vel_y
momentum_component = 'y'
pressure = pressure
porosity = porosity
[]
[energy_advection]
type = PINSFVEnergyAdvection
variable = T_fluid
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = ${advected_interp_method}
[]
[energy_diffusion]
type = PINSFVEnergyDiffusion
k = ${k}
variable = T_fluid
porosity = porosity
[]
[energy_convection]
type = PINSFVEnergyAmbientConvection
variable = T_fluid
is_solid = false
T_fluid = 'T_fluid'
T_solid = 'T_solid'
h_solid_fluid = 'h_cv'
[]
[solid_energy_diffusion]
type = FVDiffusion
coeff = ${k}
variable = T_solid
[]
[solid_energy_convection]
type = PINSFVEnergyAmbientConvection
variable = T_solid
is_solid = true
T_fluid = 'T_fluid'
T_solid = 'T_solid'
h_solid_fluid = 'h_cv'
[]
[]
[FVBCs]
inactive = 'heated-side'
[inlet-u]
type = INSFVInletVelocityBC
boundary = 'left'
variable = superficial_vel_x
function = ${u_inlet}
[]
[inlet-v]
type = INSFVInletVelocityBC
boundary = 'left'
variable = superficial_vel_y
function = 0
[]
[inlet-T]
type = FVNeumannBC
variable = T_fluid
value = '${fparse u_inlet * rho * cp * T_inlet}'
boundary = 'left'
[]
[no-slip-u]
type = INSFVNoSlipWallBC
boundary = 'top'
variable = superficial_vel_x
function = 0
[]
[no-slip-v]
type = INSFVNoSlipWallBC
boundary = 'top'
variable = superficial_vel_y
function = 0
[]
[heated-side]
type = FVDirichletBC
boundary = 'top'
variable = 'T_solid'
value = 150
[]
[symmetry-u]
type = PINSFVSymmetryVelocityBC
boundary = 'bottom'
variable = superficial_vel_x
u = superficial_vel_x
v = superficial_vel_y
mu = ${mu}
momentum_component = 'x'
[]
[symmetry-v]
type = PINSFVSymmetryVelocityBC
boundary = 'bottom'
variable = superficial_vel_y
u = superficial_vel_x
v = superficial_vel_y
mu = ${mu}
momentum_component = 'y'
[]
[symmetry-p]
type = INSFVSymmetryPressureBC
boundary = 'bottom'
variable = pressure
[]
[outlet-p]
type = INSFVOutletPressureBC
boundary = 'right'
variable = pressure
function = 0.1
[]
[]
[Materials]
[constants]
type = ADGenericFunctorMaterial
prop_names = 'h_cv'
prop_values = '1'
[]
[functor_constants]
type = ADGenericFunctorMaterial
prop_names = 'cp'
prop_values = '${cp}'
[]
[ins_fv]
type = INSFVEnthalpyMaterial
rho = ${rho}
temperature = 'T_fluid'
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
nl_rel_tol = 1e-14
[]
# Some basic Postprocessors to examine the solution
[Postprocessors]
[inlet-p]
type = SideAverageValue
variable = pressure
boundary = 'left'
[]
[outlet-u]
type = SideAverageValue
variable = superficial_vel_x
boundary = 'right'
[]
[outlet-temp]
type = SideAverageValue
variable = T_fluid
boundary = 'right'
[]
[solid-temp]
type = ElementAverageValue
variable = T_solid
[]
[]
[Outputs]
exodus = true
csv = false
[]
(moose/modules/navier_stokes/test/tests/finite_volume/pins/channel-flow/heated/2d-rc-heated.i)Careful! The utilization of central difference (average
) advected interpolation may lead to oscillatory behavior in certain scenarios. Even though it is not the case for this example, if this phenomenon arises, we recommend using first order upwind
or second order TVD chemes.
The same simulation can also be set up using the NavierStokesFV action syntax:
mu = 1
rho = 1
k = 1e-3
cp = 1
u_inlet = 1
T_inlet = 200
h_cv = 1.0
[Mesh]
[mesh]
type = CartesianMeshGenerator
dim = 2
dx = '5 5'
dy = '1.0'
ix = '50 50'
iy = '20'
subdomain_id = '1 2'
[]
[]
[Variables]
[T_solid]
type = MooseVariableFVReal
[]
[]
[AuxVariables]
[porosity]
type = MooseVariableFVReal
initial_condition = 0.5
[]
[]
[Modules]
[NavierStokesFV]
compressibility = 'incompressible'
porous_medium_treatment = true
add_energy_equation = true
density = ${rho}
dynamic_viscosity = ${mu}
thermal_conductivity = ${k}
specific_heat = ${cp}
porosity = 'porosity'
initial_velocity = '${u_inlet} 1e-6 0'
initial_pressure = 0.0
initial_temperature = 0.0
inlet_boundaries = 'left'
momentum_inlet_types = 'fixed-velocity'
momentum_inlet_function = '${u_inlet} 0'
energy_inlet_types = 'heatflux'
energy_inlet_function = '${fparse u_inlet * rho * cp * T_inlet}'
wall_boundaries = 'top bottom'
momentum_wall_types = 'noslip symmetry'
energy_wall_types = 'heatflux heatflux'
energy_wall_function = '0 0'
outlet_boundaries = 'right'
momentum_outlet_types = 'fixed-pressure'
pressure_function = '0.1'
ambient_convection_alpha = ${h_cv}
ambient_temperature = 'T_solid'
mass_advection_interpolation = 'average'
momentum_advection_interpolation = 'average'
energy_advection_interpolation = 'average'
[]
[]
[FVKernels]
[solid_energy_diffusion]
type = FVDiffusion
coeff = ${k}
variable = T_solid
[]
[solid_energy_convection]
type = PINSFVEnergyAmbientConvection
variable = 'T_solid'
is_solid = true
T_fluid = 'T_fluid'
T_solid = 'T_solid'
h_solid_fluid = ${h_cv}
[]
[]
[FVBCs]
[heated-side]
type = FVDirichletBC
boundary = 'top'
variable = 'T_solid'
value = 150
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
nl_rel_tol = 1e-14
[]
# Some basic Postprocessors to examine the solution
[Postprocessors]
[inlet-p]
type = SideAverageValue
variable = pressure
boundary = 'left'
[]
[outlet-u]
type = SideAverageValue
variable = superficial_vel_x
boundary = 'right'
[]
[outlet-temp]
type = SideAverageValue
variable = T_fluid
boundary = 'right'
[]
[solid-temp]
type = ElementAverageValue
variable = T_solid
[]
[]
[Outputs]
exodus = true
csv = false
[]
(moose/modules/navier_stokes/test/tests/finite_volume/pins/channel-flow/heated/2d-rc-heated-action.i)Compared to the previous example, we see that in this case the porous medium treatment is enabled by setting "porous_medium_treatment" to true
. The corresponding porosity can be supplied through the "porosity" parameter. Furthermore, the heat exchange between the fluid and the homogenized structure is enabled using the "ambient_temperature" and "ambient_convection_alpha" parameters.
Weakly-compressible fluid flow
The last example is dedicated to demonstrating a transient flow in a channel using a weakly-compressible approximation. The following examples shows how this simulation is set up by manually defining the kernels and boundary conditions. For more information on the weakly-compressible treatment, visit the Weakly-compressible Navier Stokes page.
rho = 'rho'
l = 10
velocity_interp_method = 'rc'
advected_interp_method = 'average'
# Artificial fluid properties
# For a real case, use a GeneralFluidFunctorProperties and a viscosity rampdown
# or initialize very well!
k = 1
cp = 1000
mu = 1e2
# Operating conditions
inlet_temp = 300
outlet_pressure = 1e5
inlet_v = 0.001
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = ${l}
ymin = 0
ymax = 1
nx = 20
ny = 10
[]
[]
[GlobalParams]
rhie_chow_user_object = 'rc'
[]
[UserObjects]
[rc]
type = INSFVRhieChowInterpolator
u = vel_x
v = vel_y
pressure = pressure
[]
[]
[Variables]
[vel_x]
type = INSFVVelocityVariable
initial_condition = ${inlet_v}
[]
[vel_y]
type = INSFVVelocityVariable
initial_condition = 1e-15
[]
[pressure]
type = INSFVPressureVariable
initial_condition = ${outlet_pressure}
[]
[T_fluid]
type = INSFVEnergyVariable
initial_condition = ${inlet_temp}
[]
[]
[AuxVariables]
[mixing_length]
type = MooseVariableFVReal
[]
[power_density]
type = MooseVariableFVReal
initial_condition = 1e4
[]
[]
[FVKernels]
inactive = 'u_turb v_turb temp_turb'
[mass_time]
type = WCNSFVMassTimeDerivative
variable = pressure
drho_dt = drho_dt
[]
[mass]
type = INSFVMassAdvection
variable = pressure
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = ${rho}
[]
[u_time]
type = WCNSFVMomentumTimeDerivative
variable = vel_x
drho_dt = drho_dt
rho = rho
momentum_component = 'x'
[]
[u_advection]
type = INSFVMomentumAdvection
variable = vel_x
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = ${advected_interp_method}
rho = ${rho}
momentum_component = 'x'
[]
[u_viscosity]
type = INSFVMomentumDiffusion
variable = vel_x
mu = ${mu}
momentum_component = 'x'
[]
[u_pressure]
type = INSFVMomentumPressure
variable = vel_x
momentum_component = 'x'
pressure = pressure
[]
[u_turb]
type = INSFVMixingLengthReynoldsStress
variable = vel_x
rho = ${rho}
mixing_length = 'mixing_length'
momentum_component = 'x'
u = vel_x
v = vel_y
[]
[v_time]
type = WCNSFVMomentumTimeDerivative
variable = vel_y
drho_dt = drho_dt
rho = rho
momentum_component = 'y'
[]
[v_advection]
type = INSFVMomentumAdvection
variable = vel_y
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = ${advected_interp_method}
rho = ${rho}
momentum_component = 'y'
[]
[v_viscosity]
type = INSFVMomentumDiffusion
variable = vel_y
momentum_component = 'y'
mu = ${mu}
[]
[v_pressure]
type = INSFVMomentumPressure
variable = vel_y
momentum_component = 'y'
pressure = pressure
[]
[v_turb]
type = INSFVMixingLengthReynoldsStress
variable = vel_y
rho = ${rho}
mixing_length = 'mixing_length'
momentum_component = 'y'
u = vel_x
v = vel_y
[]
[temp_time]
type = WCNSFVEnergyTimeDerivative
variable = T_fluid
cp = cp
rho = rho
drho_dt = drho_dt
[]
[temp_conduction]
type = FVDiffusion
coeff = 'k'
variable = T_fluid
[]
[temp_advection]
type = INSFVEnergyAdvection
variable = T_fluid
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = ${advected_interp_method}
[]
[heat_source]
type = FVCoupledForce
variable = T_fluid
v = power_density
[]
[temp_turb]
type = WCNSFVMixingLengthEnergyDiffusion
variable = T_fluid
rho = rho
cp = cp
mixing_length = 'mixing_length'
schmidt_number = 1
u = vel_x
v = vel_y
[]
[]
[FVBCs]
[no_slip_x]
type = INSFVNoSlipWallBC
variable = vel_x
boundary = 'top bottom'
function = 0
[]
[no_slip_y]
type = INSFVNoSlipWallBC
variable = vel_y
boundary = 'top bottom'
function = 0
[]
# Inlet
[inlet_u]
type = INSFVInletVelocityBC
variable = vel_x
boundary = 'left'
function = ${inlet_v}
[]
[inlet_v]
type = INSFVInletVelocityBC
variable = vel_y
boundary = 'left'
function = 0
[]
[inlet_T]
type = FVDirichletBC
variable = T_fluid
boundary = 'left'
value = ${inlet_temp}
[]
[outlet_p]
type = INSFVOutletPressureBC
variable = pressure
boundary = 'right'
function = ${outlet_pressure}
[]
[]
[FluidProperties]
[fp]
type = FlibeFluidProperties
[]
[]
[Materials]
[const_functor]
type = ADGenericFunctorMaterial
prop_names = 'cp k'
prop_values = '${cp} ${k}'
[]
[rho]
type = RhoFromPTFunctorMaterial
fp = fp
temperature = T_fluid
pressure = pressure
[]
[ins_fv]
type = INSFVEnthalpyMaterial
temperature = 'T_fluid'
rho = ${rho}
[]
[]
[AuxKernels]
inactive = 'mixing_len'
[mixing_len]
type = WallDistanceMixingLengthAux
walls = 'top'
variable = mixing_length
execute_on = 'initial'
delta = 0.5
[]
[]
[Executioner]
type = Transient
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
[TimeStepper]
type = IterationAdaptiveDT
dt = 1e-3
optimal_iterations = 6
[]
end_time = 15
nl_abs_tol = 1e-9
nl_max_its = 50
line_search = 'none'
automatic_scaling = true
off_diagonals_in_auto_scaling = true
compute_scaling_once = false
[]
[Outputs]
exodus = true
[]
(moose/modules/navier_stokes/test/tests/finite_volume/wcns/channel-flow/2d-transient.i)Careful! The utilization of central difference (average
) advected interpolation may lead to oscillatory behavior in certain scenarios. Even though it is not the case for this example, if this phenomenon arises, we recommend using first order upwind
or second order TVD chemes.
The same simulation can be set up using the action syntax as folows:
l = 10
# Artificial fluid properties
# For a real case, use a GeneralFluidFunctorProperties and a viscosity rampdown
# or initialize very well!
k = 1
cp = 1000
mu = 1e2
# Operating conditions
inlet_temp = 300
outlet_pressure = 1e5
inlet_v = 0.001
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = ${l}
ymin = 0
ymax = 1
nx = 20
ny = 10
[]
[]
[Modules]
[NavierStokesFV]
compressibility = 'weakly-compressible'
add_energy_equation = true
density = 'rho'
dynamic_viscosity = 'mu'
thermal_conductivity = 'k'
specific_heat = 'cp'
initial_velocity = '${inlet_v} 1e-15 0'
initial_temperature = '${inlet_temp}'
initial_pressure = '${outlet_pressure}'
inlet_boundaries = 'left'
momentum_inlet_types = 'fixed-velocity'
momentum_inlet_function = '${inlet_v} 0'
energy_inlet_types = 'fixed-temperature'
energy_inlet_function = '${inlet_temp}'
wall_boundaries = 'top bottom'
momentum_wall_types = 'noslip noslip'
energy_wall_types = 'heatflux heatflux'
energy_wall_function = '0 0'
outlet_boundaries = 'right'
momentum_outlet_types = 'fixed-pressure'
pressure_function = '${outlet_pressure}'
external_heat_source = 'power_density'
mass_advection_interpolation = 'average'
momentum_advection_interpolation = 'average'
energy_advection_interpolation = 'average'
[]
[]
[AuxVariables]
[power_density]
type = MooseVariableFVReal
initial_condition = 1e4
[]
[]
[FluidProperties]
[fp]
type = FlibeFluidProperties
[]
[]
[Materials]
[const_functor]
type = ADGenericFunctorMaterial
prop_names = 'cp k mu'
prop_values = '${cp} ${k} ${mu}'
[]
[rho]
type = RhoFromPTFunctorMaterial
fp = fp
temperature = T_fluid
pressure = pressure
[]
[]
[Executioner]
type = Transient
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
[TimeStepper]
type = IterationAdaptiveDT
dt = 1e-3
optimal_iterations = 6
[]
end_time = 15
nl_abs_tol = 1e-9
nl_max_its = 50
line_search = 'none'
automatic_scaling = true
off_diagonals_in_auto_scaling = true
compute_scaling_once = false
[]
[Outputs]
exodus = true
[]
(moose/modules/navier_stokes/test/tests/finite_volume/wcns/channel-flow/2d-transient-action.i)We note that the weakly-compressible handling can be enabled by setting "compressibility" to weakly-compressible
. As shown in the example, an arbitrary energy source function can also be supplied to the incorporated energy equation using the "external_heat_source" parameter.
Available Actions
- Navier Stokes App
- NSFVActionThis class allows us to set up Navier-Stokes equations for porous medium or clean fluid flows using incompressible or weakly compressible approximations with a finite volume discretization.