I want to write my own Monte Carlo code for Neutron transport

Hi, i would like to write my own MC code in order to simulate the transport of Neutrons in Nuclear reactors. I know the basics of MC and i have already written a code for homogeneus reactors, my problem is the generalization to more complex geometries made of different materials, such as fuel rods and lattices.

So my question is: there exists some libraries (in Matlab or Python environments possibly) or "easy" methodologies to define different regions in the reactor? To be more precise, how can i effectively tell the neutron in which region of the reactor it is and so which nuclide is it colliding with?

My hope is to find some pre-existing functions devoted to the geometry realization (cylinders or other simple shapes) that provide me some sort of "region label" that i can use for the identification of the appropriate cross section.

Any suggestion is really appreciated.

mathman
Many years ago (I'm long since retired) my work involved problem you described. The main issue is finding 3 dimensional modeling codes. MCNP https://mcnp.lanl.gov/ seems to be widely used.

jim mcnamara
You can search the Physics Forums (and of course, you can go to google) for Monte Carlo code and find numerous examples.

hutchphd
Homework Helper
A reactor is a pretty complicated piece of geometry, and neutrons are notoriously energy dependent in their ineractiions. Perhaps a less complicated interrim problem would be useful. I don't know enough Matlab or Python stuff to suggest a good one.

You should describe your previous achievements in this modeling.

This is wildly good by Forrest Brown. https://mcnp.lanl.gov/pdf_files/la-ur-16-29043.pdf

As for the location stuff, you need more than location ID, you need distance to the next cell boundary. So in effect a ray tracer. You might be better starting with a 2D model with limited depth and square rod cross sections. When you move to 3D you'll need a set of equations for every shape.

Thank you all for your replies,

Many years ago (I'm long since retired) my work involved problem you described. The main issue is finding 3 dimensional modeling codes. MCNP https://mcnp.lanl.gov/ seems to be widely used.
Thank you, i heard about MCNP but i don't have a licence. However I have some experience with a similar code named Serpent. But my goal is to write my own.

You can search the Physics Forums (and of course, you can go to google) for Monte Carlo code and find numerous examples.
I tried to do so but unfortunately i wasn't able to find anything useful about this specific topic.

A reactor is a pretty complicated piece of geometry, and neutrons are notoriously energy dependent in their ineractiions. Perhaps a less complicated interrim problem would be useful. I don't know enough Matlab or Python stuff to suggest a good one.

You should describe your previous achievements in this modeling.
Yes, indeed it is quite complicated but i would be happy also with simpler geometries, however my code should be as general as possible. Up to now i succeeded in 3D homogeneus core and in 1D multi region slab but in these situations is easy to identify in which region of the geometry neutrons are colliding. For the slab, for example, i just need to see in which interval neutron coordinate is: x_i < x < x_{i+1}.

This is wildly good by Forrest Brown. https://mcnp.lanl.gov/pdf_files/la-ur-16-29043.pdf

As for the location stuff, you need more than location ID, you need distance to the next cell boundary. So in effect a ray tracer. You might be better starting with a 2D model with limited depth and square rod cross sections. When you move to 3D you'll need a set of equations for every shape.
Thank you for the document, it seems to be very interesting. Probably you are right, i should start with a simpler model. But even in 2D with square objects i need some functions which provide me the region ID (and for the geometry creation). Did you think that i should write these functions by my own or i can exploit some pre-existing code?
For what concern ray tracing i'm planning to use Delta-tracking method, doing so i don't need to re-sample neutron path when it cross a surface.

I had not heard of Woodcock delta-tracking. It may take a while to sink in.

For rods you can specify using if then statements for x and y. A lattice of rods is then also fairly easy. Since all surfaces are at 90 degrees to every other surface geometry math is simpler too.

https://github.com/openmc-dev/openmc/
If you are searching for a code to read try this one. Older versions were written in FORTRAN and had MCNP like problem entry. Now things are in C++ and Python, and Python is the method for problem entry. Either way it's a proper neutronics code in 3D.

mark_bose
I'm not aware of any general libraries that are available. If you want some examples, you could start looking in the source code for Serpent or OpenMC. Geometry processing is fairly complicated, and the more geometry types you support, the more complicated it becomes.

As you are probably aware from working with Serpent, the "standard" way of representing geometry in Monte Carlo codes is combinatorial geometry where you have a set of surfaces and you define volumes using the surfaces.

My suggestion is to start simple. For example, start out by supporting planes and spheres. Start with a sphere in a box with reflective boundary conditions. In this case you have 7 different surfaces (6 planes and the sphere). You have a particle at a point in space with a certain direction. The particle is in one of two volumes - inside the sphere or outside the sphere. You should calculate the distance from the point in space to each surface along the direction of travel. This takes quite a bit of geometry calculations, but you can use Google to find formulas for the distance to different surfaces. You build up a list of distances, and then find the minimum distance to find the nearest surface (the particle may not intersect all of the surfaces). After this, move the particle. It may or may not reach the surface. If it goes past the surface, you need find the next voluime, and then start your search over again. Once it leaves the surface, you also have to find out what the new volume is.

The real trick is to make this work efficiently. As you probably figured out, as your problem gets larger, it will get harder to find what volume the particle is currently in. You can start to build up internal lists of what volumes are next to each other.

Another good place to look is in video game design. Modern video games do a lot of "ray tracing" and they've come up with some very good algorithms to find the distance to surfaces and search methods.

I hope this helps. Good luck!

Alex A, hutchphd, sysprog and 3 others
gleem
My question is why? I am not a reactor physicist, but I would check to see if something is available for free that you could use instead of trying to reinvent the wheel. As @rpp suggests. For example, Open MC seems like a good place to start.

sysprog
Thank you all for the suggestions. Probably i will start with the simplifications you suggest. It seems that "Theory and Methodology" chapter in OpenMC manual could be useful.

berkeman, hutchphd and sysprog
bigfooted
Gold Member
If you want to play with something yourself, then you need multiblock mesh support. In python, you could work with pyvista and use paraview multiblock meshes to set up something. This also makes it easy to generate cell averaged means and visualize your results with pyvista

https://docs.pyvista.org/

mark_bose
There is a Russian code named Mcu with a shitton of functions. I can share it for you but documentation is in Russian. Write me in PM if you are interested.