- pin_typeHow to pin the pressure
C++ Type:MooseEnum
Controllable:No
Description:How to pin the pressure
NSPressurePin
Pins the pressure after a solve
Overview
The NSPressurePin
can pin the pressure in two modes:
by offsetting the pressure variable to make it have an average equal to the "phi0" parameter value
by offsetting the pressure variable to make its value equal to the "phi0" parameter value in the element containing the point specified by the "point" parameter.
In the NavierStokesFV Action, a NSPressurePin
can be used by setting the "pinned_pressure_type" parameter to average-uo
or point-value-uo
respectively.
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
- phi00Pressure pin value
Default:0
C++ Type:PostprocessorName
Unit:(no unit assumed)
Controllable:No
Description:Pressure pin value
- pointThe XYZ coordinates of a point inside an element where the pinned value shall be enforced.
C++ Type:libMesh::Point
Controllable:No
Description:The XYZ coordinates of a point inside an element where the pinned value shall be enforced.
- pressure_averageA postprocessor that computes the average of the pressure variable
C++ Type:PostprocessorName
Unit:(no unit assumed)
Controllable:No
Description:A postprocessor that computes the average of the pressure variable
- variablepressurePressure variable
Default:pressure
C++ Type:NonlinearVariableName
Unit:(no unit assumed)
Controllable:No
Description:Pressure variable
Optional Parameters
- allow_duplicate_execution_on_initialFalseIn the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).
Default:False
C++ Type:bool
Controllable:No
Description:In the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).
- execute_onTIMESTEP_ENDThe list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html.
Default:TIMESTEP_END
C++ Type:ExecFlagEnum
Options:XFEM_MARK, FORWARD, ADJOINT, HOMOGENEOUS_FORWARD, ADJOINT_TIMESTEP_BEGIN, ADJOINT_TIMESTEP_END, INITIAL, LINEAR_CONVERGENCE, NONLINEAR, NONLINEAR_CONVERGENCE, POSTCHECK, TIMESTEP_END, MULTIAPP_FIXED_POINT_END, MULTIAPP_FIXED_POINT_BEGIN, MULTIAPP_FIXED_POINT_CONVERGENCE, FINAL, CUSTOM
Controllable:No
Description:The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html.
- execution_order_group0Execution order groups are executed in increasing order (e.g., the lowest number is executed first). Note that negative group numbers may be used to execute groups before the default (0) group. Please refer to the user object documentation for ordering of user object execution within a group.
Default:0
C++ Type:int
Controllable:No
Description:Execution order groups are executed in increasing order (e.g., the lowest number is executed first). Note that negative group numbers may be used to execute groups before the default (0) group. Please refer to the user object documentation for ordering of user object execution within a group.
- force_postauxFalseForces the UserObject to be executed in POSTAUX
Default:False
C++ Type:bool
Controllable:No
Description:Forces the UserObject to be executed in POSTAUX
- force_preauxFalseForces the UserObject to be executed in PREAUX
Default:False
C++ Type:bool
Controllable:No
Description:Forces the UserObject to be executed in PREAUX
Execution Scheduling 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.
- 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
- prop_getter_suffixAn optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:An optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
- use_interpolated_stateFalseFor the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.
Default:False
C++ Type:bool
Controllable:No
Description:For the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.
Material Property Retrieval Parameters
Input Files
- (modules/navier_stokes/test/tests/finite_element/ins/hdg/ip/lid-driven/lid-driven-fsp.i)
- (modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp_stokes.i)
- (modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/rayleigh-bernard-two-phase.i)
- (modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp_al.i)
- (modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_fsp.i)
- (modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_fsp_diagonal_of_a_for_scaling.i)
- (modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp_elman.i)
- (modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp.i)
- (modules/navier_stokes/test/tests/finite_element/ins/lid_driven/transient_fsp.i)
phi0
Default:0
C++ Type:PostprocessorName
Unit:(no unit assumed)
Controllable:No
Description:Pressure pin value
phi0
Default:0
C++ Type:PostprocessorName
Unit:(no unit assumed)
Controllable:No
Description:Pressure pin value
point
C++ Type:libMesh::Point
Controllable:No
Description:The XYZ coordinates of a point inside an element where the pinned value shall be enforced.
pinned_pressure_type
Default:average-uo
C++ Type:MooseEnum
Options:average, point-value, average-uo, point-value-uo
Controllable:No
Description:Types for shifting (pinning) the pressure in case of incompressible simulations.
(modules/navier_stokes/test/tests/finite_element/ins/hdg/ip/lid-driven/lid-driven-fsp.i)
mu = 1
rho = 1
l = 1
U = 1e3
n = 200
gamma = 1e5
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = ${l}
ymin = 0
ymax = ${l}
nx = ${n}
ny = ${n}
elem_type = TRI6
[]
[]
[Problem]
type = NavierStokesProblem
extra_tag_matrices = 'mass'
mass_matrix = 'mass'
use_pressure_mass_matrix = true
[]
[Variables]
[vel_x]
family = MONOMIAL
order = FIRST
[]
[vel_y]
family = MONOMIAL
order = FIRST
[]
[pressure]
family = MONOMIAL
order = CONSTANT
[]
[vel_bar_x]
family = SIDE_HIERARCHIC
order = FIRST
[]
[vel_bar_y]
family = SIDE_HIERARCHIC
order = FIRST
[]
[pressure_bar]
family = SIDE_HIERARCHIC
order = FIRST
[]
[]
[HDGKernels]
[momentum_x_convection]
type = AdvectionIPHDGKernel
variable = vel_x
face_variable = vel_bar_x
velocity = 'velocity'
coeff = ${rho}
[]
[momentum_x_diffusion]
type = NavierStokesStressIPHDGKernel
variable = vel_x
face_variable = vel_bar_x
diffusivity = 'mu'
alpha = 6
pressure_variable = pressure
pressure_face_variable = pressure_bar
component = 0
[]
[momentum_y_convection]
type = AdvectionIPHDGKernel
variable = vel_y
face_variable = vel_bar_y
velocity = 'velocity'
coeff = ${rho}
[]
[momentum_y_diffusion]
type = NavierStokesStressIPHDGKernel
variable = vel_y
face_variable = vel_bar_y
diffusivity = 'mu'
alpha = 6
pressure_variable = pressure
pressure_face_variable = pressure_bar
component = 1
[]
[pressure_convection]
type = AdvectionIPHDGKernel
variable = pressure
face_variable = pressure_bar
velocity = 'velocity'
coeff = '${fparse -rho}'
self_advection = false
[]
[]
[Kernels]
[mass_matrix_pressure]
type = MassMatrix
variable = pressure
matrix_tags = 'mass'
density = '${fparse -1/gamma}'
[]
[grad_div_x]
type = GradDiv
variable = vel_x
u = vel_x
v = vel_y
gamma = ${gamma}
component = 0
[]
[grad_div_y]
type = GradDiv
variable = vel_y
u = vel_x
v = vel_y
gamma = ${gamma}
component = 1
[]
[]
[DGKernels]
[pb_mass]
type = MassMatrixDGKernel
variable = pressure_bar
matrix_tags = 'mass'
density = '${fparse -1/gamma}'
[]
[u_jump]
type = MassFluxPenalty
variable = vel_x
u = vel_x
v = vel_y
component = 0
gamma = ${gamma}
[]
[v_jump]
type = MassFluxPenalty
variable = vel_y
u = vel_x
v = vel_y
component = 1
gamma = ${gamma}
[]
[]
[BCs]
[momentum_x_diffusion_walls]
type = NavierStokesStressIPHDGDirichletBC
boundary = 'left bottom right'
variable = vel_x
face_variable = vel_bar_x
pressure_variable = pressure
pressure_face_variable = pressure_bar
alpha = 6
functor = '0'
diffusivity = 'mu'
component = 0
[]
[momentum_x_diffusion_top]
type = NavierStokesStressIPHDGDirichletBC
boundary = 'top'
variable = vel_x
face_variable = vel_bar_x
pressure_variable = pressure
pressure_face_variable = pressure_bar
alpha = 6
functor = '${U}'
diffusivity = 'mu'
component = 0
[]
[momentum_y_diffusion_all]
type = NavierStokesStressIPHDGDirichletBC
boundary = 'left bottom right top'
variable = vel_y
face_variable = vel_bar_y
pressure_variable = pressure
pressure_face_variable = pressure_bar
alpha = 6
functor = '0'
diffusivity = 'mu'
component = 1
[]
[mass_convection]
type = AdvectionIPHDGPrescribedFluxBC
face_variable = pressure_bar
variable = pressure
velocity = 'velocity'
coeff = '${fparse -rho}'
self_advection = false
boundary = 'left bottom top right'
prescribed_normal_flux = 0
[]
[pb_mass]
type = MassMatrixIntegratedBC
variable = pressure_bar
matrix_tags = 'mass'
boundary = 'left right bottom top'
density = '${fparse -1/gamma}'
[]
[u_jump_walls]
type = MassFluxPenaltyBC
variable = vel_x
u = vel_x
v = vel_y
component = 0
boundary = 'left right bottom'
gamma = ${gamma}
dirichlet_value = 'walls'
[]
[v_jump_walls]
type = MassFluxPenaltyBC
variable = vel_y
u = vel_x
v = vel_y
component = 1
boundary = 'left right bottom'
gamma = ${gamma}
dirichlet_value = 'walls'
[]
[u_jump_top]
type = MassFluxPenaltyBC
variable = vel_x
u = vel_x
v = vel_y
component = 0
boundary = 'top'
gamma = ${gamma}
dirichlet_value = 'top'
[]
[v_jump_top]
type = MassFluxPenaltyBC
variable = vel_y
u = vel_x
v = vel_y
component = 1
boundary = 'top'
gamma = ${gamma}
dirichlet_value = 'top'
[]
[]
[Functions]
[top]
type = ParsedVectorFunction
value_x = ${U}
value_y = 0
[]
[walls]
type = ParsedVectorFunction
value_x = 0
value_y = 0
[]
[]
[Materials]
[const]
type = ADGenericConstantMaterial
prop_names = 'rho mu'
prop_values = '${rho} ${mu}'
[]
[vel]
type = ADVectorFromComponentVariablesMaterial
vector_prop_name = 'velocity'
u = vel_x
v = vel_y
[]
[]
[Preconditioning]
[FSP]
type = FSP
topsplit = 'up'
[up]
splitting = 'u p'
splitting_type = schur
petsc_options_iname = '-pc_fieldsplit_schur_fact_type -pc_fieldsplit_schur_precondition -ksp_gmres_restart -ksp_type -ksp_pc_side -ksp_rtol'
petsc_options_value = 'full self 300 fgmres right 1e-4'
[]
[u]
vars = 'vel_x vel_y vel_bar_x vel_bar_y'
petsc_options = '-ksp_converged_reason'
petsc_options_iname = '-pc_type -ksp_type -ksp_rtol -ksp_gmres_restart -ksp_pc_side -pc_factor_mat_solver_type'
petsc_options_value = 'ilu gmres 1e-2 300 right strumpack'
[]
[p]
vars = 'pressure pressure_bar'
petsc_options = '-ksp_converged_reason'
petsc_options_iname = '-ksp_type -ksp_gmres_restart -ksp_rtol -pc_type -ksp_pc_side'
petsc_options_value = 'gmres 300 1e-2 ilu right'
[]
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
[]
[Outputs]
active = ''
[out]
type = Exodus
hide = 'pressure_average vel_bar_x vel_bar_y pressure_bar'
[]
[]
[Postprocessors]
[Re]
type = ParsedPostprocessor
pp_names = ''
expression = '${rho} * ${U} * ${l} / ${mu}'
[]
[pressure_average]
type = ElementAverageValue
variable = pressure
[]
[]
[Correctors]
[set_pressure]
type = NSPressurePin
pin_type = 'average'
variable = pressure
pressure_average = 'pressure_average'
[]
[]
(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp_stokes.i)
rho=1
mu=1
U=1
l=1
prefactor=${fparse 1/(l/2)^2}
n=8
[Mesh]
[gen]
type = DistributedRectilinearMeshGenerator
dim = 2
xmin = 0
xmax = ${l}
ymin = 0
ymax = ${l}
nx = ${n}
ny = ${n}
elem_type = QUAD4
[]
second_order = true
parallel_type = distributed
[]
[Variables]
[vel]
order = SECOND
family = LAGRANGE_VEC
[]
[p]
order = FIRST
family = LAGRANGE
[]
[]
[Kernels]
[mass]
type = INSADMass
variable = p
[]
[mass_kernel]
type = MassMatrix
variable = p
matrix_tags = 'mass'
[]
[momentum_viscous]
type = INSADMomentumViscous
variable = vel
[]
[momentum_pressure]
type = INSADMomentumPressure
variable = vel
pressure = p
integrate_p_by_parts = true
[]
[]
[BCs]
[no_slip]
type = ADVectorFunctionDirichletBC
variable = vel
boundary = 'bottom right left'
[]
[lid]
type = ADVectorFunctionDirichletBC
variable = vel
boundary = 'top'
function_x = 'lid_function'
[]
[]
[Materials]
[const]
type = ADGenericConstantMaterial
prop_names = 'rho mu'
prop_values = '${rho} ${mu}'
[]
[insad]
type = INSADMaterial
velocity = vel
pressure = p
[]
[]
[Functions]
[lid_function]
# We pick a function that is exactly represented in the velocity
# space so that the Dirichlet conditions are the same regardless
# of the mesh spacing.
type = ParsedFunction
expression = '${prefactor}*${U}*x*(${l}-x)'
[]
[]
[Problem]
type = NavierStokesProblem
mass_matrix = 'mass'
extra_tag_matrices = 'mass'
use_pressure_mass_matrix = true
[]
[Preconditioning]
[FSP]
type = FSP
topsplit = 'up'
[up]
splitting = 'u p'
splitting_type = schur
petsc_options_iname = '-pc_fieldsplit_schur_fact_type -pc_fieldsplit_schur_precondition -ksp_gmres_restart -ksp_type -ksp_pc_side -ksp_rtol'
petsc_options_value = 'full self 300 fgmres right 1e-4'
[]
[u]
vars = 'vel'
# petsc_options = '-ksp_converged_reason'
petsc_options_iname = '-pc_type -pc_hypre_type -ksp_type -ksp_rtol -ksp_gmres_restart -ksp_pc_side'
petsc_options_value = 'hypre boomeramg gmres 1e-2 300 right'
[]
[p]
vars = 'p'
petsc_options = '-ksp_converged_reason'
petsc_options_iname = '-ksp_type -ksp_gmres_restart -ksp_rtol -pc_type -ksp_pc_side -pc_hypre_type'
petsc_options_value = 'fgmres 300 1e-2 hypre right boomeramg'
[]
[]
[]
[Postprocessors]
[pavg]
type = ElementAverageValue
variable = p
[]
[]
[Correctors]
[set_pressure]
type = NSPressurePin
pin_type = 'average'
variable = p
pressure_average = 'pavg'
[]
[]
[Executioner]
type = Steady
solve_type = NEWTON
nl_rel_tol = 1e-12
[]
[Outputs]
print_linear_residuals = false
[exo]
type = Exodus
execute_on = 'final'
hide = 'pavg'
[]
[]
(modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/rayleigh-bernard-two-phase.i)
mu = 1.0
rho = 1e3
mu_d = 0.3
rho_d = 1.0
dp = 0.01
U_lid = 0.0
g = -9.81
[GlobalParams]
velocity_interp_method = 'rc'
advected_interp_method = 'upwind'
rhie_chow_user_object = 'rc'
[]
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = .1
ymin = 0
ymax = .1
nx = 11
ny = 11
[]
[]
[Variables]
[vel_x]
type = INSFVVelocityVariable
[]
[vel_y]
type = INSFVVelocityVariable
[]
[pressure]
type = INSFVPressureVariable
[]
[phase_2]
type = INSFVScalarFieldVariable
[]
[]
[UserObjects]
[rc]
type = INSFVRhieChowInterpolator
u = vel_x
v = vel_y
pressure = pressure
[]
[]
[Correctors]
[pin_pressure]
type = NSPressurePin
variable = pressure
pin_type = point-value
point = '0 0 0'
[]
[]
[FVKernels]
[mass]
type = INSFVMassAdvection
variable = pressure
rho = 'rho_mixture'
[]
[u_time]
type = INSFVMomentumTimeDerivative
variable = vel_x
rho = 'rho_mixture'
momentum_component = 'x'
[]
[u_advection]
type = INSFVMomentumAdvection
variable = vel_x
rho = 'rho_mixture'
momentum_component = 'x'
[]
[u_viscosity]
type = INSFVMomentumDiffusion
variable = vel_x
mu = 'mu_mixture'
momentum_component = 'x'
[]
[u_pressure]
type = INSFVMomentumPressure
variable = vel_x
momentum_component = 'x'
pressure = pressure
[]
[u_buoyant]
type = INSFVMomentumGravity
variable = vel_x
rho = 'rho_mixture'
momentum_component = 'x'
gravity = '0 ${g} 0'
[]
[v_time]
type = INSFVMomentumTimeDerivative
variable = vel_y
rho = 'rho_mixture'
momentum_component = 'y'
[]
[v_advection]
type = INSFVMomentumAdvection
variable = vel_y
rho = 'rho_mixture'
momentum_component = 'y'
[]
[v_viscosity]
type = INSFVMomentumDiffusion
variable = vel_y
mu = 'mu_mixture'
momentum_component = 'y'
[]
[v_pressure]
type = INSFVMomentumPressure
variable = vel_y
momentum_component = 'y'
pressure = pressure
[]
[v_buoyant]
type = INSFVMomentumGravity
variable = vel_y
rho = 'rho_mixture'
momentum_component = 'y'
gravity = '0 ${g} 0'
[]
[phase_2_time]
type = FVFunctorTimeKernel
variable = phase_2
[]
[phase_2_advection]
type = INSFVScalarFieldAdvection
variable = phase_2
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
[]
[phase_2_diffusion]
type = FVDiffusion
variable = phase_2
coeff = 1e-3
[]
[]
[FVBCs]
[top_x]
type = INSFVNoSlipWallBC
variable = vel_x
boundary = 'top'
function = ${U_lid}
[]
[no_slip_x]
type = INSFVNoSlipWallBC
variable = vel_x
boundary = 'left right bottom'
function = 0
[]
[no_slip_y]
type = INSFVNoSlipWallBC
variable = vel_y
boundary = 'left right top bottom'
function = 0
[]
[bottom_phase_2]
type = FVDirichletBC
variable = phase_2
boundary = 'bottom'
value = 1.0
[]
[top_phase_2]
type = FVDirichletBC
variable = phase_2
boundary = 'top'
value = 0.0
[]
[]
[AuxVariables]
[U]
order = CONSTANT
family = MONOMIAL
fv = true
[]
[drag_coefficient]
type = MooseVariableFVReal
[]
[rho_mixture_var]
type = MooseVariableFVReal
[]
[mu_mixture_var]
type = MooseVariableFVReal
[]
[phase_1]
type = MooseVariableFVReal
[]
[]
[AuxKernels]
[mag]
type = VectorMagnitudeAux
variable = U
x = vel_x
y = vel_y
[]
[populate_cd]
type = FunctorAux
variable = drag_coefficient
functor = 'Darcy_coefficient'
[]
[populate_rho_mixture_var]
type = FunctorAux
variable = rho_mixture_var
functor = 'rho_mixture'
[]
[populate_mu_mixture_var]
type = FunctorAux
variable = mu_mixture_var
functor = 'mu_mixture'
[]
[compute_phase_1]
type = ParsedAux
variable = phase_1
coupled_variables = 'phase_2'
expression = '1 - phase_2'
[]
[]
[FunctorMaterials]
[CD]
type = NSFVDispersePhaseDragFunctorMaterial
rho = 'rho_mixture'
mu = mu_mixture
u = 'vel_x'
v = 'vel_y'
particle_diameter = ${dp}
[]
[mixing_material]
type = NSFVMixtureFunctorMaterial
phase_1_names = '${rho_d} ${mu_d}'
phase_2_names = '${rho} ${mu}'
prop_names = 'rho_mixture mu_mixture'
phase_1_fraction = 'phase_2'
[]
[populate_u_slip]
type = WCNSFV2PSlipVelocityFunctorMaterial
slip_velocity_name = 'vel_slip_x'
momentum_component = 'x'
u = 'vel_x'
v = 'vel_y'
rho = ${rho}
mu = 'mu_mixture'
rho_d = ${rho_d}
particle_diameter = ${dp}
linear_coef_name = 'Darcy_coefficient'
[]
[populate_v_slip]
type = WCNSFV2PSlipVelocityFunctorMaterial
slip_velocity_name = 'vel_slip_y'
momentum_component = 'y'
u = 'vel_x'
v = 'vel_y'
rho = ${rho}
mu = 'mu_mixture'
rho_d = ${rho_d}
particle_diameter = ${dp}
linear_coef_name = 'Darcy_coefficient'
[]
[]
[Postprocessors]
[average_void]
type = ElementAverageValue
variable = 'phase_2'
[]
[max_y_velocity]
type = ElementExtremeValue
variable = 'vel_y'
value_type = max
[]
[min_y_velocity]
type = ElementExtremeValue
variable = 'vel_y'
value_type = min
[]
[max_x_velocity]
type = ElementExtremeValue
variable = 'vel_x'
value_type = max
[]
[min_x_velocity]
type = ElementExtremeValue
variable = 'vel_x'
value_type = min
[]
[max_x_slip_velocity]
type = ElementExtremeFunctorValue
functor = 'vel_slip_x'
value_type = max
[]
[max_y_slip_velocity]
type = ElementExtremeFunctorValue
functor = 'vel_slip_y'
value_type = max
[]
[max_drag_coefficient]
type = ElementExtremeFunctorValue
functor = 'drag_coefficient'
value_type = max
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
[TimeStepper]
type = IterationAdaptiveDT
optimal_iterations = 10
iteration_window = 2
growth_factor = 2
cutback_factor = 0.5
dt = 1e-3
[]
nl_max_its = 20
nl_rel_tol = 1e-03
nl_abs_tol = 1e-9
l_max_its = 5
end_time = 1e8
[]
[Outputs]
exodus = false
[CSV]
type = CSV
execute_on = 'FINAL'
[]
[]
(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp_al.i)
rho=1
mu=1e-3
U=1
l=1
prefactor=${fparse 1/(l/2)^2}
n=8
gamma=${U}
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = ${l}
ymin = 0
ymax = ${l}
nx = ${n}
ny = ${n}
elem_type = QUAD4
[]
second_order = true
[]
[Variables]
[vel]
order = SECOND
family = LAGRANGE_VEC
[]
[p]
order = FIRST
family = LAGRANGE
[]
[]
[Kernels]
[mass]
type = INSADMass
variable = p
[]
[mass_kernel]
type = MassMatrix
variable = p
matrix_tags = 'mass'
density = ${fparse -1/(gamma + mu)}
[]
[momentum_advection]
type = INSADMomentumAdvection
variable = vel
[]
[momentum_viscous]
type = INSADMomentumViscous
variable = vel
[]
[momentum_pressure]
type = INSADMomentumPressure
variable = vel
pressure = p
integrate_p_by_parts = true
[]
[momentum_graddiv]
type = INSADMomentumGradDiv
variable = vel
gamma = ${gamma}
[]
[]
[BCs]
[no_slip]
type = ADVectorFunctionDirichletBC
variable = vel
boundary = 'bottom right left'
preset = true
[]
[lid]
type = ADVectorFunctionDirichletBC
variable = vel
boundary = 'top'
function_x = 'lid_function'
preset = true
[]
[]
[Materials]
[const]
type = ADGenericConstantMaterial
prop_names = 'rho mu'
prop_values = '${rho} ${mu}'
[]
[insad]
type = INSADTauMaterial
velocity = vel
pressure = p
[]
[]
[Functions]
[lid_function]
# We pick a function that is exactly represented in the velocity
# space so that the Dirichlet conditions are the same regardless
# of the mesh spacing.
type = ParsedFunction
expression = '${prefactor}*${U}*x*(${l}-x)'
[]
[]
[Problem]
type = NavierStokesProblem
mass_matrix = 'mass'
extra_tag_matrices = 'mass'
use_pressure_mass_matrix = true
[]
[Preconditioning]
[FSP]
type = FSP
topsplit = 'up'
[up]
splitting = 'u p'
splitting_type = schur
petsc_options_iname = '-pc_fieldsplit_schur_fact_type -pc_fieldsplit_schur_precondition -ksp_gmres_restart -ksp_type -ksp_pc_side -ksp_rtol'
petsc_options_value = 'full self 300 fgmres right 1e-4'
[]
[u]
vars = 'vel'
petsc_options = '-ksp_converged_reason'
petsc_options_iname = '-pc_type -ksp_type -ksp_rtol -ksp_gmres_restart -ksp_pc_side -pc_factor_mat_solver_type'
petsc_options_value = 'ilu gmres 1e-2 300 right strumpack'
[]
[p]
vars = 'p'
petsc_options = '-ksp_converged_reason'
petsc_options_iname = '-ksp_type -ksp_gmres_restart -ksp_rtol -pc_type -ksp_pc_side'
petsc_options_value = 'gmres 300 1e-2 ilu right'
[]
[]
[]
[Postprocessors]
[pavg]
type = ElementAverageValue
variable = p
[]
[]
[Correctors]
[set_pressure]
type = NSPressurePin
pin_type = 'average'
variable = p
pressure_average = 'pavg'
[]
[]
[Executioner]
type = Steady
solve_type = NEWTON
nl_rel_tol = 1e-12
[]
[Outputs]
print_linear_residuals = false
[exo]
type = Exodus
execute_on = 'final'
hide = 'pavg'
[]
[]
(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_fsp.i)
rho=1
mu=2e-3
U=1
l=1
prefactor=${fparse 1/(l/2)^2}
n=64
[GlobalParams]
gravity = '0 0 0'
[]
[Mesh]
[gen]
type = DistributedRectilinearMeshGenerator
dim = 2
xmin = 0
xmax = ${l}
ymin = 0
ymax = ${l}
nx = ${n}
ny = ${n}
elem_type = QUAD4
[]
second_order = true
parallel_type = distributed
[]
[Variables]
[vel_x]
order = SECOND
family = LAGRANGE
[]
[vel_y]
order = SECOND
family = LAGRANGE
[]
[p]
order = FIRST
family = LAGRANGE
[]
[]
[Kernels]
[mass]
type = INSMass
variable = p
u = vel_x
v = vel_y
pressure = p
[]
[x_momentum_space]
type = INSMomentumLaplaceForm
variable = vel_x
u = vel_x
v = vel_y
pressure = p
component = 0
[]
[momentum_x_mass]
type = MassMatrix
variable = vel_x
density = ${rho}
matrix_tags = 'mass'
[]
[y_momentum_space]
type = INSMomentumLaplaceForm
variable = vel_y
u = vel_x
v = vel_y
pressure = p
component = 1
[]
[momentum_y_mass]
type = MassMatrix
variable = vel_y
density = ${rho}
matrix_tags = 'mass'
[]
[]
[BCs]
[x_no_slip]
type = DirichletBC
variable = vel_x
boundary = 'bottom right left'
value = 0.0
[]
[lid]
type = FunctionDirichletBC
variable = vel_x
boundary = 'top'
function = 'lid_function'
[]
[y_no_slip]
type = DirichletBC
variable = vel_y
boundary = 'bottom right top left'
value = 0.0
[]
[]
[Materials]
[const]
type = GenericConstantMaterial
block = 0
prop_names = 'rho mu'
prop_values = '${rho} ${mu}'
[]
[]
[Functions]
[lid_function]
# We pick a function that is exactly represented in the velocity
# space so that the Dirichlet conditions are the same regardless
# of the mesh spacing.
type = ParsedFunction
expression = '${prefactor}*${U}*x*(${l}-x)'
[]
[]
[Problem]
type = NavierStokesProblem
mass_matrix = 'mass'
extra_tag_matrices = 'mass'
[]
[Preconditioning]
[FSP]
type = FSP
topsplit = 'up'
[up]
splitting = 'u p'
splitting_type = schur
petsc_options_iname = '-pc_fieldsplit_schur_fact_type -pc_fieldsplit_schur_precondition -ksp_gmres_restart -ksp_type -ksp_pc_side -ksp_rtol'
petsc_options_value = 'full self 300 fgmres right 1e-4'
[]
[u]
vars = 'vel_x vel_y'
# petsc_options = '-ksp_converged_reason'
petsc_options_iname = '-pc_type -pc_hypre_type -ksp_type -ksp_rtol -ksp_gmres_restart -ksp_pc_side'
petsc_options_value = 'hypre boomeramg gmres 1e-2 300 right'
[]
[p]
vars = 'p'
petsc_options = '-pc_lsc_scale_diag -ksp_converged_reason'# -lsc_ksp_converged_reason -lsc_ksp_monitor_true_residual
petsc_options_iname = '-ksp_type -ksp_gmres_restart -ksp_rtol -pc_type -ksp_pc_side -lsc_pc_type -lsc_pc_hypre_type -lsc_ksp_type -lsc_ksp_rtol -lsc_ksp_pc_side -lsc_ksp_gmres_restart'
petsc_options_value = 'fgmres 300 1e-2 lsc right hypre boomeramg gmres 1e-1 right 300'
[]
[]
[]
[Postprocessors]
[pavg]
type = ElementAverageValue
variable = p
[]
[]
[Correctors]
[set_pressure]
type = NSPressurePin
pin_type = 'average'
variable = p
pressure_average = 'pavg'
[]
[]
[Executioner]
type = Steady
solve_type = NEWTON
[]
[Outputs]
[exo]
type = Exodus
execute_on = 'final'
hide = 'pavg'
[]
[]
(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_fsp_diagonal_of_a_for_scaling.i)
rho=1
mu=2e-3
U=1
l=1
prefactor=${fparse 1/(l/2)^2}
n=64
[GlobalParams]
gravity = '0 0 0'
[]
[Mesh]
[gen]
type = DistributedRectilinearMeshGenerator
dim = 2
xmin = 0
xmax = ${l}
ymin = 0
ymax = ${l}
nx = ${n}
ny = ${n}
elem_type = QUAD4
[]
second_order = true
parallel_type = distributed
[]
[Variables]
[vel_x]
order = SECOND
family = LAGRANGE
[]
[vel_y]
order = SECOND
family = LAGRANGE
[]
[p]
order = FIRST
family = LAGRANGE
[]
[]
[Kernels]
[mass]
type = INSMass
variable = p
u = vel_x
v = vel_y
pressure = p
[]
[x_momentum_space]
type = INSMomentumLaplaceForm
variable = vel_x
u = vel_x
v = vel_y
pressure = p
component = 0
[]
[y_momentum_space]
type = INSMomentumLaplaceForm
variable = vel_y
u = vel_x
v = vel_y
pressure = p
component = 1
[]
[]
[BCs]
[x_no_slip]
type = DirichletBC
variable = vel_x
boundary = 'bottom right left'
value = 0.0
[]
[lid]
type = FunctionDirichletBC
variable = vel_x
boundary = 'top'
function = 'lid_function'
[]
[y_no_slip]
type = DirichletBC
variable = vel_y
boundary = 'bottom right top left'
value = 0.0
[]
[]
[Materials]
[const]
type = GenericConstantMaterial
block = 0
prop_names = 'rho mu'
prop_values = '${rho} ${mu}'
[]
[]
[Functions]
[lid_function]
# We pick a function that is exactly represented in the velocity
# space so that the Dirichlet conditions are the same regardless
# of the mesh spacing.
type = ParsedFunction
expression = '${prefactor}*${U}*x*(${l}-x)'
[]
[]
[Problem]
type = NavierStokesProblem
[]
[Preconditioning]
[FSP]
type = FSP
topsplit = 'up'
[up]
splitting = 'u p'
splitting_type = schur
petsc_options_iname = '-pc_fieldsplit_schur_fact_type -pc_fieldsplit_schur_precondition -ksp_gmres_restart -ksp_type -ksp_pc_side -ksp_rtol'
petsc_options_value = 'full self 300 fgmres right 1e-4'
[]
[u]
vars = 'vel_x vel_y'
# petsc_options = '-ksp_converged_reason'
petsc_options_iname = '-pc_type -pc_hypre_type -ksp_type -ksp_rtol -ksp_gmres_restart -ksp_pc_side'
petsc_options_value = 'hypre boomeramg gmres 1e-2 300 right'
[]
[p]
vars = 'p'
petsc_options = '-pc_lsc_scale_diag -ksp_converged_reason'# -lsc_ksp_converged_reason -lsc_ksp_monitor_true_residual
petsc_options_iname = '-ksp_type -ksp_gmres_restart -ksp_rtol -pc_type -ksp_pc_side -lsc_pc_type -lsc_pc_hypre_type -lsc_ksp_type -lsc_ksp_rtol -lsc_ksp_pc_side -lsc_ksp_gmres_restart'
petsc_options_value = 'fgmres 300 1e-2 lsc right hypre boomeramg gmres 1e-1 right 300'
[]
[]
[]
[Postprocessors]
[pavg]
type = ElementAverageValue
variable = p
[]
[]
[Correctors]
[set_pressure]
type = NSPressurePin
pin_type = 'average'
variable = p
pressure_average = 'pavg'
[]
[]
[Executioner]
type = Steady
solve_type = NEWTON
[]
[Outputs]
[exo]
type = Exodus
execute_on = 'final'
hide = 'pavg'
[]
[]
(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp_elman.i)
rho=1
mu=1
U=1
l=1
prefactor=${fparse 1/(l/2)^2}
n=8
[Mesh]
[gen]
type = DistributedRectilinearMeshGenerator
dim = 2
xmin = 0
xmax = ${l}
ymin = 0
ymax = ${l}
nx = ${n}
ny = ${n}
elem_type = QUAD4
[]
second_order = true
parallel_type = distributed
[]
[Variables]
[vel]
order = SECOND
family = LAGRANGE_VEC
[]
[p]
order = FIRST
family = LAGRANGE
[]
[]
[Kernels]
[mass]
type = INSADMass
variable = p
[]
[velocity_mass_kernel]
type = VectorMassMatrix
variable = vel
matrix_tags = 'mass'
[]
[momentum_convection]
type = INSADMomentumAdvection
variable = vel
[]
[momentum_viscous]
type = INSADMomentumViscous
variable = vel
[]
[momentum_pressure]
type = INSADMomentumPressure
variable = vel
pressure = p
integrate_p_by_parts = true
[]
[]
[BCs]
[no_slip]
type = ADVectorFunctionDirichletBC
variable = vel
boundary = 'bottom right left'
[]
[lid]
type = ADVectorFunctionDirichletBC
variable = vel
boundary = 'top'
function_x = 'lid_function'
[]
[]
[Materials]
[const]
type = ADGenericConstantMaterial
prop_names = 'rho mu'
prop_values = '${rho} ${mu}'
[]
[insad]
type = INSADMaterial
velocity = vel
pressure = p
[]
[]
[Functions]
[lid_function]
# We pick a function that is exactly represented in the velocity
# space so that the Dirichlet conditions are the same regardless
# of the mesh spacing.
type = ParsedFunction
expression = '${prefactor}*${U}*x*(${l}-x)'
[]
[]
[Problem]
type = NavierStokesProblem
mass_matrix = 'mass'
extra_tag_matrices = 'mass'
[]
[Preconditioning]
[FSP]
type = FSP
topsplit = 'up'
[up]
splitting = 'u p'
splitting_type = schur
petsc_options_iname = '-pc_fieldsplit_schur_fact_type -pc_fieldsplit_schur_precondition -ksp_gmres_restart -ksp_type -ksp_pc_side -ksp_rtol'
petsc_options_value = 'full self 300 fgmres right 1e-4'
[]
[u]
vars = 'vel'
# petsc_options = '-ksp_converged_reason'
petsc_options_iname = '-pc_type -pc_hypre_type -ksp_type -ksp_rtol -ksp_gmres_restart -ksp_pc_side'
petsc_options_value = 'hypre boomeramg gmres 1e-2 300 right'
[]
[p]
vars = 'p'
petsc_options = '-ksp_converged_reason -pc_lsc_scale_diag'
petsc_options_iname = '-ksp_type -ksp_gmres_restart -ksp_rtol -pc_type -ksp_pc_side -lsc_pc_type -lsc_pc_hypre_type -lsc_ksp_type -lsc_ksp_rtol -lsc_ksp_pc_side -lsc_ksp_gmres_restart'
petsc_options_value = 'fgmres 300 1e-2 lsc right hypre boomeramg gmres 1e-1 right 300'
[]
[]
[]
[Postprocessors]
[pavg]
type = ElementAverageValue
variable = p
[]
[]
[Correctors]
[set_pressure]
type = NSPressurePin
pin_type = 'average'
variable = p
pressure_average = 'pavg'
[]
[]
[Executioner]
type = Steady
solve_type = NEWTON
[]
[Outputs]
print_linear_residuals = false
[exo]
type = Exodus
execute_on = 'final'
hide = 'pavg'
[]
[]
(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/steady_vector_fsp.i)
rho=1
mu=1
U=1
l=1
prefactor=${fparse 1/(l/2)^2}
n=8
[Mesh]
[gen]
type = DistributedRectilinearMeshGenerator
dim = 2
xmin = 0
xmax = ${l}
ymin = 0
ymax = ${l}
nx = ${n}
ny = ${n}
elem_type = QUAD4
[]
second_order = true
parallel_type = distributed
[]
[Variables]
[vel]
order = SECOND
family = LAGRANGE_VEC
[]
[p]
order = FIRST
family = LAGRANGE
[]
[]
[Kernels]
[mass]
type = INSADMass
variable = p
[]
[mass_kernel]
type = MassMatrix
variable = p
matrix_tags = 'mass'
[]
[momentum_convection]
type = INSADMomentumAdvection
variable = vel
[]
[momentum_viscous]
type = INSADMomentumViscous
variable = vel
extra_matrix_tags = 'L'
[]
[momentum_pressure]
type = INSADMomentumPressure
variable = vel
pressure = p
integrate_p_by_parts = true
[]
[]
[BCs]
[no_slip]
type = ADVectorFunctionDirichletBC
variable = vel
boundary = 'bottom right left'
extra_matrix_tags = 'L'
[]
[lid]
type = ADVectorFunctionDirichletBC
variable = vel
boundary = 'top'
function_x = 'lid_function'
extra_matrix_tags = 'L'
[]
[]
[Materials]
[const]
type = ADGenericConstantMaterial
prop_names = 'rho mu'
prop_values = '${rho} ${mu}'
[]
[insad]
type = INSADMaterial
velocity = vel
pressure = p
[]
[]
[Functions]
[lid_function]
# We pick a function that is exactly represented in the velocity
# space so that the Dirichlet conditions are the same regardless
# of the mesh spacing.
type = ParsedFunction
expression = '${prefactor}*${U}*x*(${l}-x)'
[]
[]
[Problem]
type = NavierStokesProblem
mass_matrix = 'mass'
extra_tag_matrices = 'mass L'
L_matrix = 'L'
commute_lsc = true
[]
[Preconditioning]
[FSP]
type = FSP
topsplit = 'up'
[up]
splitting = 'u p'
splitting_type = schur
petsc_options_iname = '-pc_fieldsplit_schur_fact_type -pc_fieldsplit_schur_precondition -ksp_gmres_restart -ksp_type -ksp_pc_side -ksp_rtol'
petsc_options_value = 'full self 300 fgmres right 1e-4'
[]
[u]
vars = 'vel'
# petsc_options = '-ksp_converged_reason'
petsc_options_iname = '-pc_type -pc_hypre_type -ksp_type -ksp_rtol -ksp_gmres_restart -ksp_pc_side'
petsc_options_value = 'hypre boomeramg gmres 1e-2 300 right'
[]
[p]
vars = 'p'
petsc_options = '-ksp_converged_reason -pc_lsc_commute'
petsc_options_iname = '-ksp_type -ksp_gmres_restart -ksp_rtol -pc_type -ksp_pc_side -lsc_pc_type -lsc_pc_hypre_type -lsc_ksp_type -lsc_ksp_rtol -lsc_ksp_pc_side -lsc_ksp_gmres_restart -lsc_mass_pc_type -lsc_mass_pc_hypre_type -lsc_mass_ksp_rtol -lsc_mass_ksp_type'
petsc_options_value = 'fgmres 300 1e-2 lsc right hypre boomeramg fgmres 1e-1 right 300 hypre boomeramg 1e-1 gmres'
[]
[]
[]
[Postprocessors]
[pavg]
type = ElementAverageValue
variable = p
[]
[]
[Correctors]
[set_pressure]
type = NSPressurePin
pin_type = 'average'
variable = p
pressure_average = 'pavg'
[]
[]
[Executioner]
type = Steady
solve_type = NEWTON
[]
[Outputs]
print_linear_residuals = false
[exo]
type = Exodus
execute_on = 'final'
hide = 'pavg'
file_base = 'fsp_steady_low_Re_olshanskii'
[]
[]
(modules/navier_stokes/test/tests/finite_element/ins/lid_driven/transient_fsp.i)
n=64
mu=2e-3
[GlobalParams]
gravity = '0 0 0'
preset = true
supg = false
[]
[Problem]
extra_tag_matrices = 'mass'
previous_nl_solution_required = true
type = NavierStokesProblem
mass_matrix = 'mass'
schur_fs_index = '1'
[]
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = 1.0
ymin = 0
ymax = 1.0
nx = ${n}
ny = ${n}
elem_type = QUAD9
[]
[]
[Variables]
[vel_x]
order = SECOND
family = LAGRANGE
[]
[vel_y]
order = SECOND
family = LAGRANGE
[]
[p]
order = FIRST
family = LAGRANGE
[]
[]
[Kernels]
# mass
[mass]
type = INSMass
variable = p
u = vel_x
v = vel_y
pressure = p
[]
[x_time]
type = INSMomentumTimeDerivative
variable = vel_x
[]
[x_momentum_space]
type = INSMomentumLaplaceForm
variable = vel_x
u = vel_x
v = vel_y
pressure = p
component = 0
[]
[x_mass]
type = MassMatrix
variable = vel_x
matrix_tags = 'mass'
[]
[y_time]
type = INSMomentumTimeDerivative
variable = vel_y
[]
[y_momentum_space]
type = INSMomentumLaplaceForm
variable = vel_y
u = vel_x
v = vel_y
pressure = p
component = 1
[]
[y_mass]
type = MassMatrix
variable = vel_y
matrix_tags = 'mass'
[]
[]
[BCs]
[x_no_slip]
type = DirichletBC
variable = vel_x
boundary = 'bottom right left'
value = 0.0
[]
[lid]
type = FunctionDirichletBC
variable = vel_x
boundary = 'top'
function = 'lid_function'
[]
[y_no_slip]
type = DirichletBC
variable = vel_y
boundary = 'bottom right top left'
value = 0.0
[]
[]
[Materials]
[const]
type = GenericConstantMaterial
block = 0
prop_names = 'rho mu'
prop_values = '1 ${mu}'
[]
[]
[Functions]
[lid_function]
# We pick a function that is exactly represented in the velocity
# space so that the Dirichlet conditions are the same regardless
# of the mesh spacing.
type = ParsedFunction
expression = '4*x*(1-x)'
[]
[]
[Preconditioning]
[FSP]
type = FSP
topsplit = 'by_diri_others'
[by_diri_others]
splitting = 'diri others'
splitting_type = additive
petsc_options_iname = '-ksp_type'
petsc_options_value = 'preonly'
[]
[diri]
sides = 'left right top bottom'
vars = 'vel_x vel_y'
petsc_options_iname = '-pc_type'
petsc_options_value = 'jacobi'
[]
[others]
splitting = 'u p'
splitting_type = schur
petsc_options_iname = '-pc_fieldsplit_schur_fact_type -pc_fieldsplit_schur_precondition -ksp_gmres_restart -ksp_rtol -ksp_type -ksp_atol'
petsc_options_value = 'full self 300 1e-5 fgmres 1e-9'
unside_by_var_boundary_name = 'left top right bottom left top right bottom'
unside_by_var_var_name = 'vel_x vel_x vel_x vel_x vel_y vel_y vel_y vel_y'
[]
[u]
vars = 'vel_x vel_y'
unside_by_var_boundary_name = 'left top right bottom left top right bottom'
unside_by_var_var_name = 'vel_x vel_x vel_x vel_x vel_y vel_y vel_y vel_y'
# petsc_options = '-ksp_converged_reason'
petsc_options_iname = '-pc_type -ksp_pc_side -ksp_type -ksp_rtol -pc_hypre_type -ksp_gmres_restart'
petsc_options_value = 'hypre right gmres 1e-2 boomeramg 300'
[]
[p]
vars = 'p'
petsc_options = '-pc_lsc_scale_diag -ksp_converged_reason'# -lsc_ksp_converged_reason -lsc_ksp_monitor_true_residual
petsc_options_iname = '-ksp_type -ksp_gmres_restart -ksp_rtol -pc_type -ksp_pc_side -lsc_pc_type -lsc_pc_hypre_type -lsc_ksp_type -lsc_ksp_rtol -lsc_ksp_pc_side -lsc_ksp_gmres_restart'
petsc_options_value = 'fgmres 300 1e-2 lsc right hypre boomeramg gmres 1e-1 right 300'
[]
[]
[]
[Postprocessors]
[pavg]
type = ElementAverageValue
variable = p
[]
[]
[Correctors]
[set_pressure]
type = NSPressurePin
pin_type = 'average'
variable = p
pressure_average = 'pavg'
[]
[]
[Executioner]
solve_type = NEWTON
type = Transient
petsc_options_iname = '-snes_max_it'
petsc_options_value = '100'
line_search = 'none'
nl_rel_tol = 1e-8
nl_abs_tol = 1e-9
abort_on_solve_fail = true
normalize_solution_diff_norm_by_dt = false
[TimeStepper]
type = IterationAdaptiveDT
optimal_iterations = 6
dt = 1e-2
[]
steady_state_detection = true
[]
[Outputs]
[exo]
type = Exodus
execute_on = 'final'
hide = 'pavg'
[]
[]