Problem of Negative Feedback in Eigenvalue Calculation

[js-jixu]

Hi, developers.

In the steady state calculation of the neutron diffusion equation, I think the calculated neutron flux distribution and power distribution should be different when given different temperature distributions. But in my input file. I have tried several different temperature distributions in the fuel region, but the calculated neutron flux distributions are all exactly the same. Temperature does not have any negative feedback effect on the neutron flux distribution. Why does this happen?

The temperature distributions in the following two graphs are different, but the power distributions are exactly the same:

Here are input file and mesh file:

[GlobalParams]
  num_groups = 6
  num_precursor_groups = 6
  use_exp_form = false
  group_fluxes = 'group1 group2 group3 group4 group5 group6'
  pre_concs = 'pre1 pre2 pre3 pre4 pre5 pre6'
  temperature = T
  sss2_input = true
  account_delayed = true
[../]

[Mesh]
  file = 3layers_3d_4parts.msh
[]

[Problem]
  type = FEProblem
[]

[Nt]
  var_name_base = group
  vacuum_boundaries = 'cool_bottom cool_top clad_bottom clad_top fuel_bottom fuel_top'
  create_temperature_var = false
  eigen = true
  #scaling = 1e3
  block = 'fuel cool clad'
  pre_blocks = 'fuel'
[]

[Precursors]
  [./pres]
    var_name_base = pre
    outlet_boundaries = ''
    constant_velocity_values = true
    u_def = 0
    v_def = 0
    w_def = 0
    nt_exp_form = false
    family = MONOMIAL
    order = CONSTANT
    loop_precursors = false
    transient = false
    eigen = true
    #scaling = 1e3
    block = 'fuel'
  [../]
[]

[AuxVariables]
  [./T]
    family = LAGRANGE
    order = FIRST
  [../]
  [./power]
    family = MONOMIAL
    order = FIRST
	block = fuel
  [../]
[]

[ICs]
  [./fuel_ini_T]
    type = FunctionIC
	variable = T
	function = fuel_temp_func
	block = fuel
  [../]
  [./clad_ini_T]
    type = ConstantIC
	variable = T
	value = 600
	block = clad
  [../]
  [./cool_ini_T]
    type = ConstantIC
	variable = T
	value = 600
	block = cool
  [../]
[]

[Functions]
  [./fuel_temp_func]
    type = ParsedFunction
	value = '600 + 750 * (z / 90)'
  [../]
[]

[AuxKernels]
  [./heat_source]
    type = FissionHeatSourceAux
    variable = power
    tot_fission_heat = powernorm
    power = 21230
  [../]
[]

[Materials]
  [./fuel]
    type = GenericMoltresMaterial
    property_tables_root = './xs/fuel_'
    interp_type = 'linear'
    block = 'fuel'
  [../]
  [./clad]
    type = GenericMoltresMaterial
    property_tables_root = './xs/clad_'
    interp_type = 'linear'
    block = 'clad'
  [../]
  [./cool]
    type = GenericMoltresMaterial
    property_tables_root = './xs/cool_'
    interp_type = 'linear'
    block = 'cool'
  [../]
[]

[Executioner]
  type = InversePowerMethod
  max_power_iterations = 50

  # fission power normalization settings
  normalization = 'powernorm'
  normal_factor = 21230           # Watts

  xdiff = 'group1diff'
  bx_norm = 'bnorm'
  k0 = 1
  l_max_its = 200
  eig_check_tol = 1e-5

  solve_type = 'NEWTON'
  petsc_options = '-snes_converged_reason -ksp_converged_reason -snes_linesearch_monitor'
#  petsc_options_iname = '-pc_type -sub_pc_type -ksp_gmres_restart -pc_gasm_overlap -sub_pc_factor_shift_type -pc_gasm_blocks -sub_pc_factor_mat_solver_type'
#  petsc_options_value = 'gasm     lu           200                1                NONZERO                   16              superlu_dist'


#  petsc_options_iname = '-pc_type -pc_factor_shift_type -pc_factor_mat_solver_package'
#  petsc_options_value = 'lu       NONZERO               superlu_dist'

## Use the settings below instead if running on a desktop/small cluster
  petsc_options_iname = '-pc_type -sub_pc_type -ksp_gmres_restart -pc_asm_overlap -sub_pc_factor_shift_type'
  petsc_options_value = 'asm      lu           200                1               NONZERO'

  line_search = none
  automatic_scaling = true
[]

[Postprocessors]
  [./bnorm]
    type = ElmIntegTotFissNtsPostprocessor
    execute_on = linear
  [../]
  [./tot_fissions]
    type = ElmIntegTotFissPostprocessor
    execute_on = linear
  [../]
  [./powernorm]
    type = ElmIntegTotFissHeatPostprocessor
    execute_on = linear
  [../]
  [./group1norm]
    type = ElementIntegralVariablePostprocessor
    variable = group1
    execute_on = linear
  [../]
  [./group1max]
    type = NodalMaxValue
    variable = group1
    execute_on = timestep_end
  [../]
  [./group1diff]
    type = ElementL2Diff
    variable = group1
    execute_on = 'linear timestep_end'
    use_displaced_mesh = false
  [../]
  [./group2norm]
    type = ElementIntegralVariablePostprocessor
    variable = group2
    execute_on = linear
  [../]
  [./group2max]
    type = NodalMaxValue
    variable = group2
    execute_on = timestep_end
  [../]
  [./group2diff]
    type = ElementL2Diff
    variable = group2
    execute_on = 'linear timestep_end'
    use_displaced_mesh = false
  [../]
[]

[Outputs]
  perf_graph = true
  print_linear_residuals = true
  [./exodus]
    type = Exodus
    file_base = 'nts-steady'
    execute_on = 'TIMESTEP_END'
  [../]
  [./csv]
    type = CSV
    file_base = 'nts-steady'
    execute_on = 'TIMESTEP_END'
  [../]
[]

[Debug]
  show_var_residual_norms = true
[]

3layers_3d_4parts.zip

[smpark7]

As far as I can tell, your input file seems correct. Did you provide cross section data ranging from 600K to 1350K?

[js-jixu]

Yep, I calculated the cross section data with openmc.

xs.zip

[js-jixu]

I have plotted the axial distribution of neutron flux and power. There is a slight difference between them, albeit a very small one. But their temperature distributions are quite different. Should I use OpenMC to verify the power distribution under different temperature distributions?

[smpark7]

I'm not sure if you can get comparable results from OpenMC due to the statistical uncertainty, but you can try. I noticed that your fission XS in the fuel region does not vary much with temperature. Going from 850K to 1350K, it decreases by about 1.5%. For comparison, the fission XS of molten salt in the MSFR decreases by 10% (from Doppler effect and density change). Perhaps that's why we see such a small change.

[js-jixu]

Yep, the doppler fuel temperature coefficient of the lead-bismuth fast reactor is very small, about -0.14pcm/K. I will discuss this with my tutor. The discussion results will be added after the discussion.

[js-jixu]

In order to check whether it is a XS problem, I plan to do a nuclear-thermal coupled simulation of a single rod of a pressurized water reactor. I gave the fuel region an axial sinusoidal temperature distribution, and it turned out that the negative feedback effect was significant. More coupling results will be provided in the coming days.