- fpFluid properties userobject
C++ Type:UserObjectName
Controllable:No
Description:Fluid properties userobject
- momentum_componentThe component of the momentum equation that this kernel applies to.
C++ Type:MooseEnum
Controllable:No
Description:The component of the momentum equation that this kernel applies to.
- variableThe name of the finite volume variable this kernel applies to
C++ Type:NonlinearVariableName
Controllable:No
Description:The name of the finite volume variable this kernel applies to
CNSFVMomentumHLLC
Implements the momentum flux portion of the free-flow HLLC discretization.
Overview
This object implements the momentum equation inter-cell fluxes for the Harten-Lax-Van Leer-Contact (HLLC) formulation described in CNSFVHLLCBase.
Input Parameters
- blockThe list of blocks (ids or names) that this object will be applied
C++ Type:std::vector<SubdomainName>
Controllable:No
Description:The list of blocks (ids or names) that this object will be applied
- boundaries_to_forceThe set of boundaries to force execution of this FVFluxKernel on.
C++ Type:std::vector<BoundaryName>
Controllable:No
Description:The set of boundaries to force execution of this FVFluxKernel on.
- boundaries_to_not_forceThe set of boundaries to not force execution of this FVFluxKernel on.
C++ Type:std::vector<BoundaryName>
Controllable:No
Description:The set of boundaries to not force execution of this FVFluxKernel on.
- force_boundary_executionFalseWhether to force execution of this object on the boundary.
Default:False
C++ Type:bool
Controllable:No
Description:Whether to force execution of this object on the boundary.
- ghost_layers1The number of layers of elements to ghost.
Default:1
C++ Type:unsigned short
Controllable:No
Description:The number of layers of elements to ghost.
- prop_getter_suffixAn optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
C++ Type:MaterialPropertyName
Controllable:No
Description:An optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
- use_point_neighborsFalseWhether to use point neighbors, which introduces additional ghosting to that used for simple face neighbors.
Default:False
C++ Type:bool
Controllable:No
Description:Whether to use point neighbors, which introduces additional ghosting to that used for simple face neighbors.
Optional Parameters
- control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Controllable:No
Description:Adds user-defined labels for accessing object parameters via control logic.
- enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:Yes
Description:Set the enabled status of the MooseObject.
- implicitTrueDetermines whether this object is calculated using an implicit or explicit form
Default:True
C++ Type:bool
Controllable:No
Description:Determines whether this object is calculated using an implicit or explicit form
- use_displaced_meshFalseWhether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Default:False
C++ Type:bool
Controllable:No
Description:Whether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Advanced Parameters
- extra_matrix_tagsThe extra tags for the matrices this Kernel should fill
C++ Type:std::vector<TagName>
Controllable:No
Description:The extra tags for the matrices this Kernel should fill
- extra_vector_tagsThe extra tags for the vectors this Kernel should fill
C++ Type:std::vector<TagName>
Controllable:No
Description:The extra tags for the vectors this Kernel should fill
- matrix_tagssystemThe tag for the matrices this Kernel should fill
Default:system
C++ Type:MultiMooseEnum
Options:nontime, system
Controllable:No
Description:The tag for the matrices this Kernel should fill
- vector_tagsnontimeThe tag for the vectors this Kernel should fill
Default:nontime
C++ Type:MultiMooseEnum
Options:nontime, time
Controllable:No
Description:The tag for the vectors this Kernel should fill
Tagging Parameters
Input Files
- (modules/navier_stokes/test/tests/finite_volume/cns/mms/1d/hllc-mms.i)
- (modules/navier_stokes/test/tests/finite_volume/cns/symmetry_test/2D_symmetry.i)
- (modules/navier_stokes/test/tests/finite_volume/cns/mms/1d-with-bcs/free-flow-hllc.i)
- (modules/navier_stokes/test/tests/finite_volume/cns/stagnation_inlet/supersonic_nozzle_hllc.i)
- (modules/navier_stokes/test/tests/finite_volume/cns/mms/1d-with-bcs/regular-straight-channel.i)
- (modules/navier_stokes/test/tests/finite_volume/cns/implicit_bcs/hllc_sod_shocktube.i)
- (modules/navier_stokes/test/tests/finite_volume/cns/benchmark_shock_tube_1D/hllc_sod_shocktube.i)
- (modules/navier_stokes/test/tests/finite_volume/cns/shock_tube_2D_cavity/hllc_sod_shocktube_2D.i)
- (modules/navier_stokes/test/tests/finite_volume/cns/pressure_outlet/subsonic_nozzle_fixed_inflow_hllc.i)
- (modules/navier_stokes/test/tests/finite_volume/cns/mms/1d-with-bcs/straight-channel-hllc.i)
(modules/navier_stokes/test/tests/finite_volume/cns/mms/1d/hllc-mms.i)
[GlobalParams]
fp = fp
[]
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 1
xmin = .1
xmax = 1
nx = 2
[]
[]
[Variables]
[rho]
type = MooseVariableFVReal
[]
[rho_u]
type = MooseVariableFVReal
[]
[rho_et]
type = MooseVariableFVReal
[]
[]
[ICs]
[rho]
type = FunctionIC
variable = rho
function = 'cos(1.1*x)'
[]
[rho_u]
type = FunctionIC
variable = rho_u
function = '2*sin(1.1*x)'
[]
[rho_et]
type = FunctionIC
variable = rho_et
function = '3*cos(1.1*x)'
[]
[]
[FVKernels]
[mass_advection]
type = CNSFVMassHLLC
variable = rho
[]
[fake_diffusivity]
type = FVDiffusion
variable = rho
coeff = 1
[]
[mass_fn]
type = FVBodyForce
variable = rho
function = 'forcing_rho'
[]
[momentum_advection]
type = CNSFVMomentumHLLC
variable = rho_u
momentum_component = x
[]
[viscosity]
type = FVDiffusion
variable = rho_u
coeff = 1
[]
[momentum_fn]
type = FVBodyForce
variable = rho_u
function = 'forcing_rho_u'
[]
[fluid_energy_advection]
type = CNSFVFluidEnergyHLLC
variable = rho_et
[]
[fake_conduction]
type = FVDiffusion
variable = rho_et
coeff = 1
[]
[energy_fn]
type = FVBodyForce
variable = rho_et
function = 'forcing_rho_et'
[]
[]
[FVBCs]
[rho]
type = FVFunctionDirichletBC
boundary = 'left right'
variable = rho
function = 'exact_rho'
[]
[rho_u]
type = FVFunctionDirichletBC
boundary = 'left right'
variable = rho_u
function = 'exact_rho_u'
[]
[rho_et]
type = FVFunctionDirichletBC
boundary = 'left right'
variable = rho_et
function = 'exact_rho_et'
[]
[]
[Materials]
[var_mat]
type = ConservedVarValuesMaterial
rho = rho
rhou = rho_u
rho_et = rho_et
[]
[]
[Modules]
[FluidProperties]
[fp]
type = TestConservedVarFluidProperties
[]
[]
[]
[Functions]
[exact_rho]
type = ParsedFunction
value = 'cos(x)'
[]
[forcing_rho]
type = ParsedFunction
value = '3*cos(x)'
[]
[exact_rho_u]
type = ParsedFunction
value = '2*sin(x)'
[]
[forcing_rho_u]
type = ParsedFunction
value = '4*sin(x)^3/cos(x)^2 + 9*sin(x)'
[]
[exact_rho_et]
type = ParsedFunction
value = '3*cos(x)'
[]
[forcing_rho_et]
type = ParsedFunction
value = '11*cos(x)'
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
nl_rel_tol = 1e-12
[]
[Outputs]
exodus = true
csv = true
[]
[Postprocessors]
[h]
type = AverageElementSize
outputs = 'console csv'
execute_on = 'timestep_end'
[]
[L2rho]
type = ElementL2Error
variable = rho
function = exact_rho
outputs = 'console csv'
execute_on = 'timestep_end'
[]
[L2rho_u]
variable = rho_u
function = exact_rho_u
type = ElementL2Error
outputs = 'console csv'
execute_on = 'timestep_end'
[]
[L2rho_et]
variable = rho_et
function = exact_rho_et
type = ElementL2Error
outputs = 'console csv'
execute_on = 'timestep_end'
[]
[]
(modules/navier_stokes/test/tests/finite_volume/cns/symmetry_test/2D_symmetry.i)
rho_inside = 1
E_inside = 2.501505578
rho_outside = 0.125
E_outside = 1.999770935
radius = 0.1
angle = 45
[GlobalParams]
fp = fp
[]
[Debug]
show_material_props = true
[]
[Mesh]
[file]
type = GeneratedMeshGenerator
dim = 2
xmin = -0.5
xmax = 0.5
nx = 10
ymin = -0.5
ymax = 0.5
ny = 10
[../]
[rotate]
type = TransformGenerator
vector_value = '${angle} 0 0'
transform = ROTATE
input = file
[]
[]
[Modules]
[FluidProperties]
[fp]
type = IdealGasFluidProperties
allow_imperfect_jacobians = true
[]
[]
[]
[Variables]
[rho]
family = MONOMIAL
order = CONSTANT
fv = true
[../]
[rho_u]
family = MONOMIAL
order = CONSTANT
fv = true
initial_condition = 1e-15
outputs = none
[]
[rho_v]
family = MONOMIAL
order = CONSTANT
fv = true
initial_condition = 1e-15
outputs = none
[]
[rho_E]
family = MONOMIAL
order = CONSTANT
fv = true
[]
[]
[ICs]
[rho_ic]
type = FunctionIC
variable = rho
function = 'if (abs(x) < ${radius} & abs(y) < ${radius}, ${rho_inside}, ${rho_outside})'
[]
[rho_E_ic]
type = FunctionIC
variable = rho_E
function = 'if (abs(x) < ${radius} & abs(y) < ${radius}, ${fparse E_inside * rho_inside}, ${fparse E_outside * rho_outside})'
[]
[]
[FVKernels]
# Mass conservation
[mass_time]
type = FVTimeKernel
variable = rho
[]
[mass_advection]
type = CNSFVMassHLLC
variable = rho
fp = fp
[]
# Momentum x conservation
[momentum_x_time]
type = FVTimeKernel
variable = rho_u
[]
[momentum_x_advection]
type = CNSFVMomentumHLLC
variable = rho_u
momentum_component = x
fp = fp
[]
# Momentum y conservation
[momentum_y_time]
type = FVTimeKernel
variable = rho_v
[]
[./momentum_y_advection]
type = CNSFVMomentumHLLC
variable = rho_v
momentum_component = y
[]
# Fluid energy conservation
[./fluid_energy_time]
type = FVTimeKernel
variable = rho_E
[]
[./fluid_energy_advection]
type = CNSFVFluidEnergyHLLC
variable = rho_E
fp = fp
[]
[]
[FVBCs]
## outflow implicit conditions
[mass_outflow]
type = CNSFVHLLCMassImplicitBC
variable = rho
fp = fp
boundary = 'left right top bottom'
[]
[./momentum_x_outflow]
type = CNSFVHLLCMomentumImplicitBC
variable = rho_u
momentum_component = x
fp = fp
boundary = 'left right top bottom'
[]
[momentum_y_outflow]
type = CNSFVHLLCMomentumImplicitBC
variable = rho_v
momentum_component = y
fp = fp
boundary = 'left right top bottom'
[]
[fluid_energy_outflow]
type = CNSFVHLLCFluidEnergyImplicitBC
variable = rho_E
fp = fp
boundary = 'left right top bottom'
[]
[]
[AuxVariables]
[Ma]
family = MONOMIAL
order = CONSTANT
[]
[p]
family = MONOMIAL
order = CONSTANT
[]
[]
[AuxKernels]
[Ma_aux]
type = NSMachAux
variable = Ma
fluid_properties = fp
use_material_properties = true
[]
[p_aux]
type = ADMaterialRealAux
variable = p
property = pressure
[]
[]
[Materials]
[var_mat]
type = ConservedVarValuesMaterial
rho = rho
rhou = rho_u
rhov = rho_v
rho_et = rho_E
[]
[sound_speed]
type = SoundspeedMat
fp = fp
[]
[]
[Postprocessors]
[cfl_dt]
type = ADCFLTimeStepSize
c_names = 'sound_speed'
vel_names = 'speed'
CFL = 0.5
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
[]
[]
[Executioner]
type = Transient
end_time = 0.2
[TimeIntegrator]
type = ExplicitSSPRungeKutta
order = 2
[]
l_tol = 1e-8
[TimeStepper]
type = PostprocessorDT
postprocessor = cfl_dt
[]
[]
[Outputs]
exodus = true
[]
(modules/navier_stokes/test/tests/finite_volume/cns/mms/1d-with-bcs/free-flow-hllc.i)
[GlobalParams]
fp = fp
[]
[Mesh]
[cartesian]
type = GeneratedMeshGenerator
dim = 1
xmin = .1
xmax = 1.1
nx = 2
[]
[]
[Modules]
[FluidProperties]
[fp]
type = IdealGasFluidProperties
[]
[]
[]
[Variables]
[rho]
type = MooseVariableFVReal
[]
[rho_u]
type = MooseVariableFVReal
[]
[rho_et]
type = MooseVariableFVReal
[]
[]
[ICs]
[rho]
type = FunctionIC
variable = rho
function = 'exact_rho'
[]
[rho_u]
type = FunctionIC
variable = rho_u
function = 'exact_rho_u'
[]
[rho_et]
type = FunctionIC
variable = rho_et
function = 'exact_rho_et'
[]
[]
[FVKernels]
[mass_advection]
type = CNSFVMassHLLC
variable = rho
[]
[mass_fn]
type = FVBodyForce
variable = rho
function = 'forcing_rho'
[]
[momentum_x_advection]
type = CNSFVMomentumHLLC
variable = rho_u
momentum_component = x
[]
[momentum_fn]
type = FVBodyForce
variable = rho_u
function = 'forcing_rho_u'
[]
[fluid_energy_advection]
type = CNSFVFluidEnergyHLLC
variable = rho_et
[]
[energy_fn]
type = FVBodyForce
variable = rho_et
function = 'forcing_rho_et'
[]
[]
[FVBCs]
[mass_in]
variable = rho
type = CNSFVHLLCSpecifiedMassFluxAndTemperatureMassBC
boundary = left
temperature = 'exact_T'
rhou = 'exact_rho_u'
[]
[momentum_in]
variable = rho_u
type = CNSFVHLLCSpecifiedMassFluxAndTemperatureMomentumBC
boundary = left
temperature = 'exact_T'
rhou = 'exact_rho_u'
momentum_component = 'x'
[]
[energy_in]
variable = rho_et
type = CNSFVHLLCSpecifiedMassFluxAndTemperatureFluidEnergyBC
boundary = left
temperature = 'exact_T'
rhou = 'exact_rho_u'
[]
[mass_out]
variable = rho
type = CNSFVHLLCSpecifiedPressureMassBC
boundary = right
pressure = 'exact_p'
[]
[momentum_out]
variable = rho_u
type = CNSFVHLLCSpecifiedPressureMomentumBC
boundary = right
pressure = 'exact_p'
momentum_component = 'x'
[]
[energy_out]
variable = rho_et
type = CNSFVHLLCSpecifiedPressureFluidEnergyBC
boundary = right
pressure = 'exact_p'
[]
[]
[Materials]
[var_mat]
type = ConservedVarValuesMaterial
rho = rho
rhou = rho_u
rho_et = rho_et
[]
[]
[Functions]
[exact_rho]
type = ParsedFunction
value = '3.48788261470924*cos(x)'
[]
[forcing_rho]
type = ParsedFunction
value = '-3.83667087618017*sin(1.1*x)'
[]
[exact_rho_u]
type = ParsedFunction
value = '3.48788261470924*cos(1.1*x)'
[]
[forcing_rho_u]
type = ParsedFunction
value = '-(10.6975765229419*cos(1.2*x)/cos(x) - 0.697576522941849*cos(1.1*x)^2/cos(x)^2)*sin(x) + (10.6975765229419*sin(x)*cos(1.2*x)/cos(x)^2 - 1.3951530458837*sin(x)*cos(1.1*x)^2/cos(x)^3 + 1.53466835047207*sin(1.1*x)*cos(1.1*x)/cos(x)^2 - 12.8370918275302*sin(1.2*x)/cos(x))*cos(x) + 3.48788261470924*sin(x)*cos(1.1*x)^2/cos(x)^2 - 7.67334175236034*sin(1.1*x)*cos(1.1*x)/cos(x)'
[]
[exact_rho_et]
type = ParsedFunction
value = '26.7439413073546*cos(1.2*x)'
[]
[forcing_rho_et]
type = ParsedFunction
value = '1.0*(3.48788261470924*(3.06706896551724*cos(1.2*x)/cos(x) - 0.2*cos(1.1*x)^2/cos(x)^2)*cos(x) + 26.7439413073546*cos(1.2*x))*sin(x)*cos(1.1*x)/cos(x)^2 - 1.1*(3.48788261470924*(3.06706896551724*cos(1.2*x)/cos(x) - 0.2*cos(1.1*x)^2/cos(x)^2)*cos(x) + 26.7439413073546*cos(1.2*x))*sin(1.1*x)/cos(x) + 1.0*(-(10.6975765229419*cos(1.2*x)/cos(x) - 0.697576522941849*cos(1.1*x)^2/cos(x)^2)*sin(x) + (10.6975765229419*sin(x)*cos(1.2*x)/cos(x)^2 - 1.3951530458837*sin(x)*cos(1.1*x)^2/cos(x)^3 + 1.53466835047207*sin(1.1*x)*cos(1.1*x)/cos(x)^2 - 12.8370918275302*sin(1.2*x)/cos(x))*cos(x) - 32.0927295688256*sin(1.2*x))*cos(1.1*x)/cos(x)'
[]
[exact_T]
type = ParsedFunction
value = '0.0106975765229418*cos(1.2*x)/cos(x) - 0.000697576522941848*cos(1.1*x)^2/cos(x)^2'
[]
[exact_p]
type = ParsedFunction
value = '3.48788261470924*(3.06706896551724*cos(1.2*x)/cos(x) - 0.2*cos(1.1*x)^2/cos(x)^2)*cos(x)'
[]
[]
[Executioner]
solve_type = NEWTON
type = Steady
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
nl_max_its = 50
line_search = none
nl_rel_tol = 1e-12
nl_abs_tol = 1e-12
[]
[Outputs]
exodus = true
csv = true
[]
[Debug]
show_var_residual_norms = true
[]
[Postprocessors]
[h]
type = AverageElementSize
outputs = 'console csv'
execute_on = 'timestep_end'
[]
[L2rho]
type = ElementL2Error
variable = rho
function = exact_rho
outputs = 'console csv'
execute_on = 'timestep_end'
[]
[L2rho_u]
variable = rho_u
function = exact_rho_u
type = ElementL2Error
outputs = 'console csv'
execute_on = 'timestep_end'
[]
[L2rho_et]
variable = rho_et
function = exact_rho_et
type = ElementL2Error
outputs = 'console csv'
execute_on = 'timestep_end'
[]
[]
(modules/navier_stokes/test/tests/finite_volume/cns/stagnation_inlet/supersonic_nozzle_hllc.i)
stagnation_pressure = 1
stagnation_temperature = 1
[GlobalParams]
fp = fp
[]
[Debug]
show_material_props = true
[]
[Mesh]
[file]
type = FileMeshGenerator
file = supersonic_nozzle.e
[]
[]
[Modules]
[FluidProperties]
[fp]
type = IdealGasFluidProperties
[]
[]
[]
[Variables]
[rho]
family = MONOMIAL
order = CONSTANT
fv = true
initial_condition = 0.0034
[]
[rho_u]
family = MONOMIAL
order = CONSTANT
fv = true
initial_condition = 1e-4
outputs = none
[]
[rho_v]
family = MONOMIAL
order = CONSTANT
fv = true
outputs = none
[]
[rho_E]
family = MONOMIAL
order = CONSTANT
fv = true
initial_condition = 2.5
[]
[]
[FVKernels]
# Mass conservation
[mass_time]
type = FVTimeKernel
variable = rho
[]
[mass_advection]
type = CNSFVMassHLLC
variable = rho
[]
# Momentum x conservation
[momentum_x_time]
type = FVTimeKernel
variable = rho_u
[]
[momentum_x_advection]
type = CNSFVMomentumHLLC
variable = rho_u
momentum_component = x
[]
# Momentum y conservation
[momentum_y_time]
type = FVTimeKernel
variable = rho_v
[]
[momentum_y_advection]
type = CNSFVMomentumHLLC
variable = rho_v
momentum_component = y
[]
# Fluid energy conservation
[fluid_energy_time]
type = FVTimeKernel
variable = rho_E
[]
[fluid_energy_advection]
type = CNSFVFluidEnergyHLLC
variable = rho_E
[]
[]
[FVBCs]
## inflow stagnation boundaries
[mass_stagnation_inflow]
type = CNSFVHLLCMassStagnationInletBC
variable = rho
stagnation_pressure = ${stagnation_pressure}
stagnation_temperature = ${stagnation_temperature}
boundary = left
[]
[momentum_x_stagnation_inflow]
type = CNSFVHLLCMomentumStagnationInletBC
variable = rho_u
momentum_component = x
stagnation_pressure = ${stagnation_pressure}
stagnation_temperature = ${stagnation_temperature}
boundary = left
[]
[momentum_y_stagnation_inflow]
type = CNSFVHLLCMomentumStagnationInletBC
variable = rho_v
momentum_component = y
stagnation_pressure = ${stagnation_pressure}
stagnation_temperature = ${stagnation_temperature}
boundary = left
[../]
[fluid_energy_stagnation_inflow]
type = CNSFVHLLCFluidEnergyStagnationInletBC
variable = rho_E
stagnation_pressure = ${stagnation_pressure}
stagnation_temperature = ${stagnation_temperature}
boundary = left
[]
## outflow implicit conditions
[mass_outflow]
type = CNSFVHLLCMassImplicitBC
variable = rho
boundary = right
[]
[momentum_x_outflow]
type = CNSFVHLLCMomentumImplicitBC
variable = rho_u
momentum_component = x
boundary = right
[]
[momentum_y_outflow]
type = CNSFVHLLCMomentumImplicitBC
variable = rho_v
momentum_component = y
boundary = right
[]
[fluid_energy_outflow]
type = CNSFVHLLCFluidEnergyImplicitBC
variable = rho_E
boundary = right
[]
# wall conditions
[momentum_x_pressure_wall]
type = CNSFVMomImplicitPressureBC
variable = rho_u
momentum_component = x
boundary = wall
[]
[momentum_y_pressure_wall]
type = CNSFVMomImplicitPressureBC
variable = rho_v
momentum_component = y
boundary = wall
[]
[]
[AuxVariables]
[Ma]
family = MONOMIAL
order = CONSTANT
[]
[Ma_layered]
family = MONOMIAL
order = CONSTANT
[]
[]
[UserObjects]
[layered_Ma_UO]
type = LayeredAverage
variable = Ma
num_layers = 100
direction = x
[]
[]
[AuxKernels]
[Ma_aux]
type = NSMachAux
variable = Ma
fluid_properties = fp
use_material_properties = true
[]
[Ma_layered_aux]
type = SpatialUserObjectAux
variable = Ma_layered
user_object = layered_Ma_UO
[]
[]
[Materials]
[var_mat]
type = ConservedVarValuesMaterial
rho = rho
rhou = rho_u
rhov = rho_v
rho_et = rho_E
[]
[fluid_props]
type = GeneralFluidProps
porosity = 1
characteristic_length = 1
[]
[sound_speed]
type = SoundspeedMat
fp = fp
[]
[]
[Postprocessors]
[cfl_dt]
type = ADCFLTimeStepSize
c_names = 'sound_speed'
vel_names = 'speed'
CFL = 0.5
[]
[outflow_Ma]
type = SideAverageValue
variable = Ma
boundary = right
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
[]
[]
[Executioner]
type = Transient
end_time = 0.1
[TimeIntegrator]
type = ExplicitSSPRungeKutta
order = 2
[]
l_tol = 1e-8
[TimeStepper]
type = PostprocessorDT
postprocessor = cfl_dt
[]
[]
[VectorPostprocessors]
[Ma_layered]
type = LineValueSampler
variable = Ma_layered
start_point = '0 0 0'
end_point = '10 0 0'
num_points = 100
sort_by = x
[]
[]
[Outputs]
exodus = true
[]
(modules/navier_stokes/test/tests/finite_volume/cns/mms/1d-with-bcs/regular-straight-channel.i)
[GlobalParams]
fp = fp
[]
[Mesh]
[./gen_mesh]
type = GeneratedMeshGenerator
dim = 1
xmax = 1.5
nx = 15
[../]
[]
[Modules]
[FluidProperties]
[fp]
type = IdealGasFluidProperties
[]
[]
[]
[Variables]
[rho]
type = MooseVariableFVReal
initial_condition = 1.28969
scaling = 1e3
[]
[rho_u]
type = MooseVariableFVReal
initial_condition = 1.28969
[]
[rho_et]
type = MooseVariableFVReal
initial_condition = 2.525e5
scaling = 1e-2
[]
[]
[FVKernels]
[mass_advection]
type = CNSFVMassHLLC
variable = rho
fp = fp
[]
[momentum_x_advection]
type = CNSFVMomentumHLLC
variable = rho_u
momentum_component = x
fp = fp
[]
[drag]
type = FVReaction
variable = rho_u
rate = 1000
[]
[fluid_energy_advection]
type = CNSFVFluidEnergyHLLC
variable = rho_et
fp = fp
[]
[]
[FVBCs]
[mass_in]
variable = rho
type = CNSFVHLLCSpecifiedMassFluxAndTemperatureMassBC
boundary = left
temperature = 273.15
rhou = 1.28969
[]
[momentum_in]
variable = rho_u
type = CNSFVHLLCSpecifiedMassFluxAndTemperatureMomentumBC
boundary = left
temperature = 273.15
rhou = 1.28969
momentum_component = 'x'
[]
[energy_in]
variable = rho_et
type = CNSFVHLLCSpecifiedMassFluxAndTemperatureFluidEnergyBC
boundary = left
temperature = 273.15
rhou = 1.28969
[]
[mass_out]
variable = rho
type = CNSFVHLLCSpecifiedPressureMassBC
boundary = right
pressure = 1.01e5
[]
[momentum_out]
variable = rho_u
type = CNSFVHLLCSpecifiedPressureMomentumBC
boundary = right
pressure = 1.01e5
momentum_component = 'x'
[]
[energy_out]
variable = rho_et
type = CNSFVHLLCSpecifiedPressureFluidEnergyBC
boundary = right
pressure = 1.01e5
[]
[]
[Materials]
[var_mat]
type = ConservedVarValuesMaterial
rho = rho
rhou = rho_u
rho_et = rho_et
[]
[]
[Executioner]
solve_type = NEWTON
type = Steady
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
nl_max_its = 50
line_search = none
[]
[Outputs]
exodus = true
csv = true
[]
[Debug]
show_var_residual_norms = true
[]
(modules/navier_stokes/test/tests/finite_volume/cns/implicit_bcs/hllc_sod_shocktube.i)
rho_left = 1
E_left = 2.501505578
u_left = 1e-15
rho_right = 0.125
E_right = 1.999770935
u_right = 1e-15
middle = 0.5
[GlobalParams]
fp = fp
[]
[Mesh]
[cartesian]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = ${fparse 2 * middle}
nx = 5
ymin = 0
ymax = 1
ny = 2
[]
[]
[Modules]
[FluidProperties]
[fp]
type = IdealGasFluidProperties
allow_imperfect_jacobians = true
[]
[]
[]
[Variables]
[rho]
order = CONSTANT
family = MONOMIAL
fv = true
[]
[rho_u]
order = CONSTANT
family = MONOMIAL
fv = true
[]
[rho_v]
order = CONSTANT
family = MONOMIAL
fv = true
initial_condition = 1e-10
[]
[rho_E]
order = CONSTANT
family = MONOMIAL
fv = true
[]
[]
[FVKernels]
[mass_time]
type = FVTimeKernel
variable = rho
[]
[mass_advection]
type = CNSFVMassHLLC
variable = rho
[]
[momentum_x_time]
type = FVTimeKernel
variable = rho_u
[]
[momentum_x_advection]
type = CNSFVMomentumHLLC
variable = rho_u
momentum_component = x
[]
[momentum_y_time]
type = FVTimeKernel
variable = rho_v
[]
[momentum_y_advection]
type = CNSFVMomentumHLLC
variable = rho_v
momentum_component = y
[]
[fluid_energy_time]
type = FVTimeKernel
variable = rho_E
[]
[fluid_energy_advection]
type = CNSFVFluidEnergyHLLC
variable = rho_E
[]
[]
[FVBCs]
[mass_implicit]
type = CNSFVHLLCMassImplicitBC
variable = rho
fp = fp
boundary = 'left right'
[]
[mom_x_implicit]
type = CNSFVHLLCMomentumImplicitBC
variable = rho_u
momentum_component = x
fp = fp
boundary = 'left right'
[]
[wall]
type = CNSFVMomImplicitPressureBC
variable = rho_v
momentum_component = y
boundary = 'top bottom'
[]
[fluid_energy_implicit]
type = CNSFVHLLCFluidEnergyImplicitBC
variable = rho_E
fp = fp
boundary = 'left right'
[]
[]
[ICs]
[rho_ic]
type = FunctionIC
variable = rho
function = 'if (x < ${middle}, ${rho_left}, ${rho_right})'
[]
[rho_u_ic]
type = FunctionIC
variable = rho_u
function = 'if (x < ${middle}, ${fparse rho_left * u_left}, ${fparse rho_right * u_right})'
[]
[rho_E_ic]
type = FunctionIC
variable = rho_E
function = 'if (x < ${middle}, ${fparse E_left * rho_left}, ${fparse E_right * rho_right})'
[]
[]
[Materials]
[var_mat]
type = ConservedVarValuesMaterial
rho = rho
rhou = rho_u
rhov = rho_v
rho_et = rho_E
fp = fp
[]
[]
[Executioner]
type = Transient
[TimeIntegrator]
type = ExplicitSSPRungeKutta
order = 2
[]
l_tol = 1e-8
# run to t = 0.15
start_time = 0.0
dt = 1e-1
end_time = 10
abort_on_solve_fail = true
[]
[Outputs]
exodus = true
[]
(modules/navier_stokes/test/tests/finite_volume/cns/benchmark_shock_tube_1D/hllc_sod_shocktube.i)
rho_left = 1
E_left = 2.501505578
u_left = 1e-15
rho_right = 0.125
E_right = 1.999770935
u_right = 1e-15
middle = 50
[GlobalParams]
fp = fp
[]
[Mesh]
[cartesian]
type = GeneratedMeshGenerator
dim = 1
xmin = 0
xmax = ${fparse 2 * middle}
nx = 1000
[]
[]
[Modules]
[FluidProperties]
[fp]
type = IdealGasFluidProperties
[]
[]
[]
[Variables]
[rho]
order = CONSTANT
family = MONOMIAL
fv = true
[]
[rho_u]
order = CONSTANT
family = MONOMIAL
fv = true
[]
[rho_E]
order = CONSTANT
family = MONOMIAL
fv = true
[]
[]
[AuxVariables]
[rho_a]
order = CONSTANT
family = MONOMIAL
[]
[]
[FVKernels]
[mass_time]
type = FVTimeKernel
variable = rho
[]
[mass_advection]
type = CNSFVMassHLLC
variable = rho
[]
[momentum_time]
type = FVTimeKernel
variable = rho_u
[]
[momentum_advection]
type = CNSFVMomentumHLLC
variable = rho_u
momentum_component = x
[]
[fluid_energy_time]
type = FVTimeKernel
variable = rho_E
[../]
[fluid_energy_advection]
type = CNSFVFluidEnergyHLLC
variable = rho_E
[]
[]
[FVBCs]
[mass_implicit]
type = CNSFVHLLCMassImplicitBC
variable = rho
fp = fp
boundary = 'left right'
[]
[mom_implicit]
type = CNSFVHLLCMomentumImplicitBC
variable = rho_u
momentum_component = x
fp = fp
boundary = 'left right'
[]
[fluid_energy_implicit]
type = CNSFVHLLCFluidEnergyImplicitBC
variable = rho_E
fp = fp
boundary = 'left right'
[]
[]
[ICs]
[rho_ic]
type = FunctionIC
variable = rho
function = 'if (x < ${middle}, ${rho_left}, ${rho_right})'
[]
[rho_u_ic]
type = FunctionIC
variable = rho_u
function = 'if (x < ${middle}, ${fparse rho_left * u_left}, ${fparse rho_right * u_right})'
[]
[rho_E_ic]
type = FunctionIC
variable = rho_E
function = 'if (x < ${middle}, ${fparse E_left * rho_left}, ${fparse E_right * rho_right})'
[]
[]
[Materials]
[var_mat]
type = ConservedVarValuesMaterial
rho = rho
rhou = rho_u
rho_et = rho_E
fp = fp
[]
[]
[Preconditioning]
active = ''
[./smp]
type = SMP
full = true
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
[../]
[]
[Executioner]
type = Transient
[TimeIntegrator]
type = ExplicitSSPRungeKutta
order = 2
[]
l_tol = 1e-8
start_time = 0.0
dt = 1e-2
end_time = 20
abort_on_solve_fail = true
[]
[Outputs]
exodus = true
perf_graph = true
[]
(modules/navier_stokes/test/tests/finite_volume/cns/shock_tube_2D_cavity/hllc_sod_shocktube_2D.i)
rho_left = 1
E_left = 2.501505578
u_left = 1e-15
rho_right = 0.125
E_right = 1.999770935
u_right = 1e-15
x_sep = 35
[GlobalParams]
fp = fp
[]
[Mesh]
[./cartesian]
type = CartesianMeshGenerator
dim = 2
dx = '40 20'
ix = '200 100'
dy = '1 20 2 20 1'
iy = '4 100 10 100 4'
subdomain_id = '0 0
0 1
1 1
0 1
0 0'
[../]
[./wall]
type = SideSetsBetweenSubdomainsGenerator
input = cartesian
primary_block = 1
paired_block = 0
new_boundary = 'wall'
[../]
[./delete]
type = BlockDeletionGenerator
input = wall
block = 0
[../]
[]
[Modules]
[./FluidProperties]
[./fp]
type = IdealGasFluidProperties
allow_imperfect_jacobians = true
[../]
[../]
[]
[Variables]
[./rho]
order = CONSTANT
family = MONOMIAL
fv = true
[../]
[./rho_u]
order = CONSTANT
family = MONOMIAL
fv = true
[../]
[./rho_v]
order = CONSTANT
family = MONOMIAL
fv = true
[../]
[./rho_E]
order = CONSTANT
family = MONOMIAL
fv = true
[../]
[]
[AuxVariables]
[./Ma]
order = CONSTANT
family = MONOMIAL
[../]
[./p]
order = CONSTANT
family = MONOMIAL
[../]
[./v_norm]
order = CONSTANT
family = MONOMIAL
[../]
[./temperature]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./Ma_aux]
type = NSMachAux
variable = Ma
fluid_properties = fp
use_material_properties = true
[../]
[./p_aux]
type = ADMaterialRealAux
variable = p
property = pressure
[../]
[./v_norm_aux]
type = ADMaterialRealAux
variable = v_norm
property = speed
[../]
[./temperature_aux]
type = ADMaterialRealAux
variable = temperature
property = T_fluid
[../]
[]
[FVKernels]
[./mass_time]
type = FVTimeKernel
variable = rho
[../]
[./mass_advection]
type = CNSFVMassHLLC
variable = rho
[../]
[./momentum_x_time]
type = FVTimeKernel
variable = rho_u
[../]
[./momentum_x_advection]
type = CNSFVMomentumHLLC
variable = rho_u
momentum_component = x
[../]
[./momentum_y_time]
type = FVTimeKernel
variable = rho_v
[../]
[./momentum_y_advection]
type = CNSFVMomentumHLLC
variable = rho_v
momentum_component = y
[../]
[./fluid_energy_time]
type = FVTimeKernel
variable = rho_E
[../]
[./fluid_energy_advection]
type = CNSFVFluidEnergyHLLC
variable = rho_E
[../]
[]
[FVBCs]
[./mom_x_pressure]
type = CNSFVMomImplicitPressureBC
variable = rho_u
momentum_component = x
boundary = 'left right wall'
[../]
[./mom_y_pressure]
type = CNSFVMomImplicitPressureBC
variable = rho_v
momentum_component = y
boundary = 'wall'
[../]
[]
[ICs]
[./rho_ic]
type = FunctionIC
variable = rho
function = 'if (x < ${x_sep}, ${rho_left}, ${rho_right})'
[../]
[./rho_u_ic]
type = FunctionIC
variable = rho_u
function = 'if (x < ${x_sep}, ${fparse rho_left * u_left}, ${fparse rho_right * u_right})'
[../]
[./rho_E_ic]
type = FunctionIC
variable = rho_E
function = 'if (x < ${x_sep}, ${fparse E_left * rho_left}, ${fparse E_right * rho_right})'
[../]
[]
[Materials]
[./var_mat]
type = ConservedVarValuesMaterial
rho = rho
rhou = rho_u
rhov = rho_v
rho_et = rho_E
fp = fp
[../]
[./sound_speed]
type = SoundspeedMat
fp = fp
[../]
[]
[Preconditioning]
[./smp]
type = SMP
full = true
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
[../]
[]
[Postprocessors]
[./cfl_dt]
type = ADCFLTimeStepSize
c_names = 'sound_speed'
vel_names = 'speed'
[../]
[]
[Executioner]
type = Transient
end_time = 100
[TimeIntegrator]
type = ExplicitSSPRungeKutta
order = 2
[]
l_tol = 1e-8
[./TimeStepper]
type = PostprocessorDT
postprocessor = cfl_dt
[../]
[]
[Outputs]
exodus = true
[]
(modules/navier_stokes/test/tests/finite_volume/cns/pressure_outlet/subsonic_nozzle_fixed_inflow_hllc.i)
inlet_vel = 120
rho_in = 0.8719748696
H_in = 4.0138771448e+05
gamma = 1.4
R = 8.3145
molar_mass = 29e-3
R_specific = ${fparse R / molar_mass}
cp = ${fparse gamma * R_specific / (gamma - 1)}
cv = ${fparse cp / gamma}
T_in = ${fparse H_in / gamma / cv}
mass_flux = ${fparse inlet_vel * rho_in}
outlet_pressure = 0.9e5
[GlobalParams]
fp = fp
[]
[Debug]
show_material_props = true
[]
[Mesh]
[file]
type = FileMeshGenerator
file = subsonic_nozzle.e
[]
[]
[Modules]
[FluidProperties]
[fp]
type = IdealGasFluidProperties
[]
[]
[]
[Variables]
[rho]
family = MONOMIAL
order = CONSTANT
fv = true
initial_condition = 0.8719748696
[]
[rho_u]
family = MONOMIAL
order = CONSTANT
fv = true
initial_condition = 1e-4
[]
[rho_v]
family = MONOMIAL
order = CONSTANT
fv = true
[]
[rho_E]
family = MONOMIAL
order = CONSTANT
fv = true
initial_condition = 2.5e5
[]
[]
[FVKernels]
# Mass conservation
[mass_time]
type = FVTimeKernel
variable = rho
[]
[mass_advection]
type = CNSFVMassHLLC
variable = rho
[]
# Momentum x conservation
[momentum_x_time]
type = FVTimeKernel
variable = rho_u
[]
[momentum_x_advection]
type = CNSFVMomentumHLLC
variable = rho_u
momentum_component = x
[]
# Momentum y conservation
[momentum_y_time]
type = FVTimeKernel
variable = rho_v
[]
[momentum_y_advection]
type = CNSFVMomentumHLLC
variable = rho_v
momentum_component = y
[]
# Fluid energy conservation
[fluid_energy_time]
type = FVTimeKernel
variable = rho_E
[]
[fluid_energy_advection]
type = CNSFVFluidEnergyHLLC
variable = rho_E
[]
[]
[FVBCs]
## inflow boundaries
[mass_inflow]
type = CNSFVHLLCSpecifiedMassFluxAndTemperatureMassBC
variable = rho
boundary = left
rhou = ${mass_flux}
rhov = 0
temperature = ${T_in}
[]
[momentum_x_inflow]
type = CNSFVHLLCSpecifiedMassFluxAndTemperatureMomentumBC
variable = rho_u
boundary = left
rhou = ${mass_flux}
rhov = 0
temperature = ${T_in}
momentum_component = x
[]
[momentum_y_inflow]
type = CNSFVHLLCSpecifiedMassFluxAndTemperatureMomentumBC
variable = rho_v
boundary = left
rhou = ${mass_flux}
rhov = 0
temperature = ${T_in}
momentum_component = y
[]
[fluid_energy_inflow]
type = CNSFVHLLCSpecifiedMassFluxAndTemperatureFluidEnergyBC
variable = rho_E
boundary = left
rhou = ${mass_flux}
rhov = 0
temperature = ${T_in}
[]
## outflow conditions
[mass_outflow]
type = CNSFVHLLCSpecifiedPressureMassBC
variable = rho
boundary = right
pressure = ${outlet_pressure}
[]
[momentum_x_outflow]
type = CNSFVHLLCSpecifiedPressureMomentumBC
variable = rho_u
boundary = right
momentum_component = x
pressure = ${outlet_pressure}
[]
[momentum_y_outflow]
type = CNSFVHLLCSpecifiedPressureMomentumBC
variable = rho_v
boundary = right
momentum_component = y
pressure = ${outlet_pressure}
[]
[fluid_energy_outflow]
type = CNSFVHLLCSpecifiedPressureFluidEnergyBC
variable = rho_E
boundary = right
pressure = ${outlet_pressure}
[]
# wall conditions
[momentum_x_pressure_wall]
type = CNSFVMomImplicitPressureBC
variable = rho_u
momentum_component = x
boundary = wall
[]
[momentum_y_pressure_wall]
type = CNSFVMomImplicitPressureBC
variable = rho_v
momentum_component = y
boundary = wall
[]
[]
[AuxVariables]
[Ma]
family = MONOMIAL
order = CONSTANT
[]
[p]
family = MONOMIAL
order = CONSTANT
[]
[Ma_layered]
family = MONOMIAL
order = CONSTANT
[]
[]
[UserObjects]
[layered_Ma_UO]
type = LayeredAverage
variable = Ma
num_layers = 10
direction = x
[]
[]
[AuxKernels]
[Ma_aux]
type = NSMachAux
variable = Ma
fluid_properties = fp
use_material_properties = true
[]
[p_aux]
type = ADMaterialRealAux
variable = p
property = pressure
[]
[Ma_layered_aux]
type = SpatialUserObjectAux
variable = Ma_layered
user_object = layered_Ma_UO
[]
[]
[Materials]
[var_mat]
type = ConservedVarValuesMaterial
rho = rho
rhou = rho_u
rhov = rho_v
rho_et = rho_E
[]
[sound_speed]
type = SoundspeedMat
[]
[]
[Postprocessors]
[outflow_Ma]
type = SideAverageValue
variable = Ma
boundary = right
[]
[outflow_pressure]
type = SideAverageValue
variable = p
boundary = right
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
[]
[]
[Executioner]
type = Transient
end_time = 10
solve_type = NEWTON
nl_abs_tol = 1e-7
[TimeIntegrator]
type = ImplicitEuler
[]
[TimeStepper]
type = IterationAdaptiveDT
dt = 5e-3
optimal_iterations = 6
growth_factor = 1.5
[]
[]
[VectorPostprocessors]
[Ma_layered]
type = LineValueSampler
variable = Ma_layered
start_point = '0 0 0'
end_point = '3 0 0'
num_points = 100
sort_by = x
warn_discontinuous_face_values = false
[]
[]
[Outputs]
[out]
type = Exodus
execute_on = 'final'
[]
[]
(modules/navier_stokes/test/tests/finite_volume/cns/mms/1d-with-bcs/straight-channel-hllc.i)
[GlobalParams]
fp = fp
[]
[Mesh]
[./gen_mesh]
type = CartesianMeshGenerator
dim = 1
dx = '.1 .1 .1 .1 .1 .5 .1 .1 .1 .1 .1'
# dx = '.1 .1 .1 .1 .1 .1 .1 .1 .1 .1 .1 .1 .1 .1 .1'
[../]
[]
[Modules]
[FluidProperties]
[fp]
type = IdealGasFluidProperties
[]
[]
[]
[Variables]
[rho]
type = MooseVariableFVReal
initial_condition = 1.28969
scaling = 1e3
[]
[rho_u]
type = MooseVariableFVReal
initial_condition = 1.28969
[]
[rho_et]
type = MooseVariableFVReal
initial_condition = 2.525e5
scaling = 1e-2
[]
[]
[FVKernels]
[mass_advection]
type = CNSFVMassHLLC
variable = rho
fp = fp
[]
[momentum_x_advection]
type = CNSFVMomentumHLLC
variable = rho_u
momentum_component = x
fp = fp
[]
[drag]
type = FVReaction
variable = rho_u
rate = 1000
[]
[fluid_energy_advection]
type = CNSFVFluidEnergyHLLC
variable = rho_et
fp = fp
[]
[]
[FVBCs]
[mass_in]
variable = rho
type = CNSFVHLLCSpecifiedMassFluxAndTemperatureMassBC
boundary = left
temperature = 273.15
rhou = 1.28969
[]
[momentum_in]
variable = rho_u
type = CNSFVHLLCSpecifiedMassFluxAndTemperatureMomentumBC
boundary = left
temperature = 273.15
rhou = 1.28969
momentum_component = 'x'
[]
[energy_in]
variable = rho_et
type = CNSFVHLLCSpecifiedMassFluxAndTemperatureFluidEnergyBC
boundary = left
temperature = 273.15
rhou = 1.28969
[]
[mass_out]
variable = rho
type = CNSFVHLLCSpecifiedPressureMassBC
boundary = right
pressure = 1.01e5
[]
[momentum_out]
variable = rho_u
type = CNSFVHLLCSpecifiedPressureMomentumBC
boundary = right
pressure = 1.01e5
momentum_component = 'x'
[]
[energy_out]
variable = rho_et
type = CNSFVHLLCSpecifiedPressureFluidEnergyBC
boundary = right
pressure = 1.01e5
[]
[]
[Materials]
[var_mat]
type = ConservedVarValuesMaterial
rho = rho
rhou = rho_u
rho_et = rho_et
[]
[]
[Executioner]
solve_type = NEWTON
type = Steady
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
nl_max_its = 50
line_search = none
[]
[Outputs]
exodus = true
csv = true
[]
[Debug]
show_var_residual_norms = true
[]