- muMixture Density. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.
C++ Type:MooseFunctorName
Unit:(no unit assumed)
Controllable:No
Description:Mixture Density. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.
- rhoContinuous phase density. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.
C++ Type:MooseFunctorName
Unit:(no unit assumed)
Controllable:No
Description:Continuous phase density. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.
- uThe velocity in the x direction. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.
C++ Type:MooseFunctorName
Unit:(no unit assumed)
Controllable:No
Description:The velocity in the x direction. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.
NSFVDispersePhaseDragFunctorMaterial
This material computes the linear drag coefficient for a dispersed phase based on the particle Reynolds number . The particle Reynolds number is defined as follows:
where:
is the density of the dispersed phase particles,
is the characteristic diameter of the dispersed phase particles,
is the mixture velocity,
is the mixture viscosity.
Based on this Reynolds number, the linear drag coefficient for the dispersed phase is computed as follows Schiller (1933):
Input Parameters
- blockThe list of blocks (ids or names) that this object will be applied
C++ Type:std::vector<SubdomainName>
Unit:(no unit assumed)
Controllable:No
Description:The list of blocks (ids or names) that this object will be applied
- boundaryThe list of boundaries (ids or names) from the mesh where this object applies
C++ Type:std::vector<BoundaryName>
Unit:(no unit assumed)
Controllable:No
Description:The list of boundaries (ids or names) from the mesh where this object applies
- constant_onNONEWhen ELEMENT, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps.When SUBDOMAIN, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps. Evaluations on element qps will be skipped
Default:NONE
C++ Type:MooseEnum
Unit:(no unit assumed)
Options:NONE, ELEMENT, SUBDOMAIN
Controllable:No
Description:When ELEMENT, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps.When SUBDOMAIN, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps. Evaluations on element qps will be skipped
- declare_suffixAn optional suffix parameter that can be appended to any declared 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 declared properties. The suffix will be prepended with a '_' character.
- drag_coef_nameDarcy_coefficientName of the scalar friction coefficient defined. The vector coefficient is suffixed with _vec. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.
Default:Darcy_coefficient
C++ Type:MooseFunctorName
Unit:(no unit assumed)
Controllable:No
Description:Name of the scalar friction coefficient defined. The vector coefficient is suffixed with _vec. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.
- execute_onALWAYSThe 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:ALWAYS
C++ Type:ExecFlagEnum
Unit:(no unit assumed)
Options:XFEM_MARK, FORWARD, ADJOINT, HOMOGENEOUS_FORWARD, ADJOINT_TIMESTEP_BEGIN, ADJOINT_TIMESTEP_END, NONE, INITIAL, LINEAR, NONLINEAR_CONVERGENCE, NONLINEAR, POSTCHECK, TIMESTEP_END, TIMESTEP_BEGIN, MULTIAPP_FIXED_POINT_END, MULTIAPP_FIXED_POINT_BEGIN, FINAL, CUSTOM, ALWAYS
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.
- particle_diameter1Diameter of particles in the dispersed phase. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.
Default:1
C++ Type:MooseFunctorName
Unit:(no unit assumed)
Controllable:No
Description:Diameter of particles in the dispersed phase. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.
- 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
Unit:(no unit assumed)
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.
- vThe velocity in the y direction. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.
C++ Type:MooseFunctorName
Unit:(no unit assumed)
Controllable:No
Description:The velocity in the y direction. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.
- wThe velocity in the z direction. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.
C++ Type:MooseFunctorName
Unit:(no unit assumed)
Controllable:No
Description:The velocity in the z direction. A functor is any of the following: a variable, a functor material property, a function, a post-processor, or a number.
Optional Parameters
- control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Unit:(no unit assumed)
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
Unit:(no unit assumed)
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
Unit:(no unit assumed)
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
Unit:(no unit assumed)
Controllable:No
Description:The seed for the master random number generator
Advanced Parameters
- output_propertiesList of material properties, from this material, to output (outputs must also be defined to an output type)
C++ Type:std::vector<std::string>
Unit:(no unit assumed)
Controllable:No
Description:List of material properties, from this material, to output (outputs must also be defined to an output type)
- outputsnone Vector of output names where you would like to restrict the output of variables(s) associated with this object
Default:none
C++ Type:std::vector<OutputName>
Unit:(no unit assumed)
Controllable:No
Description:Vector of output names where you would like to restrict the output of variables(s) associated with this object
Outputs Parameters
Input Files
- (modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/channel-advection-slip.i)
- (modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/channel-drift-flux-physics.i)
- (modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_interface_area_model/turbulent_driven_growth.i)
- (modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/lid-driven-two-phase.i)
- (modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/channel-drift-flux-w-interface-area.i)
- (modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/channel-drift-flux-transient.i)
- (modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_interface_area_model/pressure_driven_growth_transient.i)
- (modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_interface_area_model/pressure_driven_growth.i)
- (modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/channel-drift-flux.i)
- (modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/rayleigh-bernard-two-phase.i)
References
- Links Schiller.
A drag coefficient correlation.
Zeit. Ver. Deutsch. Ing., 77:318–320, 1933.[BibTeX]
@article{schiller1933drag, author = "Schiller, Links", title = "A drag coefficient correlation", journal = "Zeit. Ver. Deutsch. Ing.", volume = "77", pages = "318--320", year = "1933" }
(modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/channel-advection-slip.i)
mu = 1.0
rho = 10.0
mu_d = 0.1
rho_d = 1.0
l = 2
U = 1
dp = 0.01
inlet_phase_2 = 0.1
advected_interp_method = 'average'
velocity_interp_method = 'rc'
[GlobalParams]
rhie_chow_user_object = 'rc'
mu_interp_method = 'average'
[]
[UserObjects]
[rc]
type = INSFVRhieChowInterpolator
u = vel_x
v = vel_y
pressure = pressure
[]
[]
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = '${fparse l * 5}'
ymin = '${fparse -l / 2}'
ymax = '${fparse l / 2}'
nx = 10
ny = 6
[]
uniform_refine = 0
[]
[Variables]
[vel_x]
type = INSFVVelocityVariable
initial_condition = 0
[]
[vel_y]
type = INSFVVelocityVariable
initial_condition = 0
[]
[pressure]
type = INSFVPressureVariable
[]
[phase_2]
type = INSFVScalarFieldVariable
[]
[]
[FVKernels]
[mass]
type = INSFVMassAdvection
variable = pressure
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
[]
[u_advection]
type = INSFVMomentumAdvection
variable = vel_x
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
momentum_component = 'x'
[]
[u_advection_slip]
type = WCNSFV2PMomentumAdvectionSlip
variable = vel_x
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = ${rho}
rho_d = ${rho_d}
fd = phase_2
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
momentum_component = 'x'
[]
[u_viscosity]
type = INSFVMomentumDiffusion
variable = vel_x
mu = 'mu_mixture'
limit_interpolation = true
momentum_component = 'x'
[]
[u_pressure]
type = INSFVMomentumPressure
variable = vel_x
momentum_component = 'x'
pressure = pressure
[]
[u_friction]
type = PINSFVMomentumFriction
Darcy_name = Darcy_coefficient_vec
is_porous_medium = false
momentum_component = x
mu = mu_mixture
rho = rho_mixture
variable = vel_x
[]
[v_advection]
type = INSFVMomentumAdvection
variable = vel_y
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
momentum_component = 'y'
[]
[v_advection_slip]
type = WCNSFV2PMomentumAdvectionSlip
variable = vel_y
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = ${rho}
rho_d = ${rho_d}
fd = phase_2
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
momentum_component = 'y'
[]
[v_viscosity]
type = INSFVMomentumDiffusion
variable = vel_y
mu = 'mu_mixture'
limit_interpolation = true
momentum_component = 'y'
[]
[v_pressure]
type = INSFVMomentumPressure
variable = vel_y
momentum_component = 'y'
pressure = pressure
[]
[v_friction]
type = PINSFVMomentumFriction
Darcy_name = Darcy_coefficient_vec
is_porous_medium = false
momentum_component = y
mu = mu_mixture
rho = rho_mixture
variable = vel_y
[]
[phase_2_advection]
type = INSFVScalarFieldAdvection
variable = phase_2
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = 'upwind'
[]
[phase_2_src]
type = NSFVMixturePhaseInterface
variable = phase_2
phase_coupled = phase_1
alpha = 0.1
[]
[]
[FVBCs]
[inlet-u]
type = INSFVInletVelocityBC
boundary = 'left'
variable = vel_x
functor = '${U}'
[]
[inlet-v]
type = INSFVInletVelocityBC
boundary = 'left'
variable = vel_y
functor = '0'
[]
[walls-u]
type = INSFVNoSlipWallBC
boundary = 'top bottom'
variable = vel_x
function = 0
[]
[walls-v]
type = INSFVNoSlipWallBC
boundary = 'top bottom'
variable = vel_y
function = 0
[]
[outlet_p]
type = INSFVOutletPressureBC
boundary = 'right'
variable = pressure
function = '0'
[]
[inlet_phase_2]
type = FVDirichletBC
boundary = 'left'
variable = phase_2
value = ${inlet_phase_2}
[]
[]
[AuxVariables]
[drag_coefficient]
type = MooseVariableFVReal
[]
[rho_mixture_var]
type = MooseVariableFVReal
[]
[mu_mixture_var]
type = MooseVariableFVReal
[]
[]
[AuxKernels]
[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'
[]
[]
[FunctorMaterials]
[phase_1]
property_name = 'phase_1'
type = ADParsedFunctorMaterial
functor_names = 'phase_2'
expression = '1 - phase_2'
outputs = 'out'
output_properties = 'phase_1'
[]
[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'
outputs = 'out'
output_properties = 'vel_slip_x'
[]
[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'
outputs = 'out'
output_properties = 'vel_slip_y'
[]
[CD]
type = NSFVDispersePhaseDragFunctorMaterial
rho = 'rho_mixture'
mu = mu_mixture
u = 'vel_x'
v = 'vel_y'
particle_diameter = ${dp}
[]
[mixing_material]
type = NSFVMixtureFunctorMaterial
phase_2_names = '${rho} ${mu}'
phase_1_names = '${rho_d} ${mu_d}'
prop_names = 'rho_mixture mu_mixture'
phase_1_fraction = 'phase_2'
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
nl_rel_tol = 1e-10
[]
[Preconditioning]
[SMP]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
[]
[]
[Outputs]
[out]
type = Exodus
hide = 'Re lin cum_lin'
[]
[]
[Postprocessors]
[Re]
type = ParsedPostprocessor
expression = '${rho} * ${l} * ${U}'
[]
[lin]
type = NumLinearIterations
[]
[cum_lin]
type = CumulativeValuePostprocessor
postprocessor = lin
[]
[]
(modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/channel-drift-flux-physics.i)
mu = 1.0
rho = 10.0
mu_d = 0.1
rho_d = 1.0
l = 2
U = 1
dp = 0.01
inlet_phase_2 = 0.1
advected_interp_method = 'average'
velocity_interp_method = 'rc'
# TODO remove need for those
cp = 1
k = 1
cp_d = 1
k_d = 1
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = '${fparse l * 5}'
ymin = '${fparse -l / 2}'
ymax = '${fparse l / 2}'
nx = 10
ny = 4
[]
uniform_refine = 0
[]
[Physics]
[NavierStokes]
[Flow]
[flow]
compressibility = 'incompressible'
density = 'rho_mixture'
dynamic_viscosity = 'mu_mixture'
# Initial conditions
initial_velocity = '0 0 0'
initial_pressure = 0
# Boundary conditions
inlet_boundaries = 'left'
momentum_inlet_types = 'fixed-velocity'
momentum_inlet_functors = '${U} 0'
wall_boundaries = 'top bottom'
momentum_wall_types = 'noslip noslip'
outlet_boundaries = 'right'
momentum_outlet_types = 'fixed-pressure'
pressure_functors = '0'
# Friction is done in drift flux term
friction_types = "Darcy"
friction_coeffs = "Darcy_vec"
standard_friction_formulation = true
mass_advection_interpolation = '${advected_interp_method}'
momentum_advection_interpolation = '${advected_interp_method}'
velocity_interpolation = '${velocity_interp_method}'
mu_interp_method = 'average'
[]
[]
[TwoPhaseMixture]
[mixture]
phase_1_fraction_name = 'phase_1'
phase_2_fraction_name = 'phase_2'
# Phase transport equation
add_phase_transport_equation = true
alpha_exchange = 0.1
phase_advection_interpolation = 'upwind'
# see flow for inlet boundaries
phase_fraction_inlet_type = 'fixed-value'
phase_fraction_inlet_functors = '${inlet_phase_2}'
# Needed for some reason
ghost_layers = 5
# Drift flux parameters
add_drift_flux_momentum_terms = true
density_interp_method = 'average'
slip_linear_friction_name = 'Darcy'
# Base phase material properties
phase_1_density_name = ${rho}
phase_1_viscosity_name = ${mu}
phase_1_specific_heat_name = ${cp}
phase_1_thermal_conductivity_name = ${k}
use_dispersed_phase_drag_model = true
particle_diameter = 0.01
# Other phase material properties
phase_2_density_name = ${rho_d}
phase_2_viscosity_name = ${mu_d}
phase_2_specific_heat_name = ${cp_d}
phase_2_thermal_conductivity_name = ${k_d}
output_all_properties = true
[]
[]
[]
[]
[FunctorMaterials]
[CD]
type = NSFVDispersePhaseDragFunctorMaterial
drag_coef_name = 'Darcy'
rho = 'rho_mixture'
mu = mu_mixture
u = 'vel_x'
v = 'vel_y'
particle_diameter = ${dp}
outputs = 'all'
output_properties = 'Darcy_coefficient'
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
nl_rel_tol = 1e-10
[]
[Preconditioning]
[SMP]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
[]
[]
[Outputs]
print_linear_residuals = true
print_nonlinear_residuals = true
dofmap = true
[out]
type = Exodus
hide = 'Re lin cum_lin'
[]
[perf]
type = PerfGraphOutput
[]
[]
[Postprocessors]
[Re]
type = ParsedPostprocessor
expression = '${rho} * ${l} * ${U}'
[]
[lin]
type = NumLinearIterations
[]
[cum_lin]
type = CumulativeValuePostprocessor
postprocessor = lin
[]
[]
(modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_interface_area_model/turbulent_driven_growth.i)
###############################################################################
# Validation test based on Hibiki and Ishii experiment [1] reported in Figure 5
# [1] Hibiki, T., & Ishii, M. (2000). One-group interfacial area transport of
# bubbly flows in vertical round tubes.
# International Journal of Heat and Mass Transfer, 43(15), 2711-2726.
###############################################################################
mu = 1.0
rho = 1000.0
mu_d = 1.0
rho_d = 1.0
l = ${fparse 50.8/1000.0}
U = 5.031429
dp = 0.005
inlet_phase_2 = 0.442
advected_interp_method = 'upwind'
velocity_interp_method = 'rc'
mass_exchange_coeff = 0.0
inlet_interface_area = ${fparse 6.0*inlet_phase_2/dp}
outlet_pressure = 1e5
[GlobalParams]
rhie_chow_user_object = 'rc'
density_interp_method = 'average'
mu_interp_method = 'average'
[]
[Problem]
identify_variable_groups_in_nl = false
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = INSFVRhieChowInterpolator
u = vel_x
v = vel_y
pressure = pressure
[]
[]
[Mesh]
coord_type = 'RZ'
rz_coord_axis = 'X'
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = '${fparse l * 60}'
ymin = 0
ymax = '${fparse l / 2}'
nx = 20
ny = 5
[]
uniform_refine = 0
[]
[Variables]
[vel_x]
type = INSFVVelocityVariable
initial_condition = 0
[]
[vel_y]
type = INSFVVelocityVariable
initial_condition = 0
[]
[pressure]
type = INSFVPressureVariable
[]
[phase_2]
type = INSFVScalarFieldVariable
initial_condition = ${inlet_phase_2}
[]
[interface_area]
type = INSFVScalarFieldVariable
initial_condition = ${inlet_interface_area}
[]
[]
[FVKernels]
[mass]
type = INSFVMassAdvection
variable = pressure
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = ${rho}
[]
[u_advection]
type = INSFVMomentumAdvection
variable = vel_x
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
momentum_component = 'x'
[]
[u_drift]
type = WCNSFV2PMomentumDriftFlux
variable = vel_x
rho_d = ${rho_d}
fd = 'rho_mixture_var'
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
momentum_component = 'x'
[]
[u_viscosity]
type = INSFVMomentumDiffusion
variable = vel_x
mu = 'mu_mixture'
limit_interpolation = true
momentum_component = 'x'
[]
[u_pressure]
type = INSFVMomentumPressure
variable = vel_x
momentum_component = 'x'
pressure = pressure
[]
[v_advection]
type = INSFVMomentumAdvection
variable = vel_y
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
momentum_component = 'y'
[]
[v_drift]
type = WCNSFV2PMomentumDriftFlux
variable = vel_y
rho_d = ${rho_d}
fd = 'rho_mixture_var'
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
momentum_component = 'y'
[]
[v_viscosity]
type = INSFVMomentumDiffusion
variable = vel_y
mu = 'mu_mixture'
limit_interpolation = true
momentum_component = 'y'
[]
[v_pressure]
type = INSFVMomentumPressure
variable = vel_y
momentum_component = 'y'
pressure = pressure
[]
[phase_2_advection]
type = INSFVScalarFieldAdvection
variable = phase_2
u_slip = 'vel_x'
v_slip = 'vel_y'
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = 'upwind'
[]
[phase_2_diffusion]
type = FVDiffusion
variable = phase_2
coeff = 1.0
[]
[phase_2_src]
type = NSFVMixturePhaseInterface
variable = phase_2
phase_coupled = phase_1
alpha = ${mass_exchange_coeff}
[]
[interface_area_advection]
type = INSFVScalarFieldAdvection
variable = interface_area
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = 'upwind'
[]
[interface_area_diffusion]
type = FVDiffusion
variable = interface_area
coeff = 0.1
[]
[interface_area_source_sink]
type = WCNSFV2PInterfaceAreaSourceSink
variable = interface_area
u = 'vel_x'
v = 'vel_y'
L = ${fparse l/2}
rho = 'rho_mixture'
rho_d = 'rho'
pressure = 'pressure'
k_c = '${fparse mass_exchange_coeff}'
fd = 'phase_2'
sigma = 1e-3
[]
[]
[FVBCs]
[inlet-u]
type = INSFVInletVelocityBC
boundary = 'left'
variable = vel_x
functor = '${U}'
[]
[inlet-v]
type = INSFVInletVelocityBC
boundary = 'left'
variable = vel_y
functor = '0'
[]
[walls-u]
type = INSFVNoSlipWallBC
boundary = 'top'
variable = vel_x
function = 0
[]
[walls-v]
type = INSFVNoSlipWallBC
boundary = 'top'
variable = vel_y
function = 0
[]
[outlet_p]
type = INSFVOutletPressureBC
boundary = 'right'
variable = pressure
function = '${outlet_pressure}'
[]
[inlet_phase_2]
type = FVDirichletBC
boundary = 'left'
variable = phase_2
value = ${inlet_phase_2}
[]
[inlet_interface_area]
type = FVDirichletBC
boundary = 'left'
variable = interface_area
value = ${inlet_interface_area}
[]
[symmetry-u]
type = PINSFVSymmetryVelocityBC
boundary = 'bottom'
variable = vel_x
u = vel_x
v = vel_y
mu = 'mu_mixture'
momentum_component = 'x'
[]
[symmetry-v]
type = PINSFVSymmetryVelocityBC
boundary = 'bottom'
variable = vel_y
u = vel_x
v = vel_y
mu = 'mu_mixture'
momentum_component = 'y'
[]
[symmetry-p]
type = INSFVSymmetryPressureBC
boundary = 'bottom'
variable = pressure
[]
[symmetry-phase-2]
type = INSFVSymmetryScalarBC
boundary = 'bottom'
variable = phase_2
[]
[symmetry-interface-area]
type = INSFVSymmetryScalarBC
boundary = 'bottom'
variable = interface_area
[]
[]
[AuxVariables]
[drag_coefficient]
type = MooseVariableFVReal
[]
[rho_mixture_var]
type = MooseVariableFVReal
[]
[mu_mixture_var]
type = MooseVariableFVReal
[]
[]
[AuxKernels]
[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'
[]
[]
[FluidProperties]
[fp]
type = IdealGasFluidProperties
[]
[]
[FunctorMaterials]
[bubble_properties]
type = GeneralFunctorFluidProps
fp = 'fp'
pressure = 'pressure'
T_fluid = 300.0
speed = 1.0
characteristic_length = 1.0
porosity = 1.0
output_properties = 'rho'
outputs = 'out'
[]
[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'
[]
[compute_phase_1]
type = ADParsedFunctorMaterial
property_name = phase_1
functor_names = 'phase_2'
expression = '1 - phase_2'
[]
[CD]
type = NSFVDispersePhaseDragFunctorMaterial
rho = 'rho_mixture'
mu = mu_mixture
u = 'vel_x'
v = 'vel_y'
particle_diameter = ${dp}
[]
[mixing_material]
type = NSFVMixtureFunctorMaterial
phase_2_names = '${rho} ${mu}'
phase_1_names = 'rho ${mu_d}'
prop_names = 'rho_mixture mu_mixture'
phase_1_fraction = 'phase_2'
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
nl_rel_tol = 1e-10
line_search = 'none'
[]
[Debug]
show_var_residual_norms = true
[]
[Preconditioning]
[SMP]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
[]
[]
[Outputs]
[out]
type = Exodus
[]
[]
[Postprocessors]
[Re]
type = ParsedPostprocessor
expression = '${rho} * ${l} * ${U}'
pp_names = ''
[]
[rho_outlet]
type = SideAverageValue
boundary = 'right'
variable = 'rho_mixture_var'
[]
[]
(modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/lid-driven-two-phase.i)
mu = 1.0
rho = 1.0e3
mu_d = 0.3
rho_d = 1.0
dp = 0.01
U_lid = 0.1
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 = 10
ny = 10
[]
[]
[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
[]
[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'
[]
# NOTE: the friction terms for u and v are missing
[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
[]
[]
[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'
[]
[]
[FunctorMaterials]
[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'
gravity = '0 ${g} 0'
[]
[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'
gravity = '0 ${g} 0'
[]
[compute_phase_1]
type = ADParsedFunctorMaterial
property_name = phase_1
functor_names = 'phase_2'
expression = '1 - phase_2'
[]
[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'
[]
[]
[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 = 7
iteration_window = 2
growth_factor = 2.0
cutback_factor = 0.5
dt = 1e-3
[]
nl_max_its = 10
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_volume/two_phase/mixture_model/channel-drift-flux-w-interface-area.i)
mu = 10.0
rho = 100.0
mu_d = 1.0
rho_d = 1.0
l = 2
U = 1
dp = 0.01
inlet_phase_2 = 0.0
advected_interp_method = 'upwind'
velocity_interp_method = 'rc'
mass_exchange_coeff = 0.01
[GlobalParams]
rhie_chow_user_object = 'rc'
density_interp_method = 'average'
mu_interp_method = 'average'
[]
[Problem]
identify_variable_groups_in_nl = false
[]
[UserObjects]
[rc]
type = INSFVRhieChowInterpolator
u = vel_x
v = vel_y
pressure = pressure
[]
[]
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = '${fparse l * 5}'
ymin = '${fparse -l / 2}'
ymax = '${fparse l / 2}'
nx = 20
ny = 5
[]
uniform_refine = 0
[]
[Variables]
[vel_x]
type = INSFVVelocityVariable
initial_condition = 0
[]
[vel_y]
type = INSFVVelocityVariable
initial_condition = 0
[]
[pressure]
type = INSFVPressureVariable
[]
[phase_2]
type = INSFVScalarFieldVariable
[]
[interface_area]
type = INSFVScalarFieldVariable
[]
[]
[FVKernels]
inactive = 'u_time v_time phase_2_time interface_area_time'
[mass]
type = INSFVMassAdvection
variable = pressure
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = ${rho}
[]
[u_time]
type = INSFVMomentumTimeDerivative
variable = vel_x
rho = 'rho_mixture'
momentum_component = 'x'
[]
[u_advection]
type = INSFVMomentumAdvection
variable = vel_x
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
momentum_component = 'x'
[]
[u_drift]
type = WCNSFV2PMomentumDriftFlux
variable = vel_x
rho_d = ${rho_d}
fd = 'rho_mixture_var'
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
momentum_component = 'x'
[]
[u_viscosity]
type = INSFVMomentumDiffusion
variable = vel_x
mu = 'mu_mixture'
limit_interpolation = true
momentum_component = 'x'
[]
[u_pressure]
type = INSFVMomentumPressure
variable = vel_x
momentum_component = 'x'
pressure = pressure
[]
[v_time]
type = INSFVMomentumTimeDerivative
variable = vel_y
rho = 'rho_mixture'
momentum_component = 'y'
[]
[v_advection]
type = INSFVMomentumAdvection
variable = vel_y
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
momentum_component = 'y'
[]
[v_drift]
type = WCNSFV2PMomentumDriftFlux
variable = vel_y
rho_d = ${rho_d}
fd = 'rho_mixture_var'
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
momentum_component = 'y'
[]
[v_viscosity]
type = INSFVMomentumDiffusion
variable = vel_y
mu = 'mu_mixture'
limit_interpolation = true
momentum_component = 'y'
[]
[v_pressure]
type = INSFVMomentumPressure
variable = vel_y
momentum_component = 'y'
pressure = pressure
[]
[phase_2_time]
type = FVFunctorTimeKernel
variable = phase_2
functor = phase_2
[]
[phase_2_advection]
type = INSFVScalarFieldAdvection
variable = phase_2
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = 'upwind'
[]
[phase_2_diffusion]
type = FVDiffusion
variable = phase_2
coeff = 1.0
[]
[phase_2_src]
type = NSFVMixturePhaseInterface
variable = phase_2
phase_coupled = phase_1
alpha = ${mass_exchange_coeff}
[]
[interface_area_time]
type = FVFunctorTimeKernel
variable = interface_area
functor = interface_area
[]
[interface_area_advection]
type = INSFVScalarFieldAdvection
variable = interface_area
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = 'upwind'
[]
[interface_area_diffusion]
type = FVDiffusion
variable = interface_area
coeff = 0.1
[]
[interface_area_source_sink]
type = WCNSFV2PInterfaceAreaSourceSink
variable = interface_area
u = 'vel_x'
v = 'vel_y'
L = 1.0
rho = 'rho_mixture'
rho_d = ${rho_d}
pressure = 'pressure'
k_c = ${fparse mass_exchange_coeff * 100.0}
fd = 'phase_2'
sigma = 1e-3
[]
[]
[FVBCs]
[inlet-u]
type = INSFVInletVelocityBC
boundary = 'left'
variable = vel_x
functor = '${U}'
[]
[inlet-v]
type = INSFVInletVelocityBC
boundary = 'left'
variable = vel_y
functor = '0'
[]
[walls-u]
type = INSFVNoSlipWallBC
boundary = 'top bottom'
variable = vel_x
function = 0
[]
[walls-v]
type = INSFVNoSlipWallBC
boundary = 'top bottom'
variable = vel_y
function = 0
[]
[outlet_p]
type = INSFVOutletPressureBC
boundary = 'right'
variable = pressure
function = '0'
[]
[inlet_phase_2]
type = FVDirichletBC
boundary = 'left'
variable = phase_2
value = ${inlet_phase_2}
[]
[inlet_interface_area]
type = FVDirichletBC
boundary = 'left'
variable = interface_area
value = 0.0
[]
[]
[AuxVariables]
[drag_coefficient]
type = MooseVariableFVReal
[]
[rho_mixture_var]
type = MooseVariableFVReal
[]
[mu_mixture_var]
type = MooseVariableFVReal
[]
[]
[AuxKernels]
[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'
[]
[]
[FunctorMaterials]
[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'
[]
[compute_phase_1]
type = ADParsedFunctorMaterial
property_name = phase_1
functor_names = 'phase_2'
expression = '1 - phase_2'
[]
[CD]
type = NSFVDispersePhaseDragFunctorMaterial
rho = 'rho_mixture'
mu = mu_mixture
u = 'vel_x'
v = 'vel_y'
particle_diameter = ${dp}
[]
[mixing_material]
type = NSFVMixtureFunctorMaterial
phase_2_names = '${rho} ${mu}'
phase_1_names = '${rho_d} ${mu_d}'
prop_names = 'rho_mixture mu_mixture'
phase_1_fraction = 'phase_2'
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
nl_rel_tol = 1e-10
# dt = 0.1
# end_time = 1.0
# nl_max_its = 10
[]
[Debug]
show_var_residual_norms = true
[]
[Preconditioning]
[SMP]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
[]
[]
[Outputs]
exodus = true
[]
[Postprocessors]
[Re]
type = ParsedPostprocessor
function = '${rho} * ${l} * ${U}'
pp_names = ''
[]
[rho_outlet]
type = SideAverageValue
boundary = 'right'
variable = 'rho_mixture_var'
[]
[]
(modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/channel-drift-flux-transient.i)
mu = 1.0
rho = 10.0
mu_d = 0.1
rho_d = 1.0
l = 2
U = 1
dp = 0.01
inlet_phase_2 = 0.1
advected_interp_method = 'average'
velocity_interp_method = 'rc'
[GlobalParams]
rhie_chow_user_object = 'rc'
density_interp_method = 'average'
mu_interp_method = 'average'
[]
[UserObjects]
[rc]
type = INSFVRhieChowInterpolator
u = vel_x
v = vel_y
pressure = pressure
[]
[]
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = '${fparse l * 5}'
ymin = '${fparse -l / 2}'
ymax = '${fparse l / 2}'
nx = 10
ny = 4
[]
uniform_refine = 0
[]
[Variables]
[vel_x]
type = INSFVVelocityVariable
initial_condition = 0
[]
[vel_y]
type = INSFVVelocityVariable
initial_condition = 0
[]
[pressure]
type = INSFVPressureVariable
[]
[phase_2]
type = INSFVScalarFieldVariable
[]
[]
[FVKernels]
[mass]
type = INSFVMassAdvection
variable = pressure
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
[]
[u_time]
type = INSFVMomentumTimeDerivative
variable = vel_x
rho = 'rho_mixture'
momentum_component = 'x'
[]
[u_advection]
type = INSFVMomentumAdvection
variable = vel_x
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
momentum_component = 'x'
[]
[u_drift]
type = WCNSFV2PMomentumDriftFlux
variable = vel_x
rho_d = ${rho_d}
fd = 'phase_2'
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
momentum_component = 'x'
[]
[u_viscosity]
type = INSFVMomentumDiffusion
variable = vel_x
mu = 'mu_mixture'
limit_interpolation = true
momentum_component = 'x'
[]
[u_pressure]
type = INSFVMomentumPressure
variable = vel_x
momentum_component = 'x'
pressure = pressure
[]
[u_friction]
type = PINSFVMomentumFriction
Darcy_name = Darcy_coefficient_vec
is_porous_medium = false
momentum_component = x
mu = mu_mixture
rho = rho_mixture
variable = vel_x
[]
[v_time]
type = INSFVMomentumTimeDerivative
variable = vel_y
rho = 'rho_mixture'
momentum_component = 'y'
[]
[v_advection]
type = INSFVMomentumAdvection
variable = vel_y
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
momentum_component = 'y'
[]
[v_drift]
type = WCNSFV2PMomentumDriftFlux
variable = vel_y
rho_d = ${rho_d}
fd = 'phase_2'
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
momentum_component = 'y'
[]
[v_viscosity]
type = INSFVMomentumDiffusion
variable = vel_y
mu = 'mu_mixture'
limit_interpolation = true
momentum_component = 'y'
[]
[v_pressure]
type = INSFVMomentumPressure
variable = vel_y
momentum_component = 'y'
pressure = pressure
[]
[v_friction]
type = PINSFVMomentumFriction
Darcy_name = Darcy_coefficient_vec
is_porous_medium = false
momentum_component = y
mu = mu_mixture
rho = rho_mixture
variable = vel_y
[]
[phase_2_time]
type = FVFunctorTimeKernel
variable = phase_2
functor = phase_2
[]
[phase_2_advection]
type = INSFVScalarFieldAdvection
variable = phase_2
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = 'upwind'
[]
[phase_2_src]
type = NSFVMixturePhaseInterface
variable = phase_2
phase_coupled = phase_1
alpha = 0.1
[]
[]
[FVBCs]
[inlet-u]
type = INSFVInletVelocityBC
boundary = 'left'
variable = vel_x
functor = '${U}'
[]
[inlet-v]
type = INSFVInletVelocityBC
boundary = 'left'
variable = vel_y
functor = '0'
[]
[walls-u]
type = INSFVNoSlipWallBC
boundary = 'top bottom'
variable = vel_x
function = 0
[]
[walls-v]
type = INSFVNoSlipWallBC
boundary = 'top bottom'
variable = vel_y
function = 0
[]
[outlet_p]
type = INSFVOutletPressureBC
boundary = 'right'
variable = pressure
function = '0'
[]
[inlet_phase_2]
type = FVDirichletBC
boundary = 'left'
variable = phase_2
value = ${inlet_phase_2}
[]
[]
[AuxVariables]
[drag_coefficient]
type = MooseVariableFVReal
[]
[rho_mixture_var]
type = MooseVariableFVReal
[]
[mu_mixture_var]
type = MooseVariableFVReal
[]
[vel_slip_x_var]
type = MooseVariableFVReal
[]
[vel_slip_y_var]
type = MooseVariableFVReal
[]
[]
[AuxKernels]
[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'
[]
[populate_vx_slip_var]
type = FunctorAux
variable = vel_slip_x_var
functor = 'vel_slip_x'
[]
[populate_vy_slip_var]
type = FunctorAux
variable = vel_slip_y_var
functor = 'vel_slip_y'
[]
[]
[FunctorMaterials]
[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'
[]
[compute_phase_1]
type = ADParsedFunctorMaterial
property_name = phase_1
functor_names = 'phase_2'
expression = '1 - phase_2'
[]
[CD]
type = NSFVDispersePhaseDragFunctorMaterial
rho = 'rho_mixture'
mu = mu_mixture
u = 'vel_x'
v = 'vel_y'
particle_diameter = ${dp}
[]
[mixing_material]
type = NSFVMixtureFunctorMaterial
phase_2_names = '${rho} ${mu}'
phase_1_names = '${rho_d} ${mu_d}'
prop_names = 'rho_mixture mu_mixture'
phase_1_fraction = 'phase_2'
[]
[]
[Executioner]
type = Transient
solve_type = 'NEWTON'
nl_rel_tol = 1e-10
dt = 0.1
end_time = 1.0
[]
[Preconditioning]
[SMP]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
[]
[]
[Outputs]
exodus = false
[CSV]
type = CSV
execute_on = 'TIMESTEP_END'
[]
[]
[Postprocessors]
[Re]
type = ParsedPostprocessor
expression = '${rho} * ${l} * ${U}'
[]
[rho_outlet]
type = SideAverageValue
boundary = 'right'
variable = 'rho_mixture_var'
[]
[vslip_x]
type = SideExtremeValue
boundary = 'left'
variable = 'vel_slip_x_var'
[]
[vslip_y]
type = SideExtremeValue
boundary = 'left'
variable = 'vel_slip_y_var'
[]
[vslip_value]
type = ParsedPostprocessor
expression = 'sqrt(vslip_x*vslip_x + vslip_y*vslip_y)*vslip_x/abs(vslip_x)'
pp_names = 'vslip_x vslip_y'
[]
[]
(modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_interface_area_model/pressure_driven_growth_transient.i)
###############################################################################
# Validation test based on Hibiki and Ishii experiment [1] reported in Figure 3
# [1] Hibiki, T., & Ishii, M. (2000). One-group interfacial area transport of bubbly flows in vertical round tubes.
# International Journal of Heat and Mass Transfer, 43(15), 2711-2726.
###############################################################################
mu = 1.0
rho = 1000.0
mu_d = 1.0
rho_d = 1.0
l = ${fparse 50.8/1000.0}
U = 0.491230114
dp = 0.001
inlet_phase_2 = 0.049
advected_interp_method = 'upwind'
velocity_interp_method = 'rc'
mass_exchange_coeff = 0.0
inlet_interface_area = ${fparse 6.0*inlet_phase_2/dp}
outlet_pressure = 1e6
[GlobalParams]
rhie_chow_user_object = 'rc'
density_interp_method = 'average'
mu_interp_method = 'average'
[]
[Problem]
identify_variable_groups_in_nl = false
[]
[UserObjects]
[rc]
type = INSFVRhieChowInterpolator
u = vel_x
v = vel_y
pressure = pressure
[]
[]
[Mesh]
coord_type = 'RZ'
rz_coord_axis = 'X'
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = '${fparse l * 60}'
ymin = 0
ymax = '${fparse l / 2}'
nx = 20
ny = 5
[]
uniform_refine = 0
[]
[Variables]
[vel_x]
type = INSFVVelocityVariable
initial_condition = 0
[]
[vel_y]
type = INSFVVelocityVariable
initial_condition = 0
[]
[pressure]
type = INSFVPressureVariable
[]
[phase_2]
type = INSFVScalarFieldVariable
initial_condition = ${inlet_phase_2}
[]
[interface_area]
type = INSFVScalarFieldVariable
initial_condition = ${inlet_interface_area}
[]
[]
[FVKernels]
[mass]
type = INSFVMassAdvection
variable = pressure
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = ${rho}
[]
[u_time]
type = INSFVMomentumTimeDerivative
variable = vel_x
rho = 'rho_mixture'
momentum_component = 'x'
[]
[u_advection]
type = INSFVMomentumAdvection
variable = vel_x
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
momentum_component = 'x'
[]
[u_drift]
type = WCNSFV2PMomentumDriftFlux
variable = vel_x
rho_d = ${rho_d}
fd = 'rho_mixture_var'
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
momentum_component = 'x'
[]
[u_viscosity]
type = INSFVMomentumDiffusion
variable = vel_x
mu = 'mu_mixture'
limit_interpolation = true
momentum_component = 'x'
[]
[u_pressure]
type = INSFVMomentumPressure
variable = vel_x
momentum_component = 'x'
pressure = pressure
[]
[v_time]
type = INSFVMomentumTimeDerivative
variable = vel_y
rho = 'rho_mixture'
momentum_component = 'y'
[]
[v_advection]
type = INSFVMomentumAdvection
variable = vel_y
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
momentum_component = 'y'
[]
[v_drift]
type = WCNSFV2PMomentumDriftFlux
variable = vel_y
rho_d = ${rho_d}
fd = 'rho_mixture_var'
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
momentum_component = 'y'
[]
[v_viscosity]
type = INSFVMomentumDiffusion
variable = vel_y
mu = 'mu_mixture'
limit_interpolation = true
momentum_component = 'y'
[]
[v_pressure]
type = INSFVMomentumPressure
variable = vel_y
momentum_component = 'y'
pressure = pressure
[]
[phase_2_time]
type = FVFunctorTimeKernel
variable = phase_2
functor = phase_2
[]
[phase_2_advection]
type = INSFVScalarFieldAdvection
variable = phase_2
u_slip = 'vel_x'
v_slip = 'vel_y'
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = 'upwind'
[]
[phase_2_diffusion]
type = FVDiffusion
variable = phase_2
coeff = 1.0
[]
[phase_2_src]
type = NSFVMixturePhaseInterface
variable = phase_2
phase_coupled = phase_1
alpha = ${mass_exchange_coeff}
[]
[interface_area_time]
type = FVFunctorTimeKernel
variable = interface_area
functor = interface_area
[]
[interface_area_advection]
type = INSFVScalarFieldAdvection
variable = interface_area
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = 'upwind'
[]
[interface_area_diffusion]
type = FVDiffusion
variable = interface_area
coeff = 0.1
[]
[interface_area_source_sink]
type = WCNSFV2PInterfaceAreaSourceSink
variable = interface_area
u = 'vel_x'
v = 'vel_y'
L = ${fparse l/2}
rho = 'rho_mixture'
rho_d = 'rho'
pressure = 'pressure'
k_c = '${fparse mass_exchange_coeff}'
fd = 'phase_2'
sigma = 1e-3
cutoff_fraction = 0.0
[]
[]
[FVBCs]
[inlet-u]
type = INSFVInletVelocityBC
boundary = 'left'
variable = vel_x
functor = '${U}'
[]
[inlet-v]
type = INSFVInletVelocityBC
boundary = 'left'
variable = vel_y
functor = '0'
[]
[walls-u]
type = INSFVNoSlipWallBC
boundary = 'top'
variable = vel_x
function = 0
[]
[walls-v]
type = INSFVNoSlipWallBC
boundary = 'top'
variable = vel_y
function = 0
[]
[outlet_p]
type = INSFVOutletPressureBC
boundary = 'right'
variable = pressure
function = '${outlet_pressure}'
[]
[inlet_phase_2]
type = FVDirichletBC
boundary = 'left'
variable = phase_2
value = ${inlet_phase_2}
[]
[inlet_interface_area]
type = FVDirichletBC
boundary = 'left'
variable = interface_area
value = ${inlet_interface_area}
[]
[symmetry-u]
type = PINSFVSymmetryVelocityBC
boundary = 'bottom'
variable = vel_x
u = vel_x
v = vel_y
mu = 'mu_mixture'
momentum_component = 'x'
[]
[symmetry-v]
type = PINSFVSymmetryVelocityBC
boundary = 'bottom'
variable = vel_y
u = vel_x
v = vel_y
mu = 'mu_mixture'
momentum_component = 'y'
[]
[symmetry-p]
type = INSFVSymmetryPressureBC
boundary = 'bottom'
variable = pressure
[]
[symmetry-phase-2]
type = INSFVSymmetryScalarBC
boundary = 'bottom'
variable = phase_2
[]
[symmetry-interface-area]
type = INSFVSymmetryScalarBC
boundary = 'bottom'
variable = interface_area
[]
[]
[AuxVariables]
[drag_coefficient]
type = MooseVariableFVReal
[]
[rho_mixture_var]
type = MooseVariableFVReal
[]
[mu_mixture_var]
type = MooseVariableFVReal
[]
[]
[AuxKernels]
[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'
[]
[]
[FluidProperties]
[fp]
type = IdealGasFluidProperties
[]
[]
[FunctorMaterials]
[bubble_properties]
type = GeneralFunctorFluidProps
fp = 'fp'
pressure = 'pressure'
T_fluid = 300.0
speed = 1.0
characteristic_length = 1.0
porosity = 1.0
output_properties = 'rho'
outputs = 'out'
[]
[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'
[]
[compute_phase_1]
type = ADParsedFunctorMaterial
property_name = phase_1
functor_names = 'phase_2'
expression = '1 - phase_2'
[]
[CD]
type = NSFVDispersePhaseDragFunctorMaterial
rho = 'rho_mixture'
mu = mu_mixture
u = 'vel_x'
v = 'vel_y'
particle_diameter = ${dp}
[]
[mixing_material]
type = NSFVMixtureFunctorMaterial
phase_2_names = '${rho} ${mu}'
phase_1_names = 'rho ${mu_d}'
prop_names = 'rho_mixture mu_mixture'
phase_1_fraction = 'phase_2'
[]
[]
[Executioner]
type = Transient
solve_type = 'NEWTON'
nl_abs_tol = 1e-7
dt = 0.1
end_time = 1.0
nl_max_its = 10
line_search = 'none'
[]
[Debug]
show_var_residual_norms = true
[]
[Preconditioning]
[SMP]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
[]
[]
[Outputs]
[out]
type = Exodus
[]
[]
[Postprocessors]
[Re]
type = ParsedPostprocessor
expression = '${rho} * ${l} * ${U}'
pp_names = ''
[]
[rho_outlet]
type = SideAverageValue
boundary = 'right'
variable = 'rho_mixture_var'
[]
[]
(modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_interface_area_model/pressure_driven_growth.i)
###############################################################################
# Validation test based on Hibiki and Ishii experiment [1] reported in Figure 3
# [1] Hibiki, T., & Ishii, M. (2000). One-group interfacial area transport of bubbly flows in vertical round tubes.
# International Journal of Heat and Mass Transfer, 43(15), 2711-2726.
###############################################################################
mu = 1.0
rho = 1000.0
mu_d = 1.0
rho_d = 1.0
l = ${fparse 50.8/1000.0}
U = 0.491230114
dp = 0.001
inlet_phase_2 = 0.049
advected_interp_method = 'upwind'
velocity_interp_method = 'rc'
mass_exchange_coeff = 0.0
inlet_interface_area = ${fparse 6.0*inlet_phase_2/dp}
outlet_pressure = 1e5
[GlobalParams]
rhie_chow_user_object = 'rc'
density_interp_method = 'average'
mu_interp_method = 'average'
[]
[Problem]
identify_variable_groups_in_nl = false
previous_nl_solution_required = true
[]
[UserObjects]
[rc]
type = INSFVRhieChowInterpolator
u = vel_x
v = vel_y
pressure = pressure
[]
[]
[Mesh]
coord_type = 'RZ'
rz_coord_axis = 'X'
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = '${fparse l * 60}'
ymin = 0
ymax = '${fparse l / 2}'
nx = 20
ny = 5
[]
uniform_refine = 0
[]
[Variables]
[vel_x]
type = INSFVVelocityVariable
initial_condition = 0
[]
[vel_y]
type = INSFVVelocityVariable
initial_condition = 0
[]
[pressure]
type = INSFVPressureVariable
[]
[phase_2]
type = INSFVScalarFieldVariable
initial_condition = ${inlet_phase_2}
[]
[interface_area]
type = INSFVScalarFieldVariable
initial_condition = ${inlet_interface_area}
[]
[]
[FVKernels]
[mass]
type = INSFVMassAdvection
variable = pressure
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = ${rho}
[]
[u_advection]
type = INSFVMomentumAdvection
variable = vel_x
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
momentum_component = 'x'
[]
[u_drift]
type = WCNSFV2PMomentumDriftFlux
variable = vel_x
rho_d = ${rho_d}
fd = 'rho_mixture_var'
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
momentum_component = 'x'
[]
[u_viscosity]
type = INSFVMomentumDiffusion
variable = vel_x
mu = 'mu_mixture'
limit_interpolation = true
momentum_component = 'x'
[]
[u_pressure]
type = INSFVMomentumPressure
variable = vel_x
momentum_component = 'x'
pressure = pressure
[]
[v_advection]
type = INSFVMomentumAdvection
variable = vel_y
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
momentum_component = 'y'
[]
[v_drift]
type = WCNSFV2PMomentumDriftFlux
variable = vel_y
rho_d = ${rho_d}
fd = 'rho_mixture_var'
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
momentum_component = 'y'
[]
[v_viscosity]
type = INSFVMomentumDiffusion
variable = vel_y
mu = 'mu_mixture'
limit_interpolation = true
momentum_component = 'y'
[]
[v_pressure]
type = INSFVMomentumPressure
variable = vel_y
momentum_component = 'y'
pressure = pressure
[]
[phase_2_advection]
type = INSFVScalarFieldAdvection
variable = phase_2
u_slip = 'vel_x'
v_slip = 'vel_y'
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = 'upwind'
[]
[phase_2_diffusion]
type = FVDiffusion
variable = phase_2
coeff = 1.0
[]
[phase_2_src]
type = NSFVMixturePhaseInterface
variable = phase_2
phase_coupled = phase_1
alpha = ${mass_exchange_coeff}
[]
[interface_area_advection]
type = INSFVScalarFieldAdvection
variable = interface_area
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = 'upwind'
[]
[interface_area_diffusion]
type = FVDiffusion
variable = interface_area
coeff = 0.1
[]
[interface_area_source_sink]
type = WCNSFV2PInterfaceAreaSourceSink
variable = interface_area
u = 'vel_x'
v = 'vel_y'
L = ${fparse l/2}
rho = 'rho_mixture'
rho_d = 'rho'
pressure = 'pressure'
k_c = '${fparse mass_exchange_coeff}'
fd = 'phase_2'
sigma = 1e-3
cutoff_fraction = 0.0
[]
[]
[FVBCs]
[inlet-u]
type = INSFVInletVelocityBC
boundary = 'left'
variable = vel_x
functor = '${U}'
[]
[inlet-v]
type = INSFVInletVelocityBC
boundary = 'left'
variable = vel_y
functor = '0'
[]
[walls-u]
type = INSFVNoSlipWallBC
boundary = 'top'
variable = vel_x
function = 0
[]
[walls-v]
type = INSFVNoSlipWallBC
boundary = 'top'
variable = vel_y
function = 0
[]
[outlet_p]
type = INSFVOutletPressureBC
boundary = 'right'
variable = pressure
function = '${outlet_pressure}'
[]
[inlet_phase_2]
type = FVDirichletBC
boundary = 'left'
variable = phase_2
value = ${inlet_phase_2}
[]
[inlet_interface_area]
type = FVDirichletBC
boundary = 'left'
variable = interface_area
value = ${inlet_interface_area}
[]
[symmetry-u]
type = PINSFVSymmetryVelocityBC
boundary = 'bottom'
variable = vel_x
u = vel_x
v = vel_y
mu = 'mu_mixture'
momentum_component = 'x'
[]
[symmetry-v]
type = PINSFVSymmetryVelocityBC
boundary = 'bottom'
variable = vel_y
u = vel_x
v = vel_y
mu = 'mu_mixture'
momentum_component = 'y'
[]
[symmetry-p]
type = INSFVSymmetryPressureBC
boundary = 'bottom'
variable = pressure
[]
[symmetry-phase-2]
type = INSFVSymmetryScalarBC
boundary = 'bottom'
variable = phase_2
[]
[symmetry-interface-area]
type = INSFVSymmetryScalarBC
boundary = 'bottom'
variable = interface_area
[]
[]
[AuxVariables]
[drag_coefficient]
type = MooseVariableFVReal
[]
[rho_mixture_var]
type = MooseVariableFVReal
[]
[mu_mixture_var]
type = MooseVariableFVReal
[]
[]
[AuxKernels]
[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'
[]
[]
[FluidProperties]
[fp]
type = IdealGasFluidProperties
[]
[]
[FunctorMaterials]
[bubble_properties]
type = GeneralFunctorFluidProps
fp = 'fp'
pressure = 'pressure'
T_fluid = 300.0
speed = 1.0
characteristic_length = 1.0
porosity = 1.0
output_properties = 'rho'
outputs = 'out'
[]
[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'
[]
[compute_phase_1]
type = ADParsedFunctorMaterial
property_name = phase_1
functor_names = 'phase_2'
expression = '1 - phase_2'
[]
[CD]
type = NSFVDispersePhaseDragFunctorMaterial
rho = 'rho_mixture'
mu = mu_mixture
u = 'vel_x'
v = 'vel_y'
particle_diameter = ${dp}
[]
[mixing_material]
type = NSFVMixtureFunctorMaterial
phase_2_names = '${rho} ${mu}'
phase_1_names = 'rho ${mu_d}'
prop_names = 'rho_mixture mu_mixture'
phase_1_fraction = 'phase_2'
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
nl_rel_tol = 1e-10
line_search = 'none'
[]
[Debug]
show_var_residual_norms = true
[]
[Preconditioning]
[SMP]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
[]
[]
[Outputs]
[out]
type = Exodus
[]
[]
[Postprocessors]
[Re]
type = ParsedPostprocessor
expression = '${rho} * ${l} * ${U}'
pp_names = ''
[]
[rho_outlet]
type = SideAverageValue
boundary = 'right'
variable = 'rho_mixture_var'
[]
[]
(modules/navier_stokes/test/tests/finite_volume/two_phase/mixture_model/channel-drift-flux.i)
mu = 1.0
rho = 10.0
mu_d = 0.1
rho_d = 1.0
l = 2
U = 1
dp = 0.01
inlet_phase_2 = 0.1
advected_interp_method = 'average'
velocity_interp_method = 'rc'
[GlobalParams]
rhie_chow_user_object = 'rc'
density_interp_method = 'average'
mu_interp_method = 'average'
[]
[UserObjects]
[rc]
type = INSFVRhieChowInterpolator
u = vel_x
v = vel_y
pressure = pressure
[]
[]
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = '${fparse l * 5}'
ymin = '${fparse -l / 2}'
ymax = '${fparse l / 2}'
nx = 10
ny = 4
[]
uniform_refine = 0
[]
[Variables]
[vel_x]
type = INSFVVelocityVariable
initial_condition = 0
[]
[vel_y]
type = INSFVVelocityVariable
initial_condition = 0
[]
[pressure]
type = INSFVPressureVariable
[]
[phase_2]
type = INSFVScalarFieldVariable
[]
[]
[FVKernels]
[mass]
type = INSFVMassAdvection
variable = pressure
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
[]
[u_advection]
type = INSFVMomentumAdvection
variable = vel_x
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
momentum_component = 'x'
[]
[u_drift]
type = WCNSFV2PMomentumDriftFlux
variable = vel_x
rho_d = ${rho_d}
fd = 'phase_2'
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
momentum_component = 'x'
[]
[u_viscosity]
type = INSFVMomentumDiffusion
variable = vel_x
mu = 'mu_mixture'
limit_interpolation = true
momentum_component = 'x'
[]
[u_pressure]
type = INSFVMomentumPressure
variable = vel_x
momentum_component = 'x'
pressure = pressure
[]
[u_friction]
type = PINSFVMomentumFriction
Darcy_name = Darcy_coefficient_vec
is_porous_medium = false
momentum_component = x
mu = mu_mixture
rho = rho_mixture
variable = vel_x
[]
[v_advection]
type = INSFVMomentumAdvection
variable = vel_y
advected_interp_method = ${advected_interp_method}
velocity_interp_method = ${velocity_interp_method}
rho = 'rho_mixture'
momentum_component = 'y'
[]
[v_drift]
type = WCNSFV2PMomentumDriftFlux
variable = vel_y
rho_d = ${rho_d}
fd = 'phase_2'
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
momentum_component = 'x'
[]
[v_viscosity]
type = INSFVMomentumDiffusion
variable = vel_y
mu = 'mu_mixture'
limit_interpolation = true
momentum_component = 'y'
[]
[v_pressure]
type = INSFVMomentumPressure
variable = vel_y
momentum_component = 'y'
pressure = pressure
[]
[v_friction]
type = PINSFVMomentumFriction
Darcy_name = Darcy_coefficient_vec
is_porous_medium = false
momentum_component = y
mu = mu_mixture
rho = rho_mixture
variable = vel_y
[]
[phase_2_advection]
type = INSFVScalarFieldAdvection
variable = phase_2
u_slip = 'vel_slip_x'
v_slip = 'vel_slip_y'
velocity_interp_method = ${velocity_interp_method}
advected_interp_method = 'upwind'
[]
[phase_2_src]
type = NSFVMixturePhaseInterface
variable = phase_2
phase_coupled = phase_1
alpha = 0.1
[]
[]
[FVBCs]
[inlet-u]
type = INSFVInletVelocityBC
boundary = 'left'
variable = vel_x
functor = '${U}'
[]
[inlet-v]
type = INSFVInletVelocityBC
boundary = 'left'
variable = vel_y
functor = '0'
[]
[walls-u]
type = INSFVNoSlipWallBC
boundary = 'top bottom'
variable = vel_x
function = 0
[]
[walls-v]
type = INSFVNoSlipWallBC
boundary = 'top bottom'
variable = vel_y
function = 0
[]
[outlet_p]
type = INSFVOutletPressureBC
boundary = 'right'
variable = pressure
function = '0'
[]
[inlet_phase_2]
type = FVDirichletBC
boundary = 'left'
variable = phase_2
value = ${inlet_phase_2}
[]
[]
[AuxVariables]
[drag_coefficient]
type = MooseVariableFVReal
[]
[rho_mixture_var]
type = MooseVariableFVReal
[]
[mu_mixture_var]
type = MooseVariableFVReal
[]
[]
[AuxKernels]
[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'
[]
[]
[FunctorMaterials]
[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'
outputs = 'out'
output_properties = 'vel_slip_x'
ghost_layers = 5
[]
[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'
outputs = 'out'
output_properties = 'vel_slip_y'
ghost_layers = 5
[]
[compute_phase_1]
type = ADParsedFunctorMaterial
property_name = phase_1
functor_names = 'phase_2'
expression = '1 - phase_2'
outputs = 'out'
output_properties = 'phase_1'
[]
[CD]
type = NSFVDispersePhaseDragFunctorMaterial
rho = 'rho_mixture'
mu = mu_mixture
u = 'vel_x'
v = 'vel_y'
particle_diameter = ${dp}
[]
[mixing_material]
type = NSFVMixtureFunctorMaterial
phase_2_names = '${rho} ${mu}'
phase_1_names = '${rho_d} ${mu_d}'
prop_names = 'rho_mixture mu_mixture'
phase_1_fraction = 'phase_2'
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
nl_rel_tol = 1e-10
[]
[Preconditioning]
[SMP]
type = SMP
full = true
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
[]
[]
[Outputs]
print_linear_residuals = true
print_nonlinear_residuals = true
[out]
type = Exodus
hide = 'Re lin cum_lin'
[]
[perf]
type = PerfGraphOutput
[]
[]
[Postprocessors]
[Re]
type = ParsedPostprocessor
expression = '${rho} * ${l} * ${U}'
[]
[lin]
type = NumLinearIterations
[]
[cum_lin]
type = CumulativeValuePostprocessor
postprocessor = lin
[]
[]
(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
[]
[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'
[]
[]