Hi,
When I try to run an example which compute the depletion of a reactor,I found that keff has not changed throughout its entire lifespan.I don’t know what caused this problem, so I tried changing the environment, changing the cross-sectional file, and setting ’depletable=True‘, but nothing changed. Finally, I plotted the changes in the number of particles for U235 and found that there was no change at all. I am certain that I have set ’depletable=True‘ of fuel. Who knows the reason for this and can you help me solve this problem?
Thanks.
Could you share your depletion code?
sure，there is my code:
import openmc
import os
from math import log10
import numpy as np
import math
import matplotlib.pyplot as pyplot
openmc.config['cross_sections'] ="/home/shuoxing/cs/mcnp_endfb71/cross_sections.xml"
fuel20=openmc.Material(name='fuel20')
fuel20.add_element('U',1.0,enrichment=20)
fuel20.add_nuclide('O16',2.0)
fuel20.set_density('g/cc',10.3)
fuel20.depletable=True
water=openmc.Material(name='water')
water.add_nuclide('H1', 2.0)
water.add_nuclide('O16',1.0)
water.set_density('g/cc', 1.0)
water.add_s_alpha_beta('c_H_in_H2O')
Er2O3=openmc.Material(name='Er2O3')
Er2O3.add_nuclide('Er167',2.0)
Er2O3.add_nuclide('O16',3.0)
Er2O3.set_density('g/cm3',8.64)
zirconium = openmc.Material(name="zirconium")
zirconium.add_element('Zr', 1.0)
zirconium.set_density('g/cm3', 6.6)
niobium=openmc.Material(name='niobium')
niobium.add_element('Nb',1.0)
niobium.set_density('g/cm3',8.57)
sn=openmc.Material(name='sn')
sn.add_element('Sn',1.0)
sn.set_density('g/cm3',7.28)
cr=openmc.Material(name='cr')
cr.add_element('Cr',1.0)
cr.set_density('g/cm3',7.2)
Fe=openmc.Material(name='Fe')
Fe.add_element('Fe',1.0)
Fe.set_density('g/cm3',7.86)
ag=openmc.Material(name='ag')
ag.add_element('Ag',1.0)
ag.set_density('g/cm3',10.49)
In=openmc.Material(name='In')
In.add_element('In',1.0)
In.set_density('g/cm3',7.30)
He=openmc.Material(name='He')
He.add_element('He',1.0)
He.set_density('g/cc',1.29)
B4C=openmc.Material(name='keran')
B4C.add_nuclide('B10',4.0)
B4C.add_nuclide('C0',1.0)
B4C.set_density('g/cc',2.52)
Al2O3=openmc.Material(name='Al2O3')
Al2O3.add_element('Al',2.0)
Al2O3.add_nuclide('O16',3.0)
Al2O3.set_density('g/cm3',3.5)
BeO=openmc.Material(name='BeO')
BeO.add_element('Be',1.0)
BeO.add_nuclide('O16',1.0)
BeO.set_density('g/cm3',3.02)
WABA_48=openmc.Material.mix_materials([Al2O3,B4C],[0.93,0.07])
WABA_48.depletable=True
zr2=openmc.Material.mix_materials([sn,Fe,cr,niobium,zirconium],[0.0145,0.00135,0.001,0.00055,0.9826])
zr4=openmc.Material.mix_materials([sn,Fe,cr,zirconium],[0.0145,0.0021,0.001,0.9824])
keran=openmc.Material.mix_materials([zr2,B4C],[0.9672,0.0328],'wo')
control=openmc.Material.mix_materials([ag,In,cr],[0.8,0.15,0.05])
fuel_Er10=openmc.Material.mix_materials([fuel20,Er2O3],[0.90,0.10])
fuel_Er10.depletable=True
materials_file = openmc.Materials([fuel20,water,zr4,keran,zr2,control,He,Er2O3,fuel_Er10,B4C,Al2O3,WABA_48,BeO])
materials_file.export_to_xml()
# Boudries and outer universe
all_water_out=openmc.Cell(cell_id=200, fill=water)
num_pitch_pin = 1.131656934
num_pitch_lat=num_pitch_pin*17
num_pitch_core=num_pitch_lat*7
left = openmc.XPlane(-num_pitch_lat/2)
right=openmc.XPlane(num_pitch_lat/2)
up=openmc.YPlane(num_pitch_lat/2)
down=openmc.YPlane(-num_pitch_lat/2)
# boundary of the assembly
assembly_z0 = openmc.ZPlane(surface_id=300, z0=-50)
assembly_z1 = openmc.ZPlane(surface_id=301, z0=50)
assembly_out=openmc.ZCylinder(r=80)
core_top=openmc.ZPlane(z0=55)
core_bot=openmc.ZPlane(z0=-55)
# boundary & bottom of the reflector
reflector_out=openmc.ZCylinder(r=85)
# boundary of shell
shell_in=openmc.ZCylinder(r=88)
shell_out=openmc.ZCylinder(r=95,boundary_type='vacuum')
shell_top=openmc.ZPlane(z0=62,boundary_type='vacuum')
shell_bot=openmc.ZPlane(z0=-62,boundary_type='vacuum')
'''reflect_cell.region = ~assembly_cell.region & +r_left&-r_right&+r_down&-r_up & -assembly_z1 & +assembly_z0
top_reflect_cell.region = +r_left&-r_right&+r_down&-r_up & +assembly_z1 & -reflector_z1
bot_reflect_cell.region = +r_left&-r_right&+r_down&-r_up & -assembly_z0 & +reflector_z0'''
# Create universes
u_root = openmc.Universe()
r_fuel=0.65/2
fuel20out=openmc.ZCylinder(r=r_fuel)
clad20out=openmc.ZCylinder(r=0.75/2)
clad20in=openmc.ZCylinder(r=0.68/2)
controlout=openmc.ZCylinder(r=0.62/2)
cladcontrolin=openmc.ZCylinder(r=0.65/2)
cladcontrolout=openmc.ZCylinder(r=0.80/2)
fuel20_region=-fuel20out&+assembly_z0&-assembly_z1
clad20_region=+clad20in&-clad20out&+assembly_z0&-assembly_z1
he_region=+fuel20out & -clad20in&+assembly_z0&-assembly_z1
reflector_region=+assembly_out & -reflector_out &+assembly_z0&-assembly_z1
shell_region = +shell_in & -shell_out&+assembly_z0&-assembly_z1
ref_water_shell_region=+reflector_out&-shell_in&+assembly_z0&-assembly_z1
fuel20_cell=openmc.Cell(cell_id=101,fill=fuel20,region=fuel20_region)
clad20_cell=openmc.Cell(cell_id=102,fill=zr4,region=clad20_region)
he_cell=openmc.Cell(cell_id=103,fill=He,region=he_region)
water_cell=openmc.Cell(fill=water,region=+clad20out)
fuel_Gd6_cell=openmc.Cell(cell_id=104,fill=fuel_Er10,region=fuel20_region)
clad20_Gd6_cell=openmc.Cell(cell_id=105,fill=zr4,region=clad20_region)
he_cell_Gd6=openmc.Cell(cell_id=106,fill=He,region=he_region)
water_cell_Gd6=openmc.Cell(fill=water,region=+clad20out)
control_region=-controlout&+assembly_z0&-assembly_z1
cladcontrol_region=+cladcontrolin & -cladcontrolout&+assembly_z0&-assembly_z1
gap_region=+controlout&-cladcontrolin&+assembly_z0&-assembly_z1
control_water_20_cell=openmc.Cell(cell_id=1001,fill=water,region=control_region)
control_clad20_cell=openmc.Cell(cell_id=1002,fill=zr4,region=cladcontrol_region)
control_water20_cell=openmc.Cell(cell_id=1003,fill=water,region=gap_region)
control_water_cell=openmc.Cell(fill=water,region=+clad20out)
'''
control_bang_20_cell=openmc.Cell(cell_id=1004,fill=control,region=control_region)
control_bang_clad20_cell=openmc.Cell(cell_id=1005,fill=zr4,region=cladcontrol_region)
control_bang_water20_cell=openmc.Cell(cell_id=1006,fill=water,region=gap_region)
control__bang_water_cell=openmc.Cell(fill=water,region=+clad20out)'''
WABA_waterout=openmc.ZCylinder(r=0.338346/2)
WABA_zr1out=openmc.ZCylinder(r=0.401504/2)
WABA_in=openmc.ZCylinder(r=0.418045/2)
WABA_out=openmc.ZCylinder(r=0.478195/2)
WABA_zr2in=openmc.ZCylinder(r=0.494737/2)
WABA_zr2out=openmc.ZCylinder(r=0.572932/2)
WABA_zr3in=openmc.ZCylinder(r=0.748872/2)
WABA_zr3rout=openmc.ZCylinder(r=0.8/2)
WABA_water_region=-WABA_waterout&+assembly_z0&-assembly_z1
WABA_zr1_region=+WABA_waterout&-WABA_zr1out&+assembly_z0&-assembly_z1
WABA_He1_region=+WABA_zr1out&-WABA_in&+assembly_z0&-assembly_z1
WABA_region=+WABA_in&-WABA_out&+assembly_z0&-assembly_z1
WABA_He2_region=+WABA_out&-WABA_zr2in&+assembly_z0&-assembly_z1
WABA_zr2_region=+WABA_zr2in&-WABA_zr2out&+assembly_z0&-assembly_z1
WABA_water2_region=+WABA_zr2out&-WABA_zr3in&+assembly_z0&-assembly_z1
WABA_zr3_region=+WABA_zr3in&-WABA_zr3rout&+assembly_z0&-assembly_z1
WABA_water_cell=openmc.Cell(fill=water,region=WABA_water_region)
WABA_zr1_cell=openmc.Cell(fill=zr4,region=WABA_zr1_region)
WABA_He1_cell=openmc.Cell(fill=He,region=WABA_He1_region)
WABA_cell=openmc.Cell(fill=WABA_48,region=WABA_region)
WABA_He2_cell=openmc.Cell(fill=He,region=WABA_He2_region)
WABA_zr2_cell=openmc.Cell(fill=zr4,region=WABA_zr2_region)
WABA_water2_cell=openmc.Cell(fill=water,region=WABA_water2_region)
WABA_zr3_cell=openmc.Cell(fill=zr4,region=WABA_zr3_region)
WABA_waterout_cell=openmc.Cell(fill=water,region=+WABA_zr3rout)
u20=openmc.Universe(cells=(fuel20_cell,clad20_cell,water_cell,he_cell))
u20_Gd6=openmc.Universe(cells=(fuel_Gd6_cell,clad20_Gd6_cell,he_cell_Gd6,water_cell_Gd6))
control_water=openmc.Universe(cells=(control_water_cell,control_water20_cell,control_clad20_cell,control_water_20_cell))
#control_bang=openmc.Universe(cells=(control_bang_20_cell,control_bang_clad20_cell,control_bang_water20_cell,control__bang_water_cell))
WABA_u=openmc.Universe(cells=(WABA_cell,WABA_zr3_cell,WABA_zr2_cell,WABA_He2_cell,WABA_water_cell,WABA_zr1_cell,WABA_He1_cell,WABA_water2_cell,WABA_waterout_cell))
all_water_cell=openmc.Cell(cell_id=301,fill=water)
outer_universe=openmc.Universe(cells=(all_water_cell,))
lattice = openmc.RectLattice()
lattice.lower_left = (-num_pitch_lat/2, -num_pitch_lat/2)
lattice.pitch = (num_pitch_pin, num_pitch_pin)
lattice.universes = [[u20, u20_Gd6, u20, u20,u20_Gd6,u20,u20,u20,u20,u20,u20,u20,u20_Gd6,u20,u20,u20_Gd6,u20],
[u20_Gd6, u20, u20, u20_Gd6,u20,u20,u20_Gd6,u20,u20,u20,u20_Gd6,u20,u20,u20_Gd6,u20,u20,u20_Gd6],
[u20, u20, u20, u20,u20,WABA_u,u20,u20,control_water,u20,u20,WABA_u,u20,u20,u20,u20,u20],
[u20, u20_Gd6, u20, control_water,u20_Gd6,u20,u20,u20_Gd6,u20,u20_Gd6,u20,u20,u20_Gd6,control_water,u20,u20_Gd6,u20],
[u20_Gd6, u20, u20, u20_Gd6,u20,u20,u20_Gd6,u20,u20,u20,u20_Gd6,u20,u20,u20_Gd6,u20,u20,u20_Gd6],
[u20, u20, WABA_u, u20,u20,control_water,u20,u20,control_water,u20,u20,control_water,u20,u20,WABA_u,u20,u20],
[u20, u20_Gd6, u20, u20,u20_Gd6,u20,u20,u20_Gd6,u20,u20_Gd6,u20,u20,u20_Gd6,u20,u20,u20_Gd6,u20],
[u20, u20, u20, u20_Gd6,u20,u20,u20_Gd6,u20,u20,u20,u20_Gd6,u20,u20,u20_Gd6,u20,u20,u20],
[u20, u20, control_water, u20,u20,control_water,u20,u20,control_water,u20,u20,control_water,u20,u20,control_water,u20,u20],
[u20, u20, u20, u20_Gd6,u20,u20,u20_Gd6,u20,u20,u20,u20_Gd6,u20,u20,u20_Gd6,u20,u20,u20],
[u20, u20_Gd6, u20, u20,u20_Gd6,u20,u20,u20_Gd6,u20,u20_Gd6,u20,u20,u20_Gd6,u20,u20,u20_Gd6,u20],
[u20, u20, WABA_u, u20,u20,control_water,u20,u20,control_water,u20,u20,control_water,u20,u20,WABA_u,u20,u20],
[u20_Gd6, u20, u20, u20_Gd6,u20,u20,u20_Gd6,u20,u20,u20,u20_Gd6,u20,u20,u20_Gd6,u20,u20,u20_Gd6],
[u20, u20_Gd6, u20, control_water,u20_Gd6,u20,u20,u20_Gd6,u20,u20_Gd6,u20,u20,u20_Gd6,control_water,u20,u20_Gd6,u20],
[u20, u20, u20, u20,u20,WABA_u,u20,u20,control_water,u20,u20,WABA_u,u20,u20,u20,u20,u20],
[u20_Gd6, u20, u20, u20_Gd6,u20,u20,u20_Gd6,u20,u20,u20,u20_Gd6,u20,u20,u20_Gd6,u20,u20,u20_Gd6],
[u20, u20_Gd6, u20, u20,u20_Gd6,u20,u20,u20,u20,u20,u20,u20,u20_Gd6,u20,u20,u20_Gd6,u20]]
fuel_lattice = openmc.Cell()
fuel_lattice.fill = lattice
Fuel20 = openmc.Universe(cells=[fuel_lattice])
lat_left=openmc.XPlane(x0=-num_pitch_lat/2)
lat_right=openmc.XPlane(x0=num_pitch_lat/2)
lat_up=openmc.YPlane(y0=num_pitch_lat/2)
lat_down=openmc.YPlane(y0=-num_pitch_lat/2)
Pri_lat = +lat_left & -lat_right & +lat_down &-lat_up
C_empty = openmc.Cell(fill=water, region=Pri_lat)
WATER = openmc.Universe(cells=[C_empty])
L_core = openmc.RectLattice()
L_core.lower_left = (-num_pitch_core / 2, -num_pitch_core / 2)
L_core.pitch = (num_pitch_lat, num_pitch_lat)
L_core.universes = [
[WATER , WATER , Fuel20, Fuel20, Fuel20, WATER , WATER ],
[WATER , Fuel20, Fuel20, Fuel20, Fuel20, Fuel20, WATER ],
[Fuel20 , Fuel20, Fuel20, Fuel20, Fuel20, Fuel20, Fuel20],
[Fuel20 , Fuel20, Fuel20, Fuel20, Fuel20, Fuel20, Fuel20],
[Fuel20 , Fuel20, Fuel20, Fuel20, Fuel20, Fuel20, Fuel20],
[WATER , Fuel20, Fuel20, Fuel20, Fuel20, Fuel20, WATER ],
[WATER , WATER , Fuel20, Fuel20, Fuel20, WATER , WATER ]]
all_water_out=openmc.Cell(fill=water)
all_water_out_u=openmc.Universe(cells=[all_water_out])
L_core.outer=all_water_out_u
C_core = openmc.Cell()
C_core.region = -assembly_out & +assembly_z0 & -assembly_z1
C_core.fill = L_core
R_flan = +assembly_out & -reflector_out & +core_bot & -core_top
C_flan = openmc.Cell(fill=BeO, region=R_flan)
R_downChannel = +reflector_out & -shell_in & +core_bot & -core_top
C_downChannel = openmc.Cell(fill=water, region=R_downChannel)
R_vess = +shell_in & -shell_out & +core_bot & -core_top
C_vess = openmc.Cell(fill=zr4, region=R_vess)
R_bottom=-assembly_out & -assembly_z0 & +core_bot
C_bot = openmc.Cell(fill=water, region=R_bottom)
R_top=-assembly_out& + assembly_z1 & -core_top
C_top = openmc.Cell(fill=water, region=R_top)
R_bottom_vess=+shell_bot&-core_bot&-shell_out
C_bot_vess = openmc.Cell(fill=zr4, region=R_bottom_vess)
R_top_vess=-shell_top&+core_top&-shell_out
C_top_vess = openmc.Cell(fill=zr4, region=R_top_vess)
u_root.add_cells([C_core, C_flan, C_downChannel, C_vess,C_bot,C_top,C_bot_vess,C_top_vess])
geometry=openmc.Geometry(u_root)
geometry.export_to_xml()
plot=openmc.Plot()
plot.origin=(0,0,0)
plot.width=(200,200)
plot.pixels=(500,500)
plot.color_by='material'
plot.colors={water:'blue',
zr4:'black',
BeO:'green'}
plot2=openmc.Plot()
plot2.basis='yz'
plot2.origin=(0,0,0)
plot2.width=(200,130)
plot2.pixels=(600,600)
plot2.color_by='material'
plot2.colors={water:'blue',
zr4: 'black',
BeO: 'green'}
plots=openmc.Plots([plot,plot2])
plots.export_to_xml()
'''openmc.plot_geometry()
lower_left, upper_right = geometry.bounding_box
cell_vol = openmc.VolumeCalculation([fuel20_cell, fuel_Gd6_cell,clad20_Gd6_cell,clad20_cell,he_cell,water_cell,all_water_out,he_cell_Gd6,water_cell_Gd6,control_water_cell,
control_water20_cell,control_water_20_cell,control_clad20_cell,WABA_water_cell,WABA_water2_cell,WABA_zr2_cell,WABA_zr3_cell,WABA_zr1_cell,
WABA_He1_cell,WABA_He2_cell,WABA_cell,WABA_waterout_cell,fuel_lattice,C_empty,C_top,C_bot_vess,C_top_vess,C_bot,C_core,C_downChannel,C_vess,C_flan] , 1000000000, lower_left, upper_right)
settings = openmc.Settings()
settings.volume_calculations = [cell_vol]
settings.run_mode = 'volume'
settings.export_to_xml()
openmc.calculate_volumes()'''
# OpenMC simulation parameters
batches = 200
particles = 20000
settings = openmc.Settings()
settings.batches = batches
settings.particles = particles
uniform_dist = openmc.stats.Box([-num_pitch_lat/2, -num_pitch_lat/2, -num_pitch_lat/2],[num_pitch_lat/2, num_pitch_lat/2, num_pitch_lat/2],only_fissionable=True)
settings.source = openmc.source.Source(space=uniform_dist)
settings.output = {'tallies': True}
settings.export_to_xml()
# 网格计数
mesh = openmc.RegularMesh()
mesh.dimension = (100, 100, 100)
mesh.lower_left = (-num_pitch_lat/2, -num_pitch_lat/2, -50)
mesh.upper_right = (num_pitch_lat/2, num_pitch_lat/2, 50)
mesh_filter = openmc.MeshFilter(mesh)
mesh_tally = openmc.Tally(name="Mesh tally")
mesh_tally.filters = [mesh_filter]
mesh_tally.scores = ['flux', 'fission', 'nu-fission']
# 能谱
e_min, e_max = 1e-5, 20.0e6
groups = 500
energies = np.logspace(log10(e_min), log10(e_max), groups + 1)
energy_filter = openmc.EnergyFilter(energies)
spectrum_tally = openmc.Tally(name="Flux spectrum")
spectrum_tally.filters = [energy_filter]
spectrum_tally.scores = ['flux']
fuel20.volume=math.pi * r_fuel ** 2 *200*37*100
fuel_Er10.volume=math.pi*r_fuel **2 *64*37*100
WABA_48.volume=math.pi*((0.478195/2)**2-(0.418045/2)**2)*8*37*100
'''cell_vol = openmc.VolumeCalculation([fuel20_cell, fuel_Gd6_cell,clad20_Gd6_cell,clad20_cell,he_cell,water_cell,all_water_out,he_cell_Gd6,water_cell_Gd6,control_water_cell,
control_water20_cell,control_water_20_cell,control_clad20_cell,WABA_water_cell,WABA_water2_cell,WABA_zr2_cell,WABA_zr3_cell,WABA_zr1_cell,
WABA_He1_cell,WABA_He2_cell,WABA_cell,WABA_waterout_cell,fuel_lattice,C_empty,C_top,C_bot_vess,C_top_vess,C_bot,C_core,C_downChannel,C_vess,C_flan] , 1000000000, lower_left, upper_right)
fuel20_cell.volume=245559.5425956
fuel_Gd6_cell.volume=78587.53963160001
clad20_Gd6_cell.volume=18613.632188
clad20_cell.volume=58188.334153200005
he_cell.volume=23176.4804248
water_cell.volume=620702.21231
all_water_out.volume=296828.06066719996
he_cell_Gd6.volume=7408.262944
water_cell_Gd6.volume=198607.75581240002
control_water_cell.volume=52785.14925
control_water20_cell.volume=1880.7997475999998
control_water_20_cell.volume=18992.2729584
control_clad20_cell.volume=6914.896518
WABA_water_cell.volume=2660.6020568
WABA_water2_cell.volume=5409.3399532
WABA_zr2_cell.volume=1937.3322031999999
WABA_zr3_cell.volume=1842.8398756
WABA_zr1_cell.volume=1088.4859004
WABA_He1_cell.volume=314.98636239999996
WABA_He2_cell.volume=375.9817888
WABA_cell.volume=1255.6749639999998
WABA_waterout_cell.volume=23038.9430348
fuel_lattice.volume=1369341.0646728
C_empty.volume=344483.219706
C_top.volume=100516.4742348
C_bot_vess.volume=198516.81327
C_top_vess.volume=198427.27184079998
C_bot.volume=100512.212702
C_core.volume= 2010652.3450459999
C_downChannel.volume=179337.2611648
C_vess.volume=442703.62751799996
C_flan.volume=285101.8757124
'''
import openmc.deplete
chain = openmc.deplete.Chain.from_xml('/home/shuoxing/chain/chain_endfb71_pwr.xml')
print(chain.nuclide_dict)
model = openmc.Model(materials=materials_file, geometry=geometry, settings=settings)
#model.differentiate_depletable_mats(diff_volume_method='match cell')
operator = openmc.deplete.CoupledOperator(model, "/home/shuoxing/chain/chain_endfb71_pwr.xml")
power = 200000000
time_step=[60]*5
keff_list=[]
integrator = openmc.deplete.PredictorIntegrator(operator, time_step, power, timestep_units='d')
integrator.integrate()
results = openmc.deplete.Results('./depletion_results.h5')
time, k = results.get_keff()
print(k)
_, u235 = results.get_atoms("1", "U235")
_, xe135 = results.get_atoms("1", "Xe135")
pyplot.plot(time, u235, label="U235")
pyplot.xlabel("Time [d]")
pyplot.ylabel("Number of atoms - U235");
I have calculated the change in nuclear density, and there has been no change in nuclear density. I suspect that the nuclear density has not been updated properly. At the same time, I conducted calculations on another computer, and in version 0.14.0, the nuclear density still cannot be updated. However, in version 0.13.3 of another computer, calculations can be performed. It is worth noting that on the computer I frequently use, version 0.13.2 still cannot update the nuclear density.
It seems to work fine for me using my depletion chain file and my cross sections library on version 0.14. Whatever the problem is it doesn’t seem to be in the code itself.
Yes, that’s what’s bothering me the most. I don’t know how to solve this problem