- variableThe name of the variable that this residual object operates on
C++ Type:NonlinearVariableName
Unit:(no unit assumed)
Controllable:No
Description:The name of the variable that this residual object operates on
Decay
Computes changing concentration due to radioactive decay.
Description
This kernel is intended for the following PDE: where is the radioactive decay constant and is the concentration. This kernel computes the right hand side portion of the PDE.
Example Input Syntax
[Kernels<<<{"href": "../../syntax/Kernels/index.html"}>>>]
[decay]
type = Decay<<<{"description": "Computes changing concentration due to radioactive decay.", "href": "Decay.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = c
radioactive_decay_constant<<<{"description": "Radioactive decay constant"}>>> = 1.0
[]
[]Input Parameters
- blockThe list of blocks (ids or names) that this object will be applied
C++ Type:std::vector<SubdomainName>
Controllable:No
Description:The list of blocks (ids or names) that this object will be applied
- displacementsThe displacements
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:The displacements
- radioactive_decay_constantRadioactive decay constant
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Radioactive decay constant
Optional Parameters
- absolute_value_vector_tagsThe tags for the vectors this residual object should fill with the absolute value of the residual contribution
C++ Type:std::vector<TagName>
Controllable:No
Description:The tags for the vectors this residual object should fill with the absolute value of the residual contribution
- extra_matrix_tagsThe extra tags for the matrices this Kernel should fill
C++ Type:std::vector<TagName>
Controllable:No
Description:The extra tags for the matrices this Kernel should fill
- extra_vector_tagsThe extra tags for the vectors this Kernel should fill
C++ Type:std::vector<TagName>
Controllable:No
Description:The extra tags for the vectors this Kernel should fill
- matrix_tagssystemThe tag for the matrices this Kernel should fill
Default:system
C++ Type:MultiMooseEnum
Options:nontime, system
Controllable:No
Description:The tag for the matrices this Kernel should fill
- vector_tagsnontimeThe tag for the vectors this Kernel should fill
Default:nontime
C++ Type:MultiMooseEnum
Options:nontime, time
Controllable:No
Description:The tag for the vectors this Kernel should fill
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.
- diag_save_inThe name of auxiliary variables to save this Kernel's diagonal Jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
C++ Type:std::vector<AuxVariableName>
Unit:(no unit assumed)
Controllable:No
Description:The name of auxiliary variables to save this Kernel's diagonal Jacobian contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
- enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:Yes
Description:Set the enabled status of the MooseObject.
- implicitTrueDetermines whether this object is calculated using an implicit or explicit form
Default:True
C++ Type:bool
Controllable:No
Description:Determines whether this object is calculated using an implicit or explicit form
- save_inThe name of auxiliary variables to save this Kernel's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
C++ Type:std::vector<AuxVariableName>
Unit:(no unit assumed)
Controllable:No
Description:The name of auxiliary variables to save this Kernel's residual contributions to. Everything about that variable must match everything about this variable (the type, what blocks it's on, etc.)
- seed0The seed for the master random number generator
Default:0
C++ Type:unsigned int
Controllable:No
Description:The seed for the master random number generator
- use_displaced_meshFalseWhether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Default:False
C++ Type:bool
Controllable:No
Description:Whether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Advanced Parameters
- prop_getter_suffixAn optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:An optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
- use_interpolated_stateFalseFor the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.
Default:False
C++ Type:bool
Controllable:No
Description:For the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.
Material Property Retrieval Parameters
Input Files
- (assessment/verification/TRISO_diffusion/inpile_stable/inpile_stable.i)
- (examples/TRISO/parfume/parfume.i)
- (assessment/verification/TRISO_diffusion/outofpile/outofpile.i)
- (assessment/verification/TRISO_diffusion/inpile/inpile.i)
- (examples/metal_fuel/mechanistic_fcci/fcci.i)
- (examples/TRISO/accident_simulation/triso1D_accident.i)
- (examples/TRISO/parfume/parfume_un.i)
- (examples/TRISO/accident_simulation/triso2D_accident.i)
- (examples/TRISO/accident_simulation/triso2D_accident_mortar.i)
- (examples/TRISO/accident_simulation/triso2D_accident_ad.i)
- (test/tests/radioactive_decay/radioactive_decay.i)
- (test/tests/triso/base_irradiation/triso1D_accident.i)
(test/tests/radioactive_decay/radioactive_decay.i)
# Tests radioactive decay kernel Decay.C against an analytical solution. The ODE is
#
# dc/dt = - lambda*c
#
# which has the solution
#
# c = c0 * exp(-lambda*t)
#
# Numerical solution is for a single linear element on a unit cube having an initial
# concentration of 1. No BC's are provided, which defaults to zero mass flux on all
# faces. To minimize time error, the solution is run for 1000 steps out to a time of 1.
#
# Analytical result: c(1) = 0.367879
# Numerical result: c(1) = 0.368063
#
[Mesh]
[mesh]
type = FileMeshGenerator
file = 1x1x1_cube.e
[]
[]
[Variables]
[c]
order = FIRST
family = LAGRANGE
initial_condition = 1.0
[]
[]
[Kernels]
[decay]
type = Decay
variable = c
radioactive_decay_constant = 1.0
[]
[ie]
type = TimeDerivative
variable = c
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
start_time = 0.0
num_steps = 1000
dt = .001
[]
[Outputs]
exodus = true
[]
(assessment/verification/TRISO_diffusion/inpile_stable/inpile_stable.i)
#
# TRISO kernel in 1D spherical
#
# In-pile, a stable fission product
# -------------------------------
# Non-zero source, p/=0
# Stable fission product, decay_const=0
# Zero initial concentration, C0=0
#
timestep = 100
p = 1.0
C0 = 0.0
decay_const = 0
kernel_radius = 212.5e-6
Dprime = 1e-4
[GlobalParams]
order = SECOND
family = LAGRANGE
[]
[Mesh]
coord_type = RSPHERICAL
[gen]
type = TRISO1DMeshGenerator
elem_type = EDGE3
coordinates = '0 ${kernel_radius}'
mesh_density = 100
block_names = 'fuel'
bias = 1
dual_bias = 0.8
[]
[]
[Problem]
type = ReferenceResidualProblem
reference_vector = 'ref'
extra_tag_vectors = 'ref'
[]
[Variables]
[conc]
initial_condition = '${C0}'
[]
[]
[Kernels]
[mass_dt]
type = TimeDerivative
variable = conc
[]
[mass]
type = MatDiffusion
variable = conc
diffusivity = diffusion_coef
extra_vector_tags = 'ref'
[]
[mass_source]
type = BodyForce
variable = conc
value = '${p}'
extra_vector_tags = 'ref'
[]
[decay]
type = Decay
variable = conc
radioactive_decay_constant = ${decay_const}
extra_vector_tags = 'ref'
[]
[]
[BCs]
[freesurf_conc]
type = DirichletBC
variable = conc
boundary = exterior
value = 0.0
[]
[]
[Materials]
[fuel_conc]
type = GenericConstantMaterial
block = fuel
prop_names = diffusion_coef
prop_values = '${fparse Dprime*kernel_radius^2}'
[]
[]
[Debug]
show_var_residual_norms = true
show_var_residual = 'conc'
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
line_search = 'none'
nl_rel_tol = 1e-10
nl_abs_tol = 1e-20
nl_max_its = 20
l_tol = 1e-8
l_max_its = 30
start_time = 0.0
end_time = '${fparse 1/Dprime}'
num_steps = 100
dtmax = '${timestep}'
dtmin = '${timestep}'
[Predictor]
type = SimplePredictor
scale = 1.0
[]
[]
[Functions]
[analytic]
type = ParsedFunction
expression = '1 - x*x'
[]
[]
[Postprocessors]
[dt]
type = TimestepSize
execute_on = timestep_end
[]
# Release rate
[release_inc]
type = SideIntegralMassFlux
variable = conc
boundary = exterior
arrhenius_prpty_name = diffusion_coef
execute_on = 'initial timestep_end'
[]
[released]
type = TimeIntegratedPostprocessor
value = release_inc
execute_on = 'initial timestep_end'
[]
[total]
type = ElementIntegralVariablePostprocessor
variable = conc
execute_on = 'initial timestep_end'
[]
[x_released]
type = FractionalRelease
released = released
retained = retained
[]
[retained]
type = ElementIntegralVariablePostprocessor
variable = conc
[]
[conc_point_25]
type = PointValue
variable = conc
point = '${fparse 0.25*kernel_radius} 0 0'
execute_on = 'initial timestep_end'
[]
[conc_point_50]
type = PointValue
variable = conc
point = '${fparse 0.50*kernel_radius} 0 0'
execute_on = 'initial timestep_end'
[]
[conc_point_75]
type = PointValue
variable = conc
point = '${fparse 0.75*kernel_radius} 0 0'
execute_on = 'initial timestep_end'
[]
[conc_point_95]
type = PointValue
variable = conc
point = '${fparse 0.95*kernel_radius} 0 0'
execute_on = 'initial timestep_end'
[]
[bison_dimensionless_25]
type = ParsedPostprocessor
expression = 'conc_point_25 / (p / (6 * Dprime))'
pp_names = conc_point_25
constant_names = 'p Dprime'
constant_expressions = '${p} ${Dprime}'
execute_on = 'initial timestep_end'
[]
[bison_dimensionless_50]
type = ParsedPostprocessor
expression = 'conc_point_50 / (p / (6 * Dprime))'
pp_names = conc_point_50
constant_names = 'p Dprime'
constant_expressions = '${p} ${Dprime}'
execute_on = 'initial timestep_end'
[]
[bison_dimensionless_75]
type = ParsedPostprocessor
expression = 'conc_point_75 / (p / (6 * Dprime))'
pp_names = conc_point_75
constant_names = 'p Dprime'
constant_expressions = '${p} ${Dprime}'
execute_on = 'initial timestep_end'
[]
[bison_dimensionless_95]
type = ParsedPostprocessor
expression = 'conc_point_95 / (p / (6 * Dprime))'
pp_names = conc_point_95
constant_names = 'p Dprime'
constant_expressions = '${p} ${Dprime}'
execute_on = 'initial timestep_end'
[]
[analytic_dimensionless_25]
type = FunctionValuePostprocessor
function = analytic
point = '0.25 0 0'
execute_on = 'initial timestep_end'
[]
[analytic_dimensionless_50]
type = FunctionValuePostprocessor
function = analytic
point = '0.50 0 0'
execute_on = 'initial timestep_end'
[]
[analytic_dimensionless_75]
type = FunctionValuePostprocessor
function = analytic
point = '0.75 0 0'
execute_on = 'initial timestep_end'
[]
[analytic_dimensionless_95]
type = FunctionValuePostprocessor
function = analytic
point = '0.95 0 0'
execute_on = 'initial timestep_end'
[]
[error_25]
type = ParsedPostprocessor
pp_names = 'bison_dimensionless_25 analytic_dimensionless_25'
expression = '(bison_dimensionless_25 - analytic_dimensionless_25)/analytic_dimensionless_25*100'
execute_on = 'initial timestep_end'
[]
[error_50]
type = ParsedPostprocessor
pp_names = 'bison_dimensionless_50 analytic_dimensionless_50'
expression = '(bison_dimensionless_50 - analytic_dimensionless_50)/analytic_dimensionless_50*100'
execute_on = 'initial timestep_end'
[]
[error_75]
type = ParsedPostprocessor
pp_names = 'bison_dimensionless_75 analytic_dimensionless_75'
expression = '(bison_dimensionless_75 - analytic_dimensionless_75)/analytic_dimensionless_75*100'
execute_on = 'initial timestep_end'
[]
[error_95]
type = ParsedPostprocessor
pp_names = 'bison_dimensionless_95 analytic_dimensionless_95'
expression = '(bison_dimensionless_95 - analytic_dimensionless_95)/analytic_dimensionless_95*100'
execute_on = 'initial timestep_end'
[]
[]
[Outputs]
time_step_interval = 1
exodus = false
print_linear_converged_reason = false
print_linear_residuals = false
print_nonlinear_converged_reason = false
[console]
type = Console
show = 'analytic_dimensionless_25 analytic_dimensionless_50 analytic_dimensionless_75 analytic_dimensionless_95
bison_dimensionless_25 bison_dimensionless_50 bison_dimensionless_75 bison_dimensionless_95
error_25 error_50 error_75 error_95 released x_released'
[]
[out]
type = CSV
execute_on = final
[]
[]
(examples/TRISO/parfume/parfume.i)
# UCO TRISO particle using several PARFUME models
initial_fuel_density = 10400
[GlobalParams]
order = SECOND
family = LAGRANGE
displacements = 'disp_x'
initial_enrichment = 0.1955 # [wt-]
flux_conversion_factor = 1.0 # convert E>0.10 to E>0.18 MeV
stress_free_temperature = 923.15 # used for thermal expansion
energy_per_fission = 3.204e-11 # [J/fission]
O_U = 1.5 # Initial Oxygen to Uranium atom ratio
C_U = 0.4 # Initial Carbon to Uranium atom ratio
[]
[Mesh]
coord_type = RSPHERICAL
[gen]
type = TRISO1DMeshGenerator
elem_type = EDGE3
coordinates = '0 2.125e-4 3.125e-4 3.125e-4 3.525e-4 3.875e-4 4.275e-4'
mesh_density = '6 6 0 6 8 6'
block_names = 'fuel buffer IPyC SiC OPyC'
[]
[]
[Problem]
type = ReferenceResidualProblem
reference_vector = 'ref'
extra_tag_vectors = 'ref'
[]
[UserObjects]
[particle_geometry]
type = TRISOGeometry
[]
[]
[Variables]
[temperature]
initial_condition = 923.15
[]
[conc]
initial_condition = 0.0
scaling = 1e18
[]
[]
[AuxVariables]
[disp_y]
[]
[disp_z]
[]
[fission_rate]
block = fuel
order = CONSTANT
family = MONOMIAL
[]
[burnup]
block = fuel
order = CONSTANT
family = MONOMIAL
[]
[fast_neutron_flux]
order = CONSTANT
family = MONOMIAL
[]
[fast_neutron_fluence]
order = CONSTANT
family = MONOMIAL
[]
[gap_condSlave]
order = CONSTANT
family = MONOMIAL
[]
[density]
order = CONSTANT
family = MONOMIAL
[]
[thermal_conductivity]
order = CONSTANT
family = MONOMIAL
[]
[swelling]
order = CONSTANT
family = MONOMIAL
[]
[specific_heat]
order = CONSTANT
family = MONOMIAL
[]
[volumetric_IIDC_strain]
order = CONSTANT
family = MONOMIAL
[]
[radial_IIDC_strain]
order = CONSTANT
family = MONOMIAL
[]
[tangential_IIDC_strain]
order = CONSTANT
family = MONOMIAL
[]
[BAF]
order = CONSTANT
family = MONOMIAL
[]
[fis_gas_produced]
order = CONSTANT
family = MONOMIAL
[]
[fis_gas_released]
order = CONSTANT
family = MONOMIAL
[]
[gap_HTC]
order = CONSTANT
family = MONOMIAL
[]
[gap_distance]
order = CONSTANT
family = MONOMIAL
[]
[]
[Functions]
[power_history]
type = PiecewiseLinear
x = '0 76e6 76.001e6'
y = '1 1 0'
[]
[fission_rate]
type = LinearCombinationFunction
functions = power_history
w = 3.89e19
[]
[temp_bc]
type = PiecewiseLinear
x = '0 76e6 76.001e6 84.641e6 84.6482e6'
y = '1500 1500 300 300 2073'
[]
[k_function]
type = PiecewiseLinear
x = '0 200e6'
y = '4e-37 4e-37'
[]
[d1_function]
type = ParsedFunction
expression = 'exp(t/4.5e25)'
[]
[d_gap]
type = PiecewiseLinear
x = '1500 2100'
y = '1e-14 1e-12'
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[fuel]
block = fuel
add_variables = true
strain = FINITE
incremental = true
generate_output = 'hydrostatic_stress stress_xx stress_yy stress_zz strain_xx strain_yy strain_zz'
eigenstrain_names = 'UCO_swelling_eigenstrain UCO_thermal_strain'
extra_vector_tags = 'ref'
[]
[buffer]
block = buffer
add_variables = true
strain = FINITE
incremental = true
generate_output = 'hydrostatic_stress stress_xx stress_yy stress_zz strain_xx strain_yy strain_zz'
eigenstrain_names = 'buffer_IIDC_strain buffer_thermal_strain'
extra_vector_tags = 'ref'
[]
[IPyC]
block = IPyC
add_variables = true
strain = FINITE
incremental = true
generate_output = 'hydrostatic_stress stress_xx stress_yy stress_zz strain_xx strain_yy strain_zz'
eigenstrain_names = 'IPyC_IIDC_strain IPyC_thermal_strain'
extra_vector_tags = 'ref'
[]
[SiC]
block = SiC
add_variables = true
strain = FINITE
incremental = true
generate_output = 'hydrostatic_stress stress_xx stress_yy stress_zz strain_xx strain_yy strain_zz'
eigenstrain_names = 'SiC_thermal_eigenstrain'
extra_vector_tags = 'ref'
[]
[OPyC]
block = OPyC
add_variables = true
strain = FINITE
incremental = true
generate_output = 'hydrostatic_stress stress_xx stress_yy stress_zz strain_xx strain_yy strain_zz'
eigenstrain_names = 'OPyC_IIDC_strain OPyC_thermal_strain'
extra_vector_tags = 'ref'
[]
[]
[Kernels]
[heat_ie]
type = HeatConductionTimeDerivative
variable = temperature
extra_vector_tags = 'ref'
[]
[heat]
type = HeatConduction
variable = temperature
extra_vector_tags = 'ref'
[]
[heat_source]
type = NeutronHeatSource
variable = temperature
block = fuel
fission_rate = fission_rate
extra_vector_tags = 'ref'
[]
[mass_ie]
type = TimeDerivative
variable = conc
extra_vector_tags = 'ref'
[]
[mass]
type = ArrheniusDiffusion
variable = conc
extra_vector_tags = 'ref'
[]
[mass_source]
type = BodyForce
variable = conc
function = power_history
value = 1.22e-5 # units of mol/m**3-s
block = fuel
extra_vector_tags = 'ref'
[]
[mass_decay]
type = Decay
variable = conc
radioactive_decay_constant = 7.297e-10 # units:(1/sec) The constant for Cesium
extra_vector_tags = 'ref'
[]
[]
[AuxKernels]
[fissionrate]
type = MaterialRealAux
variable = fission_rate
property = fission_rate
block = fuel
execute_on = timestep_begin
[]
[burnup]
type = MaterialRealAux
variable = burnup
property = burnup
block = fuel
execute_on = timestep_begin
[]
[fast_neutron_flux]
type = MaterialRealAux
variable = fast_neutron_flux
property = fast_neutron_flux
block = 'fuel buffer IPyC SiC OPyC'
execute_on = timestep_begin
[]
[fast_neutron_fluence]
type = MaterialRealAux
variable = fast_neutron_fluence
property = fast_neutron_fluence
block = 'fuel buffer IPyC SiC OPyC'
execute_on = timestep_begin
[]
[conductanceSlave]
type = MaterialRealAux
property = gap_conductance
variable = gap_condSlave
boundary = buffer_outer_boundary
execute_on = 'initial timestep_end'
[]
[density]
type = MaterialRealAux
variable = density
property = density
block = 'fuel buffer IPyC SiC OPyC'
execute_on = 'initial linear'
[]
[thermal_conductivity]
type = MaterialRealAux
variable = thermal_conductivity
property = thermal_conductivity
block = 'fuel buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[specific_heat]
type = MaterialRealAux
variable = specific_heat
property = specific_heat
block = 'fuel buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[swelling]
type = MaterialRealAux
variable = swelling
property = swelling
block = fuel
execute_on = linear
[]
[volumetric_IIDC_strain]
type = MaterialRealAux
variable = volumetric_IIDC_strain
property = volumetric_IIDC_strain
block = 'IPyC OPyC'
execute_on = timestep_end
[]
[radial_IIDC_strain]
type = MaterialRealAux
variable = radial_IIDC_strain
property = radial_IIDC_strain
block = 'IPyC OPyC'
execute_on = timestep_end
[]
[tangential_IIDC_strain]
type = MaterialRealAux
variable = tangential_IIDC_strain
property = tangential_IIDC_strain
block = 'IPyC OPyC'
execute_on = timestep_end
[]
[BAF]
type = MaterialRealAux
variable = BAF
property = BAF
block = 'IPyC OPyC'
execute_on = timestep_end
[]
[fis_gas_produced]
type = MaterialRealAux
variable = fis_gas_produced
property = fis_gas_produced
block = fuel
execute_on = linear
[]
[fis_gas_released]
type = MaterialRealAux
variable = fis_gas_released
property = fis_gas_released
block = fuel
execute_on = linear
[]
[gap_HTC]
type = MaterialRealAux
property = gap_conductance
variable = gap_HTC
boundary = buffer_outer_boundary
execute_on = 'initial timestep_end'
[]
[gap_distance]
type = PenetrationAux
variable = gap_distance
boundary = buffer_outer_boundary
paired_boundary = IPyC_inner_boundary
quantity = distance
tangential_tolerance = 1e-6
execute_on = 'initial timestep_end'
[]
[]
[Contact]
[mechanical]
primary = IPyC_inner_boundary
secondary = buffer_outer_boundary
penalty = 1e5
model = frictionless
formulation = kinematic
[]
[]
[ThermalContact]
[thermal_contact]
type = GasGapHeatTransfer
variable = temperature
primary = IPyC_inner_boundary
secondary = buffer_outer_boundary
initial_moles = initial_moles # coupling to a postprocessor which supplies the initial plenum/gap gas mass
initial_gas_types = 'Kr Xe'
initial_fractions = '0.185 0.815'
gas_released = 'fis_gas_released'
released_gas_types = 'Kr Xe'
released_fractions = '0.185 0.815'
tangential_tolerance = 1e-6
roughness_primary = 0e-6
roughness_secondary = 0e-6
jumpdistance_primary = 0
jumpdistance_secondary = 0
quadrature = true
emissivity_secondary = 0.0
emissivity_primary = 0.0
min_gap = 1e-7
max_gap = 50e-6
gap_geometry_type = sphere
[]
[cesium_contact]
type = GapHeatTransfer
variable = conc
primary = IPyC_inner_boundary
secondary = buffer_outer_boundary
tangential_tolerance = 1e-6
gap_conductivity_function = d_gap
gap_conductivity_function_variable = temperature
appended_property_name = _conc
quadrature = true
gap_geometry_type = sphere
emissivity_primary = 0.0
emissivity_secondary = 0.0
min_gap = 1e-7
[]
[]
[BCs]
# pin particle along symmetry planes
[no_disp_x]
type = DirichletBC
variable = disp_x
boundary = xzero
value = 0.0
[]
# fix temperature on free surface
[freesurf_temp]
type = FunctionDirichletBC
variable = temperature
function = temp_bc
boundary = exterior
[]
# fix concentration on free surface
[freesurf_conc]
type = DirichletBC
variable = conc
boundary = exterior
value = 0.0
[]
# exterior and internal pressures
[exterior_pressure_x]
type = Pressure
variable = disp_x
boundary = exterior
factor = 0.1e6
[]
[PlenumPressure] # apply gas pressure on buffer and IPyC boundaries
[plenumPressure]
boundary = buffer_IPyC_boundary
initial_pressure = 100.0
startup_time = 0
R = 8.3145
output_initial_moles = initial_moles
temperature = ave_gas_temp
volume = 'gap_volume buffer_void_volume kernel_void_volume'
material_input = 'fis_gas_released'
output = gas_pressure
[]
[]
[]
[Materials]
[fission_rate]
type = GenericFunctionMaterial
prop_names = fission_rate
prop_values = fission_rate
[]
[fast_neutron_flux]
type = FastNeutronFlux
calculate_fluence = true
flux_function = power_history
factor = 5e17
[]
[normal_vectors_triso]
type = NormalVectorsTRISO
block = 'buffer IPyC OPyC'
[]
### UCO properties
[UCO_burnup]
type = TRISOBurnup
initial_density = ${initial_fuel_density}
[]
[UCO_thermal]
type = UCOThermal
block = fuel
temperature = temperature
[]
[UCO_elasticity_tensor]
type = UCOElasticityTensor
block = fuel
temperature = temperature
[]
[UCO_stress]
type = ComputeFiniteStrainElasticStress
block = fuel
[]
[UCO_VolumetricSwellingEigenstrain]
type = UCOVolumetricSwellingEigenstrain
block = fuel
eigenstrain_name = UCO_swelling_eigenstrain
[]
[UCO_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = fuel
thermal_expansion_coeff = 10e-6 # check this value for UCO
eigenstrain_name = UCO_thermal_strain
temperature = temperature
[]
[UCO_density]
type = StrainAdjustedDensity
block = fuel
strain_free_density = ${initial_fuel_density}
[]
[fission_gas_release]
type = UCOFGR
block = fuel
average_grain_radius = 10e-6
temperature = temperature
triso_geometry = particle_geometry
cutoff_neutron_flux = 0.0
[]
[fuel_conc]
type = ArrheniusDiffusionCoef
block = fuel
d1 = 5.6e-8 # m^2/s
q1 = 209.0e+3 # J/mol
d2 = 5.2e-4 # m^2/s
q2 = 362.0e+3 # J/mol
temperature = temperature
[]
### Buffer Properties
[buffer_elasticity_tensor]
type = BufferElasticityTensor
block = buffer
temperature = temperature
[]
[buffer_stress]
type = BufferCEGACreep
block = buffer
temperature = temperature
[]
[buffer_thermal]
type = BufferThermal
block = buffer
initial_density = 1050.0
[]
[buffer_density]
type = StrainAdjustedDensity
block = buffer
strain_free_density = 1050.0
[]
[buffer_thermal_strain]
type = BufferThermalExpansionEigenstrain
block = buffer
eigenstrain_name = buffer_thermal_strain
temperature = temperature
[]
[buffer_IIDC_strain]
type = BufferCEGAIrradiationEigenstrain
block = buffer
eigenstrain_name = buffer_IIDC_strain
temperature = temperature
[]
[buffer_conc]
type = ArrheniusDiffusionCoef
block = buffer
d1 = 1.0e-12 # m^2/s
q1 = 0.0
d2 = 0.0
q2 = 0.0
temperature = temperature
[]
### IPyC properties
[IPyC_elasticity_tensor]
type = PyCElasticityTensor
block = IPyC
temperature = temperature
initial_BAF = 1.045
[]
[IPyC_stress]
type = PyCCEGACreep
block = IPyC
temperature = temperature
[]
[IPyC_thermal]
type = HeatConductionMaterial
block = IPyC
thermal_conductivity = 4.0
specific_heat = 720.0
[]
[IPyC_density]
type = GenericConstantMaterial
block = IPyC
prop_names = 'density'
prop_values = 1900.0
[]
[IPyC_IIDC_strain]
type = PyCCEGAIrradiationEigenstrain
block = IPyC
eigenstrain_name = IPyC_IIDC_strain
temperature = temperature
[]
[BAF_IPyC]
type = BaconAnisotropyFactor
initial_BAF = 1.045
block = IPyC
[]
[BAF_OPyC]
type = BaconAnisotropyFactor
initial_BAF = 1.045
block = OPyC
[]
[IPyC_thermal_strain]
type = PyCThermalExpansionEigenstrain
block = IPyC
eigenstrain_name = IPyC_thermal_strain
temperature = temperature
[]
[IPyC_conc]
type = ArrheniusDiffusionCoef
block = IPyC
d1 = 6.3e-8 # m^2/s
q1 = 222.0e+3 # J/mol
d2 = 0.0
q2 = 0.0
temperature = temperature
[]
### SiC properties
[SiC_elasticity_tensor]
type = MonolithicSiCElasticityTensor
block = SiC
temperature = temperature
elastic_modulus_model = miller
[]
[SiC_stress]
type = ComputeFiniteStrainElasticStress
block = SiC
[]
[SiC_thermal]
type = MonolithicSiCThermal
block = SiC
temperature = temperature
thermal_conductivity_model = miller
[]
[SiC_density]
type = GenericConstantMaterial
block = SiC
prop_names = 'density'
prop_values = 3200.0
[]
[SiC_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = SiC
thermal_expansion_coeff = 4.9e-6
temperature = temperature
eigenstrain_name = SiC_thermal_eigenstrain
[]
[SiC_conc]
type = ArrheniusDiffusionCoef
block = SiC
d1 = 5.5e-14 # m^2/s
d1_function = d1_function
d1_function_variable = fast_neutron_fluence
q1 = 125.0e+3 # J/mol
d2 = 1.6e-2 # m^2/s
q2 = 514.0e+3 # J/mol
temperature = temperature
[]
### OPyC properties
[OPyC_elasticity_tensor]
type = PyCElasticityTensor
block = OPyC
temperature = temperature
initial_BAF = 1.045
[]
[OPyC_stress]
type = PyCCEGACreep
block = OPyC
temperature = temperature
[]
[OPyC_thermal_conductivity]
type = HeatConductionMaterial
block = OPyC
thermal_conductivity = 4.0
specific_heat = 720.0
[]
[OPyC_density]
type = GenericConstantMaterial
block = OPyC
prop_names = 'density'
prop_values = 1900.0
[]
[OPyC_IIDC_strain]
type = PyCCEGAIrradiationEigenstrain
block = OPyC
eigenstrain_name = OPyC_IIDC_strain
temperature = temperature
[]
[OPyC_thermal_strain]
type = PyCThermalExpansionEigenstrain
block = OPyC
eigenstrain_name = OPyC_thermal_strain
temperature = temperature
[]
[OPyC_conc]
type = ArrheniusDiffusionCoef
block = OPyC
d1 = 6.3e-8 # m^2/s
q1 = 222.0e+3 # J/mol
d2 = 0.0
q2 = 0.0
temperature = temperature
[]
[]
[Dampers]
[temp]
type = MaxIncrement
variable = temperature
max_increment = 50
[]
[disp_x]
type = MaxIncrement
variable = disp_x
max_increment = 1e-6
[]
[]
[Debug]
show_var_residual_norms = true
show_var_residual = 'disp_x temperature conc'
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
line_search = 'none'
nl_rel_tol = 1e-8
nl_abs_tol = 1e-7 #1e-12
nl_max_its = 15
l_tol = 1e-4 #1e-8
l_max_its = 50
start_time = 0.0
end_time = 85.3682e6 #5.0e7
num_steps = 1000
dtmax = 2e6
dtmin = 1
[TimeStepper]
type = IterationAdaptiveDT
dt = 20
growth_factor = 1.5
optimal_iterations = 8 #6
linear_iteration_ratio = 100
time_t = '0 76e6 76.001e6 84.641e6 84.6482e6'
time_dt = '20 20 20 20 20'
[]
[Quadrature]
order = THIRD
[]
[]
[Postprocessors]
[_dt]
type = TimestepSize
execute_on = timestep_end
[]
[cs_release]
type = SideIntegralMassFlux
variable = conc
boundary = exterior
[]
[int_cs_release]
type = TimeIntegratedPostprocessor
value = cs_release
[]
[cs_release_fuel]
type = SideIntegralMassFlux
variable = conc
boundary = fuel_outer_boundary
[]
[int_cs_release_fuel]
type = TimeIntegratedPostprocessor
value = cs_release_fuel
[]
[cs_release_PyCGapBndry]
type = SideIntegralMassFlux
variable = conc
boundary = IPyC_inner_boundary
[]
[int_cs_release_PyCGapBndry]
type = TimeIntegratedPostprocessor
value = cs_release_PyCGapBndry
[]
[ave_gas_temp]
type = ElementAverageValue
block = buffer
variable = temperature
execute_on = 'initial timestep_end'
[]
[ave_gap_temp]
type = SideAverageValue
boundary = buffer_IPyC_boundary
variable = temperature
execute_on = 'initial timestep_end'
[]
[fis_gas_produced]
type = ElementIntegralMaterialProperty
mat_prop = fis_gas_produced
block = fuel
execute_on = 'initial timestep_end'
[]
[fis_gas_released]
type = ElementIntegralMaterialProperty
mat_prop = fis_gas_released
block = fuel
execute_on = 'initial timestep_end'
[]
[gap_volume]
type = InternalVolume
boundary = buffer_IPyC_boundary
execute_on = 'initial linear'
use_displaced_mesh = true
[]
[buffer_void_volume]
type = VoidVolume
block = buffer
theoretical_density = 2250
execute_on = 'initial timestep_end'
use_displaced_mesh = true
[]
[kernel_th_density]
type = UCOTheoreticalDensity
execute_on = initial
[]
[kernel_void_volume]
type = VoidVolume
block = fuel
theoretical_density = kernel_th_density
execute_on = 'initial timestep_end'
use_displaced_mesh = true
[]
[volumeTotal]
type = InternalVolume
boundary = exterior
execute_on = 'initial timestep_end'
scale_factor = -1
[]
[volumeFuel]
type = InternalVolume
boundary = fuel_outer_boundary
execute_on = 'initial timestep_end'
scale_factor = -1
[]
[volumeGas]
type = InternalVolume
boundary = buffer_outer_boundary
execute_on = 'initial timestep_end'
scale_factor = -1
addition = 4.67e-11
[]
[volumeBufferShell]
type = InternalVolume
boundary = buffer_outer_boundary
execute_on = 'initial timestep_end'
[]
[ave_temp_interior]
type = SideAverageValue
boundary = buffer_outer_boundary
variable = temperature
execute_on = 'initial timestep_end'
[]
[gap_HTC]
type = ElementExtremeValue
variable = gap_HTC
block = buffer
value_type = 'max'
execute_on = 'initial timestep_end'
[]
### Postprocessors for CO production
[total_fission_rate]
type = ElementIntegralPower
variable = temperature
fission_rate = fission_rate
block = fuel
energy_per_fission = 1.0
outputs = exodus
execute_on = 'initial timestep_end'
[]
[total_fissions]
type = TimeIntegratedPostprocessor
value = total_fission_rate
outputs = exodus
execute_on = 'initial timestep_end'
[]
[avg_surface_temp]
type = SideAverageValue
variable = temperature
boundary = exterior
outputs = exodus
execute_on = 'initial timestep_end'
[]
[time_int_surf_temp]
type = TimeIntegratedPostprocessor
value = avg_surface_temp
outputs = exodus
execute_on = 'initial timestep_end'
[]
##### irradiation conditions
[particle_power]
type = ElementIntegralPower
variable = temperature
use_material_fission_rate = true
fission_rate_material = fission_rate
block = fuel
execute_on = 'initial timestep_end'
[]
[max_fluence]
type = ElementExtremeValue
variable = fast_neutron_fluence
block = 'fuel buffer IPyC SiC OPyC'
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[max_burnup]
type = ElementExtremeValue
variable = burnup
block = fuel
value_type = 'max'
execute_on = 'initial timestep_end'
[]
#### II strain
[OPyC_radial_IIDC_strain]
type = ElementExtremeValue
variable = radial_IIDC_strain
block = OPyC
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[OPyC_tangential_IIDC_strain]
type = ElementExtremeValue
variable = tangential_IIDC_strain
block = OPyC
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[IPyC_radial_IIDC_strain]
type = ElementExtremeValue
variable = radial_IIDC_strain
block = IPyC
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[IPyC_tangential_IIDC_strain]
type = ElementExtremeValue
variable = tangential_IIDC_strain
block = IPyC
value_type = 'max'
execute_on = 'initial timestep_end'
[]
#### temperatures
[max_T_kernel]
type = NodalExtremeValue
variable = temperature
block = fuel
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[max_T_buffer]
type = NodalExtremeValue
variable = temperature
block = buffer
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[min_T_buffer]
type = NodalExtremeValue
variable = temperature
block = buffer
value_type = 'min'
execute_on = 'initial timestep_end'
[]
[max_T_IPyC]
type = NodalExtremeValue
variable = temperature
block = IPyC
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[max_T_SiC]
type = NodalExtremeValue
variable = temperature
block = SiC
value_type = 'max'
execute_on = 'initial timestep_end'
[]
#### displacement BCs
[max_disp_kernel]
type = NodalExtremeValue
variable = disp_x
block = fuel
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[min_disp_buffer]
type = NodalExtremeValue
variable = disp_x
block = buffer
value_type = 'min'
execute_on = 'initial timestep_end'
[]
[max_disp_IPyC]
type = NodalExtremeValue
variable = disp_x
block = IPyC
value_type = 'max'
execute_on = 'initial timestep_end'
[]
#### hoop stresses
[hoop_opyc_max]
type = ElementExtremeValue
variable = stress_yy
block = OPyC
execute_on = 'initial timestep_end'
[]
[hoop_sic_max]
type = ElementExtremeValue
variable = stress_yy
block = SiC
execute_on = 'initial timestep_end'
[]
[hoop_ipyc_max]
type = ElementExtremeValue
variable = stress_yy
block = IPyC
execute_on = 'initial timestep_end'
[]
[hoop_buffer_max]
type = ElementExtremeValue
variable = stress_yy
block = buffer
execute_on = 'initial timestep_end'
[]
[hoop_opyc_min]
type = ElementExtremeValue
variable = stress_yy
block = OPyC
value_type = min
execute_on = 'initial timestep_end'
[]
[hoop_sic_min]
type = ElementExtremeValue
variable = stress_yy
block = SiC
value_type = min
execute_on = 'initial timestep_end'
[]
[hoop_ipyc_min]
type = ElementExtremeValue
variable = stress_yy
block = IPyC
value_type = min
execute_on = 'initial timestep_end'
[]
[hoop_buffer_min]
type = ElementExtremeValue
variable = stress_yy
block = buffer
value_type = min
execute_on = 'initial timestep_end'
[]
### Check warning for Density
[oPyC_density]
type = ElementExtremeValue
variable = density
block = OPyC
execute_on = 'initial timestep_end'
[]
[sic_density]
type = ElementExtremeValue
variable = density
block = SiC
execute_on = 'initial timestep_end'
[]
[IPyC_density]
type = ElementExtremeValue
variable = density
block = IPyC
execute_on = 'initial timestep_end'
[]
[buffer_density]
type = ElementExtremeValue
variable = density
block = buffer
execute_on = 'initial timestep_end'
[]
[kernel_density]
type = ElementExtremeValue
variable = density
block = fuel
execute_on = 'initial timestep_end'
[]
[pd_penetration]
type = PdPenetration
boundary = SiC_inner_boundary
variable = temperature
execute_on = 'initial timestep_end'
[]
[]
[Outputs]
print_linear_residuals = true
time_step_interval = 1
exodus = true
csv = true
perf_graph = true
[]
(assessment/verification/TRISO_diffusion/outofpile/outofpile.i)
#
# TRISO kernel in 1D spherical
#
# Out-of-pile, a decaying fission product
# -------------------------------
# Zero source, p=0
# Decaying fission product, decay_const/=0
# Non-zero initial concentration, C0/=0
#
timestep = 10
p = 0.0
C0 = 1.0
decay_const = 1e-5
Dprime = 1e-4
kernel_radius = 212.5e-6
[GlobalParams]
order = SECOND
family = LAGRANGE
[]
[Mesh]
coord_type = RSPHERICAL
[gen]
type = TRISO1DMeshGenerator
elem_type = EDGE3
coordinates = '0 ${kernel_radius}'
mesh_density = 100
block_names = 'fuel'
bias = 1
dual_bias = 0.8
[]
[]
[Problem]
type = ReferenceResidualProblem
reference_vector = 'ref'
extra_tag_vectors = 'ref'
[]
[Variables]
[conc]
initial_condition = '${C0}'
[]
[]
[Kernels]
[mass_dt]
type = TimeDerivative
variable = conc
[]
[mass]
type = MatDiffusion
variable = conc
diffusivity = diffusion_coef
extra_vector_tags = 'ref'
[]
[mass_source]
type = BodyForce
variable = conc
value = '${p}'
extra_vector_tags = 'ref'
[]
[decay]
type = Decay
variable = conc
radioactive_decay_constant = ${decay_const}
extra_vector_tags = 'ref'
[]
[]
[BCs]
[freesurf_conc]
type = DirichletBC
variable = conc
boundary = exterior
value = 0.0
[]
[]
[Materials]
[fuel_conc]
type = GenericConstantMaterial
block = fuel
prop_names = diffusion_coef
prop_values = '${fparse Dprime*kernel_radius^2}'
[]
[]
[Debug]
show_var_residual_norms = true
show_var_residual = 'conc'
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
line_search = 'none'
nl_rel_tol = 1e-10
nl_abs_tol = 1e-20
nl_max_its = 20
l_tol = 1e-8
l_max_its = 30
start_time = 0.0
end_time = '${fparse 1/Dprime}'
num_steps = 100
dtmax = '${timestep}'
dtmin = '${timestep}'
[Predictor]
type = SimplePredictor
scale = 1.0
[]
[]
[Postprocessors]
[dt]
type = TimestepSize
execute_on = timestep_end
[]
# Release rate
[release_inc]
type = SideIntegralMassFlux
variable = conc
boundary = exterior
arrhenius_prpty_name = diffusion_coef
execute_on = 'initial timestep_end'
[]
[released]
type = TimeIntegratedPostprocessor
value = release_inc
execute_on = 'initial timestep_end'
[]
[total]
type = ElementIntegralVariablePostprocessor
variable = conc
execute_on = 'initial timestep_end'
[]
[x_released]
type = FractionalRelease
released = released
retained = retained
[]
[retained]
type = ElementIntegralVariablePostprocessor
variable = conc
[]
[conc_point_25]
type = PointValue
variable = conc
point = '${fparse 0.25*kernel_radius} 0 0'
execute_on = 'initial timestep_end'
[]
[conc_point_50]
type = PointValue
variable = conc
point = '${fparse 0.50*kernel_radius} 0 0'
execute_on = 'initial timestep_end'
[]
[conc_point_75]
type = PointValue
variable = conc
point = '${fparse 0.75*kernel_radius} 0 0'
execute_on = 'initial timestep_end'
[]
[conc_point_95]
type = PointValue
variable = conc
point = '${fparse 0.95*kernel_radius} 0 0'
execute_on = 'initial timestep_end'
[]
[bison_dimensionless_25]
type = ParsedPostprocessor
expression = 'conc_point_25 / C0'
pp_names = conc_point_25
constant_names = 'C0'
constant_expressions = '${C0}'
execute_on = 'initial timestep_end'
[]
[bison_dimensionless_50]
type = ParsedPostprocessor
expression = 'conc_point_50 / C0'
pp_names = conc_point_50
constant_names = 'C0'
constant_expressions = '${C0}'
execute_on = 'initial timestep_end'
[]
[bison_dimensionless_75]
type = ParsedPostprocessor
expression = 'conc_point_75 / C0'
pp_names = conc_point_75
constant_names = 'C0'
constant_expressions = '${C0}'
execute_on = 'initial timestep_end'
[]
[bison_dimensionless_95]
type = ParsedPostprocessor
expression = 'conc_point_95 / C0'
pp_names = conc_point_95
constant_names = 'C0'
constant_expressions = '${C0}'
execute_on = 'initial timestep_end'
[]
[analytic_dimensionless_25]
type = ParsedPostprocessor
# The analytic expression involves a sum we can't compute here
expression = '0.640190356'
pp_names = ''
execute_on = 'initial timestep_end'
[]
[analytic_dimensionless_50]
type = ParsedPostprocessor
# The analytic expression involves a sum we can't compute here
expression = '0.469766231'
pp_names = ''
execute_on = 'initial timestep_end'
[]
[analytic_dimensionless_75]
type = ParsedPostprocessor
# The analytic expression involves a sum we can't compute here
expression = '0.229613020'
pp_names = ''
execute_on = 'initial timestep_end'
[]
[analytic_dimensionless_95]
type = ParsedPostprocessor
# The analytic expression involves a sum we can't compute here
expression = '0.040674510'
pp_names = ''
execute_on = 'initial timestep_end'
[]
[error_25]
type = ParsedPostprocessor
pp_names = 'bison_dimensionless_25 analytic_dimensionless_25'
expression = '(bison_dimensionless_25 - analytic_dimensionless_25)/analytic_dimensionless_25*100'
execute_on = 'initial timestep_end'
[]
[error_50]
type = ParsedPostprocessor
pp_names = 'bison_dimensionless_50 analytic_dimensionless_50'
expression = '(bison_dimensionless_50 - analytic_dimensionless_50)/analytic_dimensionless_50*100'
execute_on = 'initial timestep_end'
[]
[error_75]
type = ParsedPostprocessor
pp_names = 'bison_dimensionless_75 analytic_dimensionless_75'
expression = '(bison_dimensionless_75 - analytic_dimensionless_75)/analytic_dimensionless_75*100'
execute_on = 'initial timestep_end'
[]
[error_95]
type = ParsedPostprocessor
pp_names = 'bison_dimensionless_95 analytic_dimensionless_95'
expression = '(bison_dimensionless_95 - analytic_dimensionless_95)/analytic_dimensionless_95*100'
execute_on = 'initial timestep_end'
[]
[]
[Outputs]
time_step_interval = 1
exodus = false
print_linear_converged_reason = false
print_linear_residuals = false
print_nonlinear_converged_reason = false
[console]
type = Console
show = 'analytic_dimensionless_25 analytic_dimensionless_50 analytic_dimensionless_75 analytic_dimensionless_95
bison_dimensionless_25 bison_dimensionless_50 bison_dimensionless_75 bison_dimensionless_95
error_25 error_50 error_75 error_95 released x_released'
[]
[out]
type = CSV
execute_on = final
[]
[]
(assessment/verification/TRISO_diffusion/inpile/inpile.i)
#
# TRISO kernel in 1D spherical
#
# In-pile, a decaying fission product
# -------------------------------
# Non-zero source, p/=0
# Decaying fission product, decay_const/=0
# Zero initial concentration, C0=0
#
timestep = 100
p = 1.0
C0 = 0.0
decay_const = 1e-5
Dprime = 1e-4
kernel_radius = 212.5e-6
[GlobalParams]
order = SECOND
family = LAGRANGE
[]
[Mesh]
coord_type = RSPHERICAL
[gen]
type = TRISO1DMeshGenerator
elem_type = EDGE3
coordinates = '0 ${kernel_radius}'
mesh_density = 100
block_names = 'fuel'
bias = 1
dual_bias = 0.8
[]
[]
[Problem]
type = ReferenceResidualProblem
reference_vector = 'ref'
extra_tag_vectors = 'ref'
[]
[Variables]
[conc]
initial_condition = '${C0}'
[]
[]
[Kernels]
[mass_dt]
type = TimeDerivative
variable = conc
[]
[mass]
type = MatDiffusion
variable = conc
diffusivity = diffusion_coef
extra_vector_tags = 'ref'
[]
[mass_source]
type = BodyForce
variable = conc
value = '${p}'
extra_vector_tags = 'ref'
[]
[decay]
type = Decay
variable = conc
radioactive_decay_constant = ${decay_const}
extra_vector_tags = 'ref'
[]
[]
[BCs]
[freesurf_conc]
type = DirichletBC
variable = conc
boundary = exterior
value = 0.0
[]
[]
[Materials]
[fuel_conc]
type = GenericConstantMaterial
block = fuel
prop_names = diffusion_coef
prop_values = '${fparse Dprime*kernel_radius^2}'
[]
[]
[Debug]
show_var_residual_norms = true
show_var_residual = 'conc'
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
line_search = 'none'
nl_rel_tol = 1e-10
nl_abs_tol = 1e-20
nl_max_its = 20
l_tol = 1e-8
l_max_its = 30
start_time = 0.0
end_time = '${fparse 1/Dprime}'
num_steps = 100
dtmax = '${timestep}'
dtmin = '${timestep}'
[Predictor]
type = SimplePredictor
scale = 1.0
[]
[]
[Functions]
[analytic]
type = ParsedFunction
symbol_names = 'mu'
symbol_values = '${fparse decay_const/Dprime}'
expression = '1-sinh(x*sqrt(mu))/(x*sinh(sqrt(mu)))'
[]
[]
[Postprocessors]
[dt]
type = TimestepSize
execute_on = timestep_end
[]
# Release rate
[release_inc]
type = SideIntegralMassFlux
variable = conc
boundary = exterior
arrhenius_prpty_name = diffusion_coef
execute_on = 'initial timestep_end'
[]
[released]
type = TimeIntegratedPostprocessor
value = release_inc
execute_on = 'initial timestep_end'
[]
[total]
type = ElementIntegralVariablePostprocessor
variable = conc
execute_on = 'initial timestep_end'
[]
[x_released]
type = FractionalRelease
released = released
retained = retained
[]
[retained]
type = ElementIntegralVariablePostprocessor
variable = conc
[]
[conc_point_25]
type = PointValue
variable = conc
point = '${fparse 0.25*kernel_radius} 0 0'
execute_on = 'initial timestep_end'
[]
[conc_point_50]
type = PointValue
variable = conc
point = '${fparse 0.50*kernel_radius} 0 0'
execute_on = 'initial timestep_end'
[]
[conc_point_75]
type = PointValue
variable = conc
point = '${fparse 0.75*kernel_radius} 0 0'
execute_on = 'initial timestep_end'
[]
[conc_point_95]
type = PointValue
variable = conc
point = '${fparse 0.95*kernel_radius} 0 0'
execute_on = 'initial timestep_end'
[]
[bison_dimensionless_25]
type = ParsedPostprocessor
expression = 'conc_point_25 / (p / lambda)'
pp_names = conc_point_25
constant_names = 'p lambda'
constant_expressions = '${p} ${decay_const}'
execute_on = 'initial timestep_end'
[]
[bison_dimensionless_50]
type = ParsedPostprocessor
expression = 'conc_point_50 / (p / lambda)'
pp_names = conc_point_50
constant_names = 'p lambda'
constant_expressions = '${p} ${decay_const}'
execute_on = 'initial timestep_end'
[]
[bison_dimensionless_75]
type = ParsedPostprocessor
expression = 'conc_point_75 / (p / lambda)'
pp_names = conc_point_75
constant_names = 'p lambda'
constant_expressions = '${p} ${decay_const}'
execute_on = 'initial timestep_end'
[]
[bison_dimensionless_95]
type = ParsedPostprocessor
expression = 'conc_point_95 / (p / lambda)'
pp_names = conc_point_95
constant_names = 'p lambda'
constant_expressions = '${p} ${decay_const}'
execute_on = 'initial timestep_end'
[]
[analytic_dimensionless_25]
type = FunctionValuePostprocessor
function = analytic
point = '0.25 0 0'
execute_on = 'initial timestep_end'
[]
[analytic_dimensionless_50]
type = FunctionValuePostprocessor
function = analytic
point = '0.50 0 0'
execute_on = 'initial timestep_end'
[]
[analytic_dimensionless_75]
type = FunctionValuePostprocessor
function = analytic
point = '0.75 0 0'
execute_on = 'initial timestep_end'
[]
[analytic_dimensionless_95]
type = FunctionValuePostprocessor
function = analytic
point = '0.95 0 0'
execute_on = 'initial timestep_end'
[]
[error_25]
type = ParsedPostprocessor
pp_names = 'bison_dimensionless_25 analytic_dimensionless_25'
expression = '(bison_dimensionless_25 - analytic_dimensionless_25)/analytic_dimensionless_25*100'
execute_on = 'initial timestep_end'
[]
[error_50]
type = ParsedPostprocessor
pp_names = 'bison_dimensionless_50 analytic_dimensionless_50'
expression = '(bison_dimensionless_50 - analytic_dimensionless_50)/analytic_dimensionless_50*100'
execute_on = 'initial timestep_end'
[]
[error_75]
type = ParsedPostprocessor
pp_names = 'bison_dimensionless_75 analytic_dimensionless_75'
expression = '(bison_dimensionless_75 - analytic_dimensionless_75)/analytic_dimensionless_75*100'
execute_on = 'initial timestep_end'
[]
[error_95]
type = ParsedPostprocessor
pp_names = 'bison_dimensionless_95 analytic_dimensionless_95'
expression = '(bison_dimensionless_95 - analytic_dimensionless_95)/analytic_dimensionless_95*100'
execute_on = 'initial timestep_end'
[]
[]
[Outputs]
time_step_interval = 1
exodus = false
print_linear_converged_reason = false
print_linear_residuals = false
print_nonlinear_converged_reason = false
[console]
type = Console
show = 'analytic_dimensionless_25 analytic_dimensionless_50 analytic_dimensionless_75 analytic_dimensionless_95
bison_dimensionless_25 bison_dimensionless_50 bison_dimensionless_75 bison_dimensionless_95
error_25 error_50 error_75 error_95 released x_released'
[]
[out]
type = CSV
execute_on = final
[]
[]
(examples/metal_fuel/mechanistic_fcci/fcci.i)
# This example exercises the mechanistic Fuel-Cladding Chemical Interaction
# (FCCI) modeling approach for U-Zr and U-Pu-Zr metallic nuclear fuels. It is
# run as a SubApp, which is called by the combination of an existing assessment
# case and XXXX_multiapp.i.
# It is a collection of BISON classes and objects needed to apply a mechanistic
# FCCI model for U-Zr and U-Pu-Zr metallic nuclear fuels. It is designed to
# complement the modeling approaches adopted by the FIPD-powered X447 assessment
# cases and manually-developed X430 and IFR1 assessment cases. Thermomechanical
# behaviors are modeled in the existing assessment cases whereas lanthanide (Ln)
# transport, FCCI, and wastage growth are modeled in the SubApp. XXXX_multiapp.i
# then overrides the empirical FCCI model currently used in the assessments and
# replaces it with the mechanistic FCCI model. Examples of how to do this are
# included in the test and example specifications.
# Note that, in this model, c indicates the stable concentration of Lns. The
# model relies on an effective diffusion coefficient and flux terms informed by
# lower length scale modeling and simulations. Decay terms are included for
# illustration purposes but do not affect the solution.
# combined effective Ln production constant in mol/fission
ln_yield = 0.3 # atom/fission
avogadro = 6.0221e23 # atom/mol
production_constant = '${fparse ln_yield / avogadro}' # mol/fission
# combined effective Ln decay constant in 1/s
decay_constant = 0
# variable order, block and boundary names, and end time used in the assessment
order = placeholder
fuel_block_name = placeholder
fuel_surface_name = placeholder
end_time = placeholder
# fuel rod geometry and material properties used in the assessment
cladding_density = placeholder
fuel_outer_radius = placeholder
gap_thickness = placeholder
[GlobalParams]
order = ${order}
family = LAGRANGE
[]
[Problem]
# nothing is being solved for on the cladding
kernel_coverage_check = false
material_coverage_check = false
[]
[Mesh]
[]
[Variables]
[c] # stable Ln concentration in mol/m^3
block = ${fuel_block_name}
[]
[]
[AuxVariables]
[T]
initial_condition = 298
[]
[fission_rate]
block = ${fuel_block_name}
family = MONOMIAL
order = CONSTANT
[]
[porosity_aux]
block = ${fuel_block_name}
family = MONOMIAL
order = CONSTANT
[]
[sodium_logged_porosity_aux]
block = ${fuel_block_name}
family = MONOMIAL
order = CONSTANT
[]
[diffusivity]
block = ${fuel_block_name}
family = MONOMIAL
order = CONSTANT
[]
[dcdx]
block = ${fuel_block_name}
family = MONOMIAL
order = CONSTANT
[]
[c_flux]
block = ${fuel_block_name}
family = MONOMIAL
order = CONSTANT
[]
[]
[AuxKernels]
[diffusivity]
type = MaterialRealAux
boundary = ${fuel_surface_name}
variable = diffusivity
property = D
[]
[dcdx]
type = VariableGradientComponent
boundary = ${fuel_surface_name}
variable = dcdx
gradient_variable = c
component = x
[]
[c_flux]
type = ParsedAux
boundary = ${fuel_surface_name}
variable = c_flux
coupled_variables = 'diffusivity dcdx'
expression = '-diffusivity * dcdx'
[]
[]
[Kernels]
[c_dt]
type = TimeDerivative
variable = c
[]
[c_diffusion]
type = MatDiffusion
block = ${fuel_block_name}
variable = c
diffusivity = D
args = 'c T'
[]
[c_source]
type = NeutronHeatSource
block = ${fuel_block_name}
variable = c
fission_rate = fission_rate # fission/m^3-s
energy_per_fission = ${production_constant} # mol/fission giving mol/m^3-s
[]
[c_decay]
type = Decay
block = ${fuel_block_name}
variable = c
radioactive_decay_constant = ${decay_constant} # 1/s
[]
[]
[Materials]
[porosity]
type = ParsedMaterial
block = ${fuel_block_name}
property_name = porosity
coupled_variables = porosity_aux
expression = porosity_aux
[]
[sodium_logged_porosity]
type = ParsedMaterial
block = ${fuel_block_name}
property_name = sodium_logged_porosity
coupled_variables = sodium_logged_porosity_aux
expression = sodium_logged_porosity_aux
[]
[diffusivity]
# informed by lower length scale simulations
type = UPuZrLanthanideDiffusivity
block = ${fuel_block_name}
property_name = D
temperature = T
porosity_name = porosity
sodium_logged_porosity_name = sodium_logged_porosity
outputs = exodus
output_properties = D
[]
[flux]
# informed by lower length scale simulations
type = UPuZrLanthanideFlux
boundary = ${fuel_surface_name}
property_name = F
temperature = T
concentration = c
outputs = exodus
output_properties = F
[]
[wastage_thickness]
# informed by lower length scale simulations
type = UPuZrLanthanideWastage
boundary = ${fuel_surface_name}
cladding_density = ${cladding_density}
fuel_outer_radius = ${fuel_outer_radius}
gap_thickness =${gap_thickness}
lanthanide_flux = F
outputs = exodus
[]
[]
[BCs]
[c_right]
type = MatNeumannBC
boundary = ${fuel_surface_name}
variable = c
boundary_material = F
value = -1
[]
[]
[Executioner]
type = Transient
solve_type = PJFNK
automatic_scaling = true
compute_scaling_once = false
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package -pc_factor_shift_type -pc_factor_shift_amount'
petsc_options_value = 'lu superlu_dist NONZERO 1e-15'
line_search = NONE
l_max_its = 30
l_tol = 1e-3
nl_max_its = 30
nl_rel_tol = 1e-10
nl_abs_tol = 1e-6
start_time = 0
dt = 1e7
end_time = ${end_time}
verbose = true
[]
[Postprocessors]
[c_in_fuel] # mol
type = ElementIntegralVariablePostprocessor
block = ${fuel_block_name}
variable = c
execute_on = 'initial timestep_end'
outputs = 'console csv'
[]
[fission_rate_total] # fission/s
type = ElementIntegralVariablePostprocessor
block = ${fuel_block_name}
variable = fission_rate
execute_on = 'initial timestep_end'
outputs = csv
[]
[fission_total] # fission
type = TimeIntegratedPostprocessor
value = fission_rate_total
execute_on = 'initial timestep_end'
outputs = csv
[]
[c_produced] # mol
type = LinearCombinationPostprocessor
pp_names = fission_total
pp_coefs = ${production_constant} # mol/fission giving mol
execute_on = 'initial timestep_end'
outputs = 'console csv'
[]
[c_decay_rate] # mol/s
type = LinearCombinationPostprocessor
pp_names = c_in_fuel # mol
pp_coefs = ${decay_constant} # 1/s
execute_on = 'initial timestep_end'
outputs = csv
[]
[c_decayed] # mol
type = TimeIntegratedPostprocessor
value = c_decay_rate
execute_on = 'initial timestep_end'
outputs = csv
[]
[c_reaction_rate] # mol/s
type = SideDiffusiveFluxIntegral
boundary = ${fuel_surface_name}
variable = c
diffusivity = D
execute_on = 'initial timestep_end'
outputs = csv
[]
[c_reacted] # mol
type = TimeIntegratedPostprocessor
value = c_reaction_rate
execute_on = 'initial timestep_end'
outputs = 'console csv'
[]
[c_flux_max]
type = SideExtremeValue
boundary = ${fuel_surface_name}
variable = c_flux
value_type = max
execute_on = 'initial timestep_end'
outputs = csv
[]
[c_flux_min]
type = SideExtremeValue
boundary = ${fuel_surface_name}
variable = c_flux
value_type = min
execute_on = 'initial timestep_end'
outputs = csv
[]
[wastage_thickness_max] # m
type = ElementExtremeValue
variable = wastage_thickness
value_type = max
execute_on = 'initial timestep_end'
outputs = 'console csv'
[]
[c_surface_max] # mol
type = SideExtremeValue
boundary = ${fuel_surface_name}
variable = c
value_type = max
execute_on = 'initial timestep_end'
outputs = csv
[]
[c_surface_min] # mol
type = SideExtremeValue
boundary = ${fuel_surface_name}
variable = c
value_type = min
execute_on = 'initial timestep_end'
outputs = csv
[]
[]
# [VectorPostprocessors]
# [surface_conditions]
# type = SideValueSampler
# variable = 'T D F c c_flux wastage_thickness'
# boundary = ${fuel_surface_name}
# sort_by = y
# execute_on = 'initial timestep_end'
# outputs = csv
# []
# []
[Outputs]
csv = true
exodus = true
[]
(examples/TRISO/accident_simulation/triso1D_accident.i)
# This example is 1D spherical analysis of a TRISO fuel particle. Fully coupled
# heat transfer and solid mechanics, plus diffusion of the fission product
# species cesium (Cs) are simulated. The mesh includes contact surfaces
# between the buffer and IPyC layers to facilitate a gap opening between
# these layers. These surfaces are initially in mechanical contact but
# are assumed to have no strength in tension. A coarse mesh is used to
# provide a short run time.
# The calculation simulates fuel-life in three steps. The first step is an
# irradiation period, where constant power and a fixed particle surface
# temperature (1500 K) are assumed over a lifetime of 76 Ms (2.4 yrs).
# For the second step, fuel removal and storage are simulated by setting
# the reactor power and Cs source terms to zero, reducing the particle
# surface temperature to ambient (300 K), and then holding it
# for 100 days. A third and final step simulates accident
# behavior by increasing the particle surface temperature from ambient
# to 2073 K over 2 hrs, and then holding it at this elevated temperature
# for an additional 200 hrs. At the particle outer boundary, the Cs
# concentration is held at zero and the pressure at ambient during the
# entire simulation. The particle is assumed to be stress-free at an
# initial temperature of 1500 K.
#
# Details about this simulation are given in Section 4 of the following
# article: J. D. Hales, R. L. Williamson, S. R. Novascone, D. M. Perez,
# B. W. Spencer and G. Pastore, "Multidimensional multiphysics simulation
# of TRISO particle fuel", Journal of Nuclear Materials, Vol. 443, p. 531,
# 2013.
# This is a version using an interface kernel to model gap mass transfer.
# Sorption constants are given in Table 1 of the following article: A.
# Londono-Hurtado, I. Szlufarska, R. Bratton and D. Morgan, "A review of
# fission product sorption in carbon structures", Journal of Nuclear
# Materials, Vol. 426, p. 254, 2012.
initial_fuel_density = 11000
[GlobalParams]
order = SECOND
family = LAGRANGE
displacements = disp_x
flux_conversion_factor = 0.85
[]
[Mesh]
coord_type = RSPHERICAL
[gen] # exclude gap to establish buffer-IPyC neighbor relationships for the sorption interface kernel
type = TRISO1DMeshGenerator
elem_type = EDGE3
coordinates = '0 2.125e-4 3.125e-4 3.525e-4 3.875e-4 4.275e-4'
mesh_density = '10 5 2 2 2'
block_names = 'fuel buffer IPyC SiC OPyC'
[]
[break] # create gap between buffer and IPyC to model mechanical and thermal contact
type = BreakMeshByBlockGenerator
input = gen
block_pairs = 'buffer IPyC'
split_interface = true
add_interface_on_two_sides = true
[]
[]
[Problem]
type = ReferenceResidualProblem
reference_vector = 'ref'
extra_tag_vectors = 'ref'
[]
[Variables]
[disp_x]
[]
[temp]
initial_condition = 1500.0
[]
[conc]
initial_condition = 0.0
[]
[]
[AuxVariables]
[fission_rate]
block = fuel
order = CONSTANT
family = MONOMIAL
[]
[fluence]
order = CONSTANT
family = MONOMIAL
[]
[burnup]
block = fuel
order = CONSTANT
family = MONOMIAL
[]
[gap_cond]
order = CONSTANT
family = MONOMIAL
[]
[creep_xx]
order = CONSTANT
family = MONOMIAL
[]
[creep_yy]
order = CONSTANT
family = MONOMIAL
[]
[creep_zz]
order = CONSTANT
family = MONOMIAL
[]
[gap_HTC]
order = CONSTANT
family = MONOMIAL
[]
[gap_distance]
order = CONSTANT
family = MONOMIAL
[]
[]
[Functions]
[power_history]
type = PiecewiseLinear
x = '0 76e6 76.001e6'
y = '1 1 0'
[]
[temp_bc]
type = PiecewiseLinear
x = '0 76e6 76.001e6 84.641e6 84.6482e6'
y = '1500 1500 300 300 2073'
[]
[k_function]
type = PiecewiseLinear
x = '0 200e6'
y = '4e-37 4e-37'
[]
[d1_function]
type = ParsedFunction
expression = 'exp(t/4.5e25)'
[]
[integral_flux_error]
type = ParsedFunction
symbol_names = 'buffer_integral_flux IPyC_integral_flux'
symbol_values = 'buffer_integral_flux IPyC_integral_flux'
expression = 'IPyC_integral_flux + buffer_integral_flux'
[]
[partial_pressure_error]
type = ParsedFunction
symbol_names = 'buffer_partial_pressure IPyC_partial_pressure'
symbol_values = 'buffer_partial_pressure IPyC_partial_pressure'
expression = 'IPyC_partial_pressure - buffer_partial_pressure'
[]
[]
[Physics/SolidMechanics/QuasiStatic]
generate_output = 'stress_xx stress_yy stress_zz stress_xy stress_yz stress_zx hydrostatic_stress'
strain = FINITE
incremental = true
add_variables = false
[default]
block = 'fuel buffer IPyC OPyC'
eigenstrain_names = 'thermal_strain swelling_strain'
extra_vector_tags = 'ref'
[]
[SiC]
block = 'SiC'
eigenstrain_names = 'thermal_strain'
extra_vector_tags = 'ref'
[]
[]
[Kernels]
[heat_ie]
type = HeatConductionTimeDerivative
variable = temp
extra_vector_tags = 'ref'
[]
[heat]
type = HeatConduction
variable = temp
extra_vector_tags = 'ref'
[]
[heat_source]
type = NeutronHeatSource
variable = temp
block = fuel
energy_per_fission = 3.2e-11 # units of J/fission
fission_rate = fission_rate
extra_vector_tags = 'ref'
[]
[mass_ie]
type = TimeDerivative
variable = conc
extra_vector_tags = 'ref'
[]
[mass]
type = ArrheniusDiffusion
variable = conc
extra_vector_tags = 'ref'
[]
[mass_source]
type = BodyForce
variable = conc
function = power_history
value = 1.22e-5 # units of moles/m**3-s
block = fuel
extra_vector_tags = 'ref'
[]
[mass_decay]
type = Decay
variable = conc
radioactive_decay_constant = 7.297e-10 # units:(1/sec) The constant for Cesium
extra_vector_tags = 'ref'
[]
[]
[AuxKernels]
[fission_rate]
type = FissionRateGeneral
fission_rate_formulation = GENERIC
variable = fission_rate
block = fuel
fission_rate_function = power_history
value = 3.89e19
execute_on = timestep_begin
[]
[fluence]
type = MaterialRealAux
property = fast_neutron_fluence
variable = fluence
[]
[burnup]
type = BurnupAux
variable = burnup
block = fuel
fission_rate = fission_rate
molecular_weight = 0.270 # units of kg/mole
execute_on = timestep_begin
density = ${initial_fuel_density}
[]
[creep_xx]
type = RankTwoAux
rank_two_tensor = creep_strain
variable = creep_xx
index_i = 0
index_j = 0
block = 'buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[creep_yy]
type = RankTwoAux
rank_two_tensor = creep_strain
variable = creep_yy
index_i = 1
index_j = 1
block = 'buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[creep_zz]
type = RankTwoAux
rank_two_tensor = creep_strain
variable = creep_zz
index_i = 2
index_j = 2
block = 'buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[conductance]
type = MaterialRealAux
property = gap_conductance
variable = gap_cond
boundary = buffer_IPyC
execute_on = linear
[]
[]
[Contact]
[pellet_clad_mechanical]
primary = IPyC_buffer
secondary = buffer_IPyC
penalty = 1e5
model = frictionless
formulation = penalty
[]
[]
[ThermalContact]
[thermal_contact]
type = GasGapHeatTransfer
variable = temp
primary = IPyC_buffer
secondary = buffer_IPyC
initial_moles = initial_moles # coupling to a postprocessor which supplies the initial plenum/gap gas mass
gas_released = 'fis_gas_released co_production' # coupling to postprocessors which supply the fission gas addition, co addition
released_gas_types = 'Kr Xe;
CO'
released_fractions = '0.153 0.847;
1'
gap_geometry_type = SPHERE
tangential_tolerance = 1e-6
roughness_coef = 0.0
quadrature = true
[]
[]
[InterfaceKernels]
[cesium_gap]
type = SorptionIsothermGapInterface
variable = conc
neighbor_var = conc
partial_pressure_name = partial_pressure
sorption_penalty = 1e5
diffusivity = arrhenius_diffusion_coef
use_flux_penalty = true
flux_penalty = 1e3
boundary = buffer_IPyC
extra_vector_tags = 'ref'
[]
[]
[BCs]
# pin particle along symmetry planes
[no_disp_x]
type = DirichletBC
variable = disp_x
boundary = xzero
value = 0.0
extra_vector_tags = 'ref'
[]
# fix temperature on free surface
[freesurf_temp]
type = FunctionDirichletBC
variable = temp
boundary = exterior
function = temp_bc
extra_vector_tags = 'ref'
[]
# fix concentration on free surface
[freesurf_conc]
type = DirichletBC
variable = conc
boundary = exterior
value = 0.0
extra_vector_tags = 'ref'
[]
[PlenumPressure] # apply plenum pressure on clad inner walls and pellet surfaces
[plenumPressure]
boundary = 'buffer_IPyC IPyC_buffer'
initial_pressure = 0
startup_time = 1.0e4
R = 8.3145
output_initial_moles = initial_moles # coupling to post processor to get initial fill gas mass
temperature = ave_temp_interior # coupling to post processor to get gas temperature approximation
volume = volumeGas # coupling to post processor to get gas volume
material_input = 'fis_gas_released co_production' # coupling to post processor to get fission gas added, co added
output = plenum_pressure # coupling to post processor to output plenum/gap pressure
[]
[]
[]
[Materials]
[flux]
type = FastNeutronFlux
calculate_fluence = true
factor = 5e17
[]
[fission_gas_release] # Sifgrs fission gas release mode
type = UO2Sifgrs
block = fuel
temperature = temp
fission_rate = fission_rate # coupling to fission_rate aux variable
grain_radius_const = 5.0e-6
[]
[fuel_thermal]
type = UO2Thermal
thermal_conductivity_model = FINK_LUCUTA
block = fuel
temperature = temp
burnup = burnup
initial_porosity = 0.0
[]
[fuel_swelling]
type = UO2VolumetricSwellingEigenstrain
gas_swelling_model_type = MATPRO
block = fuel
temperature = temp
burnup = burnup
eigenstrain_name = 'swelling_strain'
initial_fuel_density = ${initial_fuel_density}
[]
[fuel_stress]
type = ComputeFiniteStrainElasticStress
block = 'fuel'
[]
[fuel_elasticity]
type = ComputeIsotropicElasticityTensor
block = fuel
youngs_modulus = 2.2e11
poissons_ratio = .345
[]
[fuel_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = fuel
thermal_expansion_coeff = 10e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[fuel_den]
type = StrainAdjustedDensity
block = fuel
strain_free_density = ${initial_fuel_density} # kg/m^3
[]
[fuel_conc]
type = ArrheniusDiffusionCoef
block = fuel
d1 = 5.6e-8 # m^2/s
q1 = 209.0e+3 # J/mol
d2 = 5.2e-4 # m^2/s
q2 = 362.0e+3 # J/mol
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[buffer_eigenstrain]
type = PyCIrradiationEigenstrain
block = buffer
pyc_type = buffer
eigenstrain_name = 'swelling_strain'
[]
[buffer_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = buffer
thermal_expansion_coeff = 5.65e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[buffer_elasticity]
type = ComputeIsotropicElasticityTensor
block = buffer
youngs_modulus = 2e10
poissons_ratio = .23
[]
[buffer_stress]
type = PyCCreep
block = buffer
temperature = temp
[]
[buffer_temp]
type = HeatConductionMaterial
block = buffer
thermal_conductivity = 0.5 # J/m-s-K
specific_heat = 720.0 # J/kg-K
[]
[buffer_den]
type = StrainAdjustedDensity
strain_free_density = 1000.0 #kg/m^3
block = buffer
[]
[buffer_conc]
type = ArrheniusDiffusionCoef
block = buffer
d1 = 1.0e-12 # m^2/s
q1 = 0.0
d2 = 0.0
q2 = 0.0
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[buffer_partial_pressure]
type = SorptionPartialPressure
A = 19.33
B = -47290
D = 1.518
E = 4338
d1 = 3.397
d2 = 6.15e-4
unit_scale = 1e3 # convert from mol to mmol
density = density # convert from mmol/m^3 to mmol/kg
concentration = conc
temperature = temp
block = 'buffer IPyC'
outputs = 'all'
output_properties = partial_pressure
[]
[normal_vectors_triso]
type = NormalVectorsTRISO
block = 'IPyC OPyC buffer'
[]
[IPyC_eigenstrain]
type = PyCIrradiationEigenstrain
block = IPyC
pyc_type = dense
eigenstrain_name = 'swelling_strain'
[]
[IPyC_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = IPyC
thermal_expansion_coeff = 5.65e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[IPyC_elasticity]
type = ComputeIsotropicElasticityTensor
block = IPyC
youngs_modulus = 4.74e10
poissons_ratio = .23
[]
[IPyC_disp]
type = PyCCreep
block = 'IPyC OPyC'
temperature = temp
[]
[IPyC_temp]
type = HeatConductionMaterial
block = 'IPyC OPyC'
thermal_conductivity = 4.0
specific_heat = 720.0
[]
[IPyC_den]
type = StrainAdjustedDensity
block = 'IPyC OPyC'
strain_free_density = 1900.0
[]
[IPyC_conc]
type = ArrheniusDiffusionCoef
block = IPyC
d1 = 6.3e-8
q1 = 222.0e+3
d2 = 0.0
q2 = 0.0
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[SiC_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = SiC
thermal_expansion_coeff = 4.9e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[SiC_elasticity]
type = ComputeIsotropicElasticityTensor
block = SiC
youngs_modulus = 3.4e11
poissons_ratio = .13
[]
[SiC_creep]
type = MonolithicSiCCreepUpdate
block = SiC
temperature = temp
k_function = k_function
[]
[SiC_stress]
type = ComputeMultipleInelasticStress
block = SiC
tangent_operator = elastic
inelastic_models = 'SiC_creep'
[]
[SiC_temp]
type = HeatConductionMaterial
block = SiC
thermal_conductivity = 13.9 # J/m-s-K
specific_heat = 620.0 # J/kg-K
[]
[SiC_den]
type = StrainAdjustedDensity
strain_free_density = 3180.0 # kg/m^3
block = SiC
[]
[SiC_conc]
type = ArrheniusDiffusionCoef
block = SiC
d1 = 5.5e-14 # m^2/s
d1_function = d1_function
d1_function_variable = fluence
q1 = 125.0e+3 # J/mol
d2 = 1.6e-2 # m^2/s
q2 = 514.0e+3 # J/mol
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[OPyC_eigenstrain]
type = PyCIrradiationEigenstrain
block = OPyC
pyc_type = dense
eigenstrain_name = 'swelling_strain'
[]
[OPyC_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = OPyC
thermal_expansion_coeff = 5.65e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[OPyC_elasticity]
type = ComputeIsotropicElasticityTensor
block = OPyC
youngs_modulus = 4.74e10
poissons_ratio = .23
[]
[OPyC_conc]
type = ArrheniusDiffusionCoef
block = OPyC
d1 = 6.3e-8 # m^2/s
q1 = 222.0e+3 # J/mol
d2 = 0.0
q2 = 0.0
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[]
[Dampers]
[temp]
type = MaxIncrement
variable = temp
max_increment = 50
[]
[]
[Debug]
show_var_residual_norms = true
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
line_search = 'none'
nl_rel_tol = 5e-4
nl_abs_tol = 1e-9
nl_max_its = 15
l_tol = 1e-3
l_max_its = 50
start_time = 0.0
end_time = 85.3682e6
dt = 100
dtmax = 2e6
dtmin = 1
[TimeStepper]
type = IterationAdaptiveDT
dt = 100
optimal_iterations = 6
growth_factor = 1.5
linear_iteration_ratio = 100
time_t = '0 76e6 76.001e6 84.641e6 84.6482e6'
time_dt = '20 20 20 20 20'
[]
[Predictor]
type = SimplePredictor
scale = 1
skip_times_old = '0 76e6 76.001e6 84.641e6 84.6482e6'
[]
[]
[Outputs]
perf_graph = true
exodus = true
[console]
type = Console
max_rows = 25
[]
[csv]
type = CSV
sync_times = '100 6308007 75696087'
sync_only = true
[]
[]
[Postprocessors]
[Cs_release]
type = SideIntegralMassFlux
variable = conc
boundary = exterior
execute_on = timestep_end
[]
[dt]
type = TimestepSize
execute_on = timestep_end
[]
[fis_gas_produced] # fission gas produced (moles)
type = ElementIntegralFisGasGeneratedSifgrs
block = fuel
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[fis_gas_released] # fission gas released to plenum (moles)
type = ElementIntegralFisGasReleasedSifgrs
block = fuel
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[volumeTotal]
type = InternalVolume
boundary = exterior
scale_factor = -1
execute_on = 'initial timestep_end'
[]
[volumeFuel]
type = InternalVolume
boundary = fuel_outer_boundary
scale_factor = -1
execute_on = 'initial timestep_end'
[]
[volumeGas]
type = InternalVolume
boundary = 'buffer_IPyC IPyC_buffer'
# ro = 3.125e-4
# ri = 2.125e-4
# vb = 4/3*pi*(ro^3-ri^3) = 8.76e-11
# buffer density = 1000
# PyC density = 1900
# fill ratio = 10/19
# vb*10/19 = 4.6e-11
# Must remove 4.6e-11 m^3 from the volume
addition = -4.6e-11
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[volumeBufferShell]
type = InternalVolume
boundary = 'buffer_IPyC IPyC_buffer'
execute_on = 'initial timestep_end'
[]
[ave_temp_interior]
type = SideAverageValue
boundary = 'buffer_IPyC IPyC_buffer'
variable = temp
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
# Postprocessors for CO production
[total_fission_rate]
type = ElementIntegralPower
variable = temp
fission_rate = fission_rate
block = fuel
energy_per_fission = 1.0
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[total_fissions]
type = TimeIntegratedPostprocessor
value = total_fission_rate
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[avg_surface_temp]
type = SideAverageValue
variable = temp
boundary = exterior
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[time_int_surf_temp]
type = TimeIntegratedPostprocessor
value = avg_surface_temp
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[co_production]
type = CarbonMonoxideProduction
total_fissions = total_fissions
time_integrated_triso_temperature = time_int_surf_temp
initial_enrichment = 0.14029
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[num_lin_it]
type = NumLinearIterations
[]
[num_nonlin_it]
type = NumNonlinearIterations
[]
[tot_lin_it]
type = CumulativeValuePostprocessor
postprocessor = num_lin_it
[]
[tot_nonlin_it]
type = CumulativeValuePostprocessor
postprocessor = num_nonlin_it
[]
[alive_time]
type = PerfGraphData
section_name = Root
data_type = TOTAL
[]
[buffer_integral_flux]
type = SideDiffusiveFluxIntegral
variable = conc
boundary = buffer_IPyC
diffusivity = arrhenius_diffusion_coef
[]
[IPyC_integral_flux]
type = SideDiffusiveFluxIntegral
variable = conc
boundary = IPyC_buffer
diffusivity = arrhenius_diffusion_coef
[]
[buffer_partial_pressure]
type = SideAverageMaterialProperty
property = partial_pressure
boundary = buffer_IPyC
[]
[IPyC_partial_pressure]
type = SideAverageMaterialProperty
property = partial_pressure
boundary = IPyC_buffer
[]
[integral_flux_error]
type = FunctionValuePostprocessor
function = integral_flux_error
[]
[partial_pressure_error]
type = FunctionValuePostprocessor
function = partial_pressure_error
[]
[integral_Cs_release]
type = TimeIntegratedPostprocessor
value = Cs_release
[]
[Cs_production]
type = FunctionValuePostprocessor
function = power_history
scale_factor = 1.22e-5 # units of moles/m**3-s
[]
[time_integral_Cs_production]
type = TimeIntegratedPostprocessor
value = Cs_production
[]
[volumeFuel_initial]
type = InternalVolume
boundary = fuel_outer_boundary
scale_factor = -1
execute_on = initial
[]
[integral_Cs_production]
type = ParsedPostprocessor
pp_names = 'time_integral_Cs_production volumeFuel_initial'
expression = 'time_integral_Cs_production * volumeFuel_initial'
[]
[Cs_release_fraction]
type = ParsedPostprocessor
pp_names = 'integral_Cs_release integral_Cs_production'
expression = 'integral_Cs_release / integral_Cs_production'
[]
[]
[VectorPostprocessors]
[temperaturevpp]
type = SideValueSampler
boundary = 11
variable = temp
sort_by = x
outputs = 'csv'
use_displaced_mesh = true
[]
[]
(examples/TRISO/parfume/parfume_un.i)
# UN TRISO particle using several PARFUME models
[GlobalParams]
order = SECOND
family = LAGRANGE
displacements = 'disp_x'
initial_enrichment = 0.1955 # [wt-]
flux_conversion_factor = 1.0 # convert E>0.10 to E>0.18 MeV
stress_free_temperature = 923.15 # used for thermal expansion
energy_per_fission = 3.204e-11 # [J/fission]
O_U = 1.5 # Initial Oxygen to Uranium atom ratio
C_U = 0.4 # Initial Carbon to Uranium atom ratio
[]
[Mesh]
coord_type = RSPHERICAL
[gen]
type = TRISO1DMeshGenerator
elem_type = EDGE3
coordinates = '0 2.125e-4 3.125e-4 3.125e-4 3.525e-4 3.875e-4 4.275e-4'
mesh_density = '6 6 0 6 8 6'
block_names = 'fuel buffer IPyC SiC OPyC'
[]
[]
[Problem]
type = ReferenceResidualProblem
reference_vector = 'ref'
extra_tag_vectors = 'ref'
[]
[UserObjects]
[particle_geometry]
type = TRISOGeometry
[]
[]
[Variables]
[temperature]
initial_condition = 923.15
[]
[conc]
initial_condition = 0.0
scaling = 1e18
[]
[]
[AuxVariables]
[disp_y]
[]
[disp_z]
[]
[fission_rate]
block = fuel
order = CONSTANT
family = MONOMIAL
[]
[burnup]
block = fuel
order = CONSTANT
family = MONOMIAL
[]
[fast_neutron_flux]
order = CONSTANT
family = MONOMIAL
[]
[fast_neutron_fluence]
order = CONSTANT
family = MONOMIAL
[]
[gap_condSlave]
order = CONSTANT
family = MONOMIAL
[]
[density]
order = CONSTANT
family = MONOMIAL
[]
[thermal_conductivity]
order = CONSTANT
family = MONOMIAL
[]
[swelling]
order = CONSTANT
family = MONOMIAL
[]
[specific_heat]
order = CONSTANT
family = MONOMIAL
[]
[volumetric_IIDC_strain]
order = CONSTANT
family = MONOMIAL
[]
[radial_IIDC_strain]
order = CONSTANT
family = MONOMIAL
[]
[tangential_IIDC_strain]
order = CONSTANT
family = MONOMIAL
[]
[BAF]
order = CONSTANT
family = MONOMIAL
[]
[fis_gas_produced]
order = CONSTANT
family = MONOMIAL
[]
[fis_gas_released]
order = CONSTANT
family = MONOMIAL
[]
[gap_HTC]
order = CONSTANT
family = MONOMIAL
[]
[gap_distance]
order = CONSTANT
family = MONOMIAL
[]
[]
[Functions]
[power_history]
type = PiecewiseLinear
x = '0 76e6 76.001e6'
y = '1 1 0'
[]
[fission_rate]
type = LinearCombinationFunction
functions = power_history
w = 3.89e19
[]
[temp_bc]
type = PiecewiseLinear
# A final temperature ramp is not possible with the UNThermal model since
# its range of applicability ends at 1800 K
# To use the model beyond its limit but get a warning, add
# value_range_behavior = WARN in the GlobalParams block.
x = '0 76e6 76.001e6 84.641e6 84.6482e6'
y = '1500 1500 300 300 2073'
[]
[k_function]
type = PiecewiseLinear
x = '0 200e6'
y = '4e-37 4e-37'
[]
[d1_function]
type = ParsedFunction
expression = 'exp(t/4.5e25)'
[]
[d_gap]
type = PiecewiseLinear
x = '1500 2100'
y = '1e-14 1e-12'
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[fuel]
block = fuel
add_variables = true
strain = FINITE
incremental = true
generate_output = 'hydrostatic_stress stress_xx stress_yy stress_zz strain_xx strain_yy strain_zz'
eigenstrain_names = 'UN_swelling_eigenstrain UN_thermal_strain'
extra_vector_tags = 'ref'
[]
[buffer]
block = buffer
add_variables = true
strain = FINITE
incremental = true
generate_output = 'hydrostatic_stress stress_xx stress_yy stress_zz strain_xx strain_yy strain_zz'
eigenstrain_names = 'buffer_IIDC_strain buffer_thermal_strain'
extra_vector_tags = 'ref'
[]
[IPyC]
block = IPyC
add_variables = true
strain = FINITE
incremental = true
generate_output = 'hydrostatic_stress stress_xx stress_yy stress_zz strain_xx strain_yy strain_zz'
eigenstrain_names = 'IPyC_IIDC_strain IPyC_thermal_strain'
extra_vector_tags = 'ref'
[]
[SiC]
block = SiC
add_variables = true
strain = FINITE
incremental = true
generate_output = 'hydrostatic_stress stress_xx stress_yy stress_zz strain_xx strain_yy strain_zz'
eigenstrain_names = 'SiC_thermal_eigenstrain'
extra_vector_tags = 'ref'
[]
[OPyC]
block = OPyC
add_variables = true
strain = FINITE
incremental = true
generate_output = 'hydrostatic_stress stress_xx stress_yy stress_zz strain_xx strain_yy strain_zz'
eigenstrain_names = 'OPyC_IIDC_strain OPyC_thermal_strain'
extra_vector_tags = 'ref'
[]
[]
[Kernels]
[heat_ie]
type = HeatConductionTimeDerivative
variable = temperature
extra_vector_tags = 'ref'
[]
[heat]
type = HeatConduction
variable = temperature
extra_vector_tags = 'ref'
[]
[heat_source]
type = NeutronHeatSource
variable = temperature
block = fuel
fission_rate = fission_rate
extra_vector_tags = 'ref'
[]
[mass_ie]
type = TimeDerivative
variable = conc
extra_vector_tags = 'ref'
[]
[mass]
type = ArrheniusDiffusion
variable = conc
extra_vector_tags = 'ref'
[]
[mass_source]
type = BodyForce
variable = conc
function = power_history
value = 1.22e-5 # units of mol/m**3-s
block = fuel
extra_vector_tags = 'ref'
[]
[mass_decay]
type = Decay
variable = conc
radioactive_decay_constant = 7.297e-10 # units:(1/sec) The constant for Cesium
extra_vector_tags = 'ref'
[]
[]
[AuxKernels]
[fissionrate]
type = MaterialRealAux
variable = fission_rate
property = fission_rate
block = fuel
execute_on = timestep_begin
[]
[burnup]
type = MaterialRealAux
variable = burnup
property = burnup
block = fuel
execute_on = timestep_begin
[]
[fast_neutron_flux]
type = MaterialRealAux
variable = fast_neutron_flux
property = fast_neutron_flux
block = 'fuel buffer IPyC SiC OPyC'
execute_on = timestep_begin
[]
[fast_neutron_fluence]
type = MaterialRealAux
variable = fast_neutron_fluence
property = fast_neutron_fluence
block = 'fuel buffer IPyC SiC OPyC'
execute_on = timestep_begin
[]
[conductanceSlave]
type = MaterialRealAux
property = gap_conductance
variable = gap_condSlave
boundary = buffer_outer_boundary
execute_on = 'initial timestep_end'
[]
[density]
type = MaterialRealAux
variable = density
property = density
block = 'fuel buffer IPyC SiC OPyC'
execute_on = 'initial linear'
[]
[thermal_conductivity]
type = MaterialRealAux
variable = thermal_conductivity
property = thermal_conductivity
block = 'fuel buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[specific_heat]
type = MaterialRealAux
variable = specific_heat
property = specific_heat
block = 'fuel buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[swelling]
type = MaterialRealAux
variable = swelling
property = swelling
block = fuel
execute_on = linear
[]
[volumetric_IIDC_strain]
type = MaterialRealAux
variable = volumetric_IIDC_strain
property = volumetric_IIDC_strain
block = 'IPyC OPyC'
execute_on = timestep_end
[]
[radial_IIDC_strain]
type = MaterialRealAux
variable = radial_IIDC_strain
property = radial_IIDC_strain
block = 'IPyC OPyC'
execute_on = timestep_end
[]
[tangential_IIDC_strain]
type = MaterialRealAux
variable = tangential_IIDC_strain
property = tangential_IIDC_strain
block = 'IPyC OPyC'
execute_on = timestep_end
[]
[BAF]
type = MaterialRealAux
variable = BAF
property = BAF
block = 'IPyC OPyC'
execute_on = timestep_end
[]
[fis_gas_produced]
type = MaterialRealAux
variable = fis_gas_produced
property = fis_gas_produced
block = fuel
execute_on = linear
[]
[fis_gas_released]
type = MaterialRealAux
variable = fis_gas_released
property = fis_gas_released
block = fuel
execute_on = linear
[]
[gap_HTC]
type = MaterialRealAux
property = gap_conductance
variable = gap_HTC
boundary = buffer_outer_boundary
execute_on = 'initial timestep_end'
[]
[gap_distance]
type = PenetrationAux
variable = gap_distance
boundary = buffer_outer_boundary
paired_boundary = IPyC_inner_boundary
quantity = distance
tangential_tolerance = 1e-6
execute_on = 'initial timestep_end'
[]
[]
[Contact]
[mechanical]
primary = IPyC_inner_boundary
secondary = buffer_outer_boundary
penalty = 1e5
model = frictionless
formulation = kinematic
[]
[]
[ThermalContact]
[thermal_contact]
type = GasGapHeatTransfer
variable = temperature
primary = IPyC_inner_boundary
secondary = buffer_outer_boundary
initial_moles = initial_moles # coupling to a postprocessor which supplies the initial plenum/gap gas mass
initial_gas_types = 'Kr Xe'
initial_fractions = '0.185 0.815'
gas_released = 'fis_gas_released'
released_gas_types = 'Kr Xe'
released_fractions = '0.185 0.815'
tangential_tolerance = 1e-6
roughness_primary = 0e-6
roughness_secondary = 0e-6
jumpdistance_primary = 0
jumpdistance_secondary = 0
quadrature = true
emissivity_secondary = 0.0
emissivity_primary = 0.0
min_gap = 1e-7
max_gap = 50e-6
gap_geometry_type = sphere
[]
[cesium_contact]
type = GapHeatTransfer
variable = conc
primary = IPyC_inner_boundary
secondary = buffer_outer_boundary
tangential_tolerance = 1e-6
gap_conductivity_function = d_gap
gap_conductivity_function_variable = temperature
appended_property_name = _conc
quadrature = true
gap_geometry_type = sphere
emissivity_primary = 0.0
emissivity_secondary = 0.0
min_gap = 1e-7
[]
[]
[BCs]
# pin particle along symmetry planes
[no_disp_x]
type = DirichletBC
variable = disp_x
boundary = xzero
value = 0.0
[]
# fix temperature on free surface
[freesurf_temp]
type = FunctionDirichletBC
variable = temperature
function = temp_bc
boundary = exterior
[]
# fix concentration on free surface
[freesurf_conc]
type = DirichletBC
variable = conc
boundary = exterior
value = 0.0
[]
# exterior and internal pressures
[exterior_pressure_x]
type = Pressure
variable = disp_x
boundary = exterior
factor = 0.1e6
[]
[PlenumPressure] # apply gas pressure on buffer and IPyC boundaries
[plenumPressure]
boundary = buffer_IPyC_boundary
initial_pressure = 100.0
startup_time = 0
R = 8.3145
output_initial_moles = initial_moles
temperature = ave_gas_temp
volume = 'gap_volume buffer_void_volume kernel_void_volume'
material_input = 'fis_gas_released'
output = gas_pressure
[]
[]
[]
[Materials]
[fission_rate]
type = GenericFunctionMaterial
prop_names = fission_rate
prop_values = fission_rate
[]
[fast_neutron_flux]
type = FastNeutronFlux
calculate_fluence = true
flux_function = power_history
factor = 5e17
[]
[normal_vectors_triso]
type = NormalVectorsTRISO
block = 'buffer IPyC OPyC'
[]
### UN properties
[UN_burnup]
type = TRISOBurnup
initial_density = 13760.0
kernel_type = UN
[]
[UN_thermal]
type = MNThermal
block = fuel
temperature = temperature
formulation = COLLIN_BAUER
[]
[UN_elasticity_tensor]
type = UNElasticityTensor
block = fuel
temperature = temperature
[]
[UN_stress]
type = ComputeFiniteStrainElasticStress
block = fuel
[]
[UN_VolumetricSwellingEigenstrain]
type = BurnupDependentEigenstrain
block = fuel
swelling_name = swelling
eigenstrain_name = UN_swelling_eigenstrain
swelling_factor = 0.8
[]
[UN_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = fuel
thermal_expansion_coeff = 10e-6 # check this value for UN
eigenstrain_name = UN_thermal_strain
temperature = temperature
[]
[UN_density]
type = StrainAdjustedDensity
block = fuel
strain_free_density = 13760.0
[]
[fission_gas_release]
type = UNFGR
block = fuel
average_grain_radius = 10e-6
temperature = temperature
triso_geometry = particle_geometry
[]
[fuel_conc]
type = ArrheniusDiffusionCoef
block = fuel
d1 = 5.6e-8 # m^2/s
q1 = 209.0e+3 # J/mol
d2 = 5.2e-4 # m^2/s
q2 = 362.0e+3 # J/mol
temperature = temperature
[]
### Buffer Properties
[buffer_elasticity_tensor]
type = BufferElasticityTensor
block = buffer
temperature = temperature
[]
[buffer_stress]
type = BufferCEGACreep
block = buffer
temperature = temperature
[]
[buffer_thermal]
type = BufferThermal
block = buffer
initial_density = 1050.0
[]
[buffer_density]
type = StrainAdjustedDensity
block = buffer
strain_free_density = 1050.0
[]
[buffer_thermal_strain]
type = BufferThermalExpansionEigenstrain
block = buffer
eigenstrain_name = buffer_thermal_strain
temperature = temperature
[]
[buffer_IIDC_strain]
type = BufferCEGAIrradiationEigenstrain
block = buffer
eigenstrain_name = buffer_IIDC_strain
temperature = temperature
[]
[buffer_conc]
type = ArrheniusDiffusionCoef
block = buffer
d1 = 1.0e-12 # m^2/s
q1 = 0.0
d2 = 0.0
q2 = 0.0
temperature = temperature
[]
### IPyC properties
[IPyC_elasticity_tensor]
type = PyCElasticityTensor
block = IPyC
temperature = temperature
initial_BAF = 1.045
poissons_ratio = 0.23
[]
[IPyC_stress]
type = PyCCEGACreep
block = IPyC
temperature = temperature
[]
[IPyC_thermal]
type = HeatConductionMaterial
block = IPyC
thermal_conductivity = 4.0
specific_heat = 720.0
[]
[IPyC_density]
type = GenericConstantMaterial
block = IPyC
prop_names = 'density'
prop_values = 1900.0
[]
[IPyC_IIDC_strain]
type = PyCCEGAIrradiationEigenstrain
block = IPyC
eigenstrain_name = IPyC_IIDC_strain
temperature = temperature
[]
[BAF_IPyC]
type = BaconAnisotropyFactor
initial_BAF = 1.045
block = IPyC
[]
[BAF_OPyC]
type = BaconAnisotropyFactor
initial_BAF = 1.045
block = OPyC
[]
[IPyC_thermal_strain]
type = PyCThermalExpansionEigenstrain
block = IPyC
eigenstrain_name = IPyC_thermal_strain
temperature = temperature
[]
[IPyC_conc]
type = ArrheniusDiffusionCoef
block = IPyC
d1 = 6.3e-8 # m^2/s
q1 = 222.0e+3 # J/mol
d2 = 0.0
q2 = 0.0
temperature = temperature
[]
### SiC properties
[SiC_elasticity_tensor]
type = MonolithicSiCElasticityTensor
block = SiC
temperature = temperature
elastic_modulus_model = miller
[]
[SiC_stress]
type = ComputeFiniteStrainElasticStress
block = SiC
[]
[SiC_thermal]
type = MonolithicSiCThermal
block = SiC
temperature = temperature
thermal_conductivity_model = miller
[]
[SiC_density]
type = GenericConstantMaterial
block = SiC
prop_names = 'density'
prop_values = 3200.0
[]
[SiC_thermal_strain]
type = MonolithicSiCThermalExpansionEigenstrain
block = SiC
temperature = temperature
eigenstrain_name = SiC_thermal_eigenstrain
[]
[SiC_conc]
type = ArrheniusDiffusionCoef
block = SiC
d1 = 5.5e-14 # m^2/s
d1_function = d1_function
d1_function_variable = fast_neutron_fluence
q1 = 125.0e+3 # J/mol
d2 = 1.6e-2 # m^2/s
q2 = 514.0e+3 # J/mol
temperature = temperature
[]
### OPyC properties
[OPyC_elasticity_tensor]
type = PyCElasticityTensor
block = OPyC
temperature = temperature
initial_BAF = 1.045
poissons_ratio = 0.23
[]
[OPyC_stress]
type = PyCCEGACreep
block = OPyC
temperature = temperature
[]
[OPyC_thermal_conductivity]
type = HeatConductionMaterial
block = OPyC
thermal_conductivity = 4.0
specific_heat = 720.0
[]
[OPyC_density]
type = GenericConstantMaterial
block = OPyC
prop_names = 'density'
prop_values = 1900.0
[]
[OPyC_IIDC_strain]
type = PyCCEGAIrradiationEigenstrain
block = OPyC
eigenstrain_name = OPyC_IIDC_strain
temperature = temperature
[]
[OPyC_thermal_strain]
type = PyCThermalExpansionEigenstrain
block = OPyC
eigenstrain_name = OPyC_thermal_strain
temperature = temperature
[]
[OPyC_conc]
type = ArrheniusDiffusionCoef
block = OPyC
d1 = 6.3e-8 # m^2/s
q1 = 222.0e+3 # J/mol
d2 = 0.0
q2 = 0.0
temperature = temperature
[]
[]
[Dampers]
[temp]
type = MaxIncrement
variable = temperature
max_increment = 50
[]
[disp_x]
type = MaxIncrement
variable = disp_x
max_increment = 1e-6
[]
[]
[Debug]
show_var_residual_norms = true
show_var_residual = 'disp_x temperature conc'
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
line_search = 'none'
nl_rel_tol = 1e-8
nl_abs_tol = 1e-7
nl_max_its = 15
l_tol = 1e-4
l_max_its = 50
start_time = 0.0
end_time = 84.641e6 #85.3682e6
num_steps = 1000
dtmax = 2e6
dtmin = 1
[TimeStepper]
type = IterationAdaptiveDT
dt = 20
growth_factor = 1.5
optimal_iterations = 8 #6
linear_iteration_ratio = 100
time_t = '0 76e6 76.001e6 84.641e6 84.6482e6'
time_dt = '20 20 20 20 20'
[]
[Quadrature]
order = THIRD
[]
[]
[Postprocessors]
[_dt]
type = TimestepSize
execute_on = timestep_end
[]
[cs_release]
type = SideIntegralMassFlux
variable = conc
boundary = exterior
[]
[int_cs_release]
type = TimeIntegratedPostprocessor
value = cs_release
[]
[cs_release_fuel]
type = SideIntegralMassFlux
variable = conc
boundary = fuel_outer_boundary
[]
[int_cs_release_fuel]
type = TimeIntegratedPostprocessor
value = cs_release_fuel
[]
[cs_release_PyCGapBndry]
type = SideIntegralMassFlux
variable = conc
boundary = IPyC_inner_boundary
[]
[int_cs_release_PyCGapBndry]
type = TimeIntegratedPostprocessor
value = cs_release_PyCGapBndry
[]
[ave_gas_temp]
type = ElementAverageValue
block = buffer
variable = temperature
execute_on = 'initial timestep_end'
[]
[ave_gap_temp]
type = SideAverageValue
boundary = buffer_IPyC_boundary
variable = temperature
execute_on = 'initial timestep_end'
[]
[fis_gas_produced]
type = ElementIntegralMaterialProperty
mat_prop = fis_gas_produced
block = fuel
execute_on = 'initial timestep_end'
[]
[fis_gas_released]
type = ElementIntegralMaterialProperty
mat_prop = fis_gas_released
block = fuel
execute_on = 'initial timestep_end'
[]
[gap_volume]
type = InternalVolume
boundary = buffer_IPyC_boundary
execute_on = 'initial linear'
use_displaced_mesh = true
[]
[buffer_void_volume]
type = VoidVolume
block = buffer
theoretical_density = 2250
execute_on = 'initial timestep_end'
use_displaced_mesh = true
[]
[kernel_void_volume]
type = VoidVolume
block = fuel
theoretical_density = 14330
execute_on = 'initial timestep_end'
use_displaced_mesh = true
[]
[volumeTotal]
type = InternalVolume
boundary = exterior
execute_on = 'initial timestep_end'
scale_factor = -1
[]
[volumeFuel]
type = InternalVolume
boundary = fuel_outer_boundary
execute_on = 'initial timestep_end'
scale_factor = -1
[]
[volumeGas]
type = InternalVolume
boundary = buffer_outer_boundary
execute_on = 'initial timestep_end'
scale_factor = -1
addition = 4.67e-11
[]
[volumeBufferShell]
type = InternalVolume
boundary = buffer_outer_boundary
execute_on = 'initial timestep_end'
[]
[ave_temp_interior]
type = SideAverageValue
boundary = buffer_outer_boundary
variable = temperature
execute_on = 'initial timestep_end'
[]
[gap_HTC]
type = ElementExtremeValue
variable = gap_HTC
block = buffer
value_type = 'max'
execute_on = 'initial timestep_end'
[]
### Postprocessors for CO production
[total_fission_rate]
type = ElementIntegralPower
variable = temperature
fission_rate = fission_rate
block = fuel
energy_per_fission = 1.0
outputs = exodus
execute_on = 'initial timestep_end'
[]
[total_fissions]
type = TimeIntegratedPostprocessor
value = total_fission_rate
outputs = exodus
execute_on = 'initial timestep_end'
[]
[avg_surface_temp]
type = SideAverageValue
variable = temperature
boundary = exterior
outputs = exodus
execute_on = 'initial timestep_end'
[]
[time_int_surf_temp]
type = TimeIntegratedPostprocessor
value = avg_surface_temp
outputs = exodus
execute_on = 'initial timestep_end'
[]
##### irradiation conditions
[particle_power]
type = ElementIntegralPower
variable = temperature
use_material_fission_rate = true
fission_rate_material = fission_rate
block = fuel
execute_on = 'initial timestep_end'
[]
[max_fluence]
type = ElementExtremeValue
variable = fast_neutron_fluence
block = 'fuel buffer IPyC SiC OPyC'
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[max_burnup]
type = ElementExtremeValue
variable = burnup
block = fuel
value_type = 'max'
execute_on = 'initial timestep_end'
[]
#### II strain
[OPyC_radial_IIDC_strain]
type = ElementExtremeValue
variable = radial_IIDC_strain
block = OPyC
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[OPyC_tangential_IIDC_strain]
type = ElementExtremeValue
variable = tangential_IIDC_strain
block = OPyC
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[IPyC_radial_IIDC_strain]
type = ElementExtremeValue
variable = radial_IIDC_strain
block = IPyC
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[IPyC_tangential_IIDC_strain]
type = ElementExtremeValue
variable = tangential_IIDC_strain
block = IPyC
value_type = 'max'
execute_on = 'initial timestep_end'
[]
#### temperatures
[max_T_kernel]
type = NodalExtremeValue
variable = temperature
block = fuel
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[max_T_buffer]
type = NodalExtremeValue
variable = temperature
block = buffer
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[min_T_buffer]
type = NodalExtremeValue
variable = temperature
block = buffer
value_type = 'min'
execute_on = 'initial timestep_end'
[]
[max_T_IPyC]
type = NodalExtremeValue
variable = temperature
block = IPyC
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[max_T_SiC]
type = NodalExtremeValue
variable = temperature
block = SiC
value_type = 'max'
execute_on = 'initial timestep_end'
[]
#### displacement BCs
[max_disp_kernel]
type = NodalExtremeValue
variable = disp_x
block = fuel
value_type = 'max'
execute_on = 'initial timestep_end'
[]
[min_disp_buffer]
type = NodalExtremeValue
variable = disp_x
block = buffer
value_type = 'min'
execute_on = 'initial timestep_end'
[]
[max_disp_IPyC]
type = NodalExtremeValue
variable = disp_x
block = IPyC
value_type = 'max'
execute_on = 'initial timestep_end'
[]
#### hoop stresses
[hoop_opyc_max]
type = ElementExtremeValue
variable = stress_yy
block = OPyC
execute_on = 'initial timestep_end'
[]
[hoop_sic_max]
type = ElementExtremeValue
variable = stress_yy
block = SiC
execute_on = 'initial timestep_end'
[]
[hoop_ipyc_max]
type = ElementExtremeValue
variable = stress_yy
block = IPyC
execute_on = 'initial timestep_end'
[]
[hoop_buffer_max]
type = ElementExtremeValue
variable = stress_yy
block = buffer
execute_on = 'initial timestep_end'
[]
[hoop_opyc_min]
type = ElementExtremeValue
variable = stress_yy
block = OPyC
value_type = min
execute_on = 'initial timestep_end'
[]
[hoop_sic_min]
type = ElementExtremeValue
variable = stress_yy
block = SiC
value_type = min
execute_on = 'initial timestep_end'
[]
[hoop_ipyc_min]
type = ElementExtremeValue
variable = stress_yy
block = IPyC
value_type = min
execute_on = 'initial timestep_end'
[]
[hoop_buffer_min]
type = ElementExtremeValue
variable = stress_yy
block = buffer
value_type = min
execute_on = 'initial timestep_end'
[]
### Check warning for Density
[oPyC_density]
type = ElementExtremeValue
variable = density
block = OPyC
execute_on = 'initial timestep_end'
[]
[sic_density]
type = ElementExtremeValue
variable = density
block = SiC
execute_on = 'initial timestep_end'
[]
[IPyC_density]
type = ElementExtremeValue
variable = density
block = IPyC
execute_on = 'initial timestep_end'
[]
[buffer_density]
type = ElementExtremeValue
variable = density
block = buffer
execute_on = 'initial timestep_end'
[]
[kernel_density]
type = ElementExtremeValue
variable = density
block = fuel
execute_on = 'initial timestep_end'
[]
[pd_penetration]
type = PdPenetration
boundary = SiC_inner_boundary
variable = temperature
execute_on = 'initial timestep_end'
[]
[]
[Outputs]
print_linear_residuals = true
time_step_interval = 1
exodus = true
csv = true
perf_graph = true
[]
(examples/TRISO/accident_simulation/triso2D_accident.i)
# This example is 2D-RZ analysis of a TRISO fuel particle. Fully coupled
# heat transfer and solid mechanics, plus diffusion of the fission product
# species cesium (Cs) are simulated. The mesh includes contact surfaces
# between the buffer and IPyC layers to facilitate a gap opening between
# these layers. These surfaces are initially in mechanical contact but
# are assumed to have no strength in tension. A coarse mesh is used to
# provide a short run time.
# The calculation simulates fuel-life in three steps. The first step is an
# irradiation period, where constant power and a fixed particle surface
# temperature (1500 K) are assumed over a lifetime of 76 Ms (2.4 yrs).
# For the second step, fuel removal and storage are simulated by setting
# the reactor power and Cs source terms to zero, reducing the particle
# surface temperature to ambient (300 K), and then holding it
# for 100 days. A third and final step simulates accident
# behavior by increasing the particle surface temperature from ambient
# to 2073 K over 2 hrs, and then holding it at this elevated temperature
# for an additional 200 hrs. At the particle outer boundary, the Cs
# concentration is held at zero and the pressure at ambient during the
# entire simulation. The particle is assumed to be stress-free at an
# initial temperature of 1500 K.
#
# Details about this simulation are given in Section 4 of the following
# article: J. D. Hales, R. L. Williamson, S. R. Novascone, D. M. Perez,
# B. W. Spencer and G. Pastore, "Multidimensional multiphysics simulation
# of TRISO particle fuel", Journal of Nuclear Materials, Vol. 443, p. 531,
# 2013.
initial_fuel_density = 11000.0
[GlobalParams]
order = SECOND
family = LAGRANGE
displacements = 'disp_x disp_y'
flux_conversion_factor = 0.85
[]
[Mesh]
coord_type = RZ
[mesh]
type = FileMeshGenerator
file = triso2Dmed.e
[]
[]
[Problem]
type = ReferenceResidualProblem
reference_vector = 'ref'
extra_tag_vectors = 'ref'
[]
[Variables]
[disp_x]
[]
[disp_y]
[]
[temp]
initial_condition = 1500.0
[]
[conc]
initial_condition = 0.0
[]
[]
[AuxVariables]
[fission_rate]
block = fuel
order = CONSTANT
family = MONOMIAL
[]
[fluence]
order = CONSTANT
family = MONOMIAL
[]
[burnup]
block = fuel
order = CONSTANT
family = MONOMIAL
[]
[gap_condSlave]
order = CONSTANT
family = MONOMIAL
[]
[creep_xx]
order = CONSTANT
family = MONOMIAL
[]
[creep_yy]
order = CONSTANT
family = MONOMIAL
[]
[creep_zz]
order = CONSTANT
family = MONOMIAL
[]
[]
[Functions]
[power_history]
type = PiecewiseLinear
x = '0 76e6 76.001e6'
y = '1 1 0'
[]
[temp_bc]
type = PiecewiseLinear
x = '0 76e6 76.001e6 84.641e6 84.6482e6'
y = '1500 1500 300 300 2073'
[]
[k_function]
type = PiecewiseLinear
x = '0 200e6'
y = '4e-37 4e-37'
[]
[d1_function]
type = ParsedFunction
expression = 'exp(t/4.5e25)'
[]
[d_gap]
type = PiecewiseLinear
x = '1500 2100'
y = '1e-14 1e-12'
[]
[integral_flux_error]
type = ParsedFunction
symbol_names = 'buffer_integral_flux IPyC_integral_flux'
symbol_values = 'buffer_integral_flux IPyC_integral_flux'
expression = 'IPyC_integral_flux + buffer_integral_flux'
[]
[]
[Physics/SolidMechanics/QuasiStatic]
generate_output = 'stress_xx stress_yy stress_zz stress_xy stress_yz stress_zx hydrostatic_stress'
strain = FINITE
incremental = true
add_variables = false
[default]
block = 'fuel buffer IPyC OPyC'
eigenstrain_names = 'thermal_strain swelling_strain'
extra_vector_tags = 'ref'
[]
[SiC]
block = 'SiC'
eigenstrain_names = 'thermal_strain'
extra_vector_tags = 'ref'
[]
[]
[Kernels]
[heat_ie]
type = HeatConductionTimeDerivative
variable = temp
extra_vector_tags = 'ref'
[]
[heat]
type = HeatConduction
variable = temp
extra_vector_tags = 'ref'
[]
[heat_source]
type = NeutronHeatSource
variable = temp
block = fuel
energy_per_fission = 3.2e-11 # units of J/fission
fission_rate = fission_rate
extra_vector_tags = 'ref'
[]
[mass_ie]
type = TimeDerivative
variable = conc
extra_vector_tags = 'ref'
[]
[mass]
type = ArrheniusDiffusion
variable = conc
extra_vector_tags = 'ref'
[]
[mass_source]
type = BodyForce
variable = conc
function = power_history
value = 1.22e-5 # units of moles/m**3-s
block = fuel
extra_vector_tags = 'ref'
[]
[mass_decay]
type = Decay
variable = conc
radioactive_decay_constant = 7.297e-10 # units:(1/sec) The constant for Cesium
extra_vector_tags = 'ref'
[]
[]
[AuxKernels]
[fission_rate]
type = FissionRateGeneral
fission_rate_formulation = GENERIC
variable = fission_rate
block = fuel
fission_rate_function = power_history
value = 3.89e19
execute_on = timestep_begin
[]
[fluence]
type = MaterialRealAux
property = fast_neutron_fluence
variable = fluence
[]
[burnup]
type = BurnupAux
variable = burnup
block = fuel
fission_rate = fission_rate
molecular_weight = 0.270 # units of kg/mole
execute_on = timestep_begin
density = ${initial_fuel_density}
[]
[creep_xx]
type = RankTwoAux
rank_two_tensor = creep_strain
variable = creep_xx
index_i = 0
index_j = 0
block = 'buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[creep_yy]
type = RankTwoAux
rank_two_tensor = creep_strain
variable = creep_yy
index_i = 1
index_j = 1
block = 'buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[creep_zz]
type = RankTwoAux
rank_two_tensor = creep_strain
variable = creep_zz
index_i = 2
index_j = 2
block = 'buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[conductanceSlave]
type = MaterialRealAux
property = gap_conductance
variable = gap_condSlave
boundary = BufferGapBndry
execute_on = linear
[]
[]
[Contact]
[pellet_clad_mechanical]
primary = 15
secondary = 17
penalty = 1e5
model = frictionless
formulation = penalty
[]
[]
[ThermalContact]
[thermal_contact]
type = GasGapHeatTransfer
variable = temp
primary = 15
secondary = 17
initial_moles = initial_moles # coupling to a postprocessor which supplies the initial plenum/gap gas mass
gas_released = 'fis_gas_released co_production' # coupling to postprocessors which supply the fission gas addition, co addition
released_gas_types = 'Kr Xe;
CO'
released_fractions = '0.153 0.847;
1'
gap_geometry_type = CYLINDER
tangential_tolerance = 1e-6
roughness_coef = 0.0
quadrature = true
[]
[cesium_contact]
type = GapHeatTransfer
variable = conc
primary = 15
secondary = 17
tangential_tolerance = 1e-6
gap_conductivity_function = d_gap
gap_conductivity_function_variable = temp
appended_property_name = _conc
emissivity_primary = 0
emissivity_secondary = 0
quadrature = true
[]
[]
[BCs]
# pin particle along symmetry planes
[no_disp_x]
type = DirichletBC
variable = disp_x
boundary = xzero
value = 0.0
extra_vector_tags = 'ref'
[]
[no_disp_y]
type = DirichletBC
variable = disp_y
boundary = yzero
value = 0.0
extra_vector_tags = 'ref'
[]
# fix temperature on free surface
[freesurf_temp]
type = FunctionDirichletBC
variable = temp
boundary = exterior
function = temp_bc
extra_vector_tags = 'ref'
[]
# fix concentration on free surface
[freesurf_conc]
type = DirichletBC
variable = conc
boundary = exterior
value = 0.0
extra_vector_tags = 'ref'
[]
[PlenumPressure] # apply plenum pressure on clad inner walls and pellet surfaces
[plenumPressure]
boundary = BufferGapVol
initial_pressure = 0
startup_time = 1.0e4
R = 8.3145
output_initial_moles = initial_moles # coupling to post processor to get initial fill gas mass
temperature = ave_temp_interior # coupling to post processor to get gas temperature approximation
volume = volumeGas # coupling to post processor to get gas volume
material_input = 'fis_gas_released co_production' # coupling to post processor to get fission gas added, co added
output = plenum_pressure # coupling to post processor to output plenum/gap pressure
[]
[]
[]
[Materials]
[flux]
type = FastNeutronFlux
calculate_fluence = true
factor = 5e17
[]
[fission_gas_release] # Sifgrs fission gas release mode
type = UO2Sifgrs
block = fuel
temperature = temp
fission_rate = fission_rate # coupling to fission_rate aux variable
grain_radius_const = 5.0e-6
[]
[fuel_thermal]
type = UO2Thermal
thermal_conductivity_model = FINK_LUCUTA
block = fuel
temperature = temp
burnup = burnup
initial_porosity = 0.0
[]
[fuel_swelling]
type = UO2VolumetricSwellingEigenstrain
gas_swelling_model_type = MATPRO
block = fuel
temperature = temp
burnup = burnup
eigenstrain_name = 'swelling_strain'
initial_fuel_density = ${initial_fuel_density}
[]
[fuel_stress]
type = ComputeFiniteStrainElasticStress
block = 'fuel'
[]
[fuel_elasticity]
type = ComputeIsotropicElasticityTensor
block = fuel
youngs_modulus = 2.2e11
poissons_ratio = .345
[]
[fuel_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = fuel
thermal_expansion_coeff = 10e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[fuel_den]
type = StrainAdjustedDensity
block = fuel
strain_free_density = ${initial_fuel_density} # kg/m^3
[]
[fuel_conc]
type = ArrheniusDiffusionCoef
block = fuel
d1 = 5.6e-8 # m^2/s
q1 = 209.0e+3 # J/mol
d2 = 5.2e-4 # m^2/s
q2 = 362.0e+3 # J/mol
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[buffer_eigenstrain]
type = PyCIrradiationEigenstrain
block = buffer
pyc_type = buffer
eigenstrain_name = 'swelling_strain'
[]
[buffer_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = buffer
thermal_expansion_coeff = 5.65e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[buffer_elasticity]
type = ComputeIsotropicElasticityTensor
block = buffer
youngs_modulus = 2e10
poissons_ratio = .23
[]
[buffer_stress]
type = PyCCreep
block = buffer
temperature = temp
[]
[buffer_temp]
type = HeatConductionMaterial
block = buffer
thermal_conductivity = 0.5 # J/m-s-K
specific_heat = 720.0 # J/kg-K
[]
[buffer_den]
type = StrainAdjustedDensity
strain_free_density = 1000.0 #kg/m^3
block = buffer
[]
[buffer_conc]
type = ArrheniusDiffusionCoef
block = buffer
d1 = 1.0e-12 # m^2/s
q1 = 0.0
d2 = 0.0
q2 = 0.0
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[normal_vectors_triso]
type = NormalVectorsTRISO
block = 'IPyC OPyC buffer'
[]
[IPyC_eigenstrain]
type = PyCIrradiationEigenstrain
block = IPyC
pyc_type = dense
eigenstrain_name = 'swelling_strain'
[]
[IPyC_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = IPyC
thermal_expansion_coeff = 5.65e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[IPyC_elasticity]
type = ComputeIsotropicElasticityTensor
block = IPyC
youngs_modulus = 4.74e10
poissons_ratio = .23
[]
[IPyC_disp]
type = PyCCreep
block = 'IPyC OPyC'
temperature = temp
[]
[IPyC_temp]
type = HeatConductionMaterial
block = 'IPyC OPyC'
thermal_conductivity = 4.0
specific_heat = 720.0
[]
[IPyC_den]
type = StrainAdjustedDensity
block = 'IPyC OPyC'
strain_free_density = 1900.0
[]
[IPyC_conc]
type = ArrheniusDiffusionCoef
block = IPyC
d1 = 6.3e-8
q1 = 222.0e+3
d2 = 0.0
q2 = 0.0
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[SiC_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = SiC
thermal_expansion_coeff = 4.9e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[SiC_elasticity]
type = ComputeIsotropicElasticityTensor
block = SiC
youngs_modulus = 3.4e11
poissons_ratio = .13
[]
[SiC_creep]
type = MonolithicSiCCreepUpdate
block = SiC
temperature = temp
k_function = k_function
[]
[SiC_stress]
type = ComputeMultipleInelasticStress
block = SiC
tangent_operator = elastic
inelastic_models = 'SiC_creep'
[]
[SiC_temp]
type = HeatConductionMaterial
block = SiC
thermal_conductivity = 13.9 # J/m-s-K
specific_heat = 620.0 # J/kg-K
[]
[SiC_den]
type = StrainAdjustedDensity
strain_free_density = 3180.0 # kg/m^3
block = SiC
[]
[SiC_conc]
type = ArrheniusDiffusionCoef
block = SiC
d1 = 5.5e-14 # m^2/s
d1_function = d1_function
d1_function_variable = fluence
q1 = 125.0e+3 # J/mol
d2 = 1.6e-2 # m^2/s
q2 = 514.0e+3 # J/mol
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[OPyC_eigenstrain]
type = PyCIrradiationEigenstrain
block = OPyC
pyc_type = dense
eigenstrain_name = 'swelling_strain'
[]
[OPyC_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = OPyC
thermal_expansion_coeff = 5.65e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[OPyC_elasticity]
type = ComputeIsotropicElasticityTensor
block = OPyC
youngs_modulus = 4.74e10
poissons_ratio = .23
[]
[OPyC_conc]
type = ArrheniusDiffusionCoef
block = OPyC
d1 = 6.3e-8 # m^2/s
q1 = 222.0e+3 # J/mol
d2 = 0.0
q2 = 0.0
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[]
[Dampers]
[temp]
type = MaxIncrement
variable = temp
max_increment = 50
[]
[]
[Debug]
show_var_residual_norms = true
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
line_search = 'none'
nl_rel_tol = 5e-4
nl_abs_tol = 1e-9
nl_max_its = 15
l_tol = 1e-3
l_max_its = 50
start_time = 0.0
end_time = 85.3682e6
dt = 100
dtmax = 2e6
dtmin = 1
automatic_scaling = true
compute_scaling_once = false
scaling_group_variables = 'conc; disp_x disp_y; temp'
[TimeStepper]
type = IterationAdaptiveDT
dt = 100
optimal_iterations = 6
growth_factor = 1.5
linear_iteration_ratio = 100
time_t = '0 76e6 76.001e6 84.641e6 84.6482e6'
time_dt = '20 20 20 20 20'
[]
[Predictor]
type = SimplePredictor
scale = 1
skip_times_old = '0 76e6 76.001e6 84.641e6 84.6482e6'
[]
[Quadrature]
order = THIRD
side_order = FIFTH
[]
[]
[Outputs]
perf_graph = true
exodus = true
[console]
type = Console
max_rows = 25
[]
[csv]
type = CSV
sync_times = '100 6308007 75696087'
sync_only = true
[]
[]
[Postprocessors]
[Cs_release]
type = SideIntegralMassFlux
variable = conc
boundary = exterior
execute_on = timestep_end
[]
[dt]
type = TimestepSize
execute_on = timestep_end
[]
[fis_gas_produced] # fission gas produced (moles)
type = ElementIntegralFisGasGeneratedSifgrs
block = fuel
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[fis_gas_released] # fission gas released to plenum (moles)
type = ElementIntegralFisGasReleasedSifgrs
block = fuel
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[volumeTotal]
type = InternalVolume
boundary = exterior
execute_on = 'initial timestep_end'
[]
[volumeFuel]
type = InternalVolume
boundary = fuel
execute_on = 'initial timestep_end'
[]
[volumeGas]
type = InternalVolume
boundary = BufferGapVol
# ro = 3.125e-4
# ri = 2.125e-4
# vb = 4/3*pi*(ro^3-ri^3) = 8.76e-11
# buffer density = 1000
# PyC density = 1900
# fill ratio = 10/19
# vb*10/19 = 4.6e-11
# Must remove 4.6e-11 m^3 from the volume
addition = -4.6e-11
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[volumeBufferShell]
type = InternalVolume
boundary = BufferGapVol
execute_on = 'initial timestep_end'
[]
[ave_temp_interior]
type = SideAverageValue
boundary = BufferGapVol
variable = temp
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
# Postprocessors for CO production
[total_fission_rate]
type = ElementIntegralPower
variable = temp
fission_rate = fission_rate
block = fuel
energy_per_fission = 1.0
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[total_fissions]
type = TimeIntegratedPostprocessor
value = total_fission_rate
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[avg_surface_temp]
type = SideAverageValue
variable = temp
boundary = exterior
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[time_int_surf_temp]
type = TimeIntegratedPostprocessor
value = avg_surface_temp
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[co_production]
type = CarbonMonoxideProduction
total_fissions = total_fissions
time_integrated_triso_temperature = time_int_surf_temp
initial_enrichment = 0.14029
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[num_lin_it]
type = NumLinearIterations
[]
[num_nonlin_it]
type = NumNonlinearIterations
[]
[tot_lin_it]
type = CumulativeValuePostprocessor
postprocessor = num_lin_it
[]
[tot_nonlin_it]
type = CumulativeValuePostprocessor
postprocessor = num_nonlin_it
[]
[alive_time]
type = PerfGraphData
section_name = Root
data_type = TOTAL
[]
[buffer_avg_conc]
type = SideAverageValue
variable = conc
boundary = 17
[]
[IPyC_avg_conc]
type = SideAverageValue
variable = conc
boundary = 15
[]
[buffer_integral_flux]
type = SideIntegralMassFlux
variable = conc
boundary = 17
[]
[IPyC_integral_flux]
type = SideIntegralMassFlux
variable = conc
boundary = 15
[]
[integral_flux_error]
type = FunctionValuePostprocessor
function = integral_flux_error
[]
[integral_Cs_release]
type = TimeIntegratedPostprocessor
value = Cs_release
[]
[Cs_production]
type = FunctionValuePostprocessor
function = power_history
scale_factor = 1.22e-5 # units of moles/m**3-s
[]
[time_integral_Cs_production]
type = TimeIntegratedPostprocessor
value = Cs_production
[]
[volumeFuel_initial]
type = InternalVolume
boundary = fuel
execute_on = initial
[]
[integral_Cs_production]
type = ParsedPostprocessor
pp_names = 'time_integral_Cs_production volumeFuel_initial'
expression = 'time_integral_Cs_production * volumeFuel_initial'
[]
[Cs_release_fraction]
type = ParsedPostprocessor
pp_names = 'integral_Cs_release integral_Cs_production'
expression = 'integral_Cs_release / integral_Cs_production'
[]
[]
[VectorPostprocessors]
[temperaturevpp]
type = SideValueSampler
boundary = 11
variable = temp
sort_by = x
outputs = 'csv'
use_displaced_mesh = true
[]
[]
(examples/TRISO/accident_simulation/triso2D_accident_mortar.i)
# This example is 2D-RZ analysis of a TRISO fuel particle. Fully coupled
# heat transfer and solid mechanics, plus diffusion of the fission product
# species cesium (Cs) are simulated. The mesh includes contact surfaces
# between the buffer and IPyC layers to facilitate a gap opening between
# these layers. These surfaces are initially in mechanical contact but
# are assumed to have no strength in tension. A coarse mesh is used to
# provide a short run time.
# The calculation simulates fuel-life in three steps. The first step is an
# irradiation period, where constant power and a fixed particle surface
# temperature (1500 K) are assumed over a lifetime of 76 Ms (2.4 yrs).
# For the second step, fuel removal and storage are simulated by setting
# the reactor power and Cs source terms to zero, reducing the particle
# surface temperature to ambient (300 K), and then holding it
# for 100 days. A third and final step simulates accident
# behavior by increasing the particle surface temperature from ambient
# to 2073 K over 2 hrs, and then holding it at this elevated temperature
# for an additional 200 hrs. At the particle outer boundary, the Cs
# concentration is held at zero and the pressure at ambient during the
# entire simulation. The particle is assumed to be stress-free at an
# initial temperature of 1500 K.
#
# Details about this simulation are given in Section 4 of the following
# article: J. D. Hales, R. L. Williamson, S. R. Novascone, D. M. Perez,
# B. W. Spencer and G. Pastore, "Multidimensional multiphysics simulation
# of TRISO particle fuel", Journal of Nuclear Materials, Vol. 443, p. 531,
# 2013.
# This is a version using a thermomechanical mortar approach.
initial_fuel_density = 11000.0
[GlobalParams]
order = SECOND
family = LAGRANGE
displacements = 'disp_x disp_y'
flux_conversion_factor = 0.85
[]
[Mesh]
coord_type = RZ
[file]
type = FileMeshGenerator
file = triso2Dmed.e
[]
[]
[Problem]
type = ReferenceResidualProblem
reference_vector = 'ref'
extra_tag_vectors = 'ref'
converge_on = 'disp_x disp_y temp conc'
[]
[Variables]
[disp_x]
[]
[disp_y]
[]
[temp]
initial_condition = 1500.0
[]
[conc]
initial_condition = 0.0
[]
[]
[AuxVariables]
[fission_rate]
block = fuel
order = CONSTANT
family = MONOMIAL
[]
[fluence]
order = CONSTANT
family = MONOMIAL
[]
[burnup]
block = fuel
order = CONSTANT
family = MONOMIAL
[]
[creep_xx]
order = CONSTANT
family = MONOMIAL
[]
[creep_yy]
order = CONSTANT
family = MONOMIAL
[]
[creep_zz]
order = CONSTANT
family = MONOMIAL
[]
[]
[Functions]
[power_history]
type = PiecewiseLinear
x = '0 76e6 76.001e6'
y = '1 1 0'
[]
[temp_bc]
type = PiecewiseLinear
x = '0 76e6 76.001e6 84.641e6 84.6482e6'
y = '1500 1500 300 300 2073'
[]
[k_function]
type = PiecewiseLinear
x = '0 200e6'
y = '4e-37 4e-37'
[]
[d1_function]
type = ParsedFunction
expression = 'exp(t/4.5e25)'
[]
[d_gap]
type = PiecewiseLinear
x = '1500 2100'
y = '1e-14 1e-12'
[]
[integral_flux_error]
type = ParsedFunction
symbol_names = 'buffer_integral_flux IPyC_integral_flux'
symbol_values = 'buffer_integral_flux IPyC_integral_flux'
expression = 'IPyC_integral_flux + buffer_integral_flux'
[]
[]
[Physics/SolidMechanics/QuasiStatic]
generate_output = 'stress_xx stress_yy stress_zz stress_xy stress_yz stress_zx hydrostatic_stress'
strain = FINITE
incremental = true
add_variables = false
[default]
block = 'fuel buffer IPyC OPyC'
eigenstrain_names = 'thermal_strain swelling_strain'
extra_vector_tags = 'ref'
[]
[SiC]
block = 'SiC'
eigenstrain_names = 'thermal_strain'
extra_vector_tags = 'ref'
[]
[]
[Kernels]
[heat_ie]
type = HeatConductionTimeDerivative
variable = temp
extra_vector_tags = 'ref'
block = 'fuel buffer IPyC SiC OPyC'
[]
[heat]
type = HeatConduction
variable = temp
extra_vector_tags = 'ref'
block = 'fuel buffer IPyC SiC OPyC'
[]
[heat_source]
type = NeutronHeatSource
variable = temp
block = fuel
energy_per_fission = 3.2e-11 # units of J/fission
fission_rate = fission_rate
extra_vector_tags = 'ref'
[]
[mass_ie]
type = TimeDerivative
variable = conc
extra_vector_tags = 'ref'
block = 'fuel buffer IPyC SiC OPyC'
[]
[mass]
type = ArrheniusDiffusion
variable = conc
extra_vector_tags = 'ref'
block = 'fuel buffer IPyC SiC OPyC'
[]
[mass_source]
type = BodyForce
variable = conc
function = power_history
value = 1.22e-5 # units of moles/m**3-s
block = fuel
extra_vector_tags = 'ref'
[]
[mass_decay]
type = Decay
variable = conc
radioactive_decay_constant = 7.297e-10 # units:(1/sec) The constant for Cesium
block = 'fuel buffer IPyC SiC OPyC'
extra_vector_tags = 'ref'
[]
[]
[AuxKernels]
[fission_rate]
type = FissionRateGeneral
fission_rate_formulation = GENERIC
variable = fission_rate
block = fuel
fission_rate_function = power_history
value = 3.89e19
execute_on = timestep_begin
[]
[fluence]
type = MaterialRealAux
property = fast_neutron_fluence
variable = fluence
[]
[burnup]
type = BurnupAux
variable = burnup
block = fuel
fission_rate = fission_rate
molecular_weight = 0.270 # units of kg/mole
execute_on = timestep_begin
density = ${initial_fuel_density}
[]
[creep_xx]
type = RankTwoAux
rank_two_tensor = creep_strain
variable = creep_xx
index_i = 0
index_j = 0
block = 'buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[creep_yy]
type = RankTwoAux
rank_two_tensor = creep_strain
variable = creep_yy
index_i = 1
index_j = 1
block = 'buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[creep_zz]
type = RankTwoAux
rank_two_tensor = creep_strain
variable = creep_zz
index_i = 2
index_j = 2
block = 'buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[]
[ThermalContactMortar]
[thermal]
secondary_variable = temp
primary_boundary = 15
secondary_boundary = 17
initial_moles = initial_moles # coupling to a postprocessor which supplies the initial plenum/gap gas mass
gas_released = 'fis_gas_released co_production' # coupling to postprocessors which supply the fission gas addition, co addition
released_gas_types = 'Kr Xe;
CO'
released_fractions = '0.153 0.847;
1'
gap_geometry_type = CYLINDER
min_gap = 1e-7
max_gap = 50e-6
roughness_coef = 0.0
correct_edge_dropping = true
[]
[]
[Contact]
[pellet_clad_mechanical]
primary = 15
secondary = 17
model = frictionless
formulation = mortar
c_normal = 1.0e8
correct_edge_dropping = true
[]
[]
[ThermalContact]
[cesium_contact]
type = GapHeatTransfer
variable = conc
primary = 15
secondary = 17
tangential_tolerance = 1e-6
gap_conductivity_function = d_gap
gap_conductivity_function_variable = temp
appended_property_name = _conc
emissivity_primary = 0
emissivity_secondary = 0
quadrature = true
[]
[]
[BCs]
# pin particle along symmetry planes
[no_disp_x]
type = DirichletBC
variable = disp_x
boundary = xzero
value = 0.0
extra_vector_tags = 'ref'
[]
[no_disp_y]
type = DirichletBC
variable = disp_y
boundary = yzero
value = 0.0
extra_vector_tags = 'ref'
[]
# fix temperature on free surface
[freesurf_temp]
type = FunctionDirichletBC
variable = temp
boundary = exterior
function = temp_bc
extra_vector_tags = 'ref'
[]
# fix concentration on free surface
[freesurf_conc]
type = DirichletBC
variable = conc
boundary = exterior
value = 0.0
extra_vector_tags = 'ref'
[]
[PlenumPressure] # apply plenum pressure on clad inner walls and pellet surfaces
[plenumPressure]
boundary = BufferGapVol
initial_pressure = 0
startup_time = 1.0e4
R = 8.3145
output_initial_moles = initial_moles # coupling to post processor to get initial fill gas mass
temperature = ave_temp_interior # coupling to post processor to get gas temperature approximation
volume = volumeGas # coupling to post processor to get gas volume
material_input = 'fis_gas_released co_production' # coupling to post processor to get fission gas added, co added
output = plenum_pressure # coupling to post processor to output plenum/gap pressure
[]
[]
[]
[Materials]
[flux]
type = FastNeutronFlux
calculate_fluence = true
factor = 5e17
[]
[fission_gas_release] # Sifgrs fission gas release mode
type = UO2Sifgrs
block = fuel
temperature = temp
fission_rate = fission_rate # coupling to fission_rate aux variable
grain_radius_const = 5.0e-6
[]
[fuel_thermal]
type = UO2Thermal
thermal_conductivity_model = FINK_LUCUTA
block = fuel
temperature = temp
burnup = burnup
initial_porosity = 0.0
[]
[fuel_swelling]
type = UO2VolumetricSwellingEigenstrain
gas_swelling_model_type = MATPRO
block = fuel
temperature = temp
burnup = burnup
eigenstrain_name = 'swelling_strain'
initial_fuel_density = ${initial_fuel_density}
[]
[fuel_stress]
type = ComputeFiniteStrainElasticStress
block = 'fuel'
[]
[fuel_elasticity]
type = ComputeIsotropicElasticityTensor
block = fuel
youngs_modulus = 2.2e11
poissons_ratio = .345
[]
[fuel_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = fuel
thermal_expansion_coeff = 10e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[fuel_den]
type = StrainAdjustedDensity
block = fuel
strain_free_density = ${initial_fuel_density} # kg/m^3
[]
[fuel_conc]
type = ArrheniusDiffusionCoef
block = fuel
d1 = 5.6e-8 # m^2/s
q1 = 209.0e+3 # J/mol
d2 = 5.2e-4 # m^2/s
q2 = 362.0e+3 # J/mol
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[buffer_eigenstrain]
type = PyCIrradiationEigenstrain
block = buffer
pyc_type = buffer
eigenstrain_name = 'swelling_strain'
[]
[buffer_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = buffer
thermal_expansion_coeff = 5.65e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[buffer_elasticity]
type = ComputeIsotropicElasticityTensor
block = buffer
youngs_modulus = 2e10
poissons_ratio = .23
[]
[buffer_stress]
type = PyCCreep
block = buffer
temperature = temp
[]
[buffer_temp]
type = HeatConductionMaterial
block = buffer
thermal_conductivity = 0.5 # J/m-s-K
specific_heat = 720.0 # J/kg-K
[]
[buffer_den]
type = StrainAdjustedDensity
strain_free_density = 1000.0 #kg/m^3
block = buffer
[]
[buffer_conc]
type = ArrheniusDiffusionCoef
block = buffer
d1 = 1.0e-12 # m^2/s
q1 = 0.0
d2 = 0.0
q2 = 0.0
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[normal_vectors_triso]
type = NormalVectorsTRISO
block = 'IPyC OPyC buffer'
[]
[IPyC_eigenstrain]
type = PyCIrradiationEigenstrain
block = IPyC
pyc_type = dense
eigenstrain_name = 'swelling_strain'
[]
[IPyC_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = IPyC
thermal_expansion_coeff = 5.65e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[IPyC_elasticity]
type = ComputeIsotropicElasticityTensor
block = IPyC
youngs_modulus = 4.74e10
poissons_ratio = .23
[]
[IPyC_disp]
type = PyCCreep
block = 'IPyC OPyC'
temperature = temp
[]
[IPyC_temp]
type = HeatConductionMaterial
block = 'IPyC OPyC'
thermal_conductivity = 4.0
specific_heat = 720.0
[]
[IPyC_den]
type = StrainAdjustedDensity
block = 'IPyC OPyC'
strain_free_density = 1900.0
[]
[IPyC_conc]
type = ArrheniusDiffusionCoef
block = IPyC
d1 = 6.3e-8
q1 = 222.0e+3
d2 = 0.0
q2 = 0.0
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[SiC_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = SiC
thermal_expansion_coeff = 4.9e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[SiC_elasticity]
type = ComputeIsotropicElasticityTensor
block = SiC
youngs_modulus = 3.4e11
poissons_ratio = .13
[]
[SiC_creep]
type = MonolithicSiCCreepUpdate
block = SiC
temperature = temp
k_function = k_function
[]
[SiC_stress]
type = ComputeMultipleInelasticStress
block = SiC
tangent_operator = elastic
inelastic_models = 'SiC_creep'
[]
[SiC_temp]
type = HeatConductionMaterial
block = SiC
thermal_conductivity = 13.9 # J/m-s-K
specific_heat = 620.0 # J/kg-K
[]
[SiC_den]
type = StrainAdjustedDensity
strain_free_density = 3180.0 # kg/m^3
block = SiC
[]
[SiC_conc]
type = ArrheniusDiffusionCoef
block = SiC
d1 = 5.5e-14 # m^2/s
d1_function = d1_function
d1_function_variable = fluence
q1 = 125.0e+3 # J/mol
d2 = 1.6e-2 # m^2/s
q2 = 514.0e+3 # J/mol
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[OPyC_eigenstrain]
type = PyCIrradiationEigenstrain
block = OPyC
pyc_type = dense
eigenstrain_name = 'swelling_strain'
[]
[OPyC_thermal_strain]
type = ComputeThermalExpansionEigenstrain
block = OPyC
thermal_expansion_coeff = 5.65e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[OPyC_elasticity]
type = ComputeIsotropicElasticityTensor
block = OPyC
youngs_modulus = 4.74e10
poissons_ratio = .23
[]
[OPyC_conc]
type = ArrheniusDiffusionCoef
block = OPyC
d1 = 6.3e-8 # m^2/s
q1 = 222.0e+3 # J/mol
d2 = 0.0
q2 = 0.0
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[]
[Dampers]
[temp]
type = MaxIncrement
variable = temp
max_increment = 50
[]
[]
[Debug]
show_var_residual_norms = true
[]
[Preconditioning]
[SMP]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'NEWTON'
petsc_options = '-snes_converged_reason -ksp_converged_reason -snes_ksp_ew'
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package -mat_mffd_err -pc_factor_shift_type -pc_factor_shift_amount'
petsc_options_value = 'lu superlu_dist 1e-5 NONZERO 1e-14'
snesmf_reuse_base = false
line_search = 'none'
nl_rel_tol = 5e-4
nl_abs_tol = 1e-10
nl_max_its = 20
l_max_its = 8
start_time = 0.0
end_time = 85.3682e6
dt = 100
dtmax = 2e6
dtmin = 1
[TimeStepper]
type = IterationAdaptiveDT
dt = 100
optimal_iterations = 10
growth_factor = 1.5
linear_iteration_ratio = 100
time_t = '0 76e6 76.001e6 84.641e6 84.6482e6'
time_dt = '20 20 20 20 20'
[]
[Predictor]
type = SimplePredictor
scale = 0.5
skip_times_old = '0 76e6 76.001e6 84.641e6 84.6482e6'
[]
[]
[Outputs]
perf_graph = true
exodus = true
[console]
type = Console
max_rows = 25
[]
[csv]
type = CSV
sync_times = '100 6308007 75696087'
sync_only = true
[]
[]
[Postprocessors]
[Cs_release]
type = SideIntegralMassFlux
variable = conc
boundary = exterior
execute_on = timestep_end
[]
[dt]
type = TimestepSize
execute_on = timestep_end
[]
[fis_gas_produced] # fission gas produced (moles)
type = ElementIntegralFisGasGeneratedSifgrs
block = fuel
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[fis_gas_released] # fission gas released to plenum (moles)
type = ElementIntegralFisGasReleasedSifgrs
block = fuel
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[volumeTotal]
type = InternalVolume
boundary = exterior
execute_on = 'initial timestep_end'
[]
[volumeFuel]
type = InternalVolume
boundary = fuel
execute_on = 'initial timestep_end'
[]
[volumeGas]
type = InternalVolume
boundary = BufferGapVol
# ro = 3.125e-4
# ri = 2.125e-4
# vb = 4/3*pi*(ro^3-ri^3) = 8.76e-11
# buffer density = 1000
# PyC density = 1900
# fill ratio = 10/19
# vb*10/19 = 4.6e-11
# Must remove 4.6e-11 m^3 from the volume
addition = -4.6e-11
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[volumeBufferShell]
type = InternalVolume
boundary = BufferGapVol
execute_on = 'initial timestep_end'
[]
[ave_temp_interior]
type = SideAverageValue
boundary = BufferGapVol
variable = temp
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
# Postprocessors for CO production
[total_fission_rate]
type = ElementIntegralPower
variable = temp
fission_rate = fission_rate
block = fuel
energy_per_fission = 1.0
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[total_fissions]
type = TimeIntegratedPostprocessor
value = total_fission_rate
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[avg_surface_temp]
type = SideAverageValue
variable = temp
boundary = exterior
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[time_int_surf_temp]
type = TimeIntegratedPostprocessor
value = avg_surface_temp
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[co_production]
type = CarbonMonoxideProduction
total_fissions = total_fissions
time_integrated_triso_temperature = time_int_surf_temp
initial_enrichment = 0.14029
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[num_lin_it]
type = NumLinearIterations
[]
[num_nonlin_it]
type = NumNonlinearIterations
[]
[tot_lin_it]
type = CumulativeValuePostprocessor
postprocessor = num_lin_it
[]
[tot_nonlin_it]
type = CumulativeValuePostprocessor
postprocessor = num_nonlin_it
[]
[alive_time]
type = PerfGraphData
section_name = Root
data_type = TOTAL
[]
[buffer_avg_conc]
type = SideAverageValue
variable = conc
boundary = 17
[]
[IPyC_avg_conc]
type = SideAverageValue
variable = conc
boundary = 15
[]
[buffer_integral_flux]
type = SideIntegralMassFlux
variable = conc
boundary = 17
[]
[IPyC_integral_flux]
type = SideIntegralMassFlux
variable = conc
boundary = 15
[]
[integral_flux_error]
type = FunctionValuePostprocessor
function = integral_flux_error
[]
[integral_Cs_release]
type = TimeIntegratedPostprocessor
value = Cs_release
[]
[Cs_production]
type = FunctionValuePostprocessor
function = power_history
scale_factor = 1.22e-5 # units of moles/m**3-s
[]
[time_integral_Cs_production]
type = TimeIntegratedPostprocessor
value = Cs_production
[]
[volumeFuel_initial]
type = InternalVolume
boundary = fuel
execute_on = initial
[]
[integral_Cs_production]
type = ParsedPostprocessor
pp_names = 'time_integral_Cs_production volumeFuel_initial'
expression = 'time_integral_Cs_production * volumeFuel_initial'
[]
[Cs_release_fraction]
type = ParsedPostprocessor
pp_names = 'integral_Cs_release integral_Cs_production'
expression = 'integral_Cs_release / integral_Cs_production'
[]
[]
[VectorPostprocessors]
[temperaturevpp]
type = SideValueSampler
boundary = 11
variable = temp
sort_by = x
outputs = 'csv'
use_displaced_mesh = true
[]
[]
(examples/TRISO/accident_simulation/triso2D_accident_ad.i)
# This example is 2D-RZ analysis of a TRISO fuel particle. Fully coupled
# heat transfer and solid mechanics, plus diffusion of the fission product
# species cesium (Cs) are simulated. The mesh includes contact surfaces
# between the buffer and IPyC layers to facilitate a gap opening between
# these layers. These surfaces are initially in mechanical contact but
# are assumed to have no strength in tension. A coarse mesh is used to
# provide a short run time.
# The calculation simulates fuel-life in three steps. The first step is an
# irradiation period, where constant power and a fixed particle surface
# temperature (1500 K) are assumed over a lifetime of 76 Ms (2.4 yrs).
# For the second step, fuel removal and storage are simulated by setting
# the reactor power and Cs source terms to zero, reducing the particle
# surface temperature to ambient (300 K), and then holding it
# for 100 days. A third and final step simulates accident
# behavior by increasing the particle surface temperature from ambient
# to 2073 K over 2 hrs, and then holding it at this elevated temperature
# for an additional 200 hrs. At the particle outer boundary, the Cs
# concentration is held at zero and the pressure at ambient during the
# entire simulation. The particle is assumed to be stress-free at an
# initial temperature of 1500 K.
#
# Details about this simulation are given in Section 4 of the following
# article: J. D. Hales, R. L. Williamson, S. R. Novascone, D. M. Perez,
# B. W. Spencer and G. Pastore, "Multidimensional multiphysics simulation
# of TRISO particle fuel", Journal of Nuclear Materials, Vol. 443, p. 531,
# 2013.
# This is a version using a thermomechanical mortar approach. It uses
# Automatic Differentiation classes and models gap mass transfer using
# flux preserving and sorption mortar constraints. Sorption constants are
# given in Table 1 of the following article: A. Londono-Hurtado, I.
# Szlufarska, R. Bratton and D. Morgan, "A review of fission product
# sorption in carbon structures", Journal of Nuclear Materials, Vol. 426,
# p. 254, 2012.
initial_fuel_density = 11000.0
[GlobalParams]
order = SECOND
family = LAGRANGE
displacements = 'disp_x disp_y'
flux_conversion_factor = 0.85
use_automatic_differentiation = true
[]
[Mesh]
coord_type = RZ
[file]
type = FileMeshGenerator
file = triso2Dmed.e
[]
[]
[Problem]
type = ReferenceResidualProblem
reference_vector = 'ref'
extra_tag_vectors = 'ref'
converge_on = 'disp_x disp_y temp conc'
[]
[Variables]
[disp_x]
[]
[disp_y]
[]
[temp]
initial_condition = 1500.0
[]
[conc]
initial_condition = 0.0
[]
[conc_lm]
block = pellet_clad_mechanical_secondary_subdomain
[]
[conc_dx_lm]
block = pellet_clad_mechanical_secondary_subdomain
[]
[conc_dy_lm]
block = pellet_clad_mechanical_secondary_subdomain
[]
[]
[AuxVariables]
[fission_rate]
block = fuel
order = CONSTANT
family = MONOMIAL
[]
[fluence]
order = CONSTANT
family = MONOMIAL
[]
[burnup]
block = fuel
order = CONSTANT
family = MONOMIAL
[]
[creep_xx]
order = CONSTANT
family = MONOMIAL
[]
[creep_yy]
order = CONSTANT
family = MONOMIAL
[]
[creep_zz]
order = CONSTANT
family = MONOMIAL
[]
[]
[Functions]
[power_history]
type = PiecewiseLinear
x = '0 76e6 76.001e6'
y = '1 1 0'
[]
[temp_bc]
type = PiecewiseLinear
x = '0 76e6 76.001e6 84.641e6 84.6482e6'
y = '1500 1500 300 300 2073'
[]
[k_function]
type = PiecewiseLinear
x = '0 200e6'
y = '4e-37 4e-37'
[]
[d1_function]
type = ParsedFunction
expression = 'exp(t/4.5e25)'
[]
[integral_flux_error]
type = ParsedFunction
symbol_names = 'buffer_integral_flux IPyC_integral_flux'
symbol_values = 'buffer_integral_flux IPyC_integral_flux'
expression = 'IPyC_integral_flux + buffer_integral_flux'
[]
[partial_pressure_error]
type = ParsedFunction
symbol_names = 'buffer_partial_pressure IPyC_partial_pressure'
symbol_values = 'buffer_partial_pressure IPyC_partial_pressure'
expression = 'IPyC_partial_pressure - buffer_partial_pressure'
[]
[]
[Physics/SolidMechanics/QuasiStatic]
generate_output = 'stress_xx stress_yy stress_zz stress_xy stress_yz stress_zx hydrostatic_stress'
strain = FINITE
incremental = true
add_variables = false
[default]
block = 'fuel buffer IPyC OPyC'
eigenstrain_names = 'thermal_strain swelling_strain'
extra_vector_tags = 'ref'
[]
[SiC]
block = 'SiC'
eigenstrain_names = 'thermal_strain'
extra_vector_tags = 'ref'
[]
[]
[Kernels]
[heat_ie]
type = ADHeatConductionTimeDerivative
variable = temp
extra_vector_tags = 'ref'
block = 'fuel buffer IPyC SiC OPyC'
[]
[heat]
type = ADHeatConduction
variable = temp
extra_vector_tags = 'ref'
block = 'fuel buffer IPyC SiC OPyC'
[]
[heat_source]
type = ADNeutronHeatSource
variable = temp
block = fuel
energy_per_fission = 3.2e-11 # units of J/fission
fission_rate = fission_rate
extra_vector_tags = 'ref'
[]
[mass_ie]
type = ADTimeDerivative
variable = conc
extra_vector_tags = 'ref'
block = 'fuel buffer IPyC SiC OPyC'
[]
[mass]
type = ADArrheniusDiffusion
variable = conc
extra_vector_tags = 'ref'
block = 'fuel buffer IPyC SiC OPyC'
[]
[mass_source]
type = ADBodyForce
variable = conc
function = power_history
value = 1.22e-5 # units of moles/m**3-s
block = fuel
extra_vector_tags = 'ref'
[]
[mass_decay]
type = Decay
variable = conc
radioactive_decay_constant = 7.297e-10 # units:(1/sec) The constant for Cesium
block = 'fuel buffer IPyC SiC OPyC'
extra_vector_tags = 'ref'
[]
[]
[AuxKernels]
[fission_rate]
type = FissionRateGeneral
fission_rate_formulation = GENERIC
variable = fission_rate
block = fuel
fission_rate_function = power_history
value = 3.89e19
execute_on = timestep_begin
[]
[fluence]
type = ADMaterialRealAux
property = fast_neutron_fluence
variable = fluence
[]
[burnup]
type = ADBurnupAux
variable = burnup
block = fuel
fission_rate = fission_rate
molecular_weight = 0.270 # units of kg/mole
execute_on = timestep_begin
density = ${initial_fuel_density}
[]
[creep_xx]
type = ADRankTwoAux
rank_two_tensor = creep_strain
variable = creep_xx
index_i = 0
index_j = 0
block = 'buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[creep_yy]
type = ADRankTwoAux
rank_two_tensor = creep_strain
variable = creep_yy
index_i = 1
index_j = 1
block = 'buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[creep_zz]
type = ADRankTwoAux
rank_two_tensor = creep_strain
variable = creep_zz
index_i = 2
index_j = 2
block = 'buffer IPyC SiC OPyC'
execute_on = timestep_end
[]
[]
[ThermalContactMortar]
[thermal]
secondary_variable = temp
primary_boundary = 15
secondary_boundary = 17
initial_moles = initial_moles # coupling to a postprocessor which supplies the initial plenum/gap gas mass
gas_released = 'fis_gas_released co_production' # coupling to postprocessors which supply the fission gas addition, co addition
released_gas_types = 'Kr Xe;
CO'
released_fractions = '0.153 0.847;
1'
gap_geometry_type = CYLINDER
min_gap = 1e-7
max_gap = 50e-6
roughness_coef = 0.0
correct_edge_dropping = true
[]
[]
[Contact]
[pellet_clad_mechanical]
primary = 15
secondary = 17
model = frictionless
formulation = mortar
c_normal = 1.0e8
correct_edge_dropping = true
[]
[]
[Constraints]
[cesium_gap_value]
type = MassSorptionConstraint
variable = conc_lm
primary_variable = conc
primary_boundary = 15
primary_subdomain = pellet_clad_mechanical_primary_subdomain
secondary_variable = conc
secondary_boundary = 17
secondary_subdomain = pellet_clad_mechanical_secondary_subdomain
partial_pressure_name = partial_pressure
epsilon = 1e-4
correct_edge_dropping = true
[]
[cesium_gap_flux_x]
type = MassFluxConstraint
variable = conc_dx_lm
primary_variable = conc
diffusivity_primary = arrhenius_diffusion_coef
primary_boundary = 15
primary_subdomain = pellet_clad_mechanical_primary_subdomain
secondary_variable = conc
diffusivity_secondary = arrhenius_diffusion_coef
secondary_boundary = 17
secondary_subdomain = pellet_clad_mechanical_secondary_subdomain
component = 0
epsilon = 1e-5
correct_edge_dropping = true
[]
[cesium_gap_flux_y]
type = MassFluxConstraint
variable = conc_dy_lm
primary_variable = conc
diffusivity_primary = arrhenius_diffusion_coef
primary_boundary = 15
primary_subdomain = pellet_clad_mechanical_primary_subdomain
secondary_variable = conc
diffusivity_secondary = arrhenius_diffusion_coef
secondary_boundary = 17
secondary_subdomain = pellet_clad_mechanical_secondary_subdomain
component = 1
epsilon = 1e-5
correct_edge_dropping = true
[]
[]
[BCs]
# pin particle along symmetry planes
[no_disp_x]
type = ADDirichletBC
variable = disp_x
boundary = xzero
value = 0.0
extra_vector_tags = 'ref'
[]
[no_disp_y]
type = ADDirichletBC
variable = disp_y
boundary = yzero
value = 0.0
extra_vector_tags = 'ref'
[]
# fix temperature on free surface
[freesurf_temp]
type = ADFunctionDirichletBC
variable = temp
boundary = exterior
function = temp_bc
extra_vector_tags = 'ref'
[]
# fix concentration on free surface
[freesurf_conc]
type = ADDirichletBC
variable = conc
boundary = exterior
value = 0.0
extra_vector_tags = 'ref'
[]
[PlenumPressure] # apply plenum pressure on clad inner walls and pellet surfaces
[plenumPressure]
boundary = BufferGapVol
initial_pressure = 0
startup_time = 1.0e4
R = 8.3145
output_initial_moles = initial_moles # coupling to post processor to get initial fill gas mass
temperature = ave_temp_interior # coupling to post processor to get gas temperature approximation
volume = volumeGas # coupling to post processor to get gas volume
material_input = 'fis_gas_released co_production' # coupling to post processor to get fission gas added, co added
output = plenum_pressure # coupling to post processor to output plenum/gap pressure
[]
[]
[]
[Materials]
[flux]
type = ADFastNeutronFlux
calculate_fluence = true
factor = 5e17
[]
[fission_gas_release] # Sifgrs fission gas release mode
type = ADUO2Sifgrs
block = fuel
temperature = temp
fission_rate = fission_rate # coupling to fission_rate aux variable
grain_radius_const = 5.0e-6
[]
[fuel_thermal]
type = ADUO2Thermal
thermal_conductivity_model = FINK_LUCUTA
block = fuel
temperature = temp
burnup = burnup
initial_porosity = 0.0
[]
[fuel_swelling]
type = ADUO2VolumetricSwellingEigenstrain
gas_swelling_model_type = MATPRO
block = fuel
temperature = temp
burnup = burnup
eigenstrain_name = 'swelling_strain'
initial_fuel_density = ${initial_fuel_density}
[]
[fuel_stress]
type = ADComputeFiniteStrainElasticStress
block = 'fuel'
[]
[fuel_elasticity]
type = ADComputeIsotropicElasticityTensor
block = fuel
youngs_modulus = 2.2e11
poissons_ratio = .345
[]
[fuel_thermal_strain]
type = ADComputeThermalExpansionEigenstrain
block = fuel
thermal_expansion_coeff = 10e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[fuel_den]
type = ADStrainAdjustedDensity
block = fuel
strain_free_density = ${initial_fuel_density} # kg/m^3
[]
[fuel_conc]
type = ADArrheniusDiffusionCoef
block = fuel
d1 = 5.6e-8 # m^2/s
q1 = 209.0e+3 # J/mol
d2 = 5.2e-4 # m^2/s
q2 = 362.0e+3 # J/mol
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[buffer_eigenstrain]
type = ADPyCIrradiationEigenstrain
block = buffer
pyc_type = buffer
eigenstrain_name = 'swelling_strain'
[]
[buffer_thermal_strain]
type = ADComputeThermalExpansionEigenstrain
block = buffer
thermal_expansion_coeff = 5.65e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[buffer_elasticity]
type = ADComputeIsotropicElasticityTensor
block = buffer
youngs_modulus = 2e10
poissons_ratio = .23
[]
[buffer_stress]
type = ADPyCCreep
block = buffer
temperature = temp
[]
[buffer_temp]
type = ADHeatConductionMaterial
block = buffer
thermal_conductivity = 0.5 # J/m-s-K
specific_heat = 720.0 # J/kg-K
[]
[buffer_den]
type = ADStrainAdjustedDensity
strain_free_density = 1000.0 #kg/m^3
block = buffer
[]
[buffer_conc]
type = ADArrheniusDiffusionCoef
block = buffer
d1 = 1.0e-12 # m^2/s
q1 = 0.0
d2 = 0.0
q2 = 0.0
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[buffer_partial_pressure]
type = ADSorptionPartialPressure
A = 19.33
B = -47290
D = 1.518
E = 4338
d1 = 3.397
d2 = 6.15e-4
unit_scale = 1e3 # convert from mol to mmol
density = 1000 # convert from mmol/m^3 to mmol/kg, using constant for compatibility with default AD derivative container size
concentration = conc
temperature = temp
block = buffer
outputs = 'all'
output_properties = partial_pressure
[]
[normal_vectors_triso]
type = NormalVectorsTRISO
block = 'IPyC OPyC buffer'
[]
[IPyC_eigenstrain]
type = ADPyCIrradiationEigenstrain
block = IPyC
pyc_type = dense
eigenstrain_name = 'swelling_strain'
[]
[IPyC_thermal_strain]
type = ADComputeThermalExpansionEigenstrain
block = IPyC
thermal_expansion_coeff = 5.65e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[IPyC_elasticity]
type = ADComputeIsotropicElasticityTensor
block = IPyC
youngs_modulus = 4.74e10
poissons_ratio = .23
[]
[IPyC_disp]
type = ADPyCCreep
block = 'IPyC OPyC'
temperature = temp
[]
[IPyC_temp]
type = ADHeatConductionMaterial
block = 'IPyC OPyC'
thermal_conductivity = 4.0
specific_heat = 720.0
[]
[IPyC_den]
type = ADStrainAdjustedDensity
block = 'IPyC OPyC'
strain_free_density = 1900.0
[]
[IPyC_conc]
type = ADArrheniusDiffusionCoef
block = IPyC
d1 = 6.3e-8
q1 = 222.0e+3
d2 = 0.0
q2 = 0.0
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[IPyC_partial_pressure]
type = ADSorptionPartialPressure
A = 19.33
B = -47290
D = 1.518
E = 4338
d1 = 3.397
d2 = 6.15e-4
unit_scale = 1e3 # convert from mol to mmol
density = 1900 # convert from mmol/m^3 to mmol/kg, using constant for compatibility with default AD derivative container size
concentration = conc
temperature = temp
block = IPyC
outputs = 'all'
output_properties = partial_pressure
[]
[SiC_thermal_strain]
type = ADComputeThermalExpansionEigenstrain
block = SiC
thermal_expansion_coeff = 4.9e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[SiC_elasticity]
type = ADComputeIsotropicElasticityTensor
block = SiC
youngs_modulus = 3.4e11
poissons_ratio = .13
[]
[SiC_creep]
type = ADMonolithicSiCCreepUpdate
block = SiC
temperature = temp
k_function = k_function
[]
[SiC_stress]
type = ADComputeMultipleInelasticStress
block = SiC
inelastic_models = 'SiC_creep'
[]
[SiC_temp]
type = ADHeatConductionMaterial
block = SiC
thermal_conductivity = 13.9 # J/m-s-K
specific_heat = 620.0 # J/kg-K
[]
[SiC_den]
type = ADStrainAdjustedDensity
strain_free_density = 3180.0 # kg/m^3
block = SiC
[]
[SiC_conc]
type = ADArrheniusDiffusionCoef
block = SiC
d1 = 5.5e-14 # m^2/s
d1_function = d1_function
d1_function_variable = fluence
q1 = 125.0e+3 # J/mol
d2 = 1.6e-2 # m^2/s
q2 = 514.0e+3 # J/mol
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[OPyC_eigenstrain]
type = ADPyCIrradiationEigenstrain
block = OPyC
pyc_type = dense
eigenstrain_name = 'swelling_strain'
[]
[OPyC_thermal_strain]
type = ADComputeThermalExpansionEigenstrain
block = OPyC
thermal_expansion_coeff = 5.65e-6
stress_free_temperature = 1500.0
eigenstrain_name = thermal_strain
temperature = temp
[]
[OPyC_elasticity]
type = ADComputeIsotropicElasticityTensor
block = OPyC
youngs_modulus = 4.74e10
poissons_ratio = .23
[]
[OPyC_conc]
type = ADArrheniusDiffusionCoef
block = OPyC
d1 = 6.3e-8 # m^2/s
q1 = 222.0e+3 # J/mol
d2 = 0.0
q2 = 0.0
gas_constant = 8.3143 # J/K-mol
temperature = temp
[]
[]
[Dampers]
[temp]
type = MaxIncrement
variable = temp
max_increment = 50
[]
[]
[Debug]
show_var_residual_norms = true
[]
[Preconditioning]
[SMP]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'NEWTON'
petsc_options = '-snes_converged_reason -ksp_converged_reason -snes_ksp_ew'
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package -mat_mffd_err -pc_factor_shift_type -pc_factor_shift_amount'
petsc_options_value = 'lu superlu_dist 1e-5 NONZERO 1e-14'
snesmf_reuse_base = false
line_search = 'none'
nl_rel_tol = 5e-4
nl_abs_tol = 1e-10
nl_max_its = 20
l_max_its = 8
start_time = 0.0
end_time = 85.3682e6
dt = 100
dtmax = 2e6
dtmin = 1
[TimeStepper]
type = IterationAdaptiveDT
dt = 100
optimal_iterations = 10
growth_factor = 1.5
linear_iteration_ratio = 100
time_t = '0 76e6 76.001e6 84.641e6 84.6482e6'
time_dt = '20 20 20 20 20'
[]
[Predictor]
type = SimplePredictor
scale = 0.5
skip_times_old = '0 76e6 76.001e6 84.641e6 84.6482e6'
[]
[]
[Outputs]
perf_graph = true
exodus = true
[console]
type = Console
max_rows = 25
[]
[csv]
type = CSV
sync_times = '100 6308007 75696087'
sync_only = true
[]
[]
[Postprocessors]
[Cs_release]
type = ADSideDiffusiveFluxIntegral
variable = conc
diffusivity = arrhenius_diffusion_coef
boundary = exterior
execute_on = timestep_end
[]
[dt]
type = TimestepSize
execute_on = timestep_end
[]
[fis_gas_produced] # fission gas produced (moles)
type = ADElementIntegralFisGasGeneratedSifgrs
block = fuel
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[fis_gas_released] # fission gas released to plenum (moles)
type = ADElementIntegralFisGasReleasedSifgrs
block = fuel
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[volumeTotal]
type = InternalVolume
boundary = exterior
execute_on = 'initial timestep_end'
[]
[volumeFuel]
type = InternalVolume
boundary = fuel
execute_on = 'initial timestep_end'
[]
[volumeGas]
type = InternalVolume
boundary = BufferGapVol
# ro = 3.125e-4
# ri = 2.125e-4
# vb = 4/3*pi*(ro^3-ri^3) = 8.76e-11
# buffer density = 1000
# PyC density = 1900
# fill ratio = 10/19
# vb*10/19 = 4.6e-11
# Must remove 4.6e-11 m^3 from the volume
addition = -4.6e-11
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[volumeBufferShell]
type = InternalVolume
boundary = BufferGapVol
execute_on = 'initial timestep_end'
[]
[ave_temp_interior]
type = SideAverageValue
boundary = BufferGapVol
variable = temp
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
# Postprocessors for CO production
[total_fission_rate]
type = ElementIntegralPower
variable = temp
fission_rate = fission_rate
block = fuel
energy_per_fission = 1.0
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[total_fissions]
type = TimeIntegratedPostprocessor
value = total_fission_rate
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[avg_surface_temp]
type = SideAverageValue
variable = temp
boundary = exterior
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[time_int_surf_temp]
type = TimeIntegratedPostprocessor
value = avg_surface_temp
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[co_production]
type = CarbonMonoxideProduction
total_fissions = total_fissions
time_integrated_triso_temperature = time_int_surf_temp
initial_enrichment = 0.14029
execute_on = 'initial linear nonlinear timestep_begin timestep_end'
[]
[num_lin_it]
type = NumLinearIterations
[]
[num_nonlin_it]
type = NumNonlinearIterations
[]
[tot_lin_it]
type = CumulativeValuePostprocessor
postprocessor = num_lin_it
[]
[tot_nonlin_it]
type = CumulativeValuePostprocessor
postprocessor = num_nonlin_it
[]
[alive_time]
type = PerfGraphData
section_name = Root
data_type = TOTAL
[]
[buffer_integral_flux]
type = ADSideDiffusiveFluxIntegral
variable = conc
boundary = 17
diffusivity = arrhenius_diffusion_coef
[]
[IPyC_integral_flux]
type = ADSideDiffusiveFluxIntegral
variable = conc
boundary = 15
diffusivity = arrhenius_diffusion_coef
[]
[buffer_partial_pressure]
type = ADSideAverageMaterialProperty
property = partial_pressure
boundary = 17
[]
[IPyC_partial_pressure]
type = ADSideAverageMaterialProperty
property = partial_pressure
boundary = 15
[]
[integral_flux_error]
type = FunctionValuePostprocessor
function = integral_flux_error
[]
[partial_pressure_error]
type = FunctionValuePostprocessor
function = partial_pressure_error
[]
[integral_Cs_release]
type = TimeIntegratedPostprocessor
value = Cs_release
[]
[Cs_production]
type = FunctionValuePostprocessor
function = power_history
scale_factor = 1.22e-5 # units of moles/m**3-s
[]
[time_integral_Cs_production]
type = TimeIntegratedPostprocessor
value = Cs_production
[]
[volumeFuel_initial]
type = InternalVolume
boundary = fuel
execute_on = initial
[]
[integral_Cs_production]
type = ParsedPostprocessor
pp_names = 'time_integral_Cs_production volumeFuel_initial'
expression = 'time_integral_Cs_production * volumeFuel_initial'
[]
[Cs_release_fraction]
type = ParsedPostprocessor
pp_names = 'integral_Cs_release integral_Cs_production'
expression = 'integral_Cs_release / integral_Cs_production'
[]
[]
[VectorPostprocessors]
[temperaturevpp]
type = SideValueSampler
boundary = 11
variable = temp
sort_by = x
outputs = 'csv'
use_displaced_mesh = true
[]
[]
(test/tests/radioactive_decay/radioactive_decay.i)
# Tests radioactive decay kernel Decay.C against an analytical solution. The ODE is
#
# dc/dt = - lambda*c
#
# which has the solution
#
# c = c0 * exp(-lambda*t)
#
# Numerical solution is for a single linear element on a unit cube having an initial
# concentration of 1. No BC's are provided, which defaults to zero mass flux on all
# faces. To minimize time error, the solution is run for 1000 steps out to a time of 1.
#
# Analytical result: c(1) = 0.367879
# Numerical result: c(1) = 0.368063
#
[Mesh]
[mesh]
type = FileMeshGenerator
file = 1x1x1_cube.e
[]
[]
[Variables]
[c]
order = FIRST
family = LAGRANGE
initial_condition = 1.0
[]
[]
[Kernels]
[decay]
type = Decay
variable = c
radioactive_decay_constant = 1.0
[]
[ie]
type = TimeDerivative
variable = c
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
start_time = 0.0
num_steps = 1000
dt = .001
[]
[Outputs]
exodus = true
[]
(test/tests/triso/base_irradiation/triso1D_accident.i)
initial_fuel_density = 11000.0
[GlobalParams]
density = ${initial_fuel_density} # kg/m^3
order = SECOND
displacements = 'disp_x'
[]
[Mesh]
coord_type = RSPHERICAL
[mesh]
type = FileMeshGenerator
file = triso1DFineTruss3.e
[]
[]
[Problem]
type = ReferenceResidualProblem
reference_vector = 'ref'
extra_tag_vectors = 'ref'
[]
[Variables]
[disp_x]
[]
[temperature]
initial_condition = 1500.0
[]
[conc_Cs]
initial_condition = 0.0
scaling = 1e18
[]
[]
[AuxVariables]
[fission_rate]
block = fuel
order = CONSTANT
family = MONOMIAL
[]
[burnup]
block = fuel
order = CONSTANT
family = MONOMIAL
[]
[fast_neutron_flux]
order = CONSTANT
family = MONOMIAL
[]
[fast_neutron_fluence]
order = CONSTANT
family = MONOMIAL
[]
[gap_condSlave]
order = CONSTANT
family = MONOMIAL
[]
[]
[Functions]
[power_history]
type = PiecewiseLinear
x = '0 76e6 76.001e6'
y = '1 1 0'
[]
[fission_rate]
type = LinearCombinationFunction
functions = power_history
w = 3.89e19
[]
[temp_bc]
type = PiecewiseLinear
x = '0 76e6 76.001e6 84.641e6 84.6482e6'
y = '1500 1500 300 300 2073'
[]
[k_function]
type = PiecewiseLinear
x = '0 200e6'
y = '4e-37 4e-37'
[]
[d1_function]
type = ParsedFunction
expression = 'exp(t/4.5e25)'
[]
[d_gap]
type = PiecewiseLinear
x = '1500 2100'
y = '1e-14 1e-12'
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[fuel]
block = fuel
add_variables = false
strain = FINITE
eigenstrain_names = 'fuel_thermal_strain fuel_swelling'
generate_output = 'hydrostatic_stress stress_xx stress_yy stress_zz strain_xx strain_yy strain_zz'
extra_vector_tags = 'ref'
[]
[buffer]
block = buffer
add_variables = false
strain = FINITE
eigenstrain_names = 'buffer_thermal_strain buffer_eigenstrain'
generate_output = 'hydrostatic_stress stress_xx stress_yy stress_zz strain_xx strain_yy strain_zz creep_strain_xx creep_strain_yy creep_strain_zz'
extra_vector_tags = 'ref'
[]
[IPyC]
block = IPyC
add_variables = false
strain = FINITE
eigenstrain_names = 'IPyC_eigenstrain IPyC_thermal_strain'
generate_output = 'hydrostatic_stress stress_xx stress_yy stress_zz strain_xx strain_yy strain_zz creep_strain_xx creep_strain_yy creep_strain_zz'
extra_vector_tags = 'ref'
[]
[SiC]
block = SiC
add_variables = false
strain = FINITE
eigenstrain_names = 'SiC_thermal_strain'
generate_output = 'hydrostatic_stress stress_xx stress_yy stress_zz strain_xx strain_yy strain_zz creep_strain_xx creep_strain_yy creep_strain_zz'
extra_vector_tags = 'ref'
[]
[OPyC]
block = OPyC
add_variables = false
strain = FINITE
eigenstrain_names = 'OPyC_eigenstrain OPyC_thermal_strain'
generate_output = 'hydrostatic_stress stress_xx stress_yy stress_zz strain_xx strain_yy strain_zz creep_strain_xx creep_strain_yy creep_strain_zz'
extra_vector_tags = 'ref'
[]
[]
[Kernels]
[heat_ie]
type = HeatConductionTimeDerivative
variable = temperature
extra_vector_tags = 'ref'
[]
[heat]
type = HeatConduction
variable = temperature
extra_vector_tags = 'ref'
[]
[heat_source]
type = NeutronHeatSource
variable = temperature
block = fuel
energy_per_fission = 3.2e-11 # units of J/fission
fission_rate = fission_rate
extra_vector_tags = 'ref'
[]
[mass_ie]
type = TimeDerivative
variable = conc_Cs
extra_vector_tags = 'ref'
[]
[mass]
type = ArrheniusDiffusion
variable = conc_Cs
extra_vector_tags = 'ref'
[]
[mass_source]
type = BodyForce
variable = conc_Cs
function = power_history
value = 1.22e-5 # units of mol/m**3-s
block = fuel
extra_vector_tags = 'ref'
[]
[mass_decay]
type = Decay
variable = conc_Cs
radioactive_decay_constant = 7.297e-10 # units:(1/sec) The constant for Cesium
extra_vector_tags = 'ref'
[]
[]
[AuxKernels]
[fissionrate]
type = MaterialRealAux
variable = fission_rate
property = fission_rate
block = fuel
execute_on = timestep_begin
[]
[burnup]
type = BurnupAux
variable = burnup
fission_rate = fission_rate
molecular_weight = 0.270 # units of kg/mol
block = fuel
execute_on = timestep_begin
[]
[fast_neutron_flux]
type = MaterialRealAux
variable = fast_neutron_flux
property = fast_neutron_flux
block = 'fuel buffer IPyC SiC OPyC'
execute_on = timestep_begin
[]
[fast_neutron_fluence]
type = MaterialRealAux
variable = fast_neutron_fluence
property = fast_neutron_fluence
block = 'fuel buffer IPyC SiC OPyC'
execute_on = timestep_begin
[]
[conductanceSlave]
type = MaterialRealAux
property = gap_conductance
variable = gap_condSlave
boundary = BufferGapBndry
execute_on = 'initial timestep_end'
[]
[]
[Contact]
[pellet_clad_mechanical]
primary = 15
secondary = 17
penalty = 1e5
model = frictionless
formulation = kinematic
[]
[]
[ThermalContact]
[thermal_contact]
type = GasGapHeatTransfer
variable = temperature
primary = 15
secondary = 17
initial_moles = initial_moles # coupling to a postprocessor which supplies the initial plenum/gap gas mass
gas_released = 'fis_gas_released co_production' # coupling to postprocessors which supply the fission gas addition, co addition
released_gas_types = 'Kr Xe;
CO'
released_fractions = '0.153 0.847;
1'
tangential_tolerance = 1e-6
roughness_primary = 0e-6
roughness_secondary = 0e-6
jumpdistance_primary = 0
jumpdistance_secondary = 0
quadrature = true
emissivity_secondary = 0.0
emissivity_primary = 0.0
min_gap = 1e-7
max_gap = 50e-6
gap_geometry_type = sphere
[]
[cesium_contact]
type = GapHeatTransfer
variable = conc_Cs
primary = 15
secondary = 17
tangential_tolerance = 1e-6
gap_conductivity_function = d_gap
gap_conductivity_function_variable = temperature
appended_property_name = _conc
quadrature = true
gap_geometry_type = sphere
emissivity_primary = 0.0
emissivity_secondary = 0.0
min_gap = 1e-7
[]
[]
[BCs]
# pin particle along symmetry planes
[no_disp_x]
type = DirichletBC
variable = disp_x
boundary = xzero
value = 0.0
[]
# fix temperature on free surface
[freesurf_temp]
type = FunctionDirichletBC
variable = temperature
boundary = exterior
function = temp_bc
[]
# fix concentration on free surface
[freesurf_conc]
type = DirichletBC
variable = conc_Cs
boundary = exterior
value = 0.0
[]
# exterior and internal pressures
[exterior_pressure_x]
type = Pressure
variable = disp_x
boundary = exterior
factor = 0.1e6
[]
# apply plenum pressure on clad inner walls and pellet surfaces
[PlenumPressure]
[plenumPressure]
boundary = BufferGapVol
initial_pressure = 100
startup_time = 0
R = 8.3145
output_initial_moles = initial_moles
temperature = ave_temp_interior
volume = volumeGas
material_input = 'fis_gas_released co_production'
output = plenum_pressure
[]
[]
[]
[Materials]
[normal_vectors_triso]
type = NormalVectorsTRISO
block = 'IPyC OPyC buffer'
[]
[fission_rate]
type = GenericFunctionMaterial
prop_names = fission_rate
prop_values = fission_rate
[]
[fast_neutron_flux]
type = FastNeutronFlux
calculate_fluence = true
flux_function = power_history
factor = 5e17
[]
[fuel_thermal]
type = UO2Thermal
block = fuel
thermal_conductivity_model = FINK_LUCUTA
initial_porosity = 0.0
temperature = temperature
burnup = burnup
[]
[fuel_elasticity_tensor]
type = ComputeIsotropicElasticityTensor
block = fuel
youngs_modulus = 2.2e11
poissons_ratio = 0.345
[]
[fuel_elastic_stress]
type = ComputeFiniteStrainElasticStress
block = fuel
[]
[fuel_swelling]
type = UO2VolumetricSwellingEigenstrain
gas_swelling_model_type = MATPRO
block = fuel
temperature = temperature
burnup = burnup
eigenstrain_name = fuel_swelling
initial_fuel_density = ${initial_fuel_density}
[]
[fuel_thermal_expansion]
type = ComputeThermalExpansionEigenstrain
block = fuel
thermal_expansion_coeff = 10.0e-6
temperature = temperature
stress_free_temperature = 1500.0
eigenstrain_name = fuel_thermal_strain
[]
[fuel_den]
type = StrainAdjustedDensity
block = fuel
strain_free_density = ${initial_fuel_density}
[]
[fission_gas_release]
type = UO2Sifgrs
block = fuel
temperature = temperature
fission_rate = fission_rate
[]
[fuel_conc]
type = ArrheniusDiffusionCoef
block = fuel
d1 = 5.6e-8 # m^2/s
q1 = 209.0e+3 # J/mol
d2 = 5.2e-4 # m^2/s
q2 = 362.0e+3 # J/mol
gas_constant = 8.3143 # J/K-mol
temperature = temperature
[]
[buffer_elasticity_tensor]
type = ComputeIsotropicElasticityTensor
block = buffer
youngs_modulus = 2.0e10
poissons_ratio = 0.23
[]
[buffer_stress]
type = PyCCreep
block = buffer
flux_conversion_factor = 1.0
temperature = temperature
[]
[buffer_temp]
type = HeatConductionMaterial
block = buffer
thermal_conductivity = 0.5 # J/m-s-K
specific_heat = 720.0 # J/kg-K
[]
[buffer_den]
type = StrainAdjustedDensity
strain_free_density = 1000.0 #kg/m^3
block = buffer
[]
[buffer_thermal_expansion]
type = ComputeThermalExpansionEigenstrain
block = buffer
thermal_expansion_coeff = 5.65e-6
temperature = temperature
stress_free_temperature = 1500.0
eigenstrain_name = buffer_thermal_strain
[]
[buffer_irraditation]
type = PyCIrradiationEigenstrain
block = buffer
pyc_type = buffer
eigenstrain_name = buffer_eigenstrain
[]
[buffer_conc]
type = ArrheniusDiffusionCoef
block = buffer
d1 = 1.0e-12 # m^2/s
q1 = 0.0
d2 = 0.0
q2 = 0.0
temperature = temperature
[]
[IPyC_elasticity_tensor]
type = ComputeIsotropicElasticityTensor
block = IPyC
youngs_modulus = 4.74e10
poissons_ratio = 0.23
[]
[IPyC_stress]
type = PyCCreep
block = IPyC
flux_conversion_factor = 1.0
temperature = temperature
[]
[IPyC_temp]
type = HeatConductionMaterial
block = IPyC
thermal_conductivity = 4.0 # J/m-s-K
specific_heat = 720.0 # J/kg-K
[]
[IPyC_den]
type = StrainAdjustedDensity
strain_free_density = 1900.0 # kg/m^3
block = IPyC
[]
[IPyC_densification]
type = PyCIrradiationEigenstrain
block = IPyC
pyc_type = dense
eigenstrain_name = IPyC_eigenstrain
[]
[IPyC_thermal_expansion]
type = ComputeThermalExpansionEigenstrain
block = IPyC
thermal_expansion_coeff = 5.65e-6
temperature = temperature
stress_free_temperature = 1500.0
eigenstrain_name = IPyC_thermal_strain
[]
[IPyC_conc]
type = ArrheniusDiffusionCoef
block = IPyC
d1 = 6.3e-8 # m^2/s
q1 = 222.0e+3 # J/mol
d2 = 0.0
q2 = 0.0
temperature = temperature
[]
[SiC_elasticity_tensor]
type = ComputeIsotropicElasticityTensor
block = SiC
youngs_modulus = 3.4e11
poissons_ratio = 0.13
[]
[monolithic_SiC_creep]
type = MonolithicSiCCreepUpdate
block = SiC
fast_neutron_flux = fast_neutron_flux
temperature = temperature
k_function = k_function
[]
[stress]
type = ComputeMultipleInelasticStress
inelastic_models = monolithic_SiC_creep
block = SiC
[]
[SiC_temp]
type = HeatConductionMaterial
block = SiC
thermal_conductivity = 13.9 # J/m-s-K
specific_heat = 620.0 # J/kg-K
[]
[SiC_den]
type = StrainAdjustedDensity
strain_free_density = 3180.0 # kg/m^3
block = SiC
[]
[SiC_thermal_expansion]
type = ComputeThermalExpansionEigenstrain
block = SiC
thermal_expansion_coeff = 4.9e-6
temperature = temperature
stress_free_temperature = 1500.0
eigenstrain_name = SiC_thermal_strain
[]
[SiC_conc]
type = ArrheniusDiffusionCoef
block = SiC
d1 = 5.5e-14 # m^2/s
d1_function = d1_function
d1_function_variable = fast_neutron_fluence
q1 = 125.0e+3 # J/mol
d2 = 1.6e-2 # m^2/s
q2 = 514.0e+3 # J/mol
temperature = temperature
[]
[OPyC_elasticity_tensor]
type = ComputeIsotropicElasticityTensor
block = OPyC
youngs_modulus = 4.74e10
poissons_ratio = 0.23
[]
[OPyC_stress]
type = PyCCreep
block = OPyC
flux_conversion_factor = 1.0
temperature = temperature
[]
[OPyC_temp]
type = HeatConductionMaterial
block = OPyC
thermal_conductivity = 4.0 # J/m-s-K
specific_heat = 720.0 # J/kg-K
[]
[OPyC_den]
type = StrainAdjustedDensity
strain_free_density = 1900.0 # kg/m^3
block = OPyC
[]
[OPyC_densification]
type = PyCIrradiationEigenstrain
block = OPyC
pyc_type = dense
eigenstrain_name = OPyC_eigenstrain
[]
[OPyC_thermal_expansion]
type = ComputeThermalExpansionEigenstrain
block = OPyC
thermal_expansion_coeff = 5.65e-6
temperature = temperature
stress_free_temperature = 1500.0
eigenstrain_name = OPyC_thermal_strain
[]
[OPyC_conc]
type = ArrheniusDiffusionCoef
block = OPyC
d1 = 6.3e-8 # m^2/s
q1 = 222.0e+3 # J/mol
d2 = 0.0
q2 = 0.0
temperature = temperature
[]
[]
[Dampers]
[temp]
type = MaxIncrement
variable = temperature
max_increment = 50
[]
[]
[Debug]
show_var_residual_norms = true
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
line_search = 'none'
nl_rel_tol = 5e-8
nl_abs_tol = 1e-7
nl_max_its = 15
l_tol = 1e-8
l_max_its = 50
start_time = 0.0
#end_time = 85.3682e6
end_time = 1e3
num_steps = 1000
dtmax = 2e6
dtmin = 1
[TimeStepper]
type = IterationAdaptiveDT
dt = 20
time_t = '0 76e6 76.001e6 84.641e6 84.6482e6'
time_dt = '20 20 20 20 20'
growth_factor = 1.5
optimal_iterations = 8
linear_iteration_ratio = 100
[]
[Quadrature]
order = THIRD
[]
[]
[Postprocessors]
[release_Cs_inc]
type = SideIntegralMassFlux
variable = conc_Cs
boundary = exterior
[]
[Int_Cs_release]
type = TimeIntegratedPostprocessor
value = release_Cs_inc
[]
[release_fuel_Cs]
type = SideIntegralMassFlux
variable = conc_Cs
boundary = fuel
[]
[Int_Cs_release_fuel]
type = TimeIntegratedPostprocessor
value = release_fuel_Cs
[]
[release_PyCGapBndry_Cs]
type = SideIntegralMassFlux
variable = conc_Cs
boundary = PyCGapBndry
[]
[Int_Cs_release_PyCGapBndry]
type = TimeIntegratedPostprocessor
value = release_PyCGapBndry_Cs
[]
[fis_gas_produced]
type = ElementIntegralFisGasGeneratedSifgrs
block = fuel
execute_on = linear
[]
[fis_gas_released]
type = ElementIntegralFisGasReleasedSifgrs
block = fuel
execute_on = linear
[]
[volumeTotal]
type = InternalVolume
boundary = exterior
execute_on = 'initial timestep_end'
[]
[volumeFuel]
type = InternalVolume
boundary = fuel
execute_on = 'initial timestep_end'
[]
[volumeGas]
type = InternalVolume
boundary = BufferGapVol
addition = -4.6e-11
execute_on = 'initial linear'
[]
[volumeBufferShell]
type = InternalVolume
boundary = BufferGapVol
execute_on = 'initial timestep_end'
[]
[ave_temp_interior]
type = SideAverageValue
boundary = BufferGapVol
variable = temperature
execute_on = 'initial timestep_end'
[]
# Postprocessors for CO production
[total_fission_rate]
type = ElementIntegralPower
variable = temperature
fission_rate = fission_rate
block = fuel
energy_per_fission = 1.0
outputs = exodus
execute_on = 'initial linear'
[]
[total_fissions]
type = TimeIntegratedPostprocessor
value = total_fission_rate
outputs = exodus
execute_on = 'initial timestep_end'
[]
[avg_surface_temp]
type = SideAverageValue
variable = temperature
boundary = exterior
outputs = exodus
execute_on = 'initial timestep_end'
[]
[time_int_surf_temp]
type = TimeIntegratedPostprocessor
value = avg_surface_temp
outputs = exodus
execute_on = 'initial timestep_end'
[]
[co_production]
type = CarbonMonoxideProduction
total_fissions = total_fissions
time_integrated_triso_temperature = time_int_surf_temp
initial_enrichment = 0.14029
execute_on = 'initial nonlinear'
[]
[]
[Outputs]
print_linear_residuals = false
[console]
type = Console
max_rows = 5
outlier_variable_norms = false
[]
[exodus]
type = Exodus
file_base = triso1D_accident_out
[]
[]