TrackParticles has different fusion neutron weights at low energies

[philmartin01]

Hi,

I am running simulations of DD fusion and am trying to extract spectra from tracked neutrons. However, when attempting this there seems to be an issue when summing the weights of neutrons at low energies (<2.5 MeV). I have attached a plot with 2 spectra, one extracted from a ParticleBinning diagnostic (black), and the other was created from tracked particles (red). Both diagnostics include all neutrons in the simulation box, contained within the stated angular range in the plot, although this problem is still seen even when including all neutrons in the simulation.
Is there something that I'm missing when creating the spectrum with tracked particles? I am creating it with a histogram and summing their weights in each bin, which should work, considering the rest of the spectrum above 2.5 MeV matches fine.

[mccoys]

Please provide some more details on how your diagnostics where implemented and post-processed. As the two types of diagnostics process the same particles, they should match.

[philmartin01]

ParticleBinning diagnostic:

def angles:
return np.arctan2(particles.py, particles.px)
DiagParticleBinning(
name = 'neutron angle',
deposited_quantity = "weight",
every = every,
species = ["neutron"],
axes = [
[angles , -1*math.pi, math.pi, 360],
["ekin", 0., 50, 100, "edge_inclusive"]
]
)

And then to get spectrum I extracted the data and integrated across the appropriate angular bins. With TrackParticles, I output all variables and create a similar histogram with scipy.stats.binned_statistic_dd, something like:

stats_data = stats.binned_statistic_dd([Theta, En], values = wn, statistic = 'sum', bins=[theta_bins, energy_bins])

where Theta and energy were calculated from the tracked Px and Py components, and wn is the macroparticle weights.

[mccoys]

My guess is that your method with tracked particles divides by different bin sizes than that done in happi for particle binning. As you are using a custom axis, it cannot know the geometry of your binning. It's just a guess

[philmartin01]

Bin sizes are the same for both tracking and particlebinning, so all energies would be affected equally when dividing by bin size, not just those <2.5 MeV. I cannot see anything strange with the weights, their values are broadly consistent above and below 2.5MeV, just that there are a significant number of tracked neutrons in that low energy range, perhaps that isn't being counted by the binning diag? Or their weights are calculated differently between the two diagnostics.

[mccoys]

I have no clue what is going on here. The particle binning really just does a histogram and divides by the bin size which is indeed constant in your case.

Just in case you may want to check the raw data in the hdf5 file.

[mccoys]

Just an idea that crossed my mind. You wrote "Theta and energy were calculated from the tracked Px and Py components". But there should be Pz as well for calculating the energy. Did you take it into account?

Another detail: what is the geometry, and what are boundary conditions?

[philmartin01]

Geometry is 2D, so I didn't output Pz component. Boundaries are silver-muller.

[mccoys]

2d means x and y, but there is definitely a contribution from px.

I meant boundary conditions for particles

[philmartin01]

Ah ok, I can rerun the simulation outputting Pz components as well and test it out. Boundaries were open.

[philmartin01]

Hi, including the Pz component has resolved the issue, thanks!