INSBoussinesqBodyForce

construction:Undocumented Class

The INSBoussinesqBodyForce has not been documented. The content listed below should be used as a starting point for documenting the class, which includes the typical automatic documentation associated with a MooseObject; however, what is contained is ultimately determined by what is necessary to make the documentation clear for users.

Computes a body force for natural convection buoyancy.

Overview

Example Input File Syntax

Input Parameters

component0,1,2 depending on if we are solving the x,y,z component of the momentum equation
C++ Type:unsigned int
Controllable:No
Description:0,1,2 depending on if we are solving the x,y,z component of the momentum equation
gravityDirection of the gravity vector
C++ Type:libMesh::VectorValue<double>
Controllable:No
Description:Direction of the gravity vector
temperaturetemperature variable, for off diagonal jacobian
C++ Type:std::vector<VariableName>
Controllable:No
Description:temperature variable, for off diagonal jacobian
variableThe name of the variable that this residual object operates on
C++ Type:NonlinearVariableName
Controllable:No
Description:The name of the variable that this residual object operates on

Required Parameters

alpha_namealphaThe name of the thermal expansion coefficientthis is of the form rho = rho*(1-alpha (T-T_ref))
Default:alpha
C++ Type:MaterialPropertyName
Controllable:No
Description:The name of the thermal expansion coefficientthis is of the form rho = rho*(1-alpha (T-T_ref))
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
displacementsThe displacements
C++ Type:std::vector<VariableName>
Controllable:No
Description:The displacements
prop_getter_suffixAn optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
C++ Type:MaterialPropertyName
Controllable:No
Description:An optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
ref_temptemp_refThe name of the reference temperature
Default:temp_ref
C++ Type:MaterialPropertyName
Controllable:No
Description:The name of the reference temperature
rho_namerhoThe name of the density
Default:rho
C++ Type:MaterialPropertyName
Controllable:No
Description:The name of the density
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.

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_tagsnontimeThe tag for the vectors this Kernel should fill
Default:nontime
C++ Type:MultiMooseEnum
Options:nontime, time
Controllable:No
Description:The tag for the vectors this Kernel should fill

Tagging Parameters

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.
diag_save_inThe name of auxiliary variables to save this Kernel's diagonal Jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
C++ Type:std::vector<AuxVariableName>
Controllable:No
Description:The name of auxiliary variables to save this Kernel's diagonal Jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
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
save_inThe name of auxiliary variables to save this Kernel's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
C++ Type:std::vector<AuxVariableName>
Controllable:No
Description:The name of auxiliary variables to save this Kernel's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
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

Input Files

(problems/publication_level_cases/manuBench/manubox.i)
(tests/laminar_flow/boussinesq_square.i)
(problems/2021-cnrs-benchmark/phase-2/transient.i)
(problems/2021-cnrs-benchmark/phase-1/buoyancy.i)
(problems/2021-cnrs-benchmark/phase-1/full-coupling.i)

(problems/publication_level_cases/manuBench/manubox.i)

library_path=/home/gavin/projects/moose/modules/navier_stokes/lib/

[GlobalParams]
  # neutronics
  num_groups=6
  num_precursor_groups=8
  use_exp_form = false
  group_fluxes = 'group1 group2 group3 group4 group5 group6'
  temperature = temp
  sss2_input = false
  pre_concs = 'pre1 pre2 pre3 pre4 pre5 pre6 pre7 pre8'
  account_delayed = true

  # fluids
  gravity = '0 -9.81 0'
  pspg = true
  convective_term = true
  integrate_p_by_parts = true
  laplace = true
  u = vel_x
  v = vel_y
  p = p
  alpha = 1e-6 # stabilisation parameter
  order = FIRST
  family = LAGRANGE
[]

[Mesh]
  type = GeneratedMesh
  dim = 2
  xmin = 0
  xmax = 1.0
  ymin = 0
  ymax = 1.0
  elem_type = QUAD9
  nx = 20
  ny = 20
[]

[MeshModifiers]
  [./corner_node]
    type = AddExtraNodeset
    new_boundary = 'pinned_node'
    nodes = '0'
  [../]
[]

[Variables]

  # velocities
  [./vel_x]
    family = LAGRANGE
    order = FIRST
  [../]

  [./vel_y]
    family = LAGRANGE
    order = FIRST
  [../]

  # pressure
  [./p]
    family = LAGRANGE
    order = FIRST
  [../]

  # temperature
  [./temp]
    family = LAGRANGE
    order = FIRST
    initial_condition = 900
    scaling = 1e-4
  [../]
[]

# nt -> neutron
[Nt]
  var_name_base = group
  vacuum_boundaries = 'top bottom left right'
  create_temperature_var = false
  eigen = true
[]

[Precursors]
  [./preblock]
    var_name_base = pre
    # block = 'fuel'
    constant_velocity_values = false
    nt_exp_form = false
    family = MONOMIAL
    order = CONSTANT
    outlet_boundaries = ''
    uvel = vel_x
    vvel = vel_y
    transient = false
  [../]
[]

[Kernels]
  # mass
  [./mass]
    type = INSMass
    variable = p
  [../]

  # x-momentum, space
  [./x_momentum_space]
    type = INSMomentumLaplaceForm
    variable = vel_x
    component = 0
  [../]

  # y-momentum, space
  [./y_momentum_space]
    type = INSMomentumLaplaceForm
    variable = vel_y
    component = 1
  [../]

  [./tempAdvectionDiffusion]
    type = INSTemperature
    variable = temp
    u = vel_x
    v = vel_y
  [../]
  [./buoyancy_y]
    type = INSBoussinesqBodyForce
    variable = vel_y
    dT = deltaT
    component = 1
    temperature = temp
  [../]

  [./heatRemoval]
    type = ConvectiveHeatExchanger
    variable = temp
    tref = 900.0
    htc  = 1.0
  [../]
[]

[BCs]
  [./ux_dirichlet]
    type = DirichletBC
    boundary = 'left right bottom top'
    variable = vel_x
    value = 0
  [../]
  [./uy_dirichlet]
    type = DirichletBC
    boundary = 'left right bottom top'
    variable = vel_y
    value = 0
  [../]
  [./p_pin]
    type=DirichletBC
    boundary = 'pinned_node'
    variable = p
    value = 0
  [../]
[]

[Materials]
  # alpha = coefficient of thermal expansion where rho  = rho0 -alpha * rho0 * delta T
  # Dr. Aufiero uses nu rather than mu, (kinematic rather than absolute viscosity)
  # likewise, \rho c_p is given rather than separate \rho, c_p
  # accordingly, we let \rho = 1 to obtain equivalent results
  # also, Manu gives a length expansion coefficient. Boussinesq approximations
  # typically use an expansion coefficient with units 1/K, but for mass instead.
  # The relation between the two is alpha = -3 beta, if linearized.
  # Also, had to use Pr = mu * cp / k
  [./const]
    type = GenericConstantMaterial
    prop_names = 'mu     rho alpha   k         cp'
    prop_values = '0.025 1.0 6.0e-4  0.5  6.15e6'
  [../]

  [./fuelneutronicproperties]
    type = GenericMoltresMaterial
    property_tables_root = './groupconstants/'
    interp_type = 'none' # 'spline'
  [../]
[]

[Preconditioning]
  [./SMP]
    type = SMP
    full = true
    solve_type = 'PJFNK'
  [../]
[]

[Executioner]
  type = InversePowerMethod
  max_power_iterations = 50

  # fission power normalization:
  normalization = powernorm
  normal_factor = 1e9 # watts

  # power iteration options
  xdiff = 'group1diff'
  bx_norm = 'bnorm'
  pfactor = 1e-2

  petsc_options = '-snes_converged_reason -ksp_converged_reason'
  petsc_options_iname = '-pc_type -pc_factor_shift_type'
  petsc_options_value = 'lu NONZERO'
  line_search = 'none'
  nl_rel_tol = 1e-12
  nl_abs_tol = 1e-13
  nl_max_its = 6
  l_tol = 1e-6
  l_max_its = 500
[]

[Outputs]
  [./exodus]
    type = Exodus
  [../]
  [./csv]
    type = CSV
  [../]
[]

[Postprocessors]
  [./bnorm]
    type = ElmIntegTotFissNtsPostprocessor
    execute_on = linear
  [../]
  [./group1diff]
    type = ElementL2Diff
    variable = group1
    execute_on = 'linear timestep_end'
    use_displaced_mesh = false
  [../]
  [./powernorm]
    type = ElmIntegTotFissHeatPostprocessor
    execute_on = linear
  [../]
[]

[AuxVariables]
  [./deltaT]
    type = ConstantDifferenceAux
    variable = deltaT
    comparevar = temp
    constant = 900
  [../]
[]

[AuxKernels]
  [./deltaTCalc]
    type =  ConstantDifferenceAux
    variable = deltaT
    compareVar = temp
    constant = 900
  [../]
[]

[Debug]
  show_var_residual_norms = true
[]

(tests/laminar_flow/boussinesq_square.i)

xmax = .05
ymax = .05

[GlobalParams]
  #integrate_p_by_parts = false
  gravity = '0 -9.81 0'
[]

[Mesh]
  [generate]
    type = GeneratedMeshGenerator
    dim = 2
    xmax = ${xmax}
    ymax = ${ymax}
    nx = 20
    ny = 20
    elem_type = QUAD9
  []
  [bottom_left]
    type = ExtraNodesetGenerator
    input = generate
    new_boundary = corner
    coord = '0 0'
  []
[]

[Problem]
[]

[Preconditioning]
  [Newton_SMP]
    type = SMP
    full = true
    solve_type = 'NEWTON'
  []
[]

[Executioner]
  type = Transient
  # Run for 100+ timesteps to reach steady state.
  num_steps = 5
  dt = 0.01
  petsc_options = '-snes_converged_reason -ksp_converged_reason -snes_test_display'
  petsc_options_iname = '-pc_type -pc_factor_shift_type -pc_factor_shift_amount'
  petsc_options_value = 'lu	  NONZERO		1e-10'
  # petsc_options_iname = '-snes_type'
  # petsc_options_value = 'test'
  line_search = 'none'
  nl_rel_tol = 1e-12
  nl_max_its = 20
  l_tol = 1e-6
  l_max_its = 500

  dtmin = 1e-5
  #[./TimeStepper]
  #  type = IterationAdaptiveDT
  #  dt = 1e-3
  #  cutback_factor = 0.4
  #  growth_factor = 1.2
  #  optimal_iterations = 20
  #[../]

[]

[Debug]
  show_var_residual_norms = true
[]

[Outputs]
  perf_graph = true
  print_linear_residuals = true
  [out]
    type = Exodus
    execute_on = 'final'
  []
[]

[Variables]
  [ux]
    family = LAGRANGE
    order = SECOND
  []
  [uy]
    family = LAGRANGE
    order = SECOND
  []
  [p]
    family = LAGRANGE
    order = FIRST
  []
  [temp]
    family = LAGRANGE
    order = SECOND
    initial_condition = 340
    scaling = 1e-4
  []
[]

[BCs]
  [ux_dirichlet]
    type = DirichletBC
    boundary = 'left right bottom top'
    variable = ux
    value = 0
  []
  [uy_dirichlet]
    type = DirichletBC
    boundary = 'left right bottom top'
    variable = uy
    value = 0
  []
  [p_zero]
    type = DirichletBC
    boundary = corner
    variable = p
    value = 0
  []
  [temp_insulate]
    type = NeumannBC
    variable = temp
    value = 0 # no conduction through side walls
    boundary = 'top bottom' # not top
  []
  [coldOnTop]
    type = DirichletBC
    variable = temp
    boundary = left
    value = 300
  []
  [hotOnBottom]
    type = DirichletBC
    variable = temp
    boundary = right
    value = 400
  []
[]

[Kernels]
  [mass]
    type = INSMass
    variable = p
    u = ux
    v = uy
    pressure = p
  []
  [x_time_deriv]
    type = INSMomentumTimeDerivative
    variable = ux
  []
  [y_time_deriv]
    type = INSMomentumTimeDerivative
    variable = uy
  []
  [x_momentum_space]
    type = INSMomentumLaplaceForm
    variable = ux
    u = ux
    v = uy
    pressure = p
    component = 0
  []
  [y_momentum_space]
    type = INSMomentumLaplaceForm
    variable = uy
    u = ux
    v = uy
    pressure = p
    component = 1
  []
  [tempTimeDeriv]
    type = MatINSTemperatureTimeDerivative
    variable = temp
  []
  [tempAdvectionDiffusion]
    type = INSTemperature
    variable = temp
    u = ux
    v = uy
  []
  [buoyancy_x]
    type = INSBoussinesqBodyForce
    variable = ux
    component = 0
    temperature = temp
  []
  [buoyancy_y]
    type = INSBoussinesqBodyForce
    variable = uy
    component = 1
    temperature = temp
  []
[]

[Materials]
  [const]
    type = GenericConstantMaterial
    # alpha = coefficient of thermal expansion where rho  = rho0 -alpha * rho0 * delta T
    prop_names = 'mu        rho    temp_ref alpha   k        cp'
    prop_values = '30.74e-6  .5757  900      2.9e-3  46.38e-3 1054'
  []
[]

(problems/2021-cnrs-benchmark/phase-2/transient.i)

density = .002      # kg cm-3
cp = 3075           # J kg-1 K-1, 6.15 / 2.0e-3
k = .005            # W cm-1 K-1
gamma = 1           # W cm-3 K-1, Volumetric heat transfer coefficient
viscosity = .5      # dynamic viscosity
alpha = 1           # SUPG stabilization parameter
t_alpha = 2e-4      # K-1, Thermal expansion coefficient

# Change "freq" and "dt" values for different perturbation frequencies
freq = 0.8          # Perturbation frequency = 0.8 Hz
dt = 0.00625        # Timestep size = 1 / freq / 200 = 0.00625 s

[GlobalParams]
  num_groups = 6
  num_precursor_groups = 8
  use_exp_form = false
  group_fluxes = 'group1 group2 group3 group4 group5 group6'
  pre_concs = 'pre1 pre2 pre3 pre4 pre5 pre6 pre7 pre8'
  temperature = temp
  sss2_input = true
  account_delayed = true
  integrate_p_by_parts = true
  eigenvalue_scaling = 0.9927821802
  ## Use the eigenvalue scaling factor below if running on a 40x40 mesh
  # eigenvalue_scaling = 0.9926551482
[]

[Mesh]
  [file_mesh]
    type = FileMeshGenerator
    file = '../phase-1/full-coupling_out_ntsApp0_exodus.e'
    use_for_exodus_restart = true
  []
  [corner_node]
    type = ExtraNodesetGenerator
    input = file_mesh
    new_boundary = 'pinned_node'
    coord = '200 200'
  []
[]

[Problem]
  type = FEProblem
[]

[Variables]
  [temp]
    family = LAGRANGE
    order = FIRST
  []
  [vel]
    family = LAGRANGE_VEC
    order = FIRST
  []
  [p]
    family = LAGRANGE
    order = FIRST
  []
[]

[Nt]
  var_name_base = group
  vacuum_boundaries = 'bottom left right top'
  create_temperature_var = false
  init_nts_from_file = true
[]

[Precursors]
  [pres]
    var_name_base = pre
    outlet_boundaries = ''
    constant_velocity_values = false
    uvel = velx
    vvel = vely
    nt_exp_form = false
    family = MONOMIAL
    order = CONSTANT
    loop_precursors = false
    transient = true
    init_from_file = true
  []
[]

[AuxVariables]
  [velx]
    family = LAGRANGE
    order = FIRST
  []
  [vely]
    family = LAGRANGE
    order = FIRST
  []
  [heat]
    family = MONOMIAL
    order = FIRST
  []
[]

[Kernels]
  [mass]
    type = INSADMass
    variable = p
  []
  [mass_pspg]
    type = INSADMassPSPG
    variable = p
  []
  [momentum_time]
    type = INSADMomentumTimeDerivative
    variable = vel
  []
  [momentum_advection]
    type = INSADMomentumAdvection
    variable = vel
  []
  [momentum_viscous]
    type = INSADMomentumViscous
    variable = vel
  []
  [momentum_pressure]
    type = INSADMomentumPressure
    variable = vel
    pressure = p
  []
  [momentum_supg]
    type = INSADMomentumSUPG
    variable = vel
    velocity = vel
  []
  [buoyancy]
    type = INSADBoussinesqBodyForce
    variable = vel
    # gravity vector, cm s-2
    gravity = '0 -981 0'
    alpha_name = 't_alpha'
    ref_temp = 'temp_ref'
  []
  [gravity]
    type = INSADGravityForce
    variable = vel
    # gravity vector, cm s-2
    gravity = '0 -981 0'
  []

  [temp_time]
    type = INSADHeatConductionTimeDerivative
    variable = temp
  []
  [temp_source]
    type = INSADEnergySource
    variable = temp
    source_variable = heat
  []
  [temp_advection]
    type = INSADEnergyAdvection
    variable = temp
  []
  [temp_conduction]
    type = ADHeatConduction
    variable = temp
    thermal_conductivity = 'k'
  []
  [temp_supg]
    type = INSADEnergySUPG
    variable = temp
    velocity = vel
  []
  [temp_sink]
    type = ConvectiveHeatExchanger
    variable = temp
    htc = ${gamma}
    tref = 900
  []
[]

[AuxKernels]
  [vel_x]
    type = VectorVariableComponentAux
    variable = velx
    vector_variable = vel
    component = 'x'
  []
  [vel_y]
    type = VectorVariableComponentAux
    variable = vely
    vector_variable = vel
    component = 'y'
  []
  [heat_source]
    type = FissionHeatSourceTransientAux
    variable = heat
  []
[]

[UserObjects]
  [initial_th]
    type = SolutionUserObject
    mesh = '../phase-1/full-coupling_exodus.e'
    system_variables = 'velx vely p temp heat'
    timestep = LATEST
    execute_on = INITIAL
  []
[]

[Controls]
  [func_control]
    type = RealFunctionControl
    parameter = 'Kernels/temp_sink/htc'
    function = func_alpha
    execute_on = 'initial timestep_begin'
  []
[]

[Functions]
  [func_alpha]
    type = ParsedFunction
    value = '1 + 0.1 * sin(2*pi*t*${freq})' # Perturbation frequency = 0.8Hz
  []
  [velxf]
    type = SolutionFunction
    from_variable = velx
    solution = initial_th
  []
  [velyf]
    type = SolutionFunction
    from_variable = vely
    solution = initial_th
  []
  [pf]
    type = SolutionFunction
    from_variable = p
    solution = initial_th
  []
  [tempf]
    type = SolutionFunction
    from_variable = temp
    solution = initial_th
  []
  [heatf]
    type = SolutionFunction
    from_variable = heat
    solution = initial_th
  []
[]

[ICs]
  [vel_ic]
    type = VectorFunctionIC
    variable = vel
    function_x = velxf
    function_y = velyf
  []
  [p_ic]
    type = FunctionIC
    variable = p
    function = pf
  []
  [temp_ic]
    type = FunctionIC
    variable = temp
    function = tempf
  []
  [heat_ic]
    type = FunctionIC
    variable = heat
    function = heatf
  []
[]

[BCs]
  [no_slip]
    type = VectorDirichletBC
    variable = vel
    boundary = 'bottom left right'
    values = '0 0 0'
  []
  [lid]
    type = VectorDirichletBC
    variable = vel
    boundary = 'top'
    values = '50 0 0'
  []
  [pressure_pin]
    type = DirichletBC
    variable = p
    boundary = 'pinned_node'
    value = 0
  []
[]

[Materials]
  [fuel]
    type = GenericMoltresMaterial
    property_tables_root = '../../../property_file_dir/cnrs-benchmark/benchmark_'
    interp_type = 'linear'
    prop_names = 'temp_ref'
    prop_values = '900'
  []
  [ad_mat]
    type = ADGenericConstantMaterial
    prop_names = 'k rho cp mu t_alpha'
    prop_values = '${k} ${density} ${cp} ${viscosity} ${t_alpha}'
  []
  [ins_temp]
    type = INSADStabilized3Eqn
    alpha = ${alpha}
    velocity = vel
    pressure = p
    temperature = temp
  []
[]

[Executioner]
  type = Transient
  scheme = bdf2
  num_steps = 5000

  solve_type = 'NEWTON'
  petsc_options = '-snes_converged_reason -ksp_converged_reason -snes_linesearch_monitor'
  petsc_options_iname = '-pc_type -sub_pc_type -ksp_gmres_restart -pc_asm_overlap -sub_pc_factor_shift_type'
  petsc_options_value = 'asm      lu           200                1                NONZERO'
  line_search = 'none'

  automatic_scaling = true
  compute_scaling_once = false
  resid_vs_jac_scaling_param = 0.1
  scaling_group_variables = 'group1 group2 group3 group4 group5 group6; pre1 pre2 pre3 pre4 pre5 pre6 pre7 pre8; vel; p; temp'

  nl_abs_tol = 1e-6
  nl_rel_tol = 1e-6
  nl_max_its = 25
  l_max_its = 400
  l_tol = 1e-4

  dt = ${dt}
[]

[Preconditioning]
  [SMP]
    type = SMP
    full = true
  []
[]

[Postprocessors]
  [memory]
    type = MemoryUsage
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [tot_fissions]
    type = ElmIntegTotFissPostprocessor
    execute_on = linear
  []
  [powernorm]
    type = ElmIntegTotFissHeatPostprocessor
    execute_on = linear
  []
  [group1norm]
    type = ElementIntegralVariablePostprocessor
    variable = group1
    execute_on = linear
  []
  [group1max]
    type = NodalExtremeValue
    value_type = max
    variable = group1
    execute_on = timestep_end
  []
  [group2norm]
    type = ElementIntegralVariablePostprocessor
    variable = group2
    execute_on = linear
  []
  [group2max]
    type = NodalExtremeValue
    value_type = max
    variable = group2
    execute_on = timestep_end
  []
[]

[VectorPostprocessors]
[]

[Outputs]
  perf_graph = true
  [exodus]
    type = Exodus
    interval = 20
  []
  [csv]
    type = CSV
    append_restart = true
  []
  [checkpoint]
    type = Checkpoint
    num_files = 3
    interval = 20
  []
[]

[Debug]
  show_var_residual_norms = true
[]

(problems/2021-cnrs-benchmark/phase-1/buoyancy.i)

# Step 1.3: Buoyancy input file

density = .002 # kg cm-3
cp = 3075 # J kg-1 K-1, 6.15 / 2.0e-3
k = .005 # W cm-1 K-1
gamma = 1 # W cm-3 K-1, Volumetric heat transfer coefficient
viscosity = .5 # dynamic viscosity
alpha = 1 # SUPG stabilization parameter
t_alpha = 2e-4 # K-1, Thermal expansion coefficient

[GlobalParams]
  temperature = temp
  integrate_p_by_parts = true
[]

[Mesh]
  [square]
    type = GeneratedMeshGenerator
    dim = 2

    nx = 200
    ny = 200
    ## Use a 40-by-40 mesh instead if running on a less capable computer
    #    nx = 40
    #    ny = 40

    xmin = 0
    xmax = 200
    ymin = 0
    ymax = 200
    elem_type = QUAD4
  []
  [corner_node]
    type = ExtraNodesetGenerator
    input = square
    new_boundary = 'pinned_node'
    coord = '200 200'
  []
[]

[Problem]
  type = FEProblem
[]

[Variables]
  [temp]
    family = LAGRANGE
    order = FIRST
    initial_condition = 900
  []
  [vel]
    family = LAGRANGE_VEC
    order = FIRST
  []
  [p]
    family = LAGRANGE
    order = FIRST
  []
[]

[AuxVariables]
  [velx]
    family = LAGRANGE
    order = FIRST
  []
  [vely]
    family = LAGRANGE
    order = FIRST
  []
  [heat]
    family = MONOMIAL
    order = FIRST
    initial_condition = 250
  []
[]

[Kernels]
  [mass]
    type = INSADMass
    variable = p
  []
  [mass_pspg]
    type = INSADMassPSPG
    variable = p
  []
  [momentum_time]
    type = INSADMomentumTimeDerivative
    variable = vel
  []
  [momentum_advection]
    type = INSADMomentumAdvection
    variable = vel
  []
  [momentum_viscous]
    type = INSADMomentumViscous
    variable = vel
  []
  [momentum_pressure]
    type = INSADMomentumPressure
    variable = vel
    pressure = p
  []
  [momentum_supg]
    type = INSADMomentumSUPG
    variable = vel
    velocity = vel
  []
  [buoyancy]
    type = INSADBoussinesqBodyForce
    variable = vel
    # gravity vector, cm s-2
    gravity = '0 -981 0'
    alpha_name = 't_alpha'
    ref_temp = 'temp_ref'
  []
  [gravity]
    type = INSADGravityForce
    variable = vel
    # gravity vector, cm s-2
    gravity = '0 -981 0'
  []

  [temp_time]
    type = INSADHeatConductionTimeDerivative
    variable = temp
  []
  [temp_source]
    type = INSADEnergySource
    variable = temp
    source_variable = heat
  []
  [temp_advection]
    type = INSADEnergyAdvection
    variable = temp
  []
  [temp_conduction]
    type = ADHeatConduction
    variable = temp
    thermal_conductivity = 'k'
  []
  [temp_supg]
    type = INSADEnergySUPG
    variable = temp
    velocity = vel
  []
  [temp_sink]
    type = INSADEnergyAmbientConvection
    variable = temp
    # alpha is the heat transfer coefficient.
    alpha = ${gamma}
    T_ambient = 900
  []
[]

[AuxKernels]
  [vel_x]
    type = VectorVariableComponentAux
    variable = velx
    vector_variable = vel
    component = 'x'
  []
  [vel_y]
    type = VectorVariableComponentAux
    variable = vely
    vector_variable = vel
    component = 'y'
  []
[]

[ICs]
  [vel_ic]
    type = VectorConstantIC
    x_value = 1e-15
    y_value = 1e-15
    variable = vel
  []
[]

[BCs]
  [no_slip]
    type = VectorDirichletBC
    variable = vel
    boundary = 'bottom left right top'
    values = '0 0 0'
  []
  [pressure_pin]
    type = DirichletBC
    variable = p
    boundary = 'pinned_node'
    value = 0
  []
[]

[Materials]
  [fuel]
    type = GenericConstantMaterial
    prop_names = 'temp_ref'
    prop_values = '900'
  []
  [ad_mat]
    type = ADGenericConstantMaterial
    prop_names = 'k rho cp mu t_alpha'
    prop_values = '${k} ${density} ${cp} ${viscosity} ${t_alpha}'
  []
  [ins_temp]
    type = INSADStabilized3Eqn
    alpha = ${alpha}
    velocity = vel
    pressure = p
    temperature = temp
  []
[]

[Executioner]
  type = Transient
  end_time = 2000

  solve_type = 'PJFNK'
  petsc_options = '-snes_converged_reason -ksp_converged_reason -snes_linesearch_monitor'
  petsc_options_iname = '-pc_type -sub_pc_type -ksp_gmres_restart -pc_asm_overlap -sub_pc_factor_shift_type'
  petsc_options_value = 'asm      lu           200                1               NONZERO'
  line_search = none

  automatic_scaling = true
  compute_scaling_once = false
  off_diagonals_in_auto_scaling = true
  resid_vs_jac_scaling_param = .1

  nl_abs_tol = 1e-10

  dtmin = 1e-1
  dtmax = 10
  steady_state_detection = true
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 1e-1
    cutback_factor = .5
    growth_factor = 1.5
    optimal_iterations = 20
    iteration_window = 4
    linear_iteration_ratio = 1000
  []
[]

[VectorPostprocessors]
  [vel_aa]
    type = LineValueSampler
    variable = 'velx vely'
    start_point = '0 100 0'
    end_point = '200 100 0'
    num_points = 201
    sort_by = x
    execute_on = FINAL
    outputs = 'csv'
  []
  [vel_bb]
    type = LineValueSampler
    variable = 'velx vely'
    start_point = '100 0 0'
    end_point = '100 200 0'
    num_points = 201
    sort_by = y
    execute_on = FINAL
    outputs = 'csv'
  []
  [temp_aa]
    type = LineValueSampler
    variable = 'temp'
    start_point = '0 100 0'
    end_point = '200 100 0'
    num_points = 201
    sort_by = x
    execute_on = FINAL
    outputs = 'csv'
  []
  [temp_bb]
    type = LineValueSampler
    variable = 'temp'
    start_point = '100 0 0'
    end_point = '100 200 0'
    num_points = 201
    sort_by = y
    execute_on = FINAL
    outputs = 'csv'
  []
[]

[MultiApps]
  [ntsApp]
    type = FullSolveMultiApp
    app_type = MoltresApp
    execute_on = timestep_end
    positions = '0 0 0'
    input_files = 'buoyancy-nts.i'
  []
[]

[Transfers]
  [to_sub_temp]
    type = MultiAppCopyTransfer
    to_multi_app = ntsApp
    source_variable = temp
    variable = temp
  []
  [to_sub_vel_x]
    type = MultiAppCopyTransfer
    to_multi_app = ntsApp
    source_variable = velx
    variable = vel_x
  []
  [to_sub_vel_y]
    type = MultiAppCopyTransfer
    to_multi_app = ntsApp
    source_variable = vely
    variable = vel_y
  []
  [from_sub]
    type = MultiAppCopyTransfer
    from_multi_app = ntsApp
    source_variable = heat
    variable = heat
  []
[]

[Outputs]
  perf_graph = true
  print_linear_residuals = true
  [exodus]
    type = Exodus
  []
  [csv]
    type = CSV
    execute_on = FINAL
  []
[]

(problems/2021-cnrs-benchmark/phase-1/full-coupling.i)

# Step 1.4: Full coupling input file

density = .002 # kg cm-3
cp = 3075 # J kg-1 K-1, 6.15 / 2.0e-3
k = .005 # W cm-1 K-1
gamma = 1 # W cm-3 K-1, Volumetric heat transfer coefficient
viscosity = .5 # dynamic viscosity
alpha = 1 # SUPG stabilization parameter
t_alpha = 2e-4 # K-1, Thermal expansion coefficient

[GlobalParams]
  temperature = temp
  integrate_p_by_parts = true
[]

[Mesh]
  [square]
    type = GeneratedMeshGenerator
    dim = 2

    nx = 200
    ny = 200
    ## Use a 40-by-40 mesh instead if running on a less capable computer
    #    nx = 40
    #    ny = 40

    xmin = 0
    xmax = 200
    ymin = 0
    ymax = 200
    elem_type = QUAD4
  []
  [corner_node]
    type = ExtraNodesetGenerator
    input = square
    new_boundary = 'pinned_node'
    coord = '200 200'
  []
[]

[Problem]
  type = FEProblem
[]

[Variables]
  [temp]
    family = LAGRANGE
    order = FIRST
    initial_condition = 900
  []
  [vel]
    family = LAGRANGE_VEC
    order = FIRST
  []
  [p]
    family = LAGRANGE
    order = FIRST
  []
[]

[AuxVariables]
  [velx]
    family = LAGRANGE
    order = FIRST
  []
  [vely]
    family = LAGRANGE
    order = FIRST
  []
  [heat]
    family = MONOMIAL
    order = FIRST
    initial_condition = 250
  []
[]

[Kernels]
  [mass]
    type = INSADMass
    variable = p
  []
  [mass_pspg]
    type = INSADMassPSPG
    variable = p
  []
  [momentum_time]
    type = INSADMomentumTimeDerivative
    variable = vel
  []
  [momentum_advection]
    type = INSADMomentumAdvection
    variable = vel
  []
  [momentum_viscous]
    type = INSADMomentumViscous
    variable = vel
  []
  [momentum_pressure]
    type = INSADMomentumPressure
    variable = vel
    pressure = p
  []
  [momentum_supg]
    type = INSADMomentumSUPG
    variable = vel
    velocity = vel
  []
  [buoyancy]
    type = INSADBoussinesqBodyForce
    variable = vel
    # gravity vector, cm s-2
    gravity = '0 -981 0'
    alpha_name = 't_alpha'
    ref_temp = 'temp_ref'
  []
  [gravity]
    type = INSADGravityForce
    variable = vel
    # gravity vector, cm s-2
    gravity = '0 -981 0'
  []

  [temp_time]
    type = INSADHeatConductionTimeDerivative
    variable = temp
  []
  [temp_source]
    type = INSADEnergySource
    variable = temp
    source_variable = heat
  []
  [temp_advection]
    type = INSADEnergyAdvection
    variable = temp
  []
  [temp_conduction]
    type = ADHeatConduction
    variable = temp
    thermal_conductivity = 'k'
  []
  [temp_supg]
    type = INSADEnergySUPG
    variable = temp
    velocity = vel
  []
  [temp_sink]
    type = INSADEnergyAmbientConvection
    variable = temp
    # alpha is the heat transfer coefficient.
    alpha = ${gamma}
    T_ambient = 900
  []
[]

[AuxKernels]
  [vel_x]
    type = VectorVariableComponentAux
    variable = velx
    vector_variable = vel
    component = 'x'
  []
  [vel_y]
    type = VectorVariableComponentAux
    variable = vely
    vector_variable = vel
    component = 'y'
  []
[]

[ICs]
  [vel_ic]
    type = VectorConstantIC
    x_value = 1e-15
    y_value = 1e-15
    variable = vel
  []
[]

[BCs]
  [no_slip]
    type = VectorDirichletBC
    variable = vel
    boundary = 'bottom left right'
    values = '0 0 0'
  []
  [lid]
    type = VectorDirichletBC
    variable = vel
    boundary = 'top'
    # '50 0 0' corresponds to velx = 50 cm s-1 along the top boundary.
    # Change to 10, 20, 30, 40 for different velocity boundary conditions
    # at U_lid = 0.1, 0.2, 0.3, 0.4 m s-1.
    values = '50 0 0'
  []
  [pressure_pin]
    type = DirichletBC
    variable = p
    boundary = 'pinned_node'
    value = 0
  []
[]

[Materials]
  [fuel]
    type = GenericConstantMaterial
    prop_names = 'temp_ref'
    prop_values = '900'
  []
  [ad_mat]
    type = ADGenericConstantMaterial
    prop_names = 'k rho cp mu t_alpha'
    prop_values = '${k} ${density} ${cp} ${viscosity} ${t_alpha}'
  []
  [ins_temp]
    type = INSADStabilized3Eqn
    alpha = ${alpha}
    velocity = vel
    pressure = p
    temperature = temp
  []
[]

[Executioner]
  type = Transient
  end_time = 2000

  solve_type = 'PJFNK'
  petsc_options = '-snes_converged_reason -ksp_converged_reason -snes_linesearch_monitor'
  petsc_options_iname = '-pc_type -sub_pc_type -ksp_gmres_restart -pc_asm_overlap -sub_pc_factor_shift_type'
  petsc_options_value = 'asm      lu           200                1               NONZERO'
  line_search = none

  automatic_scaling = true
  compute_scaling_once = false
  off_diagonals_in_auto_scaling = true
  resid_vs_jac_scaling_param = .1

  nl_abs_tol = 1e-10

  dtmin = 1e-1
  dtmax = 10
  steady_state_detection = true
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 1e-1
    cutback_factor = .5
    growth_factor = 1.5
    optimal_iterations = 20
    iteration_window = 4
    linear_iteration_ratio = 1000
  []
[]

[MultiApps]
  [ntsApp]
    type = FullSolveMultiApp
    app_type = MoltresApp
    execute_on = timestep_end
    positions = '0 0 0'
    input_files = 'full-coupling-nts.i'
  []
[]

[Transfers]
  [to_sub_temp]
    type = MultiAppCopyTransfer
    to_multi_app = ntsApp
    source_variable = temp
    variable = temp
  []
  [to_sub_vel_x]
    type = MultiAppCopyTransfer
    to_multi_app = ntsApp
    source_variable = velx
    variable = vel_x
  []
  [to_sub_vel_y]
    type = MultiAppCopyTransfer
    to_multi_app = ntsApp
    source_variable = vely
    variable = vel_y
  []
  [from_sub]
    type = MultiAppCopyTransfer
    from_multi_app = ntsApp
    source_variable = heat
    variable = heat
  []
[]

[Outputs]
  perf_graph = true
  print_linear_residuals = true
  [exodus]
    type = Exodus
  []
[]

Overview
Example Input File Syntax
Input Parameters
Input Files