Step 2
Step 2 modifies the input for Step 1 by adding a pressure drop in the bed. The pressure drop in the bed physically originates from the interaction of the fluid with the coexisting solid, which is essentially a transfer of momentum from the fluid to the solid. In this tutorial, we use the German Kerntechnischer Ausschuss (KTA) correlations to compute the pressure drop. The KTA pressure drop for a one-dimensional flow is given by:
where:
where is the porosity, is the hydraulic diameter, is the density, is the mass flow rate, is the area of the cylindrical flow channel, and is the Reynolds number. The Reynolds number is given by:
where is the dynamic viscosity, is the pebble diameter of m, and m. We estimate averages for and using postprocessors. These are included in the step2.i
file but are not discussed here. We find kg/m and Pa-s. With these numerical values, we obtain:
As the length of the bed is 10 m, we expect a pressure drop of kPa.
Adding pressure drop to the channel flow
We will use the KTA drag coefficient material in Pronghorn. To use this material, the characteristic length (aka hydraulic diameter) of the bed has to be set. For pebble beds, the characteristic length is the pebble diameter, which we assume to be m consistent with typical gas-cooled reactor pebbles. We create a real valued variable called pebble_diameter
for convenience at the top of the input file:
T_fluid = 300
density = 8.6545
pebble_diameter = 0.06
mass_flow_rate = 60.0
flow_area = '${fparse pi * bed_radius * bed_radius}'
(htgr/generic-pbr-tutorial/step2.i)Then we need to add the friction coefficients to the finite volume Navier-Stokes action. This is done by adding the friction_types
and friction_coeffs
parameters.
[Modules]
[NavierStokesFV]
# general control parameters
compressibility = 'weakly-compressible'
porous_medium_treatment = true
# material property parameters
density = rho
dynamic_viscosity = mu
# porous medium treatment parameters
porosity = porosity
# initial conditions
initial_velocity = '1e-6 1e-6 0'
initial_pressure = 5.4e6
# boundary conditions
inlet_boundaries = top
momentum_inlet_types = fixed-velocity
momentum_inlet_function = '0 -${flow_vel}'
wall_boundaries = 'left right'
momentum_wall_types = 'slip slip'
outlet_boundaries = bottom
momentum_outlet_types = fixed-pressure
pressure_function = ${outlet_pressure}
# friction control parameters
friction_types = 'darcy forchheimer'
friction_coeffs = 'Darcy_coefficient Forchheimer_coefficient'
[]
[]
(htgr/generic-pbr-tutorial/step2.i)Drag friction sources are often split up into a component that depends linearly on fluid speed, called the Darcy friction (Darcy coefficient does not depend on velocity) and one that depends quadratically on flow speed, called the Forchheimer friction (Forchheimer depends linearly on velocity). The same distinction is made in Pronghorn. The two parameters simply express that there is a Darcy coefficient called Darcy_coefficient
and a Forchheimer coefficient called Forchheimer_coefficient
. The Darcy and Forchheimer coefficients are internally computed as functor material properties and almost always have the names given above.
Finally, we need to add the material that computes the Darcy_coefficient
and Forchheimer_coefficient
. This is done here:
[FunctorMaterials]
[drag_pebble_bed]
type = FunctorKTADragCoefficients
fp = fluid_properties_obj
pebble_diameter = ${pebble_diameter}
porosity = porosity
T_fluid = ${T_fluid}
T_solid = ${T_fluid}
[]
[]
(htgr/generic-pbr-tutorial/step2.i)Note how we use the ${pebble_diameter} variable to set the parameter of the same name. The fluid and solid temperatures are simply set to constants here because we have not added energy equations for the fluid and solid phases. The fluid and solid temperatures are used to evaluate material properties.
Measuring Pressure Drop
We measure the pressure drop by averaging the pressure over the inlet and outlet and then taking the difference. This is performed here:
[inlet_pressure]
type = SideAverageValue
variable = pressure
boundary = top
outputs = none
[]
[outlet_pressure]
type = SideAverageValue
variable = pressure
boundary = bottom
outputs = none
[]
[pressure_drop]
type = ParsedPostprocessor
pp_names = 'inlet_pressure outlet_pressure'
expression = 'inlet_pressure - outlet_pressure'
[]
(htgr/generic-pbr-tutorial/step2.i)Note that inlet_pressure
and outlet_pressure
are neither printed to screen nor written to any possible csv file because the parameter outputs
is set to none.
Executing and comparing the pressure drop
We execute by:
./pronghorn-opt -i step2.i
and find that the pressure is , which is exactly what we estimated it should be.
(htgr/generic-pbr-tutorial/step2.i)
# ==============================================================================
# Model description:
# Step2 - Step1 plus pressure drop.
# ------------------------------------------------------------------------------
# Idaho Falls, INL, August 15, 2023 04:03 PM
# Author(s): Joseph R. Brennan, Dr. Sebastian Schunert, Dr. Mustafa K. Jaradat
# and Dr. Paolo Balestra.
# ==============================================================================
bed_height = 10.0
bed_radius = 1.2
bed_porosity = 0.39
outlet_pressure = 5.5e6
T_fluid = 300
density = 8.6545
pebble_diameter = 0.06
mass_flow_rate = 60.0
flow_area = '${fparse pi * bed_radius * bed_radius}'
flow_vel = '${fparse mass_flow_rate / flow_area / density}'
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = ${bed_radius}
ymin = 0
ymax = ${bed_height}
nx = 6
ny = 40
[]
coord_type = RZ
[]
[FluidProperties]
[fluid_properties_obj]
type = HeliumFluidProperties
[]
[]
[Modules]
[NavierStokesFV]
# general control parameters
compressibility = 'weakly-compressible'
porous_medium_treatment = true
# material property parameters
density = rho
dynamic_viscosity = mu
# porous medium treatment parameters
porosity = porosity
# initial conditions
initial_velocity = '1e-6 1e-6 0'
initial_pressure = 5.4e6
# boundary conditions
inlet_boundaries = top
momentum_inlet_types = fixed-velocity
momentum_inlet_function = '0 -${flow_vel}'
wall_boundaries = 'left right'
momentum_wall_types = 'slip slip'
outlet_boundaries = bottom
momentum_outlet_types = fixed-pressure
pressure_function = ${outlet_pressure}
# friction control parameters
friction_types = 'darcy forchheimer'
friction_coeffs = 'Darcy_coefficient Forchheimer_coefficient'
[]
[]
[FunctorMaterials]
[fluid_props_to_mat_props]
type = GeneralFunctorFluidProps
fp = fluid_properties_obj
porosity = porosity
pressure = pressure
T_fluid = ${T_fluid}
speed = speed
characteristic_length = ${pebble_diameter}
neglect_derivatives_of_density_time_derivative = false
[]
[drag_pebble_bed]
type = FunctorKTADragCoefficients
fp = fluid_properties_obj
pebble_diameter = ${pebble_diameter}
porosity = porosity
T_fluid = ${T_fluid}
T_solid = ${T_fluid}
[]
[]
[AuxVariables]
[porosity]
family = MONOMIAL
order = CONSTANT
fv = true
initial_condition = ${bed_porosity}
[]
[]
[Executioner]
type = Transient
end_time = 100
[TimeStepper]
type = IterationAdaptiveDT
iteration_window = 2
optimal_iterations = 8
cutback_factor = 0.8
growth_factor = 2
dt = 1e-3
[]
dtmax = 5
line_search = l2
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_factor_shift_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu NONZERO superlu_dist'
nl_rel_tol = 1e-6
nl_abs_tol = 1e-6
[]
[Postprocessors]
[inlet_mfr]
type = VolumetricFlowRate
advected_quantity = rho
vel_x = 'superficial_vel_x'
vel_y = 'superficial_vel_y'
boundary = 'top'
rhie_chow_user_object = pins_rhie_chow_interpolator
[]
[outlet_mfr]
type = VolumetricFlowRate
advected_quantity = rho
vel_x = 'superficial_vel_x'
vel_y = 'superficial_vel_y'
boundary = 'bottom'
rhie_chow_user_object = pins_rhie_chow_interpolator
[]
[inlet_pressure]
type = SideAverageValue
variable = pressure
boundary = top
outputs = none
[]
[outlet_pressure]
type = SideAverageValue
variable = pressure
boundary = bottom
outputs = none
[]
[pressure_drop]
type = ParsedPostprocessor
pp_names = 'inlet_pressure outlet_pressure'
expression = 'inlet_pressure - outlet_pressure'
[]
[integral_density]
type = ADElementIntegralFunctorPostprocessor
functor = rho
outputs = none
[]
[average_density]
type = ParsedPostprocessor
pp_names = 'volume integral_density'
expression = 'integral_density / volume'
[]
[integral_mu]
type = ADElementIntegralFunctorPostprocessor
functor = mu
outputs = none
[]
[average_mu]
type = ParsedPostprocessor
pp_names = 'volume integral_mu'
expression = 'integral_mu / volume'
[]
[area]
type = AreaPostprocessor
boundary = bottom
outputs = none
[]
[volume]
type = VolumePostprocessor
[]
[]
[Outputs]
exodus = true
[]
(htgr/generic-pbr-tutorial/step2.i)
# ==============================================================================
# Model description:
# Step2 - Step1 plus pressure drop.
# ------------------------------------------------------------------------------
# Idaho Falls, INL, August 15, 2023 04:03 PM
# Author(s): Joseph R. Brennan, Dr. Sebastian Schunert, Dr. Mustafa K. Jaradat
# and Dr. Paolo Balestra.
# ==============================================================================
bed_height = 10.0
bed_radius = 1.2
bed_porosity = 0.39
outlet_pressure = 5.5e6
T_fluid = 300
density = 8.6545
pebble_diameter = 0.06
mass_flow_rate = 60.0
flow_area = '${fparse pi * bed_radius * bed_radius}'
flow_vel = '${fparse mass_flow_rate / flow_area / density}'
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = ${bed_radius}
ymin = 0
ymax = ${bed_height}
nx = 6
ny = 40
[]
coord_type = RZ
[]
[FluidProperties]
[fluid_properties_obj]
type = HeliumFluidProperties
[]
[]
[Modules]
[NavierStokesFV]
# general control parameters
compressibility = 'weakly-compressible'
porous_medium_treatment = true
# material property parameters
density = rho
dynamic_viscosity = mu
# porous medium treatment parameters
porosity = porosity
# initial conditions
initial_velocity = '1e-6 1e-6 0'
initial_pressure = 5.4e6
# boundary conditions
inlet_boundaries = top
momentum_inlet_types = fixed-velocity
momentum_inlet_function = '0 -${flow_vel}'
wall_boundaries = 'left right'
momentum_wall_types = 'slip slip'
outlet_boundaries = bottom
momentum_outlet_types = fixed-pressure
pressure_function = ${outlet_pressure}
# friction control parameters
friction_types = 'darcy forchheimer'
friction_coeffs = 'Darcy_coefficient Forchheimer_coefficient'
[]
[]
[FunctorMaterials]
[fluid_props_to_mat_props]
type = GeneralFunctorFluidProps
fp = fluid_properties_obj
porosity = porosity
pressure = pressure
T_fluid = ${T_fluid}
speed = speed
characteristic_length = ${pebble_diameter}
neglect_derivatives_of_density_time_derivative = false
[]
[drag_pebble_bed]
type = FunctorKTADragCoefficients
fp = fluid_properties_obj
pebble_diameter = ${pebble_diameter}
porosity = porosity
T_fluid = ${T_fluid}
T_solid = ${T_fluid}
[]
[]
[AuxVariables]
[porosity]
family = MONOMIAL
order = CONSTANT
fv = true
initial_condition = ${bed_porosity}
[]
[]
[Executioner]
type = Transient
end_time = 100
[TimeStepper]
type = IterationAdaptiveDT
iteration_window = 2
optimal_iterations = 8
cutback_factor = 0.8
growth_factor = 2
dt = 1e-3
[]
dtmax = 5
line_search = l2
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_factor_shift_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu NONZERO superlu_dist'
nl_rel_tol = 1e-6
nl_abs_tol = 1e-6
[]
[Postprocessors]
[inlet_mfr]
type = VolumetricFlowRate
advected_quantity = rho
vel_x = 'superficial_vel_x'
vel_y = 'superficial_vel_y'
boundary = 'top'
rhie_chow_user_object = pins_rhie_chow_interpolator
[]
[outlet_mfr]
type = VolumetricFlowRate
advected_quantity = rho
vel_x = 'superficial_vel_x'
vel_y = 'superficial_vel_y'
boundary = 'bottom'
rhie_chow_user_object = pins_rhie_chow_interpolator
[]
[inlet_pressure]
type = SideAverageValue
variable = pressure
boundary = top
outputs = none
[]
[outlet_pressure]
type = SideAverageValue
variable = pressure
boundary = bottom
outputs = none
[]
[pressure_drop]
type = ParsedPostprocessor
pp_names = 'inlet_pressure outlet_pressure'
expression = 'inlet_pressure - outlet_pressure'
[]
[integral_density]
type = ADElementIntegralFunctorPostprocessor
functor = rho
outputs = none
[]
[average_density]
type = ParsedPostprocessor
pp_names = 'volume integral_density'
expression = 'integral_density / volume'
[]
[integral_mu]
type = ADElementIntegralFunctorPostprocessor
functor = mu
outputs = none
[]
[average_mu]
type = ParsedPostprocessor
pp_names = 'volume integral_mu'
expression = 'integral_mu / volume'
[]
[area]
type = AreaPostprocessor
boundary = bottom
outputs = none
[]
[volume]
type = VolumePostprocessor
[]
[]
[Outputs]
exodus = true
[]
(htgr/generic-pbr-tutorial/step2.i)
# ==============================================================================
# Model description:
# Step2 - Step1 plus pressure drop.
# ------------------------------------------------------------------------------
# Idaho Falls, INL, August 15, 2023 04:03 PM
# Author(s): Joseph R. Brennan, Dr. Sebastian Schunert, Dr. Mustafa K. Jaradat
# and Dr. Paolo Balestra.
# ==============================================================================
bed_height = 10.0
bed_radius = 1.2
bed_porosity = 0.39
outlet_pressure = 5.5e6
T_fluid = 300
density = 8.6545
pebble_diameter = 0.06
mass_flow_rate = 60.0
flow_area = '${fparse pi * bed_radius * bed_radius}'
flow_vel = '${fparse mass_flow_rate / flow_area / density}'
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = ${bed_radius}
ymin = 0
ymax = ${bed_height}
nx = 6
ny = 40
[]
coord_type = RZ
[]
[FluidProperties]
[fluid_properties_obj]
type = HeliumFluidProperties
[]
[]
[Modules]
[NavierStokesFV]
# general control parameters
compressibility = 'weakly-compressible'
porous_medium_treatment = true
# material property parameters
density = rho
dynamic_viscosity = mu
# porous medium treatment parameters
porosity = porosity
# initial conditions
initial_velocity = '1e-6 1e-6 0'
initial_pressure = 5.4e6
# boundary conditions
inlet_boundaries = top
momentum_inlet_types = fixed-velocity
momentum_inlet_function = '0 -${flow_vel}'
wall_boundaries = 'left right'
momentum_wall_types = 'slip slip'
outlet_boundaries = bottom
momentum_outlet_types = fixed-pressure
pressure_function = ${outlet_pressure}
# friction control parameters
friction_types = 'darcy forchheimer'
friction_coeffs = 'Darcy_coefficient Forchheimer_coefficient'
[]
[]
[FunctorMaterials]
[fluid_props_to_mat_props]
type = GeneralFunctorFluidProps
fp = fluid_properties_obj
porosity = porosity
pressure = pressure
T_fluid = ${T_fluid}
speed = speed
characteristic_length = ${pebble_diameter}
neglect_derivatives_of_density_time_derivative = false
[]
[drag_pebble_bed]
type = FunctorKTADragCoefficients
fp = fluid_properties_obj
pebble_diameter = ${pebble_diameter}
porosity = porosity
T_fluid = ${T_fluid}
T_solid = ${T_fluid}
[]
[]
[AuxVariables]
[porosity]
family = MONOMIAL
order = CONSTANT
fv = true
initial_condition = ${bed_porosity}
[]
[]
[Executioner]
type = Transient
end_time = 100
[TimeStepper]
type = IterationAdaptiveDT
iteration_window = 2
optimal_iterations = 8
cutback_factor = 0.8
growth_factor = 2
dt = 1e-3
[]
dtmax = 5
line_search = l2
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_factor_shift_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu NONZERO superlu_dist'
nl_rel_tol = 1e-6
nl_abs_tol = 1e-6
[]
[Postprocessors]
[inlet_mfr]
type = VolumetricFlowRate
advected_quantity = rho
vel_x = 'superficial_vel_x'
vel_y = 'superficial_vel_y'
boundary = 'top'
rhie_chow_user_object = pins_rhie_chow_interpolator
[]
[outlet_mfr]
type = VolumetricFlowRate
advected_quantity = rho
vel_x = 'superficial_vel_x'
vel_y = 'superficial_vel_y'
boundary = 'bottom'
rhie_chow_user_object = pins_rhie_chow_interpolator
[]
[inlet_pressure]
type = SideAverageValue
variable = pressure
boundary = top
outputs = none
[]
[outlet_pressure]
type = SideAverageValue
variable = pressure
boundary = bottom
outputs = none
[]
[pressure_drop]
type = ParsedPostprocessor
pp_names = 'inlet_pressure outlet_pressure'
expression = 'inlet_pressure - outlet_pressure'
[]
[integral_density]
type = ADElementIntegralFunctorPostprocessor
functor = rho
outputs = none
[]
[average_density]
type = ParsedPostprocessor
pp_names = 'volume integral_density'
expression = 'integral_density / volume'
[]
[integral_mu]
type = ADElementIntegralFunctorPostprocessor
functor = mu
outputs = none
[]
[average_mu]
type = ParsedPostprocessor
pp_names = 'volume integral_mu'
expression = 'integral_mu / volume'
[]
[area]
type = AreaPostprocessor
boundary = bottom
outputs = none
[]
[volume]
type = VolumePostprocessor
[]
[]
[Outputs]
exodus = true
[]
(htgr/generic-pbr-tutorial/step2.i)
# ==============================================================================
# Model description:
# Step2 - Step1 plus pressure drop.
# ------------------------------------------------------------------------------
# Idaho Falls, INL, August 15, 2023 04:03 PM
# Author(s): Joseph R. Brennan, Dr. Sebastian Schunert, Dr. Mustafa K. Jaradat
# and Dr. Paolo Balestra.
# ==============================================================================
bed_height = 10.0
bed_radius = 1.2
bed_porosity = 0.39
outlet_pressure = 5.5e6
T_fluid = 300
density = 8.6545
pebble_diameter = 0.06
mass_flow_rate = 60.0
flow_area = '${fparse pi * bed_radius * bed_radius}'
flow_vel = '${fparse mass_flow_rate / flow_area / density}'
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = ${bed_radius}
ymin = 0
ymax = ${bed_height}
nx = 6
ny = 40
[]
coord_type = RZ
[]
[FluidProperties]
[fluid_properties_obj]
type = HeliumFluidProperties
[]
[]
[Modules]
[NavierStokesFV]
# general control parameters
compressibility = 'weakly-compressible'
porous_medium_treatment = true
# material property parameters
density = rho
dynamic_viscosity = mu
# porous medium treatment parameters
porosity = porosity
# initial conditions
initial_velocity = '1e-6 1e-6 0'
initial_pressure = 5.4e6
# boundary conditions
inlet_boundaries = top
momentum_inlet_types = fixed-velocity
momentum_inlet_function = '0 -${flow_vel}'
wall_boundaries = 'left right'
momentum_wall_types = 'slip slip'
outlet_boundaries = bottom
momentum_outlet_types = fixed-pressure
pressure_function = ${outlet_pressure}
# friction control parameters
friction_types = 'darcy forchheimer'
friction_coeffs = 'Darcy_coefficient Forchheimer_coefficient'
[]
[]
[FunctorMaterials]
[fluid_props_to_mat_props]
type = GeneralFunctorFluidProps
fp = fluid_properties_obj
porosity = porosity
pressure = pressure
T_fluid = ${T_fluid}
speed = speed
characteristic_length = ${pebble_diameter}
neglect_derivatives_of_density_time_derivative = false
[]
[drag_pebble_bed]
type = FunctorKTADragCoefficients
fp = fluid_properties_obj
pebble_diameter = ${pebble_diameter}
porosity = porosity
T_fluid = ${T_fluid}
T_solid = ${T_fluid}
[]
[]
[AuxVariables]
[porosity]
family = MONOMIAL
order = CONSTANT
fv = true
initial_condition = ${bed_porosity}
[]
[]
[Executioner]
type = Transient
end_time = 100
[TimeStepper]
type = IterationAdaptiveDT
iteration_window = 2
optimal_iterations = 8
cutback_factor = 0.8
growth_factor = 2
dt = 1e-3
[]
dtmax = 5
line_search = l2
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_factor_shift_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu NONZERO superlu_dist'
nl_rel_tol = 1e-6
nl_abs_tol = 1e-6
[]
[Postprocessors]
[inlet_mfr]
type = VolumetricFlowRate
advected_quantity = rho
vel_x = 'superficial_vel_x'
vel_y = 'superficial_vel_y'
boundary = 'top'
rhie_chow_user_object = pins_rhie_chow_interpolator
[]
[outlet_mfr]
type = VolumetricFlowRate
advected_quantity = rho
vel_x = 'superficial_vel_x'
vel_y = 'superficial_vel_y'
boundary = 'bottom'
rhie_chow_user_object = pins_rhie_chow_interpolator
[]
[inlet_pressure]
type = SideAverageValue
variable = pressure
boundary = top
outputs = none
[]
[outlet_pressure]
type = SideAverageValue
variable = pressure
boundary = bottom
outputs = none
[]
[pressure_drop]
type = ParsedPostprocessor
pp_names = 'inlet_pressure outlet_pressure'
expression = 'inlet_pressure - outlet_pressure'
[]
[integral_density]
type = ADElementIntegralFunctorPostprocessor
functor = rho
outputs = none
[]
[average_density]
type = ParsedPostprocessor
pp_names = 'volume integral_density'
expression = 'integral_density / volume'
[]
[integral_mu]
type = ADElementIntegralFunctorPostprocessor
functor = mu
outputs = none
[]
[average_mu]
type = ParsedPostprocessor
pp_names = 'volume integral_mu'
expression = 'integral_mu / volume'
[]
[area]
type = AreaPostprocessor
boundary = bottom
outputs = none
[]
[volume]
type = VolumePostprocessor
[]
[]
[Outputs]
exodus = true
[]