Nuclear fission - chain reaction

This model simulates the process of nuclear fission by chain reaction, such as occurs in a fission reactor or weapon. It models a sample of a fissile material, such as plutonium, which has the property that when a nucleus is hit by a neutron it splits, releasing two more neutrons which are capable of splitting other nuclei. In a fission reactor this process is continuous, and regulated by the insertion of rods of a material that absorbs neutrons, in order to release a controllable amount of energy. In a nuclear weapon, two pieces of fissile material that separately have too large a surface area to volume ratio to sustain a chain reaction are brought together to create a single piece that exceeds the "critical mass". This then undergoes a runaway chain reaction releasing all the energy at once in an explosion.

The model uses a 2-D world rather than 3-D simply because it makes the process easier to visualize. The submodel on the left represents a hexagonal grid of atoms, but rather than using the hexagonal grid submodel it is a population where the individuals get their grid coordinates from another submodel. This is because we want to remove individual atoms as the simulation runs, and we want their indices to be numbers rather than vectors, so they can be made the values of discrete events. (Later versions of Simile allow vector-valued events). 

The submodel on the right represents the neutrons. Initially there are none, but the fissile material also undergoes spontaneous radioactive decay, which produces neutrons. When an atom splits, the event is passed to the immigration channel for neutrons, causing two to be added. (Note that if population channels are driven from events, the exact number of individuals are added or removed immediately). Also a state variable outside the submodels is set to the index of the atom that split, so that when the new neutron-individuals are initialized they can get the position of the atom from the coordinate submodel to initialize their own positions. Their trajectory is set at random.

To handle splitting of atoms by neutrons, there is an association submodel between atoms and neutrons which has the existence condition that the distance between them is less than 1 unit. Inside this submodel is a limit event that triggers when the distance goes down to 0.5 unit. This event represents the actual splitting. It would be possible to have the existence condition be 0.5 unit and the event always true on initialization, but that would make the model behave differently. The existence of the association is only calculated on aregular time step boundary, so if the splitting condition was the same, all the splits would occur on time step boundaries, meaning several at once which would interfere with the mechanism for initializing neutron positions. Having the association come into existence before the actual split means that extrapolation is used to calculate the exact time of the split, so all splits occur separately. Unfortunately this makes the model run more slowly, since the whole model state is evaluated each time an event occurs.

The split event in the association is used to trigger the fission event in the corresponding atom, and to remove the neutron that does the splitting. The influences for these connections have been edited to disable the 'all values' interpretation, otherwise one split event would remove all neutrons and trigger fission in all atoms. Finally, neutrons are removed from the model if they move out of the lattice to speed up execution.

Update -- base instance lookup

A new version of the model has been added, using base instance lookup to make it run faster. The problem with the original version is that, at the height of the chain reaction there will be many hundreds of neutrons and many hundreds of remaining atoms of fissile material. This means that the membership condition of the association submodel 'Hit' must be checked hundreds of thousands of times each time step, causing the model to slow down at this point. The problem is especially severe because each event counts as a time step boundary, and fission events are much more frequent than the periodic time step at this point.

The solution is base instance lookup, also known as one-sided relation enumeration. Instead of checking every neutron's position against that of every fissile atom, we make use of the fact that the atoms' positions are fixed and regular to tell which one a neutron is closest to, then we need only check for colission with that atom. To do this we add 'col' and 'row' variables to the neutron submodel, with influences from the X and Y compartments respectively. The equations are

col = round(X/3)

row = round(Y/3.464 - 0.25)

These reverse the calculations done in the 'start posns' submodel to derive the atom positions from the index number. Now we add influences from these variables to the condition in the 'Hit' submodel, and set the condition's equation to 

index(1) is 32*(row-1)+col

The use of 'is' rather than '==' in the condition means that rather than comparing the index of each individual atom with that generated from the neutron position, we should simply pick the atom with that index. Note we no longer need an influence to the condition from 'dist sqd'. 

Finally we need to check 'allow base instance lookup' in the properties of the role arrow 'is' connecting the atom model to the Hit model. This ensures that index(1) in the Hit model (which can be looked up) corresponds to the index of the atom individual. Because the atoms are a population, there is still an iteration involved in finding the one with the given index (if it still exists) but this is quicker than calculating the distance for each atom.