Creation of a continuous UCN source

User Support
1

Dear colleagues,

For my thesis, I’m working on a UCN source that will constantly generate new UCN. The actual goal of my work is to make access to openmc easier. Good engineers are trained in Germany, but also in a very classic and conservative way. little reference to programming and the use of programs like openmc or openmfoam. I am developing a didactic concept for this.

I have designed my own source for this based on the UCN source from TRIGA in Mainz Germany. It shouldn’t be about the technical feasibility of my work!This roughly corresponds to the geometrical conditions of the HFR from the Laue-Langevin Institute.

The development of this source from the outside to the inside. The source is, for example, formed by 3 Zirkaloy bottles that are inserted into one another. There should be a vacuum between Zircaloy layer 1 and Zircaloy layer 2 to reduce heat transfer. Between Zirkaloy layer 2 and Zirkaloy layer 3 there is said to be a pre-moderator of mesitylene with a thickness that has yet to be determined. There should be frozen deuterium in Zirkaloy bottle 3 with a thickness that also needs to be determined. The aim is to determine the necessary ratio between the two materials in terms of thickness.

I created my own crossection library and determined the required materials with njoy at 20K and 120K. the pre-moderator layer should have 120K. the deuterium 20K. The question here is whether I need to go further down to 15K? I checked the geometry with various plots. it appears to me without any overlaps. I’m not sure about the temperatures and the material data calculated via njoy. Is my approach to interpolating these with Njoy correct?

The first bottle must be at the temperature of the reactor’s surroundings. I define the vacuum layer between bottles 1 and 2 as vacuum for the respective geometries. The 1st Zircaloy layer should have a temperature that corresponds to that of the reactor. The 2nd Zircaloy layer should correspond to the temperature of my pre-moderator of 120 k. I would like to set the 3rd Zircaloy layer at a flat rate of 20K. I will calculate the actual temperature curves and the wall thickness of the vacuum layer in more detail later.

what do I need help with? What is missing in my code to make it executable?

I always get WARNING: After particle 8055563 crossed surface 37 it could not be located in
WARNING: After particle 9305561 crossed surface 8 it could not be located in
WARNING: any cell and it did not leak.
After particle 7222230 crossed surface 37 it could not be located in WARNING: After particle 3333337 crossed surface 8 it could not be located in any cell and it did not leak.

Getting used to OpenMC was a steep process. I’ve managed to do a lot of things now.

However, I don’t know what to do now. I need a hint as to where I should look further. or asked in a completely amateurish way. If there is nothing fissile there, then nothing will happen? The source only sends thermal neutrons, which are then supposed to generate UCN through the premoderator and SD2 layer. I am grateful for every hint. The filters for counting the ucn and similar are not included in this code snippet. I know how that works. Please let me know if you have any further questions about understanding!

This is my code:

edit: settings.run_mode = ‘fixed source’ have to be set.
import openmc
import numpy as np

Materialien definieren (Zircaloy-4)

zircaloy = openmc.Material(name=‘Zircaloy-4’)
zircaloy.set_density(‘g/cm3’, 6.55)
zircaloy.add_nuclide(‘Zr90’, 0.5145)
zircaloy.add_nuclide(‘Zr91’, 0.1139)
zircaloy.add_nuclide(‘Zr92’, 0.1715)
zircaloy.add_nuclide(‘Zr94’, 0.0177)
zircaloy.add_nuclide(‘Zr96’, 0.1824)

Mesitylen als Pre-Moderator definieren

premod = openmc.Material(name=‘Mesitylene’)
premod.set_density(‘g/cm3’, 0.866) # Dichte von Mesitylen bei Raumtemperatur
premod.add_nuclide(‘C12’, 0.989) # Kohlenstoff-12
premod.add_nuclide(‘C13’, 0.011) # Kohlenstoff-13 0.011
premod.add_nuclide(‘H1’, 1)
premod.add_s_alpha_beta(‘c_H_in_Mesitylene’)
premod.temperature = 120

Festes Deuterium als SD2-Schicht definieren

sD2 = openmc.Material(name=‘Solid Deuterium’)
sD2.set_density(‘g/cm3’, 0.167) # Dichte von solidem Deuterium
sD2.add_nuclide(‘H2’, 1.0)
sD2.add_s_alpha_beta(‘c_para_D’)
sD2.temperature = 15.0

kantenlaenge = +8.0
#boxunten = openmc.ZPlane(z0=-55.0) # Untere Grenze 1
#boxoben = openmc.ZPlane(z0=+55.0) # Obere Grenze 2

box = openmc.model.RectangularPrism(kantenlaenge, kantenlaenge, boundary_type=‘reflective’) #6
#box_cell = openmc.Cell(fill=zircaloy, region=-box & +boxunten & -boxoben)

Zirka 1 mit Halbkugel

zirka1_ir = openmc.ZCylinder(r=3.60, boundary_type=‘vacuum’)
zirka1_or = openmc.ZCylinder(r=3.70)
zirka1 = openmc.Cell(fill=zircaloy, region=-zirka1_or & +zirka1_ir & +openmc.ZPlane(z0=-50) & -openmc.ZPlane(z0=+50))
ohalb1 = openmc.Sphere(r=3.7, x0=0, y0=0, z0=-50)
ihalb1 = openmc.Sphere(r=3.6, x0=0, y0=0, z0=-50, boundary_type=‘vacuum’)
halb1 = openmc.Cell (fill=zircaloy, region=+ihalb1 & -ohalb1 & -openmc.ZPlane(z0=-50))

Zirka2 mit Halbkugel

zirka2_ir = openmc.ZCylinder(r=3.30)
zirka2_or = openmc.ZCylinder(r=3.40, boundary_type=‘vacuum’)
zirka2 = openmc.Cell(fill=zircaloy, region=-zirka2_or & +zirka2_ir & +openmc.ZPlane(z0=-49.9) & -openmc.ZPlane(z0=+50.0))
ohalb2 = openmc.Sphere(r=3.4, x0=0, y0=0, z0=-49.9, boundary_type=‘vacuum’) #, name=‘Äussere Halbkugelzirka2unten’
ihalb2 = openmc.Sphere(r=3.3, x0=0, y0=0, z0=-49.9,) #, name=‘Innere Halbkugelzirka2unten’
halb2 = openmc.Cell (fill=zircaloy, region=+ihalb2 & -ohalb2 & -openmc.ZPlane(z0=-49.9))

Vakuum zwischen Zirka 1 und Zirka 2 mit Halbkugel (unnötig da vacuum in randbedingungen)

#vac = openmc.Cell(fill=None, region=-zirka1_ir & +zirka2_or & +openmc.ZPlane(z0=-50.0) & -openmc.ZPlane(z0=50.0))
#vachk = openmc.Cell(fill=None, region=-ihalb1 & +ohalb2 &-openmc.ZPlane(z0=-50))

Zirka 3 mit Halbkugel

zirka3_ir = openmc.ZCylinder(r=2.9)
zirka3_or = openmc.ZCylinder(r=3.0)
zirka3 = openmc.Cell(fill=zircaloy, region=+zirka3_ir & -zirka3_or & +openmc.ZPlane(z0=-49.8) & -openmc.ZPlane(z0=50.0))
ohalb3 = openmc.Sphere(r=3.0, x0=0, y0=0, z0=-49.8,) #, name=‘Äussere Halbkugelzirka3unten’
ihalb3 = openmc.Sphere(r=2.9, x0=0, y0=0, z0=-49.8) #, name=‘Innere Halbkugelzirka3unten’
halb3= openmc.Cell (fill=zircaloy, region=+ihalb3 & -ohalb3 & -openmc.ZPlane(z0=-49.8))

preomd mit Halbkugel

premodz = openmc.Cell(fill=premod, region=-zirka2_ir & +zirka3_or & +openmc.ZPlane(z0=-49.9) & -openmc.ZPlane(z0=50.0))
halb_premod = openmc.Cell(fill=premod, region=-ihalb2 & +ohalb3 & -openmc.ZPlane(z0=-49.9))

#Innenraum mit Halbkugel
innen_or = openmc.ZCylinder(r=0.25, boundary_type=‘vacuum’ )
innen = openmc.Cell(fill=None, region=-innen_or & +openmc.ZPlane(z0=-49.8) & -openmc.ZPlane(z0=+50.0))
halbinnen_ir = openmc.Sphere(r=0.25, x0=0, y0=0, z0=-49.8,)
halbI = openmc.Cell(fill=None, region=-halbinnen_ir & -openmc.ZPlane(z0=-49.8))

Sd2 mit Halbkugel

sD2s = openmc.Cell(fill=sD2, region=-zirka3_ir & +innen_or & +openmc.ZPlane(z0=-49.8) & -openmc.ZPlane(z0=+50.0))
halb_sD2 = openmc.Cell(fill=sD2, region=-ihalb3 & +halbinnen_ir & -openmc.ZPlane(z0=-49.8))

geometrie = openmc.Geometry([zirka1, zirka2,zirka3, innen, premodz, halb_premod, halb_sD2, sD2s, halbI, halb3, halb2, halb1])
geometrie.export_to_xml()
materials = openmc.Materials([zircaloy, premod, sD2])
materials.export_to_xml()

Simulationssettings definieren (nach dem Exportieren von Geometrie und Materialien)

settings = openmc.Settings()
settings.batches = 1000
settings.inactive = 100
settings.particles = 10000
settings.max_particle_events = 100000

Quelle

lower_left=(-kantenlaenge / 2, -kantenlaenge / 2, -55) #Zirka1Laenge / 2
upper_right=(kantenlaenge / 2, kantenlaenge / 2, 55) #-Zirka1Laenge / 2
uniform_dist = openmc.stats.Box(lower_left, upper_right)
quelle=openmc.IndependentSource(space=uniform_dist)
quelle.angle = openmc.stats.Isotropic()
energiewerte = [0.025] # in eV
gewichte = [1.0] # Gewichtung der Energie
quelle.energy = openmc.stats.Discrete(energiewerte, gewichte)
settings.source = quelle

Einstellungen exportieren

settings.export_to_xml()

#print(“Region von zirka1:”)
#print(zirka1.region)

#print(“\Region von zirka2:”)
#print(zirka2.region)

#print(“\Region von zirka3:”)
#print(zirka3.region)

#print(“\Region von halb1:”)
#print(halb1.region)

#print(“\Region von halb2:”)
#rint(halb2.region)

#print(“\Region von halb3:”)
#print(halb3.region)

#print(“\Region von halb_sD2:”)
#print(halb_sD2.region)

#print(“\Region von sD2:”)
#print(sD2s.region)

#print(“\Region von halbI:”)
#print(halbI.region)

#print(“\Region von Innen:”)
#print(innen.region)

#print(“\Region von box_cell:”)
#print(box_cell.region)

plot = openmc.Plot()
plot.filename = ‘geometry_plot’ # Name der Ausgabedatei
plot.width = (150.0, 150.0) # Breite und Höhe des Plots
plot.origin = (0.0, 0.0, 0.0) # Ursprung des Plots
plot.basis = ‘yz’ # Ebene für den Plot (xy-Ebene)
plot.pixels = (10000, 10000)

Farben definieren

plot.color_by = ‘cell’ # Färbung nach Zelle
plot.colors = {
#box_cell: ‘blue’, # Setze die Box-Zelle auf blau
zirka1: ‘black’,
halb1: ‘black’,
zirka2: ‘red’,
halb2: ‘red’,
#zirka3: ‘blue’,
halb3: ‘blue’,
#vac: ‘yellow’,
#vachk: ‘yellow’,
premodz: ‘brown’,
halb_premod: ‘brown’,
sD2s: ‘orange’,
halb_sD2: ‘orange’,
innen: ‘purple’,

halbI: 'purple',

}

Plotsammlung erstellen und exportieren

plots = openmc.Plots([plot])
plots.export_to_xml()

Bottle
Bottle753×752 81.4 KB

2

Hi Sascha, welcome to the openmc community.
I think the openmc has told you that the problem is in your cell definition.
like when openmc said

WARNING: After particle 8055563 crossed surface 37 it could not be located in
WARNING: After particle 9305561 crossed surface 8 it could not be located in
WARNING: any cell and it did not leak.
After particle 7222230 crossed surface 37 it could not be located in WARNING: After particle 3333337 crossed surface 8 it could not be located in any cell and it did not leak.

then you need to check whether any surface or cell has an improper definition.

When you have void cells, I think you didn’t need to add vacuum BC, but you need to add the void cells. so the vacuum BC is only used for the outer surface of your model.
I also see that your cell definition is incomplete, especially at the edge of your zirlo material, which has been moved 0.01 mm above.
image

here is the notebook that I use to check your geometry, I am adding and modifying some of your cell definitions, and I hope that it can help you
UCN.ipynb (332.5 KB)

3

Hi wahidluthfi,

thanks for your advice. I will work on that. keep you updated on my progress. I was bit worried about getting no response. I worked a lot on the examples. I also played around with different geometries to get a better understanding. furthermore to get a feeling, what is a good didactic way for newbies, after that I started on this design. have read different paper. in case of needed I will increase the thickness of booth layers a lot to improve the possibility of appearance of ucn. iam sure the layers aren’t thick enough. but first I will work on your advices. sascha

1 Like
4

it worked. corrected the geometries. killed the box. understood the scource new. understood that you have to define the zplanes first. with the boundary condition. the inersection stuff for zirka1, zirka2 and zirka3 did the trick. i thought about this, but couldnt connect it intellectual to my problem. had i feeling about this as part of the solution. understood again how important clean structured comments are. many, many thanks. will work further through my task. but so far, great. it was really, really helpfull.

1 Like
5

Glad to hear that,
You can also move the ‘small but tricky region’ that I add in your zirlo cell to the cell beside it, i.e, the void and moderator if you want because I think this should be part of the void and moderator cells in real life.
I am adding that small region to the zirlo because it is a lot easier to imagine by me if that region is part of the zirlo at first. and yeah you can move it.

so, instead of changing the zirlo cells, you can focus on the void (vac) and the pre-moderator (premodz).

image
image1241×677 11.4 KB

you can see the slight curve here in the yellow and maroon cells.
image
image1209×538 6.54 KB

here is the script for this modification

# Vakuum zwischen Zirka 1 und Zirka 2 mit Halbkugel (unnötig da vacuum in randbedingungen)
vac = openmc.Cell(name='vac', fill=None, region=(-zirka1_ir & +zirka2_or & +zplane50min & -zplane50plus) | (-zirka2_or & +ohalb2 & +zplane50min & -zplane499min))
vachk = openmc.Cell(name='vachk', fill=None, region=-ihalb1 & +ohalb2 & -zplane50min)

# Zirka2 mit Halbkugel
zirka2 = openmc.Cell(name='zirka2', fill=zircaloy, region=-zirka2_or & +zirka2_ir &  +zplane499min & -zplane50plus)
halb2 = openmc.Cell (name='halb2', fill=zircaloy, region=-ohalb2 & +ihalb2 & -zplane499min)

# preomd mit Halbkugel
premodz = openmc.Cell(name='premodz', fill=premod, region=(-zirka2_ir & +zirka3_or & +zplane499min & -zplane50plus) | (-zirka3_or &  +ohalb3 & +zplane499min & -zplane498min))
halb_premod = openmc.Cell(name='halb_premod', fill=premod, region=-ihalb2 & +ohalb3 & -zplane499min)

# Zirka 3 mit Halbkugel
zirka3 = openmc.Cell(name='zirka3', fill=zircaloy, region=-zirka3_or & +zirka3_ir & +zplane498min & -zplane50plus)
halb3= openmc.Cell (name='halb3', fill=zircaloy, region=-ohalb3 & +ihalb3 & -zplane498min)
6

thanks again. that’s very helpful. I will change this tomorrow. many, many thanks. great input and support.

7

Dear @Sascha

Aside of the geometry problems you should be careful with the definitions of the materials, in particular with the use of scattering kernels for the right state at the right temperature. For solid deuterium we have developed a new scattering kernel at 5K that might be useful.

https://git.esss.dk/spallation-physics-group/solid-deuterium

Also, the interaction with Zr and other materials in the beam should probably be better computed with NCrystal materials in OpenMC:

https://docs.openmc.org/en/latest/methods/cross_sections.html#ncrystal-materials

8

hi ignacio,

yeah I have thought also about this. its good to get helpful advices like this.
its difficult for me to imagine that a delta Theta of 15 K should have such an impact, but yeah I will include this.

i will keep you guys updated. its interesting to see which parameters have which influence on the simulation. I ran 2 cycles of my code. 1000 batches a million particles. had a huge leakage 0.94 at the first run. at the second run with premed and sd2 more even its reduced to 0.7. both times no ucn. I didn’t expected that. changed the dimensions to more realistic layer thickness. all with the 20K deuterium. made a third calculation. with 8 h . leakage down to. 0.11.
still no ucn.

i also had a error warning in joy as I calculated Zr96 and I have a warning for negative values at the startup of openmc. so your second link is highly welcome. It will take a lot of tweaks and thinking to get a ucn. waiting for an appointment happen with my professor. Thank you guys. Sascha

9

Well, maybe you should start with a simplified model and try to calculate the equilibrium spectra inside of the source, which will be computationally easier. Also, I do not know what is your definition of UCN (different people have different definitions), but OpenMC does not transport below 1e-5 eV (~900 \AA). It might be possible to extend the code further, but it would also require to add proper UCN physics (for instance, to include the interaction with the walls).