- diffusion_tensorFunctor describing a diagonal diffusion tensor. A functor is any of the following: a variable, a functor material property, a function, a postprocessor or a number.
C++ Type:MooseFunctorName
Unit:(no unit assumed)
Controllable:No
Description:Functor describing a diagonal diffusion tensor. A functor is any of the following: a variable, a functor material property, a function, a postprocessor or a number.
- variableThe name of the variable whose linear system this object contributes to
C++ Type:LinearVariableName
Unit:(no unit assumed)
Controllable:No
Description:The name of the variable whose linear system this object contributes to
LinearFVAnisotropicDiffusion
Description
This kernel contributes to the system matrix and the right hand side of a system which is solved for a linear finite volume variable. The difference between this kernel and LinearFVDiffusion is that this kernel requires a vector of diffusion coefficients, where every entry describes the diffusion coefficient for a principal direction. This is equivalent to supplying a diagonal tensor in the fully anisotropic diffusion case.
The implementation in this kernel is based on the derivation in Liu et al. (2015). The contributions of the system matrix and right hand side can be derived using the divergence theorem on a volumetric integral over cell C:
VC∫∇⋅D∇udV=f∑Sf∫D∇u⋅ndS,where D denotes a space dependent diagonal diffusion tensor, while the right hand side describes the sum of the surface integrals on each side of cell C. Following Liu et al. (2015) and using the assumption that the tensor is diagonal, we can manipulate this expression to arrive to the following form:
f∑Sf∫Dn⋅∇udS,where mathbbDn can be split into two contributions:
Dn=(Dn⋅n)n+(Dn−Dn⋅nn).Plugging this expression back to the surface integral, we get the following:
f∑Sf∫(Dn⋅n)n⋅∇u+(Dn−Dn⋅nn)⋅∇udS,where we can treat the normal projection (Sf∫(Dn⋅n)n⋅∇udS) the same way as described in LinearFVDiffusion with Dn⋅n replacing the diffusion coefficient. The second term (Sf∫(Dn−Dn⋅nn)⋅∇udSdS) can be treated explicitly, similarly to the nonorthogonal correction in LinearFVDiffusion.
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
- force_boundary_executionFalseWhether to force execution of this object on all external boundaries.
Default:False
C++ Type:bool
Controllable:No
Description:Whether to force execution of this object on all external boundaries.
- matrix_onlyFalseWhether this object is only doing assembly to matrices (no vectors)
Default:False
C++ Type:bool
Controllable:No
Description:Whether this object is only doing assembly to matrices (no vectors)
- use_nonorthogonal_correctionTrueIf the nonorthogonal correction should be used when computing the normal gradient.
Default:True
C++ Type:bool
Controllable:No
Description:If the nonorthogonal correction should be used when computing the normal gradient.
- use_nonorthogonal_correction_on_boundaryFalseIf the nonorthogonal correction should be used when computing the normal gradient.
Default:False
C++ Type:bool
Controllable:No
Description:If the nonorthogonal correction should be used when computing the normal gradient.
Optional Parameters
- absolute_value_vector_tagsThe tags for the vectors this residual object should fill with the absolute value of the residual contribution
C++ Type:std::vector<TagName>
Controllable:No
Description:The tags for the vectors this residual object should fill with the absolute value of the residual contribution
- 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_tagsrhsThe tag for the vectors this Kernel should fill
Default:rhs
C++ Type:MultiMooseEnum
Options:rhs, time
Controllable:No
Description:The tag for the vectors this Kernel should fill
Contribution To Tagged Field Data 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
- seed0The seed for the master random number generator
Default:0
C++ Type:unsigned int
Controllable:No
Description:The seed for the master random number generator
- 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
- 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.
- 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.
Parallel Ghosting Parameters
Input Files
- (modules/navier_stokes/examples/flow-over-circle/flow_over_circle_linearfv.i)
- (modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/SIMPLEC/2d/2d-velocity-pressure.i)
- (modules/navier_stokes/test/tests/finite_volume/ins/mms/linear-segregated/2d-vortex/2d-vortex.i)
- (modules/navier_stokes/test/tests/finite_volume/ins/cht/flow-around-square-linear.i)
- (modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/3d/3d-velocity-pressure.i)
- (modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/channel.i)
- (modules/navier_stokes/test/tests/finite_volume/ins/cht/flow-around-square-linear-fluidonly.i)
- (modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d/2d-velocity-pressure.i)
- (modules/navier_stokes/test/tests/finite_volume/wcns/enthalpy_equation/enthalpy_equation.i)
- (modules/navier_stokes/test/tests/postprocessors/flow_rates/conservation_LinearFV.i)
- (modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d/2d-boussinesq-transient.i)
- (modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/segregated/channel-drift-flux.i)
- (modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/segregated-comparison/segregated-linear.i)
- (modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d/2d-boussinesq.i)
- (modules/navier_stokes/test/tests/finite_volume/ins/natural_convection/linear_segregated/2d/diff_heated_cavity_linear_segregated.i)
- (modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/fluid.i)
- (modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/steady-transient-compare/common-blocks.i)
- (modules/navier_stokes/test/tests/finite_volume/wcns/enthalpy_equation/1d_test_h_fp.i)
- (modules/navier_stokes/test/tests/finite_volume/wcns/enthalpy_equation/1d_test_h.i)
- (modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/SIMPLEC/3d/3d-velocity-pressure.i)
- (modules/navier_stokes/test/tests/finite_volume/ins/lid-driven/linear-segregated/lid-driven-segregated.i)
- (test/tests/linearfvkernels/anisotropic-diffusion/anisotropic-diffusion-2d.i)
References
- Xiaogang Liu, Pingjian Ming, Wenping Zhang, Lirong Fu, and Lilong Jing.
Finite-volume methods for anisotropic diffusion problems on skewed meshes.
Numerical Heat Transfer, Part B: Fundamentals, 68(3):239–256, 2015.[BibTeX]
@article{liu2015finite, author = "Liu, Xiaogang and Ming, Pingjian and Zhang, Wenping and Fu, Lirong and Jing, Lilong", title = "Finite-volume methods for anisotropic diffusion problems on skewed meshes", journal = "Numerical Heat Transfer, Part B: Fundamentals", volume = "68", number = "3", pages = "239--256", year = "2015", publisher = "Taylor \\& Francis" }
(modules/navier_stokes/examples/flow-over-circle/flow_over_circle_linearfv.i)
[Problem]
linear_sys_names = 'u_system v_system pressure_system'
previous_nl_solution_required = true
[]
[Functions]
[inlet_function]
type = ParsedFunction
expression = '4*U*(y-ymin)*(ymax-y)/(ymax-ymin)/(ymax-ymin)'
symbol_names = 'U ymax ymin'
symbol_values = '${inlet_velocity} ${y_max} ${y_min}'
[]
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
solver_sys = u_system
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
[]
[pressure]
type = MooseLinearVariableFVReal
initial_condition = 0
solver_sys = pressure_system
[]
[]
[LinearFVKernels]
[u_time]
type = LinearFVTimeDerivative
variable = vel_x
factor = ${rho}
[]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = true
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_time]
type = LinearFVTimeDerivative
variable = vel_y
factor = ${rho}
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = true
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = true
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[]
[LinearFVBCs]
[inlet_x]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = vel_x
boundary = 'left_boundary'
functor = 'inlet_function'
[]
[inlet_y]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = vel_y
boundary = 'left_boundary'
functor = 0
[]
[circle_x]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = vel_x
boundary = 'circle'
functor = 0
[]
[circle_y]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = vel_y
boundary = 'circle'
functor = 0
[]
[walls_x]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = vel_x
boundary = 'top_boundary bottom_boundary'
functor = 0
[]
[walls_y]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = vel_y
boundary = 'top_boundary bottom_boundary'
functor = 0
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right_boundary'
variable = pressure
functor = 0
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = 'right_boundary'
[]
[outlet_v]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_y
use_two_term_expansion = false
boundary = 'right_boundary'
[]
[]
[Postprocessors]
[drag_force]
type = IntegralDirectedSurfaceForce
vel_x = vel_x
vel_y = vel_y
mu = ${mu}
pressure = pressure
principal_direction = '1 0 0'
boundary = 'circle'
outputs = none
execute_on = 'INITIAL TIMESTEP_END'
[]
[drag_coeff]
type = ParsedPostprocessor
expression = '2*drag_force/rho/(avgvel*avgvel)/D'
constant_names = 'rho avgvel D'
constant_expressions = '${rho} ${fparse 2/3*inlet_velocity} ${fparse 2*circle_radius}'
pp_names = 'drag_force'
execute_on = 'INITIAL TIMESTEP_END'
[]
[lift_force]
type = IntegralDirectedSurfaceForce
vel_x = vel_x
vel_y = vel_y
mu = ${mu}
pressure = pressure
principal_direction = '0 1 0'
boundary = 'circle'
outputs = none
execute_on = 'INITIAL TIMESTEP_END'
[]
[lift_coeff]
type = ParsedPostprocessor
expression = '2*lift_force/rho/(avgvel*avgvel)/D'
constant_names = 'rho avgvel D'
constant_expressions = '${rho} ${fparse 2/3*inlet_velocity} ${fparse 2*circle_radius}'
pp_names = 'lift_force'
execute_on = 'INITIAL TIMESTEP_END'
[]
[]
[Executioner]
type = PIMPLE
momentum_l_abs_tol = 1e-12
pressure_l_abs_tol = 1e-12
momentum_l_tol = 1e-12
pressure_l_tol = 1e-12
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
momentum_equation_relaxation = 0.90
pressure_variable_relaxation = 0.4
num_iterations = 100
pressure_absolute_tolerance = 1e-10
momentum_absolute_tolerance = 1e-10
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
print_fields = false
continue_on_max_its = true
dt = 0.01
num_steps = 500
[]
[Outputs]
exodus = true
csv = true
[]
(modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/SIMPLEC/2d/2d-velocity-pressure.i)
mu = 2.6
rho = 1.0
advected_interp_method = 'average'
[Mesh]
[mesh]
type = CartesianMeshGenerator
dim = 2
dx = '1.0'
dy = '0.1'
ix = '20'
iy = '3'
[]
[]
[Problem]
linear_sys_names = 'u_system v_system pressure_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
pressure_projection_method = consistent
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
initial_condition = 0.5
solver_sys = u_system
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
initial_condition = 0.0
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = 0.2
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[]
[LinearFVBCs]
[inlet-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = '1.1'
[]
[inlet-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_y
functor = '0.0'
[]
[walls-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_x
functor = 0.0
[]
[walls-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_y
functor = 0.0
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = 1.4
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = right
[]
[outlet_v]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_y
use_two_term_expansion = false
boundary = right
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-10
pressure_l_abs_tol = 1e-10
momentum_l_tol = 0
pressure_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
momentum_equation_relaxation = 0.8
pressure_variable_relaxation = 1.0
num_iterations = 100
pressure_absolute_tolerance = 1e-10
momentum_absolute_tolerance = 1e-10
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
print_fields = false
[]
[Outputs]
exodus = true
[]
(modules/navier_stokes/test/tests/finite_volume/ins/mms/linear-segregated/2d-vortex/2d-vortex.i)
mu = 1
rho = 1
advected_interp_method = 'average'
[Problem]
linear_sys_names = 'u_system v_system pressure_system'
previous_nl_solution_required = true
[]
[Mesh]
[gmg]
type = GeneratedMeshGenerator
dim = 2
nx = 2
ny = 2
[]
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
initial_condition = 0.0
solver_sys = u_system
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
initial_condition = 0.0
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = 0
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[u_forcing]
type = LinearFVSource
variable = vel_x
source_density = forcing_u
[]
[v_forcing]
type = LinearFVSource
variable = vel_y
source_density = forcing_v
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
use_nonorthogonal_correction_on_boundary = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[]
[LinearFVBCs]
[no-slip-wall-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left right top bottom'
variable = vel_x
functor = '0'
[]
[no-slip-wall-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left right top bottom'
variable = vel_y
functor = '0'
[]
[pressure-extrapolation]
type = LinearFVExtrapolatedPressureBC
boundary = 'left right top bottom'
variable = pressure
use_two_term_expansion = true
[]
[]
[Functions]
[exact_u]
type = ParsedFunction
expression = 'x^2*(1-x)^2*(2*y-6*y^2+4*y^3)'
[]
[exact_v]
type = ParsedFunction
expression = '-y^2*(1-y)^2*(2*x-6*x^2+4*x^3)'
[]
[exact_p]
type = ParsedFunction
expression = 'x*(1-x)'
[]
[forcing_u]
type = ParsedFunction
expression = '-4*mu*(-1+2*y)*(y^2-6*x*y^2+6*x^2*y^2-y+6*x*y-6*x^2*y+3*x^2-6*x^3+3*x^4)+1-2*x+rho*4*x^3'
'*y^2*(2*y^2-2*y+1)*(y-1)^2*(-1+2*x)*(x-1)^3'
symbol_names = 'mu rho'
symbol_values = '${mu} ${rho}'
[]
[forcing_v]
type = ParsedFunction
expression = '4*mu*(-1+2*x)*(x^2-6*y*x^2+6*x^2*y^2-x+6*x*y-6*x*y^2+3*y^2-6*y^3+3*y^4)+rho*4*y^3*x^2*(2'
'*x^2-2*x+1)*(x-1)^2*(-1+2*y)*(y-1)^3'
symbol_names = 'mu rho'
symbol_values = '${mu} ${rho}'
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-8
pressure_l_abs_tol = 1e-8
momentum_l_tol = 0
pressure_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
momentum_equation_relaxation = 0.8
pressure_variable_relaxation = 0.3
num_iterations = 2000
pressure_absolute_tolerance = 1e-8
momentum_absolute_tolerance = 1e-8
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
print_fields = false
pin_pressure = true
pressure_pin_value = 0.25
pressure_pin_point = '0.5 0.5 0.0'
[]
[Outputs]
exodus = true
[csv]
type = CSV
execute_on = FINAL
[]
[]
[Postprocessors]
[h]
type = AverageElementSize
outputs = 'csv'
execute_on = FINAL
[]
[L2u]
type = ElementL2FunctorError
approximate = vel_x
exact = exact_u
outputs = 'csv'
execute_on = FINAL
[]
[L2v]
type = ElementL2FunctorError
approximate = vel_y
exact = exact_v
outputs = 'csv'
execute_on = FINAL
[]
[L2p]
approximate = pressure
exact = exact_p
type = ElementL2FunctorError
outputs = 'csv'
execute_on = FINAL
[]
[]
(modules/navier_stokes/test/tests/finite_volume/ins/cht/flow-around-square-linear.i)
mu = 0.01
rho = 1.1
k = 0.0005
cp = 10
k_s = 3.0
h_conv = 5
power_density = 10000
advected_interp_method = 'upwind'
[Mesh]
[generated_mesh]
type = GeneratedMeshGenerator
dim = 2
nx = 10
ny = 10
xmin = 0
ymin = 0
ymax = 0.1
xmax = 0.1
[]
[subdomain1]
type = SubdomainBoundingBoxGenerator
input = generated_mesh
block_name = subdomain1
bottom_left = '0.04 0.04 0'
block_id = 1
top_right = '0.06 0.06 0'
[]
[interface]
type = SideSetsBetweenSubdomainsGenerator
input = subdomain1
primary_block = 0
paired_block = 1
new_boundary = interface
[]
[]
[Problem]
linear_sys_names = 'u_system v_system pressure_system energy_system solid_energy_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
block = 0
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
initial_condition = 0.1
solver_sys = u_system
block = 0
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
initial_condition = 0.01
block = 0
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = 0.2
block = 0
[]
[T_fluid]
type = MooseLinearVariableFVReal
solver_sys = energy_system
initial_condition = 300
block = 0
[]
[T_solid]
type = MooseLinearVariableFVReal
solver_sys = solid_energy_system
initial_condition = 500
block = 1
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = true
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = true
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = true
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[h_advection]
type = LinearFVEnergyAdvection
variable = T_fluid
advected_quantity = temperature
cp = ${cp}
advected_interp_method = ${advected_interp_method}
rhie_chow_user_object = 'rc'
[]
[conduction]
type = LinearFVDiffusion
variable = T_fluid
diffusion_coeff = ${k}
use_nonorthogonal_correction = true
[]
[solid-conduction]
type = LinearFVDiffusion
variable = T_solid
diffusion_coeff = ${k_s}
use_nonorthogonal_correction = true
[]
[solid-source]
type = LinearFVSource
variable = T_solid
source_density = ${power_density}
[]
[]
[LinearFVBCs]
[inlet-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = '0.1'
[]
[inlet-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_y
functor = '0.0'
[]
[walls-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom interface'
variable = vel_x
functor = 0.0
[]
[walls-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom interface'
variable = vel_y
functor = 0.0
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = 1.4
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = right
[]
[outlet_v]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_y
use_two_term_expansion = false
boundary = right
[]
[inlet_T]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = T_fluid
functor = 300
boundary = 'left'
[]
[walls_T]
type = LinearFVAdvectionDiffusionFunctorNeumannBC
variable = T_fluid
functor = 0.0
boundary = 'top bottom'
[]
[outlet_T]
type = LinearFVAdvectionDiffusionOutflowBC
variable = T_fluid
use_two_term_expansion = false
boundary = right
[]
[fluid_solid]
type = LinearFVConvectiveHeatTransferBC
variable = T_fluid
T_solid = T_solid
T_fluid = T_fluid
boundary = interface
h = ${h_conv}
[]
[solid_fluid]
type = LinearFVConvectiveHeatTransferBC
variable = T_solid
T_solid = T_solid
T_fluid = T_fluid
boundary = interface
h = ${h_conv}
[]
[]
[FunctorMaterials]
[rhocpT]
property_name = 'rhocpT'
type = ParsedFunctorMaterial
functor_names = 'T_fluid'
expression = '${rho}*${cp}*T_fluid'
[]
[]
[Postprocessors]
[h_in]
type = VolumetricFlowRate
boundary = left
vel_x = vel_x
vel_y = vel_y
rhie_chow_user_object = rc
advected_quantity = 'rhocpT'
subtract_mesh_velocity = false
[]
[h_out]
type = VolumetricFlowRate
boundary = right
vel_x = vel_x
vel_y = vel_y
rhie_chow_user_object = rc
advected_quantity = 'rhocpT'
advected_interp_method = upwind
subtract_mesh_velocity = false
[]
[power]
type = ElementIntegralFunctorPostprocessor
functor = ${power_density}
block = 1
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-13
pressure_l_abs_tol = 1e-13
energy_l_abs_tol = 1e-13
solid_energy_l_abs_tol = 1e-13
momentum_l_tol = 0
pressure_l_tol = 0
energy_l_tol = 0
solid_energy_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
energy_system = 'energy_system'
solid_energy_system = 'solid_energy_system'
momentum_equation_relaxation = 0.8
energy_equation_relaxation = 1.0
pressure_variable_relaxation = 0.3
num_iterations = 1000
pressure_absolute_tolerance = 1e-10
momentum_absolute_tolerance = 1e-10
energy_absolute_tolerance = 1e-10
solid_energy_absolute_tolerance = 1e-10
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
energy_petsc_options_iname = '-pc_type -pc_hypre_type'
energy_petsc_options_value = 'hypre boomeramg'
solid_energy_petsc_options_iname = '-pc_type -pc_hypre_type'
solid_energy_petsc_options_value = 'hypre boomeramg'
print_fields = false
continue_on_max_its = true
[]
[Outputs]
exodus = true
execute_on = timestep_end
[]
(modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/3d/3d-velocity-pressure.i)
mu = 2.6
rho = 1.0
advected_interp_method = 'average'
[Mesh]
[mesh]
type = CartesianMeshGenerator
dim = 3
dx = '0.3'
dy = '0.3'
dz = '0.3'
ix = '3'
iy = '3'
iz = '3'
[]
[]
[Problem]
linear_sys_names = 'u_system v_system w_system pressure_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
w = vel_z
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
initial_condition = 0.5
solver_sys = u_system
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
initial_condition = 0.0
[]
[vel_z]
type = MooseLinearVariableFVReal
solver_sys = w_system
initial_condition = 0.0
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = 0.2
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
w = vel_z
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
w = vel_z
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[w_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_z
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
w = vel_z
momentum_component = 'z'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[w_pressure]
type = LinearFVMomentumPressure
variable = vel_z
pressure = pressure
momentum_component = 'z'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[]
[LinearFVBCs]
[inlet-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = '1.1'
[]
[inlet-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_y
functor = '0.0'
[]
[inlet-w]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_z
functor = '0.0'
[]
[walls-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom back front'
variable = vel_x
functor = 0.0
[]
[walls-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom back front'
variable = vel_y
functor = 0.0
[]
[walls-w]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom back front'
variable = vel_z
functor = 0.0
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = 1.4
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = right
[]
[outlet_v]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_y
use_two_term_expansion = false
boundary = right
[]
[outlet_w]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_z
use_two_term_expansion = false
boundary = right
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-10
pressure_l_abs_tol = 1e-10
momentum_l_tol = 0
pressure_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system w_system'
pressure_system = 'pressure_system'
momentum_equation_relaxation = 0.8
pressure_variable_relaxation = 0.3
num_iterations = 100
pressure_absolute_tolerance = 1e-10
momentum_absolute_tolerance = 1e-10
momentum_petsc_options_iname = '-pc_type -pc_hypre_type -pc_hypre_boomeramg_agg_nl -pc_hypre_boomeramg_agg_num_paths -pc_hypre_boomeramg_truncfactor -pc_hypre_boomeramg_strong_threshold -pc_hypre_boomeramg_coarsen_type -pc_hypre_boomeramg_interp_type'
momentum_petsc_options_value = 'hypre boomeramg 4 1 0.1 0.6 HMIS ext+i'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type -pc_hypre_boomeramg_agg_nl -pc_hypre_boomeramg_agg_num_paths -pc_hypre_boomeramg_truncfactor -pc_hypre_boomeramg_strong_threshold -pc_hypre_boomeramg_coarsen_type -pc_hypre_boomeramg_interp_type'
pressure_petsc_options_value = 'hypre boomeramg 2 1 0.1 0.6 HMIS ext+i'
print_fields = false
[]
[Outputs]
exodus = true
[]
(modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-scalar/channel.i)
mu = 2.6
rho = 1.0
advected_interp_method = 'upwind'
k1 = 0.1
k2 = 0.2
[Mesh]
[mesh]
type = CartesianMeshGenerator
dim = 2
dx = '0.25 0.25'
dy = '0.2'
ix = '5 5'
iy = '5'
subdomain_id = '0 1'
[]
[]
[Problem]
linear_sys_names = 'u_system v_system pressure_system s1_system s2_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
initial_condition = 0.5
solver_sys = u_system
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
initial_condition = 0.0
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = 0.2
[]
[scalar1]
type = MooseLinearVariableFVReal
solver_sys = s1_system
initial_condition = 1.1
[]
[scalar2]
type = MooseLinearVariableFVReal
solver_sys = s2_system
initial_condition = 3
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[s1_advection]
type = LinearFVScalarAdvection
variable = scalar1
advected_interp_method = ${advected_interp_method}
rhie_chow_user_object = 'rc'
[]
[s1_diffusion]
type = LinearFVDiffusion
variable = scalar1
diffusion_coeff = ${k1}
use_nonorthogonal_correction = false
[]
[s2_diffusion]
type = LinearFVScalarAdvection
variable = scalar2
advected_interp_method = ${advected_interp_method}
rhie_chow_user_object = 'rc'
[]
[s2_conduction]
type = LinearFVDiffusion
variable = scalar2
diffusion_coeff = ${k2}
use_nonorthogonal_correction = false
[]
[]
[LinearFVBCs]
[inlet-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = '1.1'
[]
[inlet-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_y
functor = '0.0'
[]
[walls-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_x
functor = 0.0
[]
[walls-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_y
functor = 0.0
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = 1.4
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = right
[]
[outlet_v]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_y
use_two_term_expansion = false
boundary = right
[]
[inlet_wall_scalar1]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = scalar1
functor = 1
boundary = 'left'
[]
[outlet_scalar1]
type = LinearFVAdvectionDiffusionOutflowBC
variable = scalar1
use_two_term_expansion = false
boundary = right
[]
[inlet_wall_scalar2]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = scalar2
functor = 2
boundary = 'left'
[]
[outlet_scalar2]
type = LinearFVAdvectionDiffusionOutflowBC
variable = scalar2
use_two_term_expansion = false
boundary = right
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-13
pressure_l_abs_tol = 1e-13
passive_scalar_l_abs_tol = 1e-13
momentum_l_tol = 0
pressure_l_tol = 0
passive_scalar_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
passive_scalar_systems = 's1_system s2_system'
momentum_equation_relaxation = 0.8
passive_scalar_equation_relaxation = '0.9 0.9'
pressure_variable_relaxation = 0.3
# We need to converge the problem to show conservation
num_iterations = 200
pressure_absolute_tolerance = 1e-10
momentum_absolute_tolerance = 1e-10
# The solution being flat, the normalization factor based on fluxes goes to
# 0. The convergence criteria being multiplied by said factor, we won't do any
# better than this. For a non-flat solution, use tighter tolerances
passive_scalar_absolute_tolerance = '1e-2 1e-2'
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
passive_scalar_petsc_options_iname = '-pc_type -pc_hypre_type'
passive_scalar_petsc_options_value = 'hypre boomeramg'
print_fields = false
continue_on_max_its = true
[]
[Outputs]
exodus = true
execute_on = timestep_end
hide = 'pressure vel_x vel_y'
[]
(modules/navier_stokes/test/tests/finite_volume/ins/cht/flow-around-square-linear-fluidonly.i)
mu = 0.01
rho = 1.1
k = 0.0005
cp = 10
h_conv = 5
advected_interp_method = 'upwind'
[Mesh]
[generated_mesh]
type = GeneratedMeshGenerator
dim = 2
nx = 10
ny = 10
xmin = 0
ymin = 0
ymax = 0.1
xmax = 0.1
[]
[subdomain1]
type = SubdomainBoundingBoxGenerator
input = generated_mesh
block_name = subdomain1
bottom_left = '0.04 0.04 0'
block_id = 1
top_right = '0.06 0.06 0'
[]
[interface]
type = SideSetsBetweenSubdomainsGenerator
input = subdomain1
primary_block = 0
paired_block = 1
new_boundary = interface
[]
[delete]
type = BlockDeletionGenerator
input = interface
block = 1
[]
[]
[Problem]
linear_sys_names = 'u_system v_system pressure_system energy_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
block = 0
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
initial_condition = 0.1
solver_sys = u_system
block = 0
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
initial_condition = 0.01
block = 0
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = 0.2
block = 0
[]
[T_fluid]
type = MooseLinearVariableFVReal
solver_sys = energy_system
initial_condition = 300
block = 0
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = true
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = true
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = true
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[h_advection]
type = LinearFVEnergyAdvection
variable = T_fluid
advected_quantity = temperature
cp = ${cp}
advected_interp_method = ${advected_interp_method}
rhie_chow_user_object = 'rc'
[]
[conduction]
type = LinearFVDiffusion
variable = T_fluid
diffusion_coeff = ${k}
use_nonorthogonal_correction = true
[]
[]
[LinearFVBCs]
[inlet-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = '0.1'
[]
[inlet-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_y
functor = '0.0'
[]
[walls-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom interface'
variable = vel_x
functor = 0.0
[]
[walls-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom interface'
variable = vel_y
functor = 0.0
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = 1.4
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = right
[]
[outlet_v]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_y
use_two_term_expansion = false
boundary = right
[]
[inlet_T]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = T_fluid
functor = 300
boundary = 'left'
[]
[walls_T]
type = LinearFVAdvectionDiffusionFunctorNeumannBC
variable = T_fluid
functor = 0.0
boundary = 'top bottom'
[]
[outlet_T]
type = LinearFVAdvectionDiffusionOutflowBC
variable = T_fluid
use_two_term_expansion = false
boundary = right
[]
[fluid_solid]
type = LinearFVConvectiveHeatTransferBC
variable = T_fluid
T_solid = boundary_value
T_fluid = T_fluid
boundary = interface
h = ${h_conv}
[]
[]
[FunctorMaterials]
[rhocpT]
property_name = 'rhocpT'
type = ParsedFunctorMaterial
functor_names = 'T_fluid'
expression = '${rho}*${cp}*T_fluid'
[]
[]
[Functions]
[boundary_value]
type = ConstantFunction
value = 350
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-13
pressure_l_abs_tol = 1e-13
energy_l_abs_tol = 1e-13
momentum_l_tol = 0
pressure_l_tol = 0
energy_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
energy_system = 'energy_system'
momentum_equation_relaxation = 0.8
energy_equation_relaxation = 1.0
pressure_variable_relaxation = 0.3
num_iterations = 1000
pressure_absolute_tolerance = 1e-10
momentum_absolute_tolerance = 1e-10
energy_absolute_tolerance = 1e-10
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
energy_petsc_options_iname = '-pc_type -pc_hypre_type'
energy_petsc_options_value = 'hypre boomeramg'
print_fields = false
continue_on_max_its = true
[]
[Outputs]
exodus = true
execute_on = timestep_end
[]
(modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d/2d-velocity-pressure.i)
mu = 2.6
rho = 1.0
advected_interp_method = 'average'
[Mesh]
[mesh]
type = CartesianMeshGenerator
dim = 2
dx = '0.3'
dy = '0.3'
ix = '3'
iy = '3'
[]
[]
[Problem]
linear_sys_names = 'u_system v_system pressure_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
initial_condition = 0.5
solver_sys = u_system
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
initial_condition = 0.0
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = 0.2
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[]
[LinearFVBCs]
[inlet-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = '1.1'
[]
[inlet-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_y
functor = '0.0'
[]
[walls-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_x
functor = 0.0
[]
[walls-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_y
functor = 0.0
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = 1.4
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = right
[]
[outlet_v]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_y
use_two_term_expansion = false
boundary = right
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-10
pressure_l_abs_tol = 1e-10
momentum_l_tol = 0
pressure_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
momentum_equation_relaxation = 0.8
pressure_variable_relaxation = 0.3
num_iterations = 100
pressure_absolute_tolerance = 1e-10
momentum_absolute_tolerance = 1e-10
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
print_fields = false
[]
[Outputs]
exodus = true
[]
(modules/navier_stokes/test/tests/finite_volume/wcns/enthalpy_equation/enthalpy_equation.i)
H = 0.015 #halfwidth of the channel, 10 cm of channel height
L = 1
bulk_u = 0.01
p_ref = 101325.0
advected_interp_method = 'upwind'
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = ${L}
ymin = -${H}
ymax = ${H}
nx = 30
ny = 15
[]
[]
[Problem]
linear_sys_names = 'u_system v_system pressure_system energy_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
pressure = pressure
rho = 'rho'
p_diffusion_kernel = p_diffusion
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
solver_sys = u_system
initial_condition = ${bulk_u}
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
initial_condition = 0
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = ${p_ref}
[]
[h]
type = MooseLinearVariableFVReal
solver_sys = energy_system
initial_condition = 44000 # 1900 is an approx of cp(T)
[]
[]
[AuxVariables]
[rho_var]
type = MooseLinearVariableFVReal
[]
[cp_var]
type = MooseLinearVariableFVReal
[]
[mu_var]
type = MooseLinearVariableFVReal
[]
[k_var]
type = MooseLinearVariableFVReal
[]
[T]
type = MooseLinearVariableFVReal
initial_condition = 777.
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
mu = 'mu'
momentum_component = 'x'
use_nonorthogonal_correction = false
advected_interp_method = ${advected_interp_method}
rhie_chow_user_object = 'rc'
u = vel_x
v = vel_y
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
mu = 'mu'
momentum_component = 'y'
use_nonorthogonal_correction = false
advected_interp_method = ${advected_interp_method}
rhie_chow_user_object = 'rc'
u = vel_x
v = vel_y
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[temp_conduction]
type = LinearFVDiffusion
diffusion_coeff = 'alpha'
variable = h
[]
[temp_advection]
type = LinearFVEnergyAdvection
variable = h
advected_interp_method = ${advected_interp_method}
rhie_chow_user_object = 'rc'
[]
[]
[LinearFVBCs]
[inlet_u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = ${bulk_u} #${bulk_u} #'fully_developed_velocity'
[]
[inlet-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_y
functor = 0
[]
[inlet_h]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = h
boundary = 'left'
functor = h_from_p_T # ${fparse 1900.*860.}
[]
[inlet_T]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = T
boundary = 'left'
functor = 860.
[]
[walls-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = vel_x
boundary = 'top bottom'
functor = 0.
[]
[walls-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = vel_y
boundary = 'top bottom'
functor = 0.
[]
[walls_h]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = h
boundary = 'top bottom'
functor = h_from_p_T # ${fparse 1900. * 950}
[]
[walls_T]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = T
boundary = 'top bottom'
functor = 950.
[]
[walls_p]
type = LinearFVExtrapolatedPressureBC
boundary = 'top bottom'
variable = pressure
use_two_term_expansion = false
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = ${p_ref}
[]
[outlet_h]
type = LinearFVAdvectionDiffusionOutflowBC
variable = h
use_two_term_expansion = false
boundary = 'right'
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = right
[]
[outlet_v]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_y
use_two_term_expansion = false
boundary = right
[]
[]
[FluidProperties]
[lead]
type = LeadFluidProperties
[]
[]
[FunctorMaterials]
[fluid_props_to_mat_props]
type = GeneralFunctorFluidProps
fp = lead
pressure = ${p_ref}
T_fluid = 'T'
speed = 1
porosity = 1
characteristic_length = 1
[]
[alpha]
type = ADParsedFunctorMaterial
property_name = 'alpha'
functor_names = 'k cp'
expression = 'k/cp'
[]
[enthalpy_material]
type = LinearFVEnthalpyFunctorMaterial
pressure = ${p_ref}
T_fluid = T
h = h
fp = lead
[]
[]
[AuxKernels]
[rho_out]
type = FunctorAux
functor = 'rho'
variable = 'rho_var'
execute_on = 'NONLINEAR'
[]
[cp_out]
type = FunctorAux
functor = 'cp'
variable = 'cp_var'
execute_on = 'NONLINEAR'
[]
[mu_out]
type = FunctorAux
functor = 'mu'
variable = 'mu_var'
execute_on = 'NONLINEAR'
[]
[k_out]
type = FunctorAux
functor = 'k'
variable = 'k_var'
execute_on = 'NONLINEAR'
[]
[T_from_h_functor]
type = FunctorAux
functor = 'T_from_p_h'
variable = 'T'
execute_on = 'NONLINEAR'
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-6
pressure_l_abs_tol = 1e-6
energy_l_abs_tol = 1e-8
momentum_l_tol = 0
pressure_l_tol = 0
energy_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
energy_system = 'energy_system'
momentum_equation_relaxation = 0.7
pressure_variable_relaxation = 0.3
energy_equation_relaxation = 0.9
num_iterations = 200
pressure_absolute_tolerance = 1e-6
momentum_absolute_tolerance = 1e-6
energy_absolute_tolerance = 1e-6
print_fields = false
momentum_l_max_its = 1000
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
energy_petsc_options_iname = '-pc_type -pc_hypre_type'
energy_petsc_options_value = 'hypre boomeramg'
continue_on_max_its = true
[]
[Outputs]
exodus = true
execute_on = 'TIMESTEP_BEGIN FINAL'
[]
(modules/navier_stokes/test/tests/postprocessors/flow_rates/conservation_LinearFV.i)
mu = 2.6
rho = 1.1
advected_interp_method='average'
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
[]
[]
[Mesh]
[mesh]
type = CartesianMeshGenerator
dim = 2
dx = '0.3'
dy = '0.3'
ix = '3'
iy = '3'
[]
[]
[Problem]
linear_sys_names = 'u_system v_system pressure_system'
previous_nl_solution_required = true
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
initial_condition = 0.5
solver_sys = u_system
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
initial_condition = 0.0
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = 0.2
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[]
[LinearFVBCs]
[inlet-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = '1.1'
[]
[inlet-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_y
functor = '0.0'
[]
[walls-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_x
functor = 0.0
[]
[walls-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_y
functor = 0.0
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = 1.4
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = right
[]
[outlet_v]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_y
use_two_term_expansion = false
boundary = right
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-10
pressure_l_abs_tol = 1e-10
momentum_l_tol = 0
pressure_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
momentum_equation_relaxation = 0.8
pressure_variable_relaxation = 0.3
num_iterations = 100
pressure_absolute_tolerance = 1e-10
momentum_absolute_tolerance = 1e-10
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
print_fields = false
[]
[Postprocessors]
[inlet_mass]
type = VolumetricFlowRate
boundary = left
rhie_chow_user_object = rc
advected_quantity = ${rho}
vel_x = vel_x
vel_y = vel_y
subtract_mesh_velocity = false
[]
[outlet_mass]
type = VolumetricFlowRate
boundary = right
rhie_chow_user_object = rc
advected_quantity = ${rho}
vel_x = vel_x
vel_y = vel_y
subtract_mesh_velocity = false
[]
[]
[Outputs]
csv = true
[]
(modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d/2d-boussinesq-transient.i)
mu = 2.6
rho = 1.0
advected_interp_method = 'average'
cp = 300
k = 10
alpha_b = 1e-4
[Mesh]
[mesh]
type = CartesianMeshGenerator
dim = 2
dx = '1.'
dy = '0.2'
ix = '10'
iy = '5'
[]
[]
[Problem]
linear_sys_names = 'u_system v_system pressure_system energy_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
initial_condition = 0.5
solver_sys = u_system
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
initial_condition = 0.0
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = 0.2
[]
[T]
type = MooseLinearVariableFVReal
solver_sys = energy_system
initial_condition = 300
[]
[]
[LinearFVKernels]
[u_time]
type = LinearFVTimeDerivative
variable = vel_x
factor = ${rho}
[]
[v_time]
type = LinearFVTimeDerivative
variable = vel_y
factor = ${rho}
[]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[u_boussinesq]
type = LinearFVMomentumBoussinesq
variable = vel_x
rho = ${rho}
gravity = '0 -9.81 0'
alpha_name = ${alpha_b}
ref_temperature = 300.0
T_fluid = T
momentum_component = 'x'
[]
[v_boussinesq]
type = LinearFVMomentumBoussinesq
variable = vel_y
rho = ${rho}
gravity = '0 -9.81 0'
alpha_name = ${alpha_b}
ref_temperature = 300.0
T_fluid = T
momentum_component = 'y'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[h_time]
type = LinearFVTimeDerivative
variable = T
factor = ${fparse rho*cp}
[]
[h_advection]
type = LinearFVEnergyAdvection
variable = T
advected_quantity = temperature
cp = ${cp}
advected_interp_method = ${advected_interp_method}
rhie_chow_user_object = 'rc'
[]
[conduction]
type = LinearFVDiffusion
variable = T
diffusion_coeff = ${k}
use_nonorthogonal_correction = false
[]
[]
[FunctorMaterials]
[constant_functors]
type = GenericFunctorMaterial
prop_names = 'cp alpha_b'
prop_values = '${cp} ${alpha_b}'
[]
[]
[LinearFVBCs]
[inlet-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = '1.1'
[]
[inlet-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_y
functor = '0.0'
[]
[walls-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_x
functor = 0.0
[]
[walls-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_y
functor = 0.0
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = 1.4
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = right
[]
[outlet_v]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_y
use_two_term_expansion = false
boundary = right
[]
[inlet_top_T]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = T
functor = 300.0
boundary = 'left top'
[]
[bottom_T]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = T
functor = wall-temperature
boundary = bottom
[]
[outlet_T]
type = LinearFVAdvectionDiffusionOutflowBC
variable = T
use_two_term_expansion = false
boundary = right
[]
[]
[Functions]
[wall-temperature]
type = ParsedFunction
expression = '350 + 50 * sin(6.28*t)'
[]
[]
[Executioner]
type = PIMPLE
momentum_l_abs_tol = 1e-12
pressure_l_abs_tol = 1e-12
energy_l_abs_tol = 1e-12
momentum_l_tol = 1e-12
pressure_l_tol = 1e-12
energy_l_tol = 1e-12
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
energy_system = 'energy_system'
momentum_equation_relaxation = 0.8
pressure_variable_relaxation = 0.3
energy_equation_relaxation = 0.9
num_iterations = 100
pressure_absolute_tolerance = 1e-11
momentum_absolute_tolerance = 1e-11
energy_absolute_tolerance = 1e-11
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
energy_petsc_options_iname = '-pc_type -pc_hypre_type'
energy_petsc_options_value = 'hypre boomeramg'
print_fields = false
continue_on_max_its = true
dt = 0.01
num_steps = 6
num_piso_iterations = 0
[]
[Outputs]
exodus = true
[]
(modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/segregated/channel-drift-flux.i)
mu = 1.0
rho = 10.0
mu_d = 0.1
rho_d = 1.0
l = 2
# 'average' leads to slight oscillations, upwind may be preferred
# This method is selected for consistency with the original nonlinear input
advected_interp_method = 'average'
# TODO remove need for those
cp = 1
k = 1
cp_d = 1
k_d = 1
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = '${fparse l * 5}'
ymin = '${fparse -l / 2}'
ymax = '${fparse l / 2}'
nx = 10
ny = 4
[]
uniform_refine = 0
[]
[Problem]
linear_sys_names = 'u_system v_system pressure_system phi_system'
previous_nl_solution_required = true
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
solver_sys = u_system
initial_condition = 1
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
[]
[phase_2]
type = MooseLinearVariableFVReal
solver_sys = phi_system
[]
[]
[LinearFVKernels]
[flow_p_diffusion]
type = LinearFVAnisotropicDiffusion
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
variable = pressure
[]
[flow_HbyA_divergence]
type = LinearFVDivergence
face_flux = HbyA
force_boundary_execution = true
variable = pressure
[]
[flow_ins_momentum_flux_x]
type = LinearWCNSFVMomentumFlux
advected_interp_method = ${advected_interp_method}
momentum_component = x
mu = mu_mixture
rhie_chow_user_object = ins_rhie_chow_interpolator
u = vel_x
use_deviatoric_terms = false
use_nonorthogonal_correction = false
v = vel_y
variable = vel_x
[]
[flow_ins_momentum_flux_y]
type = LinearWCNSFVMomentumFlux
advected_interp_method = ${advected_interp_method}
momentum_component = y
mu = mu_mixture
rhie_chow_user_object = ins_rhie_chow_interpolator
u = vel_x
use_deviatoric_terms = false
use_nonorthogonal_correction = false
v = vel_y
variable = vel_y
[]
[mixture_drift_flux_x]
type = LinearWCNSFV2PMomentumDriftFlux
density_interp_method = average
fraction_dispersed = phase_2
momentum_component = x
rhie_chow_user_object = ins_rhie_chow_interpolator
rho_d = ${rho_d}
u_slip = vel_slip_x
v_slip = vel_slip_y
variable = vel_x
[]
[mixture_drift_flux_y]
type = LinearWCNSFV2PMomentumDriftFlux
density_interp_method = average
fraction_dispersed = phase_2
momentum_component = y
rhie_chow_user_object = ins_rhie_chow_interpolator
rho_d = ${rho_d}
u_slip = vel_slip_x
v_slip = vel_slip_y
variable = vel_y
[]
[flow_ins_momentum_pressure_x]
type = LinearFVMomentumPressure
momentum_component = x
pressure = pressure
variable = vel_x
[]
[flow_ins_momentum_pressure_y]
type = LinearFVMomentumPressure
momentum_component = y
pressure = pressure
variable = vel_y
[]
[flow_momentum_friction_0_x]
type = LinearFVMomentumFriction
Darcy_name = Darcy_coefficient_vec
momentum_component = x
mu = mu_mixture
variable = vel_x
[]
[flow_momentum_friction_0_y]
type = LinearFVMomentumFriction
Darcy_name = Darcy_coefficient_vec
momentum_component = y
mu = mu_mixture
variable = vel_y
[]
# Mixture phase equation
[mixture_ins_phase_2_advection]
type = LinearFVScalarAdvection
advected_interp_method = upwind
rhie_chow_user_object = ins_rhie_chow_interpolator
u_slip = vel_slip_x
v_slip = vel_slip_y
variable = phase_2
[]
[mixture_phase_interface_reaction]
type = LinearFVReaction
coeff = 0.1
variable = phase_2
[]
[mixture_phase_interface_source]
type = LinearFVSource
scaling_factor = 0.1
source_density = phase_1
variable = phase_2
[]
[]
[LinearFVBCs]
[vel_x_left]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = left
functor = 1
variable = vel_x
[]
[vel_y_left]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = left
functor = 0
variable = vel_y
[]
[pressure_extrapolation_inlet_left]
type = LinearFVExtrapolatedPressureBC
boundary = left
use_two_term_expansion = true
variable = pressure
[]
[vel_x_right]
type = LinearFVAdvectionDiffusionOutflowBC
boundary = right
use_two_term_expansion = true
variable = vel_x
[]
[vel_y_right]
type = LinearFVAdvectionDiffusionOutflowBC
boundary = right
use_two_term_expansion = true
variable = vel_y
[]
[pressure_right]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = right
functor = 0
variable = pressure
[]
[vel_x_bottom]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = bottom
functor = 0
variable = vel_x
[]
[vel_y_bottom]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = bottom
functor = 0
variable = vel_y
[]
[vel_x_top]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = top
functor = 0
variable = vel_x
[]
[vel_y_top]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = top
functor = 0
variable = vel_y
[]
[pressure_extrapolation_top_bottom]
type = LinearFVExtrapolatedPressureBC
boundary = 'top bottom'
use_two_term_expansion = true
variable = pressure
[]
[phase_2_left]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = left
functor = 0.1
variable = phase_2
[]
[phase_2_right]
type = LinearFVAdvectionDiffusionOutflowBC
boundary = right
use_two_term_expansion = true
variable = phase_2
[]
[]
[FunctorMaterials]
[flow_ins_speed_material]
type = ADVectorMagnitudeFunctorMaterial
execute_on = ALWAYS
outputs = none
vector_magnitude_name = speed
x_functor = vel_x
y_functor = vel_y
[]
[mixture_phase_1_fraction]
type = ParsedFunctorMaterial
execute_on = ALWAYS
expression = '1 - phase_2'
functor_names = phase_2
output_properties = phase_1
outputs = all
property_name = phase_1
[]
[mixture_mixture_material]
type = WCNSLinearFVMixtureFunctorMaterial
execute_on = ALWAYS
limit_phase_fraction = true
outputs = all
phase_1_fraction = phase_2
phase_1_names = '${rho_d} ${mu_d} ${cp_d} ${k_d}'
phase_2_names = '${rho} ${mu} ${cp} ${k}'
prop_names = 'rho_mixture mu_mixture cp_mixture k_mixture'
[]
[mixture_slip_x]
type = WCNSFV2PSlipVelocityFunctorMaterial
execute_on = ALWAYS
gravity = '0 0 0'
linear_coef_name = Darcy_coefficient
momentum_component = x
mu = mu_mixture
outputs = all
particle_diameter = 0.01
rho = ${rho}
rho_d = ${rho_d}
slip_velocity_name = vel_slip_x
u = vel_x
v = vel_y
[]
[mixture_slip_y]
type = WCNSFV2PSlipVelocityFunctorMaterial
execute_on = ALWAYS
gravity = '0 0 0'
linear_coef_name = Darcy_coefficient
momentum_component = y
mu = mu_mixture
outputs = all
particle_diameter = 0.01
rho = ${rho}
rho_d = ${rho_d}
slip_velocity_name = vel_slip_y
u = vel_x
v = vel_y
[]
[mixture_dispersed_drag]
type = NSFVDispersePhaseDragFunctorMaterial
drag_coef_name = Darcy_coefficient
execute_on = ALWAYS
mu = mu_mixture
outputs = all
particle_diameter = 0.01
rho = rho_mixture
u = vel_x
v = vel_y
[]
[]
[UserObjects]
[ins_rhie_chow_interpolator]
type = RhieChowMassFlux
p_diffusion_kernel = flow_p_diffusion
pressure = pressure
rho = rho_mixture
u = vel_x
v = vel_y
[]
[]
[Executioner]
type = SIMPLE
rhie_chow_user_object = 'ins_rhie_chow_interpolator'
# Systems
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
active_scalar_systems = 'phi_system'
momentum_equation_relaxation = 0.8
active_scalar_equation_relaxation = '0.7'
pressure_variable_relaxation = 0.3
# We need to converge the problem to show conservation
num_iterations = 200
pressure_absolute_tolerance = 1e-10
momentum_absolute_tolerance = 1e-10
active_scalar_absolute_tolerance = '1e-10'
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
active_scalar_petsc_options_iname = '-pc_type -pc_hypre_type'
active_scalar_petsc_options_value = 'hypre boomeramg'
momentum_l_abs_tol = 1e-13
pressure_l_abs_tol = 1e-13
active_scalar_l_abs_tol = 1e-13
momentum_l_tol = 0
pressure_l_tol = 0
active_scalar_l_tol = 0
# print_fields = true
continue_on_max_its = true
[]
[Outputs]
csv = true
[]
[Postprocessors]
[Re]
type = ParsedPostprocessor
expression = '10.0 * 2 * 1'
[]
[average_phase2]
type = ElementAverageValue
variable = phase_2
[]
[dp]
type = PressureDrop
boundary = 'left right'
downstream_boundary = right
pressure = pressure
upstream_boundary = left
[]
[max_phase2]
type = ElementExtremeValue
variable = phase_2
[]
[]
(modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/segregated-comparison/segregated-linear.i)
mu = 2.6
rho = 1.0
advected_interp_method = 'average'
[Mesh]
[mesh]
type = CartesianMeshGenerator
dim = 2
dx = '0.3'
dy = '0.3'
ix = '3'
iy = '3'
[]
[]
[Problem]
linear_sys_names = 'u_system v_system pressure_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
initial_condition = 0.5
solver_sys = u_system
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
initial_condition = 0.0
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = 0.2
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[]
[LinearFVBCs]
[inlet-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = '1.1'
[]
[inlet-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_y
functor = '0.0'
[]
[walls-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_x
functor = 0.0
[]
[walls-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_y
functor = 0.0
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = 1.4
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = right
[]
[outlet_v]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_y
use_two_term_expansion = false
boundary = right
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-10
pressure_l_abs_tol = 1e-10
momentum_l_tol = 0
pressure_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
momentum_equation_relaxation = 0.8
pressure_variable_relaxation = 0.3
num_iterations = 2
pressure_absolute_tolerance = 1e-10
momentum_absolute_tolerance = 1e-10
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
print_fields = false
[]
[Outputs]
exodus = true
execute_on = TIMESTEP_END
[]
(modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d/2d-boussinesq.i)
mu = 2.6
rho = 1.0
advected_interp_method = 'average'
cp = 300
k = 20
alpha_b = 1e-4
[Mesh]
[mesh]
type = CartesianMeshGenerator
dim = 2
dx = '1.'
dy = '0.2'
ix = '10'
iy = '5'
[]
[]
[Problem]
linear_sys_names = 'u_system v_system pressure_system energy_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
initial_condition = 0.5
solver_sys = u_system
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
initial_condition = 0.0
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = 0.2
[]
[T]
type = MooseLinearVariableFVReal
solver_sys = energy_system
initial_condition = 300
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[u_boussinesq]
type = LinearFVMomentumBoussinesq
variable = vel_x
rho = ${rho}
gravity = '0 -9.81 0'
alpha_name = ${alpha_b}
ref_temperature = 300.0
T_fluid = T
momentum_component = 'x'
[]
[v_boussinesq]
type = LinearFVMomentumBoussinesq
variable = vel_y
rho = ${rho}
gravity = '0 -9.81 0'
alpha_name = ${alpha_b}
ref_temperature = 300.0
T_fluid = T
momentum_component = 'y'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[h_advection]
type = LinearFVEnergyAdvection
variable = T
advected_quantity = temperature
cp = ${cp}
advected_interp_method = ${advected_interp_method}
rhie_chow_user_object = 'rc'
[]
[conduction]
type = LinearFVDiffusion
variable = T
diffusion_coeff = ${k}
use_nonorthogonal_correction = false
[]
[]
[FunctorMaterials]
[constant_functors]
type = GenericFunctorMaterial
prop_names = 'cp alpha_b'
prop_values = '${cp} ${alpha_b}'
[]
[]
[LinearFVBCs]
[inlet-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = '1.1'
[]
[inlet-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_y
functor = '0.0'
[]
[walls-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_x
functor = 0.0
[]
[walls-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_y
functor = 0.0
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = 1.4
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = right
[]
[outlet_v]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_y
use_two_term_expansion = false
boundary = right
[]
[inlet_top_T]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = T
functor = ${fparse 300.0}
boundary = 'left top'
[]
[bottom_T]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = T
functor = ${fparse 350.0}
boundary = bottom
[]
[outlet_T]
type = LinearFVAdvectionDiffusionOutflowBC
variable = T
use_two_term_expansion = false
boundary = right
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-10
pressure_l_abs_tol = 1e-10
energy_l_abs_tol = 1e-10
momentum_l_tol = 0
pressure_l_tol = 0
energy_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
energy_system = 'energy_system'
momentum_equation_relaxation = 0.8
energy_equation_relaxation = 0.9
pressure_variable_relaxation = 0.3
num_iterations = 200
pressure_absolute_tolerance = 1e-10
momentum_absolute_tolerance = 1e-10
energy_absolute_tolerance = 1e-10
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
energy_petsc_options_iname = '-pc_type -pc_hypre_type'
energy_petsc_options_value = 'hypre boomeramg'
print_fields = false
[]
[Outputs]
exodus = true
[]
(modules/navier_stokes/test/tests/finite_volume/ins/natural_convection/linear_segregated/2d/diff_heated_cavity_linear_segregated.i)
################################################################################
# MATERIAL PROPERTIES
################################################################################
rho = 3279.
T_0 = 875.0
mu = 1.
k_cond = 38.0
cp = 640.
alpha = 3.26e-4
walls = 'right left top bottom'
[GlobalParams]
rhie_chow_user_object = 'ins_rhie_chow_interpolator'
advected_interp_method = 'upwind'
u = superficial_vel_x
v = superficial_vel_y
[]
[Problem]
linear_sys_names = 'u_system v_system pressure_system energy_system'
previous_nl_solution_required = true
[]
################################################################################
# GEOMETRY
################################################################################
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = 1
ymin = 0
ymax = 1
nx = 30
ny = 30
[]
[]
################################################################################
# EQUATIONS: VARIABLES, KERNELS & BCS
################################################################################
[UserObjects]
[ins_rhie_chow_interpolator]
type = RhieChowMassFlux
u = superficial_vel_x
v = superficial_vel_y
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
[]
[]
[Variables]
[superficial_vel_x]
type = MooseLinearVariableFVReal
solver_sys = u_system
[]
[superficial_vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
[]
[pressure]
type = MooseLinearVariableFVReal
initial_condition = 0
solver_sys = pressure_system
[]
[T_fluid]
type = MooseLinearVariableFVReal
solver_sys = energy_system
initial_condition = 875
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = superficial_vel_x
mu = ${mu}
momentum_component = 'x'
use_nonorthogonal_correction = true
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = superficial_vel_x
pressure = pressure
momentum_component = 'x'
[]
[u_buoyancy]
type = LinearFVMomentumBoussinesq
variable = superficial_vel_x
T_fluid = T_fluid
gravity = '0 -9.81 0'
rho = ${rho}
ref_temperature = ${T_0}
alpha_name = ${alpha}
momentum_component = 'x'
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = superficial_vel_y
mu = ${mu}
momentum_component = 'y'
use_nonorthogonal_correction = true
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = superficial_vel_y
pressure = pressure
momentum_component = 'y'
[]
[v_buoyancy]
type = LinearFVMomentumBoussinesq
variable = superficial_vel_y
T_fluid = T_fluid
gravity = '0 -9.81 0'
rho = ${rho}
ref_temperature = ${T_0}
alpha_name = ${alpha}
momentum_component = 'y'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = true
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
####### FUEL ENERGY EQUATION #######
[heat_advection]
type = LinearFVEnergyAdvection
variable = T_fluid
advected_quantity = temperature
cp = ${cp}
[]
[conduction]
type = LinearFVDiffusion
variable = T_fluid
diffusion_coeff = ${fparse k_cond}
[]
[]
[LinearFVBCs]
[no-slip-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = superficial_vel_x
boundary = ${walls}
functor = 0
[]
[no-slip-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = superficial_vel_y
boundary = ${walls}
functor = 0
[]
[T_cold]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = T_fluid
boundary = 'right'
functor = 870.0
[]
[T_hot]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = T_fluid
boundary = 'left'
functor = 880.0
[]
[T_all]
type = LinearFVAdvectionDiffusionExtrapolatedBC
variable = T_fluid
boundary = 'top bottom'
use_two_term_expansion = false
[]
[pressure-extrapolation]
type = LinearFVExtrapolatedPressureBC
boundary = ${walls}
variable = pressure
use_two_term_expansion = false
[]
[]
[FunctorMaterials]
[constant_functors]
type = GenericFunctorMaterial
prop_names = 'cp alpha_b'
prop_values = '${cp} ${alpha}'
[]
[]
################################################################################
# EXECUTION / SOLVE
################################################################################
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-11
pressure_l_abs_tol = 1e-11
energy_l_abs_tol = 1e-11
momentum_l_tol = 0
pressure_l_tol = 0
energy_l_tol = 0
rhie_chow_user_object = 'ins_rhie_chow_interpolator'
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
energy_system = 'energy_system'
momentum_equation_relaxation = 0.3
pressure_variable_relaxation = 0.7
energy_equation_relaxation = 0.95
num_iterations = 3000
pressure_absolute_tolerance = 1e-8
momentum_absolute_tolerance = 1e-8
energy_absolute_tolerance = 1e-8
print_fields = false
momentum_l_max_its = 300
pin_pressure = true
pressure_pin_value = 0.0
pressure_pin_point = '0.5 0.0 0.0'
# momentum_petsc_options = '-ksp_monitor'
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
energy_petsc_options_iname = '-pc_type -pc_hypre_type'
energy_petsc_options_value = 'hypre boomeramg'
[]
################################################################################
# SIMULATION OUTPUTS
################################################################################
[Outputs]
exodus = true
[]
(modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/2d-heated/fluid.i)
mu = 2.6
rho = 1.0
advected_interp_method = 'upwind'
cp = 1000
k = 1
[Mesh]
[mesh]
type = CartesianMeshGenerator
dim = 2
dx = '0.25 0.25'
dy = '0.2'
ix = '5 5'
iy = '5'
subdomain_id = '0 1'
[]
[]
[Problem]
linear_sys_names = 'u_system v_system pressure_system energy_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
initial_condition = 0.5
solver_sys = u_system
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
initial_condition = 0.0
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = 0.2
[]
[T_fluid]
type = MooseLinearVariableFVReal
solver_sys = energy_system
initial_condition = 300
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[h_advection]
type = LinearFVEnergyAdvection
variable = T_fluid
advected_quantity = temperature
cp = ${cp}
advected_interp_method = ${advected_interp_method}
rhie_chow_user_object = 'rc'
[]
[conduction]
type = LinearFVDiffusion
variable = T_fluid
diffusion_coeff = ${k}
use_nonorthogonal_correction = false
[]
[heat_exchange]
type = LinearFVVolumetricHeatTransfer
variable = T_fluid
h_solid_fluid = 100
T_fluid = T_fluid
T_solid = T_solid
is_solid = false
block = 1
[]
[]
[FunctorMaterials]
[constant_functors]
type = GenericFunctorMaterial
prop_names = 'cp'
prop_values = '${cp}'
[]
[]
[LinearFVBCs]
[inlet-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = '1.1'
[]
[inlet-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_y
functor = '0.0'
[]
[walls-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_x
functor = 0.0
[]
[walls-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_y
functor = 0.0
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = 1.4
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = right
[]
[outlet_v]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_y
use_two_term_expansion = false
boundary = right
[]
[inlet_wall_T]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = T_fluid
functor = 300
boundary = 'left top bottom'
[]
[outlet_T]
type = LinearFVAdvectionDiffusionOutflowBC
variable = T_fluid
use_two_term_expansion = false
boundary = right
[]
[]
[AuxVariables]
[T_solid]
type = MooseLinearVariableFVReal
initial_condition = 300
block = 1
[]
[]
[MultiApps]
inactive = 'solid'
[solid]
type = FullSolveMultiApp
input_files = solid.i
execute_on = timestep_begin
no_restore = true
[]
[]
[Transfers]
inactive = 'from_solid to_solid'
[from_solid]
type = MultiAppGeneralFieldShapeEvaluationTransfer
from_multi_app = solid
source_variable = 'T_solid'
variable = 'T_solid'
execute_on = timestep_begin
from_blocks = 1
[]
[to_solid]
type = MultiAppGeneralFieldShapeEvaluationTransfer
to_multi_app = solid
source_variable = 'T_fluid'
variable = 'T_fluid'
execute_on = timestep_begin
to_blocks = 1
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-13
pressure_l_abs_tol = 1e-13
energy_l_abs_tol = 1e-13
momentum_l_tol = 0
pressure_l_tol = 0
energy_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
energy_system = 'energy_system'
momentum_equation_relaxation = 0.8
energy_equation_relaxation = 0.9
pressure_variable_relaxation = 0.3
num_iterations = 20
pressure_absolute_tolerance = 1e-10
momentum_absolute_tolerance = 1e-10
energy_absolute_tolerance = 1e-10
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
energy_petsc_options_iname = '-pc_type -pc_hypre_type'
energy_petsc_options_value = 'hypre boomeramg'
print_fields = false
continue_on_max_its = true
[]
[Outputs]
exodus = true
execute_on = timestep_end
[]
(modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/linear-segregated/steady-transient-compare/common-blocks.i)
mu = 2.6
rho = 1.0
advected_interp_method = 'average'
cp = 300
k = 20
alpha_b = 1e-4
[Mesh]
[mesh]
type = CartesianMeshGenerator
dim = 2
dx = '1.'
dy = '0.2'
ix = '10'
iy = '5'
[]
[]
[Problem]
linear_sys_names = 'u_system v_system pressure_system energy_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
initial_condition = 0.5
solver_sys = u_system
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
initial_condition = 0.0
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = 0.2
[]
[T]
type = MooseLinearVariableFVReal
solver_sys = energy_system
initial_condition = 300
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[u_boussinesq]
type = LinearFVMomentumBoussinesq
variable = vel_x
rho = ${rho}
gravity = '0 -9.81 0'
alpha_name = ${alpha_b}
ref_temperature = 300.0
T_fluid = T
momentum_component = 'x'
[]
[v_boussinesq]
type = LinearFVMomentumBoussinesq
variable = vel_y
rho = ${rho}
gravity = '0 -9.81 0'
alpha_name = ${alpha_b}
ref_temperature = 300.0
T_fluid = T
momentum_component = 'y'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[h_advection]
type = LinearFVEnergyAdvection
variable = T
advected_quantity = temperature
cp = ${cp}
advected_interp_method = ${advected_interp_method}
rhie_chow_user_object = 'rc'
[]
[conduction]
type = LinearFVDiffusion
variable = T
diffusion_coeff = ${k}
use_nonorthogonal_correction = false
[]
[]
[FunctorMaterials]
[constant_functors]
type = GenericFunctorMaterial
prop_names = 'cp alpha_b'
prop_values = '${cp} ${alpha_b}'
[]
[]
[LinearFVBCs]
[inlet-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = '1.1'
[]
[inlet-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_y
functor = '0.0'
[]
[walls-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_x
functor = 0.0
[]
[walls-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom'
variable = vel_y
functor = 0.0
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = 1.4
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = right
[]
[outlet_v]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_y
use_two_term_expansion = false
boundary = right
[]
[inlet_top_T]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = T
functor = ${fparse 300.0}
boundary = 'left top'
[]
[bottom_T]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = T
functor = ${fparse 350.0}
boundary = bottom
[]
[outlet_T]
type = LinearFVAdvectionDiffusionOutflowBC
variable = T
use_two_term_expansion = false
boundary = right
[]
[]
[Outputs]
exodus = true
execute_on = FINAL
[]
(modules/navier_stokes/test/tests/finite_volume/wcns/enthalpy_equation/1d_test_h_fp.i)
L = 30
nx = 600
bulk_u = 0.01
p_ref = 101325.0
T_in = 860.
q_source = 20000.
advected_interp_method = 'upwind'
[Mesh]
[gmg]
type = GeneratedMeshGenerator
dim = 1
xmin = 0
xmax = ${L}
nx = ${nx}
[]
allow_renumbering = false
[]
[GlobalParams]
rhie_chow_user_object = 'rc'
advected_interp_method = ${advected_interp_method}
u = vel_x
[]
[Problem]
linear_sys_names = 'u_system pressure_system energy_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
pressure = pressure
rho = 'rho'
p_diffusion_kernel = p_diffusion
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
solver_sys = u_system
initial_condition = ${bulk_u}
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = ${p_ref}
[]
[h]
type = MooseLinearVariableFVReal
solver_sys = energy_system
initial_condition = ${fparse 860.*240.}
[]
[]
[AuxVariables]
[rho_var]
type = MooseLinearVariableFVReal
[]
[cp_var]
type = MooseLinearVariableFVReal
[]
[mu_var]
type = MooseLinearVariableFVReal
[]
[k_var]
type = MooseLinearVariableFVReal
[]
[alpha_var]
type = MooseLinearVariableFVReal
[]
[T]
type = MooseLinearVariableFVReal
initial_condition = ${T_in}
[]
[h_aux]
type = MooseLinearVariableFVReal
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
mu = 'mu'
momentum_component = 'x'
use_nonorthogonal_correction = false
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[temp_advection]
type = LinearFVEnergyAdvection
variable = h
[]
[source]
type = LinearFVSource
variable = h
source_density = source_func
[]
[]
[LinearFVBCs]
[inlet_u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = ${bulk_u}
[]
[inlet_h]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = h
boundary = 'left'
functor = 'h_from_p_T'
[]
[inlet_T]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = T
boundary = 'left'
functor = ${T_in}
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = ${p_ref}
[]
[outlet_h]
type = LinearFVAdvectionDiffusionOutflowBC
variable = h
use_two_term_expansion = false
boundary = 'right'
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = 'right'
[]
[]
[FluidProperties]
[lead]
type = LeadFluidProperties
[]
[]
[FunctorMaterials]
[enthalpy_material]
type = LinearFVEnthalpyFunctorMaterial
pressure = ${p_ref}
T_fluid = T
h = h
fp = lead
[]
[fluid_props_to_mat_props]
type = GeneralFunctorFluidProps
fp = lead
pressure = ${p_ref}
T_fluid = 'T'
speed = 1
porosity = 1
characteristic_length = 1
[]
[source_func]
type = ADParsedFunctorMaterial
property_name = source_func
functor_names = 'rho'
expression = ${q_source}
[]
[]
[AuxKernels]
[rho_out]
type = FunctorAux
functor = 'rho'
variable = 'rho_var'
execute_on = 'NONLINEAR'
[]
[cp_out]
type = FunctorAux
functor = 'cp'
variable = 'cp_var'
execute_on = 'NONLINEAR'
[]
[mu_out]
type = FunctorAux
functor = 'mu'
variable = 'mu_var'
execute_on = 'NONLINEAR'
[]
[k_out]
type = FunctorAux
functor = 'k'
variable = 'k_var'
execute_on = 'NONLINEAR'
[]
[T_from_h_functor_aux]
type = FunctorAux
functor = 'T_from_p_h'
variable = 'T'
execute_on = 'NONLINEAR'
[]
[h_from_T_functor_aux]
type = FunctorAux
functor = 'h_from_p_T'
variable = 'h_aux'
execute_on = 'NONLINEAR'
[]
[]
[Postprocessors]
[T_out_sim]
type = ElementalVariableValue
variable = T
elementid = ${fparse nx-1}
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-12
pressure_l_abs_tol = 1e-12
energy_l_abs_tol = 1e-12
momentum_l_tol = 0
pressure_l_tol = 0
energy_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system'
pressure_system = 'pressure_system'
energy_system = 'energy_system'
momentum_equation_relaxation = 0.7
pressure_variable_relaxation = 0.3
energy_equation_relaxation = 0.95
num_iterations = 100
pressure_absolute_tolerance = 1e-8
momentum_absolute_tolerance = 1e-8
energy_absolute_tolerance = 1e-6
print_fields = false
momentum_l_max_its = 200
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
energy_petsc_options_iname = '-pc_type -pc_hypre_type'
energy_petsc_options_value = 'hypre boomeramg'
continue_on_max_its = true
[]
[Outputs]
[out]
type = CSV
[]
[]
(modules/navier_stokes/test/tests/finite_volume/wcns/enthalpy_equation/1d_test_h.i)
L = 30
nx = 600
bulk_u = 0.01
q_source = 50000.
A_cp = 976.78
B_cp = 1.0634
T_in = 860.
p_ref = 101325.0
rho = 2000.
advected_interp_method = 'upwind'
[Mesh]
[gmg]
type = GeneratedMeshGenerator
dim = 1
xmin = 0
xmax = ${L}
nx = ${nx}
[]
allow_renumbering = false
[]
[GlobalParams]
rhie_chow_user_object = 'rc'
advected_interp_method = ${advected_interp_method}
u = vel_x
[]
[Problem]
linear_sys_names = 'u_system pressure_system energy_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
pressure = pressure
rho = 'rho'
p_diffusion_kernel = p_diffusion
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
solver_sys = u_system
initial_condition = ${bulk_u}
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = ${p_ref}
[]
[h]
type = MooseLinearVariableFVReal
solver_sys = energy_system
initial_condition = ${fparse 860.*1900.}
[]
[]
[AuxVariables]
[rho_var]
type = MooseLinearVariableFVReal
[]
[cp_var]
type = MooseLinearVariableFVReal
[]
[mu_var]
type = MooseLinearVariableFVReal
[]
[k_var]
type = MooseLinearVariableFVReal
[]
[alpha_var]
type = MooseLinearVariableFVReal
[]
[T]
type = MooseLinearVariableFVReal
initial_condition = 860.
[]
[h_aux]
type = MooseLinearVariableFVReal
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
mu = 'mu'
momentum_component = 'x'
use_nonorthogonal_correction = false
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[temp_advection]
type = LinearFVEnergyAdvection
variable = h
[]
[source]
type = LinearFVSource
variable = h
source_density = source_func
[]
[]
[LinearFVBCs]
[inlet_u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = ${bulk_u} #${bulk_u} #'fully_developed_velocity'
[]
[inlet_h]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = h
boundary = 'left'
functor = 'h_from_p_T'
[]
[inlet_T]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = T
boundary = 'left'
functor = ${T_in}
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = ${p_ref}
[]
[outlet_h]
type = LinearFVAdvectionDiffusionOutflowBC
variable = h
use_two_term_expansion = false
boundary = 'right'
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = 'right'
[]
[]
[Functions]
[source_func]
type = ParsedFunction
expression = ${q_source}
[]
[T_analytical]
type = ParsedFunction
expression = ${fparse (-A_cp+sqrt(A_cp^2-2*B_cp*(-q_source/rho/bulk_u*L-A_cp*T_in-B_cp/2*T_in*T_in)))/B_cp}
[]
[]
[FunctorMaterials]
[enthalpy_material]
type = LinearFVEnthalpyFunctorMaterial
pressure = ${p_ref}
T_fluid = T
h = h
h_from_p_T_functor = h_from_p_T_functor
T_from_p_h_functor = T_from_p_h_functor
[]
[h_from_p_T_functor]
type = ParsedFunctorMaterial
property_name = 'h_from_p_T_functor'
functor_names = 'T'
expression = '${A_cp}*T+${B_cp}/2*(T^2)'
[]
[T_from_p_h_functor]
type = ParsedFunctorMaterial
property_name = 'T_from_p_h_functor'
functor_names = 'h'
expression = '(-${A_cp}+sqrt(${A_cp}^2+2*h*${B_cp}))/${B_cp}'
[]
[rho]
type = ADParsedFunctorMaterial
property_name = 'rho'
functor_names = 'T'
expression = ${rho}
[]
[cp]
type = ADParsedFunctorMaterial
property_name = 'cp'
functor_names = 'T'
expression = '${A_cp}+${B_cp}*T'
[]
[mu]
type = ADParsedFunctorMaterial
property_name = 'mu'
functor_names = 'T'
expression = '4.5e-3'
[]
[k]
type = ADParsedFunctorMaterial
property_name = 'k'
functor_names = 'T'
expression = 0.7
[]
[]
[AuxKernels]
[rho_out]
type = FunctorAux
functor = 'rho'
variable = 'rho_var'
execute_on = 'NONLINEAR'
[]
[cp_out]
type = FunctorAux
functor = 'cp'
variable = 'cp_var'
execute_on = 'NONLINEAR'
[]
[mu_out]
type = FunctorAux
functor = 'mu'
variable = 'mu_var'
execute_on = 'NONLINEAR'
[]
[k_out]
type = FunctorAux
functor = 'k'
variable = 'k_var'
execute_on = 'NONLINEAR'
[]
[T_from_h_functor_aux]
type = FunctorAux
functor = 'T_from_p_h'
variable = 'T'
execute_on = 'NONLINEAR'
[]
[h_from_T_functor_aux]
type = FunctorAux
functor = 'h_from_p_T'
variable = 'h_aux'
execute_on = 'NONLINEAR'
[]
[]
[Postprocessors]
[T_out_sim]
type = ElementalVariableValue
variable = T
elementid = ${fparse nx-1}
[]
[T_out_analytic]
type = FunctionValuePostprocessor
function = T_analytical
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-12
pressure_l_abs_tol = 1e-12
energy_l_abs_tol = 1e-12
momentum_l_tol = 0
pressure_l_tol = 0
energy_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system'
pressure_system = 'pressure_system'
energy_system = 'energy_system'
momentum_equation_relaxation = 0.7
pressure_variable_relaxation = 0.3
energy_equation_relaxation = 0.95
num_iterations = 100
pressure_absolute_tolerance = 1e-8
momentum_absolute_tolerance = 1e-8
energy_absolute_tolerance = 1e-6
print_fields = false
momentum_l_max_its = 200
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
energy_petsc_options_iname = '-pc_type -pc_hypre_type'
energy_petsc_options_value = 'hypre boomeramg'
continue_on_max_its = true
[]
[Outputs]
[out]
type = CSV
[]
[]
(modules/navier_stokes/test/tests/finite_volume/ins/channel-flow/SIMPLEC/3d/3d-velocity-pressure.i)
mu = 2.6
rho = 1.0
advected_interp_method = 'average'
[Mesh]
[mesh]
type = CartesianMeshGenerator
dim = 3
dx = '1.0'
dy = '0.3'
dz = '0.3'
ix = '10'
iy = '3'
iz = '3'
[]
[]
[Problem]
linear_sys_names = 'u_system v_system w_system pressure_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
w = vel_z
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
pressure_projection_method = consistent
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
initial_condition = 0.5
solver_sys = u_system
[]
[vel_y]
type = MooseLinearVariableFVReal
solver_sys = v_system
initial_condition = 0.0
[]
[vel_z]
type = MooseLinearVariableFVReal
solver_sys = w_system
initial_condition = 0.0
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = 0.2
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
w = vel_z
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
w = vel_z
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[w_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_z
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
w = vel_z
momentum_component = 'z'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[w_pressure]
type = LinearFVMomentumPressure
variable = vel_z
pressure = pressure
momentum_component = 'z'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[]
[LinearFVBCs]
[inlet-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_x
functor = '1.1'
[]
[inlet-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_y
functor = '0.0'
[]
[inlet-w]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'left'
variable = vel_z
functor = '0.0'
[]
[walls-u]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom back front'
variable = vel_x
functor = 0.0
[]
[walls-v]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom back front'
variable = vel_y
functor = 0.0
[]
[walls-w]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'top bottom back front'
variable = vel_z
functor = 0.0
[]
[outlet_p]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
boundary = 'right'
variable = pressure
functor = 1.4
[]
[outlet_u]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_x
use_two_term_expansion = false
boundary = right
[]
[outlet_v]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_y
use_two_term_expansion = false
boundary = right
[]
[outlet_w]
type = LinearFVAdvectionDiffusionOutflowBC
variable = vel_z
use_two_term_expansion = false
boundary = right
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-10
pressure_l_abs_tol = 1e-10
momentum_l_tol = 0
pressure_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system w_system'
pressure_system = 'pressure_system'
momentum_equation_relaxation = 0.8
pressure_variable_relaxation = 1.0
num_iterations = 100
pressure_absolute_tolerance = 1e-10
momentum_absolute_tolerance = 1e-10
momentum_petsc_options_iname = '-pc_type -pc_hypre_type -pc_hypre_boomeramg_agg_nl -pc_hypre_boomeramg_agg_num_paths -pc_hypre_boomeramg_truncfactor -pc_hypre_boomeramg_strong_threshold -pc_hypre_boomeramg_coarsen_type -pc_hypre_boomeramg_interp_type'
momentum_petsc_options_value = 'hypre boomeramg 4 1 0.1 0.6 HMIS ext+i'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type -pc_hypre_boomeramg_agg_nl -pc_hypre_boomeramg_agg_num_paths -pc_hypre_boomeramg_truncfactor -pc_hypre_boomeramg_strong_threshold -pc_hypre_boomeramg_coarsen_type -pc_hypre_boomeramg_interp_type'
pressure_petsc_options_value = 'hypre boomeramg 2 1 0.1 0.6 HMIS ext+i'
print_fields = false
[]
[Outputs]
exodus = true
[]
(modules/navier_stokes/test/tests/finite_volume/ins/lid-driven/linear-segregated/lid-driven-segregated.i)
mu = .01
rho = 1
advected_interp_method = 'average'
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = .1
ymin = 0
ymax = .1
nx = 3
ny = 3
[]
[]
[Problem]
linear_sys_names = 'u_system v_system pressure_system'
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = RhieChowMassFlux
u = vel_x
v = vel_y
pressure = pressure
rho = ${rho}
p_diffusion_kernel = p_diffusion
[]
[]
[Variables]
[vel_x]
type = MooseLinearVariableFVReal
initial_condition = 0.0
solver_sys = u_system
[]
[vel_y]
type = MooseLinearVariableFVReal
initial_condition = 0.0
solver_sys = v_system
[]
[pressure]
type = MooseLinearVariableFVReal
solver_sys = pressure_system
initial_condition = 0.2
[]
[]
[LinearFVKernels]
[u_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_x
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'x'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[v_advection_stress]
type = LinearWCNSFVMomentumFlux
variable = vel_y
advected_interp_method = ${advected_interp_method}
mu = ${mu}
u = vel_x
v = vel_y
momentum_component = 'y'
rhie_chow_user_object = 'rc'
use_nonorthogonal_correction = false
[]
[u_pressure]
type = LinearFVMomentumPressure
variable = vel_x
pressure = pressure
momentum_component = 'x'
[]
[v_pressure]
type = LinearFVMomentumPressure
variable = vel_y
pressure = pressure
momentum_component = 'y'
[]
[p_diffusion]
type = LinearFVAnisotropicDiffusion
variable = pressure
diffusion_tensor = Ainv
use_nonorthogonal_correction = false
[]
[HbyA_divergence]
type = LinearFVDivergence
variable = pressure
face_flux = HbyA
force_boundary_execution = true
[]
[]
[LinearFVBCs]
[top_x]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = vel_x
boundary = 'top'
functor = 1
[]
[no_slip_x]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = vel_x
boundary = 'left right bottom'
functor = 0
[]
[no_slip_y]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = vel_y
boundary = 'left right top bottom'
functor = 0
[]
[pressure-extrapolation]
type = LinearFVExtrapolatedPressureBC
boundary = 'left right top bottom'
variable = pressure
use_two_term_expansion = true
[]
[]
[Executioner]
type = SIMPLE
momentum_l_abs_tol = 1e-10
pressure_l_abs_tol = 1e-10
momentum_l_tol = 0
pressure_l_tol = 0
rhie_chow_user_object = 'rc'
momentum_systems = 'u_system v_system'
pressure_system = 'pressure_system'
momentum_equation_relaxation = 0.8
pressure_variable_relaxation = 0.3
num_iterations = 500
pressure_absolute_tolerance = 1e-10
momentum_absolute_tolerance = 1e-10
momentum_petsc_options_iname = '-pc_type -pc_hypre_type'
momentum_petsc_options_value = 'hypre boomeramg'
pressure_petsc_options_iname = '-pc_type -pc_hypre_type'
pressure_petsc_options_value = 'hypre boomeramg'
print_fields = false
pin_pressure = true
pressure_pin_value = 0.0
pressure_pin_point = '0.05 0.05 0.0'
[]
[Outputs]
exodus = true
csv = false
perf_graph = false
print_nonlinear_residuals = false
print_linear_residuals = true
[]
(test/tests/linearfvkernels/anisotropic-diffusion/anisotropic-diffusion-2d.i)
[Mesh]
[gmg]
type = GeneratedMeshGenerator
dim = 2
nx = 2
ny = 1
ymax = 0.5
[]
[]
[Problem]
linear_sys_names = 'u_sys'
[]
[Variables]
[u]
type = MooseLinearVariableFVReal
solver_sys = 'u_sys'
initial_condition = 1.0
[]
[]
[LinearFVKernels]
[diffusion]
type = LinearFVAnisotropicDiffusion
variable = u
diffusion_tensor = diffusivity_tensor
use_nonorthogonal_correction = false
[]
[source]
type = LinearFVSource
variable = u
source_density = source_func
[]
[]
[LinearFVBCs]
[dir]
type = LinearFVAdvectionDiffusionFunctorDirichletBC
variable = u
boundary = "left right top bottom"
functor = analytic_solution
[]
[]
[FunctorMaterials]
[diff_tensor]
type = GenericVectorFunctorMaterial
prop_names = diffusivity_tensor
prop_values = 'coeff_func_x coeff_func_y 0.0'
[]
[]
[Functions]
[coeff_func_x]
type = ParsedFunction
expression = '1+0.5*x*y'
[]
[coeff_func_y]
type = ParsedFunction
expression = '1+x*y'
[]
[source_func]
type = ParsedFunction
expression = '(1.5-y*y)*(2+2*x*y)+(1.5-x*x)*(2+4*x*y)'
[]
[analytic_solution]
type = ParsedFunction
expression = '(1.5-x*x)*(1.5-y*y)'
[]
[]
[Postprocessors]
[h]
type = AverageElementSize
execute_on = FINAL
[]
[error]
type = ElementL2FunctorError
approximate = u
exact = analytic_solution
execute_on = FINAL
[]
[]
[Convergence]
[linear]
type = IterationCountConvergence
max_iterations = 1
converge_at_max_iterations = true
[]
[]
[Executioner]
type = Steady
system_names = u_sys
l_tol = 1e-10
multi_system_fixed_point=true
multi_system_fixed_point_convergence=linear
petsc_options_iname = '-pc_type -pc_hypre_type'
petsc_options_value = 'hypre boomeramg'
[]
[Outputs]
[csv]
type = CSV
execute_on = FINAL
[]
[]