arry size not right

tasks/task_11_CSG_shut_down_dose_tallies/3_mesh_based_shut_down_dose_rate_example.py

@@ -65,7 +65,7 @@
 box = openmc.model.RectangularParallelepiped(
  xmin=-20, xmax=20, ymin=-20, ymax=20, zmin=-20, zmax=20, boundary_type="vacuum"
 )
-sphere_cell_4 = openmc.Cell(region=-box & +sphere_surf_1)
+sphere_cell_4 = openmc.Cell(region=-box & +sphere_surf_1, fill=mat_aluminum)
 
 my_geometry = openmc.Geometry([sphere_cell_1, sphere_cell_2, sphere_cell_3, sphere_cell_4])
 
@@ -101,32 +101,9 @@
 
 model_neutron = openmc.Model(my_geometry, my_materials, my_neutron_settings)
 
-hour_in_seconds = 60*60
 
-# This section defines the neutron pulse schedule.
-# If the method made use of the CoupledOperator then there would need to be a
-# transport simulation for each timestep. However as the IndependentOperator is
-# used then just a single transport simulation is done, thus speeding up the
-# simulation considerably.
-timesteps_and_source_rates = [
- (1, 1e18), # 1 second
- (hour_in_seconds, 0), # 1 hour
- (1, 1e18), # 1 second
- (hour_in_seconds, 0), # 2 hour
- (1, 1e18), # 1 second
- (hour_in_seconds, 0), # 3 hour
- (1, 1e18), # 1 second
- (hour_in_seconds, 0), # 4 hour
- (1, 1e18), # 1 second
- (hour_in_seconds, 0), # 4 hour
- (1, 1e18), # 1 second
- (hour_in_seconds, 0), # 4 hour
- (1, 1e18), # 1 second
- (hour_in_seconds, 0), # 5 hour
-]
 
-timesteps = [item[0] for item in timesteps_and_source_rates]
-source_rates = [item[1] for item in timesteps_and_source_rates]
+
 
 model_neutron.export_to_xml(directory=statepoints_folder/ "neutrons")
 model_neutron.export_to_xml()
@@ -139,7 +116,7 @@
 # print(set(all_nuclides))
 
 # this does perform transport but just to get the flux and micro xs
-flux_in_each_group, micro_xs = openmc.deplete.get_microxs_and_flux(
+flux_in_each_group, all_micro_xs = openmc.deplete.get_microxs_and_flux(
  model=model_neutron,
  domains=regular_mesh,
  energies=[0,30e6], # different group structures see this file for all the groups available https://github.com/openmc-dev/openmc/blob/develop/openmc/mgxs/__init__.py
@@ -153,14 +130,14 @@
 openmc.lib.init()
 
 
-mesh = openmc.lib.RegularMesh()
-mesh.dimension = regular_mesh.dimension
-mesh.set_parameters(
+lib_mesh = openmc.lib.RegularMesh()
+lib_mesh.dimension = regular_mesh.dimension
+lib_mesh.set_parameters(
  lower_left=regular_mesh.lower_left,
  upper_right=regular_mesh.upper_right
 )
 
-vols = mesh.material_volumes(n_samples = 1000000)
+mesh_material_volumes = lib_mesh.material_volumes(n_samples = 1000000)
 
 mesh_voxel_material = []
 
@@ -169,12 +146,14 @@
 volume_of_material_in_voxel = []
 material_in_voxel = []
 
+selective_flux_in_each_group = []
+selective_micro_xs = []
 
-for i, entry in enumerate(vols):
- print(entry)
+for i, (mesh_material_volume, micro_xs, flux_in_group) in enumerate(zip(mesh_material_volumes, all_micro_xs, flux_in_each_group) ):
+ print(mesh_material_volume)
  materials_in_voxel = []
  volumes_in_voxel = []
- for material_volume_tuple in entry:
+ for material_volume_tuple in mesh_material_volume:
  material = material_volume_tuple[0]
  if material != None:
  volume_in_cm3 = material_volume_tuple[1]
@@ -186,22 +165,33 @@
 
  print(materials_in_voxel, volumes_in_voxel)
  if len(materials_in_voxel)==1:
+ print('2 materials found')
  voxel_mat = materials_in_voxel[0].clone() #TODO find out why cloning crashes code
+ selective_flux_in_each_group.append(flux_in_group)
+ selective_micro_xs.append(micro_xs)
+
+ voxel_mat.volume = sum(volumes_in_voxel)
+ voxel_mat.id = i+mat_number_offset
+ voxel_mat.depletable = True
+ material_in_voxel.append(voxel_mat)
 
  elif len(materials_in_voxel)>1:
  # todo check this volume fraction is correct
 
  norm_fracs = [v*(1/sum(volumes_in_voxel)) for v in volumes_in_voxel]
- print('norm_fracs',norm_fracs)
+ print('1 material found')
  voxel_mat = openmc.Material.mix_materials(materials_in_voxel, norm_fracs)
+ selective_flux_in_each_group.append(flux_in_group)
+ selective_micro_xs.append(micro_xs)
+
+ voxel_mat.volume = sum(volumes_in_voxel)
+ voxel_mat.id = i+mat_number_offset
+ voxel_mat.depletable = True
+ material_in_voxel.append(voxel_mat)
+
  else:
  print('no material in voxel')
- voxel_mat = openmc.Material()
-
- voxel_mat.volume = sum(volumes_in_voxel)
- voxel_mat.id = i+mat_number_offset
- voxel_mat.depletable = True
- material_in_voxel.append(voxel_mat)
+
 
 openmc.lib.finalize()
 
@@ -210,13 +200,39 @@
 # constructing the operator, note we pass in the flux and micro xs
 operator = openmc.deplete.IndependentOperator(
  materials=openmc.Materials(material_in_voxel),
- fluxes=[flux[0] for flux in flux_in_each_group],
- micros=micro_xs,
+ fluxes=[flux[0] for flux in selective_flux_in_each_group],
+ micros=selective_micro_xs,
  reduce_chain=True, # reduced to only the isotopes present in depletable materials and their possible progeny
  reduce_chain_level=4,
  normalization_mode="source-rate"
 )
 
+# This section defines the neutron pulse schedule.
+# If the method made use of the CoupledOperator then there would need to be a
+# transport simulation for each timestep. However as the IndependentOperator is
+# used then just a single transport simulation is done, thus speeding up the
+# simulation considerably.
+hour_in_seconds = 60*60
+timesteps_and_source_rates = [
+ (1, 1e18), # 1 second
+ (hour_in_seconds, 0), # 1 hour
+ (1, 1e18), # 1 second
+ (hour_in_seconds, 0), # 2 hour
+ (1, 1e18), # 1 second
+ (hour_in_seconds, 0), # 3 hour
+ (1, 1e18), # 1 second
+ (hour_in_seconds, 0), # 4 hour
+ (1, 1e18), # 1 second
+ (hour_in_seconds, 0), # 4 hour
+ (1, 1e18), # 1 second
+ (hour_in_seconds, 0), # 4 hour
+ (1, 1e18), # 1 second
+ (hour_in_seconds, 0), # 5 hour
+]
+
+timesteps = [item[0] for item in timesteps_and_source_rates]
+source_rates = [item[1] for item in timesteps_and_source_rates]
+
 integrator = openmc.deplete.PredictorIntegrator(
  operator=operator,
  timesteps=timesteps,