Cell search via a tree

[rzehumat]

(note I am new to OpenMC, so sorry in advance if this idea is non-sense; I also did not make much research about this yet)

Hi,

the cell search algorithm is described such that when the universe is known

loop over each cell within that universe. For each cell, we check whether the specified point is inside the cell using the algorithm described in Finding a Lattice Tile.

This is intuitive. However, it’s possible that some cells are checked multiple times.

I’ll use the BWR example’s, universe for fuel pin

# Geometry definitions for the fuel rod

fuel_or  = openmc.ZCylinder(r=1.0414/2, name='Fuel OR')
fclad_ir = openmc.ZCylinder(r=1.06426/2, name='Clad IR')
fclad_or = openmc.ZCylinder(r=1.22682/2, name='Clad OR')

fuel_region = -fuel_or 
gap_region  = +fuel_or & -fclad_ir
fclad_region  = +fclad_ir & -fclad_or
fwater_region = pin_cell_box & +fclad_or

fuel_cell = openmc.Cell(name='fuel')
fuel_cell.fill = uo2
fuel_cell.region = fuel_region 

gap_cell = openmc.Cell(name='air gap')
gap_cell.region = gap_region

clad_cell = openmc.Cell(name='clad')
clad_cell.fill = zircaloy
clad_cell.region = fclad_region

fwater_cell = openmc.Cell(name='fwater')
fwater_cell.fill = water
fwater_cell.region = fwater_region

fuel_pin_universe = openmc.Universe(cells=[fuel_cell, gap_cell, clad_cell, fwater_cell])

Not sure how are the cells and surfaces arranged, but if the order was [fuel_cell, gap_cell, clad_cell, fwater_cell], the possible scenarios in the cell search routine would be (each bullet point corresponds to possible scenario depending on where the particle is)

• test_fuel_cell(-fuel_or) … return fuel_cell
• test_fuel_cell(+fuel_or), test_gap_cell(+fuel_or, -fclad_ir) … return gap_cell
• test_fuel_cell(+fuel_or), test_gap_cell(+fuel_or, +fclad_ir), test_clad_cell(+fclad_ir -fclad_or) … return clad_cell
• test_fuel_cell(+fuel_or), test_gap_cell(+fuel_or), +fclad_ir), test_clad_cell(+fclad_ir +fclad_or), test_fwater_cell(pin_cell_box, +fclad_or) … return fwater_cell

The scenario above is rather pessimistic – thanks to the short-circuiting (PR 2178), it could be simpler:

• test_fuel_cell(-fuel_or) … return fuel_cell
• test_fuel_cell(+fuel_or), test_gap_cell(+fuel_or, -fclad_ir) … return gap_cell
• test_fuel_cell(+fuel_or), test_gap_cell(+fclad_ir), test_clad_cell(+fclad_ir, -fclad_or) … return clad_cell
• test_fuel_cell(+fuel_or), test_gap_cell(+fclad_ir), test_clad_cell(+fclad_or), test_fwater_cell(pin_cell_box, +fclad_or) … return fwater_cell

Still, especially in cases when the particle is in cladding or water, some surface tests are repetitive

• fuel_or is tested 2 times if particle is in gap_cell,
• fclad_ir is tested 2 times if particle is in clad_cell,
• fclad_or is tested 2 times if particle is in fwater_cell.

I was wondering if cell search could be done via binary tree of surfaces, as illuestrated in figure below, where branch up is + and down is -. In such case, the scenarios would look like

• -fclad_ir, -fuel_or … return fuel_cell
• -fclad_ir, +fuel_or … return gap_cell
• +fclad_ir, -fclad_or … return clad_cell
• +fclad_ir, +fclad_or … return fwater_cell

The advantage of such approach would be

• theoretically faster cell search – best scenario would be O(log_2(N_surfaces)) surface tests, while the current would be O(N_cells*N_surfaces_per_cell) or something like that

The disadvantages would be

• likely more compicated data structure to store surfaces and cells,
• information about arrangement of surfaces – e.g. tree above was done with knowledge that fclad_ir is inside fclad_or and fuel_or is inside fclad_ir, so that they don’t need to be tested if the particle is outside of the outer surface,
• the balanced tree above is best case scenario when surfaces are inside each other – but it might detoriate from tree to line and ending de facto with the same algorithm as now,
• it does not check for potentially undefined regions,
• in some cases, there are more surface checks – in the BWR example, this is true for fuel_cell, which would have just 1 surface check, but now has 2… Technically, the tree could be manually biased to check some surfaces sooner,
• it would be useful only in cases where cell search occupies most

Let me know if such idea is worth further investigation or not – or if someone has already investigated this (which I guess someone has done already). If it would be worth trying, I’ll try implementing that – but I did not want to do unless I have some confidence it’s worth it.

Thanks for any thoughts on that.

[vitmog]

Hi @rzehumat,

This theme is interesting for me, but I’m not an OpenMC developer, so I will speak just from general considerations.

You have noted absolutely right about the principal possibility of improvement of the CSG cell search scheme via the performed binary search algorithm. However, in my opinion, exactly this scheme will not be efficient due to the following reasons:

1. If the universe concept is used for creation of a tree structure (breaking of a entire model) then, as a rule, O(log_2(N_surfaces)) is quite close to O(N_cells*N_surfaces_per_cell) for practical problems (see your example).

2. In well-designed code, the computations of a current cell are required just in the beginning of history for the fixed-source problems and in the beginning of computation for the criticality problems – further computations (searches) are redundant. I don’t get yet how it currently works in OpenMC.

3. Anyway, as a rule, the spent computational resources to the calculation of testing point-surface position are not high because the calculation doesn’t require even the floating point division operations (the fp-multiplications are the heaviest ones). As a rule, the most resources of geometry processing are spent to the ray tracing when each calculation of the ray distance to a bounding surface requires the fp-divisions and, except for planar ones, at least one the square root operation – one of the most greedy to resources.

Although, it is principally possible to use specialized non-CSG geometry modules which provide limited but efficient geometry processing including optimizations like the one proposed by you. As examples, they can be the lattices implemented in OpenMC – extremely limited but very efficient. In the same manner, it can be something like a parametrized and absolutely optimized (but predominantly for ray tracing) “fuel pin in a rectangle cell” geometry block.

In my opinion, the noted by you approach with ad hoc optimizations is a potentially perspective direction of Monte Carlo particle transport geometries development, but I’m not sure that such the optimizations of geometry calculations can dramatically improve the usability of the entire code, because the expected computation cost gains are limited just by several times in an optimistic (maybe even unreal) case from my experience. However, it also may simplify the geometry building process due to the parametrization of model parts and, consequently, improve the usability.