LinearFVPhaseChangeSource

Linear FV elemental kernel that adds the apparent heat-capacity phase-change source term: rho * L * (df/dT) * T_dot, with f the liquid fraction.

Overview

This elemental kernel adds the apparent heat-capacity latent-heat contribution to a temperature-based energy equation solved with the linear finite volume discretization:

$var element = document.getElementById("moose-equation-be140c00-dd87-4629-a8b5-ac74aebf4d0c");katex.render("\\rho L \\frac{\\partial f}{\\partial T} \\frac{\\partial T}{\\partial t},", element, {displayMode:true,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$

where $var element = document.getElementById("moose-equation-ecfa7ef2-9673-47f2-b882-cbba9936ec43");katex.render("\\rho", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ is the mixture density, $var element = document.getElementById("moose-equation-44410000-9d03-4218-ab52-cb223983cb57");katex.render("L", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ the latent heat of fusion, and $var element = document.getElementById("moose-equation-0e117cc0-2059-4e52-9b59-7aa49a8cb8b1");katex.render("f", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ the liquid fraction defined over the mushy interval $var element = document.getElementById("moose-equation-c47715bd-3384-4679-ab31-63c51c0075c1");katex.render("[T_{sol}, T_{liq}]", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ using the normalized temperature

$var element = document.getElementById("moose-equation-815a1d2a-4bf6-4dd4-b521-697734592ca5");katex.render("s = \\min\\left(\\max\\left(\\frac{T - T_{sol}}{T_{liq} - T_{sol}}, 0\\right), 1\\right).", element, {displayMode:true,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$

Two liquid-fraction shapes are available through the "smoothing" parameter:

smooth (default): a smoothstep $var element = document.getElementById("moose-equation-26f929d1-dd0f-4f8f-b7de-62452529a79d");katex.render("f(s) = 3s^2 - 2s^3", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ , giving $var element = document.getElementById("moose-equation-a418c92e-3775-4f21-b4c3-ebdafedf8901");katex.render("\\partial f / \\partial T = 6 s (1 - s) / (T_{liq} - T_{sol})", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ , which vanishes is a regularization using the cubic Hermite interpolation employed for solid–liquid phase interpolation Wang et al. (1993),
sharp: a linear liquid fraction $var element = document.getElementById("moose-equation-b7fbbacd-fcd9-4603-a203-cf5f7d814ce9");katex.render("f(s) = s", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ , giving a constant $var element = document.getElementById("moose-equation-74b1d156-7a27-4070-a729-852581336277");katex.render("\\partial f / \\partial T = 1 / (T_{liq} - T_{sol})", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ inside the interval and zero outside.

Both variants release the same total latent heat across the interval. The matrix contribution treats the temperature time derivative implicitly using the time integrator of the system, while $var element = document.getElementById("moose-equation-49d4ffef-f59b-47f0-bf77-7d8f2af88925");katex.render("\\partial f / \\partial T", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ is evaluated using the current temperature iterate.

note

Ensure that the simulation has the proper spatial-temporal discretization, to converge the liquid-solid front advancement.

note

This kernel is intended for temperature-based solves only. It must not be combined with a total-enthalpy formulation, in which the latent heat is already contained in the solved variable.

Example input syntax

[LinearFVKernels<<<{"href": "../../syntax/LinearFVKernels/index.html"}>>>]

  [energy_time]
    type = LinearFVTimeDerivative<<<{"description": "Represents the matrix and right hand side contributions of a time derivative term in a partial differential equation.", "href": "LinearFVTimeDerivative.html"}>>>
    variable<<<{"description": "The name of the variable whose linear system this object contributes to"}>>> = T
    factor<<<{"description": "A multiplier on the variable within the time derivative. A functor is any of the following: a variable, a functor material property, a function, a postprocessor or a number."}>>> = '${fparse rho * cp}'
  []
  [energy_conduction]
    type = LinearFVDiffusion<<<{"description": "Represents the matrix and right hand side contributions of a diffusion term in a partial differential equation.", "href": "LinearFVDiffusion.html"}>>>
    variable<<<{"description": "The name of the variable whose linear system this object contributes to"}>>> = T
    diffusion_coeff<<<{"description": "The diffusion coefficient. A functor is any of the following: a variable, a functor material property, a function, a postprocessor or a number."}>>> = ${k}
    use_nonorthogonal_correction<<<{"description": "If the nonorthogonal correction should be used when computing the normal gradient."}>>> = false
  []
  [energy_source]
    type = LinearFVPhaseChangeSource<<<{"description": "Linear FV elemental kernel that adds the apparent heat-capacity phase-change source term: rho * L * (df/dT) * T_dot, with f the liquid fraction.", "href": "LinearFVPhaseChangeSource.html"}>>>
    variable<<<{"description": "The name of the variable whose linear system this object contributes to"}>>> = T
    L<<<{"description": "Latent heat. A functor is any of the following: a variable, a functor material property, a function, a postprocessor or a number."}>>> = ${L}
    rho<<<{"description": "The mixture density. A functor is any of the following: a variable, a functor material property, a function, a postprocessor or a number."}>>> = ${rho}
    T_solidus<<<{"description": "The solidus temperature. A functor is any of the following: a variable, a functor material property, a function, a postprocessor or a number."}>>> = ${T_solidus}
    T_liquidus<<<{"description": "The liquidus temperature. A functor is any of the following: a variable, a functor material property, a function, a postprocessor or a number."}>>> = ${T_liquidus}
    smoothing<<<{"description": "Shape of the liquid fraction over the mushy interval. 'smooth' (default) uses a smoothstep cubic function and 'sharp' uses a linear liquid fraction."}>>> = smooth
  []
[]

Input Parameters

LLatent heat. 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:Latent heat. A functor is any of the following: a variable, a functor material property, a function, a postprocessor or a number.
T_liquidusThe liquidus temperature. 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:The liquidus temperature. A functor is any of the following: a variable, a functor material property, a function, a postprocessor or a number.
T_solidusThe solidus temperature. 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:The solidus temperature. A functor is any of the following: a variable, a functor material property, a function, a postprocessor or a number.
rhoThe mixture density. 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:The mixture density. 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

Required 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
smoothingsmoothShape of the liquid fraction over the mushy interval. 'smooth' (default) uses a smoothstep cubic function and 'sharp' uses a linear liquid fraction.
Default:smooth
C++ Type:MooseEnum
Options:smooth, sharp
Controllable:No
Description:Shape of the liquid fraction over the mushy interval. 'smooth' (default) uses a smoothstep cubic function and 'sharp' uses a linear liquid fraction.

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_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)
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/solidification/1d_stefan_T.i)
(modules/navier_stokes/test/tests/finite_volume/ins/solidification/1d-stefan_linearfv.i)
(modules/navier_stokes/examples/solidification/gallium_melting_linear.i)

References

S.-L. Wang, R. F. Sekerka, A. A. Wheeler, B. T. Murray, S. R. Coriell, R. J. Braun, and G. B. McFadden. Thermodynamically-consistent phase-field models for solidification. Physica D: Nonlinear Phenomena, 69(1–2):189–200, 1993. doi:10.1016/0167-2789(93)90189-8.

BibTeX

@article{wang1993thermodynamically,
    author = "Wang, S.-L. and Sekerka, R. F. and Wheeler, A. A. and Murray, B. T. and Coriell, S. R. and Braun, R. J. and McFadden, G. B.",
    title = "Thermodynamically-consistent phase-field models for solidification",
    journal = "Physica D: Nonlinear Phenomena",
    volume = "69",
    number = "1--2",
    pages = "189--200",
    year = "1993",
    doi = "10.1016/0167-2789(93)90189-8"
}

RIS

TY  - JOUR
AU  - Wang, S.-L.
AU  - Sekerka, R. F.
AU  - Wheeler, A. A.
AU  - Murray, B. T.
AU  - Coriell, S. R.
AU  - Braun, R. J.
AU  - McFadden, G. B.
TI  - Thermodynamically-consistent phase-field models for solidification
JO  - Physica D: Nonlinear Phenomena
PY  - 1993
VL  - 69
IS  - 1–2
DO  - 10.1016/0167-2789(93)90189-8
SP  - 189
EP  - 200
ER  -

Plain Text

S.-L. Wang, R. F. Sekerka, A. A. Wheeler, B. T. Murray, S. R. Coriell, R. J. Braun, and G. B. McFadden. Thermodynamically-consistent phase-field models for solidification. Physica D: Nonlinear Phenomena, 69(1–2):189–200, 1993. doi:10.1016/0167-2789(93)90189-8.

(modules/navier_stokes/test/tests/finite_volume/ins/solidification/1d-stefan_linearfv.i)

# Regression test for LinearFVPhaseChangeSource based on the 1D two-phase
# Stefan problem: melting driven by a hot wall at x = 0.
#
# Nondimensional properties rho = cp = k = L = 1 (Stefan number = 1). The
# sharp-interface analytical solution places the front at
#
#   s(tau) = 2 * lambda * sqrt(alpha * tau),
#
# where lambda ~ 0.62 and tau = t + t0. The initial condition is the
# analytical temperature profile at tau = t0, for which the front is at
# s0 = 0.1.
#
# Compared with examples/solidification/1d_stefan_T.i, this input is
# coarsened for speed and uses a wider mushy interval so that the phase-change
# source remains resolved on the coarse mesh. This is a regression test, not
# a validation case.

rho = 1
cp = 1
k = 1
L = 1

T_solidus = 0
T_liquidus = 0.1

T_hot_bc = 1

# Time offset corresponding to an initial front position s0 = 0.1.
t0 = 0.0065023282

Lx = 0.5
nx = 100

[Problem]
  linear_sys_names = 'energy_system'

  # LinearFVPhaseChangeSource evaluates its temperature-dependent
  # coefficient using the previous nonlinear-iteration solution.
  previous_nl_solution_required = true
[]

[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 1
    nx = ${nx}
    xmin = 0.0
    xmax = ${Lx}
  []
[]

[Functions]
  # Analytical sharp-interface temperature profile at tau = t0.
  [T_exact_initial]
    type = ParsedFunction
    expression = 'if(x < 0.1, 1.0 - erf(x / (2.0 * sqrt(${t0}))) / 0.6194595791, 0.0)'
  []
[]

[Variables]
  [T]
    type = MooseLinearVariableFVReal
    solver_sys = energy_system
  []
[]

[AuxVariables]
  [fl]
    type = MooseVariableFVReal
    initial_condition = 0.0
  []
[]

[AuxKernels]
  [compute_fl]
    type = NSLiquidFractionAux
    variable = fl
    temperature = T
    T_liquidus = ${T_liquidus}
    T_solidus = ${T_solidus}
    execute_on = TIMESTEP_END
  []
[]

[ICs]
  [T_ic]
    type = FunctionIC
    variable = T
    function = T_exact_initial
  []
[]

[LinearFVKernels]

  [energy_time]
    type = LinearFVTimeDerivative
    variable = T
    factor = ${fparse rho * cp}
  []
  [energy_conduction]
    type = LinearFVDiffusion
    variable = T
    diffusion_coeff = ${k}
    use_nonorthogonal_correction = false
  []
  [energy_source]
    type = LinearFVPhaseChangeSource
    variable = T
    L = ${L}
    rho = ${rho}
    T_solidus = ${T_solidus}
    T_liquidus = ${T_liquidus}
    smoothing = smooth
  []
[]

[LinearFVBCs]
  # Hot wall at x = 0.
  [hot_left]
    type = LinearFVAdvectionDiffusionFunctorDirichletBC
    variable = T
    boundary = left
    functor = ${T_hot_bc}
  []

  # Far-field solid at the melting temperature.
  [cold_right]
    type = LinearFVAdvectionDiffusionFunctorDirichletBC
    variable = T
    boundary = right
    functor = 0.0
  []
[]

[Postprocessors]
  # Integral of the liquid fraction. In this 1D problem, it approximates
  # the position of the melting front.
  [front_position]
    type = ElementIntegralVariablePostprocessor
    variable = fl
  []

  [T_avg]
    type = ElementAverageValue
    variable = T
  []
[]

[Executioner]
  type = Transient

  dt = 1e-3
  end_time = 1e-2

  system_names = 'energy_system'

  l_abs_tol = 1e-13
  l_tol = 1e-13
  l_max_its = 1000

  petsc_options_iname = '-energy_system_pc_type -energy_system_pc_hypre_type'
  petsc_options_value = 'hypre boomeramg'

  # Repeat the temperature solve so that the lagged apparent heat-capacity
  # coefficient is updated within each time step.
  multi_system_fixed_point = true
  multi_system_fixed_point_convergence = phase_change_fixed_point
[]

[Convergence]
  # Keep the number of coefficient-update iterations deterministic for the
  # regression gold files.
  [phase_change_fixed_point]
    type = IterationCountConvergence
    max_iterations = 25
    converge_at_max_iterations = true
  []
[]

[Outputs]
  exodus = true
  csv = true
[]

(modules/navier_stokes/examples/solidification/1d_stefan_T.i)

rho_salt = 1
cp_salt = 1
k_salt = 1
L = 1

# Numerical regularization of sharp melting at T_m = 0.
# For a sharper run, reduce this to 1e-4 or 1e-5.
# Setting it exactly to zero activates your degenerate sharp branch.
T_solidus = 0
T_liquidus = 1e-2
T_hot_bc = 1

erf_lambda = 0.6194595791366345
t0 = 0.0065023282

simulation_end_time = 0.25

Lx = 2.0
nx = 1000

[Problem]
  previous_nl_solution_required = true
  linear_sys_names = 'energy_system p_system u_system'
[]

[Mesh]
  [salt_mesh]
    type = GeneratedMeshGenerator
    dim = 1
    nx = ${nx}
    xmin = 0.0
    xmax = ${Lx}
    subdomain_ids = 0
    bias_x = 1.0
  []

  [name_salt]
    type = RenameBlockGenerator
    input = salt_mesh
    old_block = 0
    new_block = 'salt'
  []
[]

[Functions]
  # Initial temperature corresponding to the sharp analytical solution.
  [T_exact_initial]
    type = ParsedFunction
    expression = 'if(x < 0.1, 1.0 - erf(x / (2.0 * sqrt(${t0}))) / ${erf_lambda}, 0.0)'
  []
[]

[UserObjects]
  [rc]
    type = RhieChowMassFlux
    u = vel_x
    v = vel_y
    w = vel_z
    pressure = pressure
    rho = ${rho_salt}
    p_diffusion_kernel = p_diffusion
    block = salt
  []
[]

[Variables]
  [T]
    type = MooseLinearVariableFVReal
    solver_sys = energy_system
    block = salt
  []

  # Dummy velocity/pressure variables kept because your PIMPLE block
  # references u_system and p_system, even though they are not solved.
  [vel_x]
    type = MooseLinearVariableFVReal
    solver_sys = 'u_system'
    initial_condition = 0
    block = salt
  []

  [pressure]
    type = MooseLinearVariableFVReal
    solver_sys = 'p_system'
    initial_condition = 0
    block = salt
  []
[]

[AuxVariables]
  [fl]
    type = MooseVariableFVReal
    initial_condition = 0.0
  []
[]

[AuxKernels]
  [compute_fl]
    type = NSLiquidFractionAux
    variable = fl
    temperature = T
    T_liquidus = '${T_liquidus}'
    T_solidus = '${T_solidus}'
    execute_on = 'TIMESTEP_END'
  []
[]

[ICs]
  [h_salt_ic]
    type = FunctionIC
    variable = T
    function = T_exact_initial
    block = salt
  []
[]

[LinearFVKernels]
  # Pressure diffusion kernel kept only because the Rhie-Chow object
  # references it. Pressure solve is disabled in the Executioner.
  [p_diffusion]
    type = LinearFVPressureCorrectionDiffusion
    variable = pressure
    diffusion_tensor = Ainv
    block = salt
    use_nonorthogonal_correction = false
  []

  # rho * dh/dt
  [h_time]
    type = LinearFVTimeDerivative
    variable = T
    factor = ${fparse rho_salt*cp_salt}
    block = salt
  []

  # div(k * dT/dh * grad(h))
  [h_conduction]
    type = LinearFVDiffusion
    variable = T
    diffusion_coeff = ${k_salt}
    use_nonorthogonal_correction = false
    block = salt
    # coeff_interp_method = harm
  []
  [energy_source]
    type = LinearFVPhaseChangeSource
    variable = T
    L = ${L}
    T_liquidus = ${T_liquidus}
    T_solidus = ${T_solidus}
    rho = ${rho_salt}
    smoothing = 'sharp'
  []
[]

[LinearFVBCs]
  # Hot wall at x = 0.
  [hot_left]
    type = LinearFVAdvectionDiffusionFunctorDirichletBC
    variable = T
    boundary = left
    functor = ${T_hot_bc}
  []

  # Far-field solid at T_m = 0.
  [cold_right]
    type = LinearFVAdvectionDiffusionFunctorDirichletBC
    variable = T
    boundary = right
    functor = 0.0
  []
[]

[FunctorMaterials]
  [ins_fv]
    type = INSFVEnthalpyFunctorMaterial
    rho = ${rho_salt}
    cp = ${cp_salt}
    temperature = 'T'
  []
[]

[VectorPostprocessors]
  [horizontal_center]
    type = LineValueSampler
    start_point = '0 0 0'
    end_point = '0.5 0 0'
    num_points = ${fparse nx*0.5/Lx}
    variable = 'fl'
    sort_by = 'x'
    execute_on = 'FINAL'
    warn_discontinuous_face_values = false
  []
[]


[Executioner]
  type = PIMPLE

  dt = 1.25e-4

  # This corresponds to analytical final time tau_end = 0.25
  # because tau = t + t0.
  end_time = ${simulation_end_time}

  num_iterations = 60
  continue_on_max_its = true
  print_fields = false

  energy_system = energy_system

  rhie_chow_user_object = 'rc'
  momentum_systems = 'u_system'
  pressure_system = 'p_system'
  should_solve_momentum = false
  should_solve_pressure = false

  energy_l_abs_tol = 1e-18
  energy_l_tol = 1e-18

  energy_absolute_tolerance = 1e-18

  energy_equation_relaxation = 0.9

  energy_petsc_options_iname = '-pc_type -pc_hypre_type'
  energy_petsc_options_value = 'hypre boomeramg'
[]

[Outputs]
  [out_test]
    type = Exodus
    time_step_interval = 400
  []
[]

(modules/navier_stokes/test/tests/finite_volume/ins/solidification/1d-stefan_linearfv.i)

# Regression test for LinearFVPhaseChangeSource based on the 1D two-phase
# Stefan problem: melting driven by a hot wall at x = 0.
#
# Nondimensional properties rho = cp = k = L = 1 (Stefan number = 1). The
# sharp-interface analytical solution places the front at
#
#   s(tau) = 2 * lambda * sqrt(alpha * tau),
#
# where lambda ~ 0.62 and tau = t + t0. The initial condition is the
# analytical temperature profile at tau = t0, for which the front is at
# s0 = 0.1.
#
# Compared with examples/solidification/1d_stefan_T.i, this input is
# coarsened for speed and uses a wider mushy interval so that the phase-change
# source remains resolved on the coarse mesh. This is a regression test, not
# a validation case.

rho = 1
cp = 1
k = 1
L = 1

T_solidus = 0
T_liquidus = 0.1

T_hot_bc = 1

# Time offset corresponding to an initial front position s0 = 0.1.
t0 = 0.0065023282

Lx = 0.5
nx = 100

[Problem]
  linear_sys_names = 'energy_system'

  # LinearFVPhaseChangeSource evaluates its temperature-dependent
  # coefficient using the previous nonlinear-iteration solution.
  previous_nl_solution_required = true
[]

[Mesh]
  [mesh]
    type = GeneratedMeshGenerator
    dim = 1
    nx = ${nx}
    xmin = 0.0
    xmax = ${Lx}
  []
[]

[Functions]
  # Analytical sharp-interface temperature profile at tau = t0.
  [T_exact_initial]
    type = ParsedFunction
    expression = 'if(x < 0.1, 1.0 - erf(x / (2.0 * sqrt(${t0}))) / 0.6194595791, 0.0)'
  []
[]

[Variables]
  [T]
    type = MooseLinearVariableFVReal
    solver_sys = energy_system
  []
[]

[AuxVariables]
  [fl]
    type = MooseVariableFVReal
    initial_condition = 0.0
  []
[]

[AuxKernels]
  [compute_fl]
    type = NSLiquidFractionAux
    variable = fl
    temperature = T
    T_liquidus = ${T_liquidus}
    T_solidus = ${T_solidus}
    execute_on = TIMESTEP_END
  []
[]

[ICs]
  [T_ic]
    type = FunctionIC
    variable = T
    function = T_exact_initial
  []
[]

[LinearFVKernels]

  [energy_time]
    type = LinearFVTimeDerivative
    variable = T
    factor = ${fparse rho * cp}
  []
  [energy_conduction]
    type = LinearFVDiffusion
    variable = T
    diffusion_coeff = ${k}
    use_nonorthogonal_correction = false
  []
  [energy_source]
    type = LinearFVPhaseChangeSource
    variable = T
    L = ${L}
    rho = ${rho}
    T_solidus = ${T_solidus}
    T_liquidus = ${T_liquidus}
    smoothing = smooth
  []
[]

[LinearFVBCs]
  # Hot wall at x = 0.
  [hot_left]
    type = LinearFVAdvectionDiffusionFunctorDirichletBC
    variable = T
    boundary = left
    functor = ${T_hot_bc}
  []

  # Far-field solid at the melting temperature.
  [cold_right]
    type = LinearFVAdvectionDiffusionFunctorDirichletBC
    variable = T
    boundary = right
    functor = 0.0
  []
[]

[Postprocessors]
  # Integral of the liquid fraction. In this 1D problem, it approximates
  # the position of the melting front.
  [front_position]
    type = ElementIntegralVariablePostprocessor
    variable = fl
  []

  [T_avg]
    type = ElementAverageValue
    variable = T
  []
[]

[Executioner]
  type = Transient

  dt = 1e-3
  end_time = 1e-2

  system_names = 'energy_system'

  l_abs_tol = 1e-13
  l_tol = 1e-13
  l_max_its = 1000

  petsc_options_iname = '-energy_system_pc_type -energy_system_pc_hypre_type'
  petsc_options_value = 'hypre boomeramg'

  # Repeat the temperature solve so that the lagged apparent heat-capacity
  # coefficient is updated within each time step.
  multi_system_fixed_point = true
  multi_system_fixed_point_convergence = phase_change_fixed_point
[]

[Convergence]
  # Keep the number of coefficient-update iterations deterministic for the
  # regression gold files.
  [phase_change_fixed_point]
    type = IterationCountConvergence
    max_iterations = 25
    converge_at_max_iterations = true
  []
[]

[Outputs]
  exodus = true
  csv = true
[]

(modules/navier_stokes/examples/solidification/gallium_melting_linear.i)

##########################################################
# Simulation of Gallium Melting Experiment
# Ref: Gau, C., & Viskanta, R. (1986). Melting and solidification of a pure metal on a vertical wall.
# Key physics: melting/solidification, convective heat transfer, natural convection
##########################################################

mu = 1.81e-3
rho_solid = 6093
rho_liquid = 6093
k_solid = 32
k_liquid = 32
cp_solid = 381.5
cp_liquid = 381.5
L = 80160
alpha_b = 1.2e-4
T_solidus = 302.93
T_liquidus = '${fparse T_solidus + 0.1}'
T_cold = 301.15
T_hot = 311.15
Nx = 200
Ny = 120

[Mesh]
  [gen]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = 0
    xmax = 88.9e-3
    ymin = 0
    ymax = 63.5e-3
    nx = ${Nx}
    ny = ${Ny}
    bias_x = 1.005
  []
[]

[UserObjects]
  [rc]
    type = RhieChowMassFlux
    u = vel_x
    v = vel_y
    pressure = pressure
    rho = ${rho_liquid}
    p_diffusion_kernel = p_diffusion
  []
[]

[Problem]
  linear_sys_names = 'u_system v_system pressure_system energy_system'
  previous_nl_solution_required = true
[]

[AuxVariables]
  [U]
    type = MooseVariableFVReal
  []
  [fl]
    type = MooseVariableFVReal
    initial_condition = 0.0
  []
  [density]
    type = MooseVariableFVReal
  []
  [th_cond]
    type = MooseVariableFVReal
  []
  [cp_var]
    type = MooseVariableFVReal
  []
  [darcy_coef]
    type = MooseVariableFVReal
  []
  [fch_coef]
    type = MooseVariableFVReal
  []
[]

[AuxKernels]
  [mag]
    type = VectorMagnitudeAux
    variable = U
    x = vel_x
    y = vel_y
  []
  [compute_fl]
    type = NSLiquidFractionAux
    variable = fl
    temperature = T
    T_liquidus = '${T_liquidus}'
    T_solidus = '${T_solidus}'
    execute_on = 'TIMESTEP_END'
  []
  [rho_out]
    type = FunctorAux
    functor = 'rho_mixture'
    variable = 'density'
  []
  [th_cond_out]
    type = FunctorAux
    functor = 'k_mixture'
    variable = 'th_cond'
  []
  [cp_out]
    type = FunctorAux
    functor = 'cp_mixture'
    variable = 'cp_var'
  []
  [darcy_out]
    type = FunctorAux
    functor = 'Darcy_coefficient'
    variable = 'darcy_coef'
  []
  [fch_out]
    type = FunctorAux
    functor = 'Forchheimer_coefficient'
    variable = 'fch_coef'
  []
[]

[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.2
  []
  [T]
    type = MooseLinearVariableFVReal
    solver_sys = energy_system
    initial_condition = ${T_cold}
  []
[]

[FVInterpolationMethods]
  [upwind]
    type = FVAdvectedUpwind
  []
[]

[LinearFVKernels]
  [u_time]
    type = LinearFVTimeDerivative
    variable = vel_x
    factor = 'rho_mixture'
  []
  [v_time]
    type = LinearFVTimeDerivative
    variable = vel_y
    factor = 'rho_mixture'
  []
  [u_advection_stress]
    type = LinearWCNSFVMomentumFlux
    variable = vel_x
    mu = ${mu}
    u = vel_x
    v = vel_y
    momentum_component = 'x'
    rhie_chow_user_object = 'rc'
    use_nonorthogonal_correction = false
    advected_interp_method_name = upwind
  []
  [v_advection_stress]
    type = LinearWCNSFVMomentumFlux
    variable = vel_y
    mu = ${mu}
    u = vel_x
    v = vel_y
    momentum_component = 'y'
    rhie_chow_user_object = 'rc'
    use_nonorthogonal_correction = false
    advected_interp_method_name = upwind
  []
  [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_liquid}'
    gravity = '0 -9.81 0'
    alpha_name = ${alpha_b}
    ref_temperature = ${T_solidus}
    T_fluid = T
    momentum_component = 'x'
  []
  [v_boussinesq]
    type = LinearFVMomentumBoussinesq
    variable = vel_y
    rho = '${rho_liquid}'
    gravity = '0 -9.81 0'
    alpha_name = ${alpha_b}
    ref_temperature = ${T_solidus}
    T_fluid = T
    momentum_component = 'y'
  []
  [u_friction_darcy]
    type = LinearFVReaction
    variable = vel_x
    coeff = 'darcy_coef_friction'
  []
  [v_friction_darcy]
    type = LinearFVReaction
    variable = vel_y
    coeff = 'darcy_coef_friction'
  []
  [u_friction_forch]
    type = LinearFVReaction
    variable = vel_x
    coeff = 'forch_coef_friction'
  []
  [v_friction_forch]
    type = LinearFVReaction
    variable = vel_y
    coeff = 'forch_coef_friction'
  []

  [p_diffusion]
    type = LinearFVPressureCorrectionDiffusion
    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 should be mixture. For simplicity and
    # because the properties are the same for liquid and solid
    # the factor is defined as below.
    factor = ${fparse rho_liquid*cp_liquid}
  []
  [h_advection]
    type = LinearFVEnergyAdvection
    variable = T
    advected_quantity = temperature
    cp = ${cp_liquid}
    rhie_chow_user_object = 'rc'
  []
  [conduction]
    type = LinearFVDiffusion
    variable = T
    diffusion_coeff = k_mixture
    use_nonorthogonal_correction = false
  []
  [energy_source]
    type = LinearFVPhaseChangeSource
    variable = T
    L = ${L}
    T_liquidus = ${T_liquidus}
    T_solidus = ${T_solidus}
    rho = 'rho_mixture'
  []
[]

[LinearFVBCs]
  [walls-u]
    type = LinearFVAdvectionDiffusionFunctorDirichletBC
    boundary = 'left right top bottom'
    variable = vel_x
    functor = 0.0
  []
  [walls-v]
    type = LinearFVAdvectionDiffusionFunctorDirichletBC
    boundary = 'left right top bottom'
    variable = vel_y
    functor = 0.0
  []
  [hot_wall]
    type = LinearFVAdvectionDiffusionFunctorDirichletBC
    variable = T
    functor = '${T_hot}'
    boundary = 'left'
  []
  [cold_wall]
    type = LinearFVAdvectionDiffusionFunctorDirichletBC
    variable = T
    functor = '${T_cold}'
    boundary = 'right'
  []
[]

[FunctorMaterials]
  [eff_cp]
    type = NSFVMixtureFunctorMaterial
    phase_2_names = '${cp_solid} ${k_solid} ${rho_solid}'
    phase_1_names = '${cp_liquid} ${k_liquid} ${rho_liquid}'
    prop_names = 'cp_mixture k_mixture rho_mixture'
    phase_1_fraction = fl
  []
  [mushy_zone_resistance]
    type = INSFVMushyPorousFrictionFunctorMaterial
    liquid_fraction = 'fl'
    mu = '${mu}'
    rho_l = '${rho_liquid}'
    dendrite_spacing_scaling = 1e-3
  []
  [darcy_coeff_friction]
    type = ParsedFunctorMaterial
    property_name = 'darcy_coef_friction'
    functor_names = 'darcy_coef'
    functor_symbols = 'darcy_coef'
    expression = 'darcy_coef'
  []
  [forch_coeff_friction]
    type = ParsedFunctorMaterial
    property_name = 'forch_coef_friction'
    functor_names = 'U fch_coef'
    functor_symbols = 'U fch_coef'
    expression = 'fch_coef * U'
  []
[]

[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.7
  pressure_variable_relaxation = 0.3
  energy_equation_relaxation = 0.9
  num_iterations = 8
  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.005

  start_time = 0.0
  end_time = 60.0
  num_piso_iterations = 0

  pin_pressure = true
  pressure_pin_value = 0.0
  pressure_pin_point = '0.05 0.05 0.0'
[]

[Postprocessors]
  [ave_p]
    type = ElementAverageValue
    variable = 'pressure'
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [ave_fl]
    type = ElementAverageValue
    variable = 'fl'
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [ave_T]
    type = ElementAverageValue
    variable = 'T'
    execute_on = 'INITIAL TIMESTEP_END'
  []
[]

[Outputs]
  [out_avg]
    type = Exodus
    time_step_interval = 2000
  []
  csv = false
[]

Overview
Example input syntax
Input Parameters
Input Files
References