Monte Carlo workflow

In order to simulate the \(\nu_\tau-\tau\) transport through a section of the Earth, Danton employs a Monte Carlo technique. In the traditional approach, primary neutrinos are injected at the top of the atmosphere and the resulting tau decays are collected.

In addition to the traditional forward method, Danton also allows for the backward sampling of tau decays. In practice, tau leptons are injected at decay points of interest, and potential primary neutrinos are collected at the top of the atmosphere. In the backward case, the candidate primary neutrinos have associated Monte Carlo weights indicating their respective likelihoods.

Backward or Forward

The default setting for Danton is a backward Monte Carlo simulation. To perform a conventional forward simulation, one must set the mode attribute of the simulation to "forward".

Backward sampling can enhance CPU efficiency by orders of magnitude in certain use cases (see, e.g., [NM18]). However, this approach is more more intricate, given that nature progresses forwards.

Decays or fluxes

The default configuration of Danton is to sample decay densities. However, it is also possible to sample fluxes of particles instead, should this be required. This can be achieved by setting the tau_decays attribute of the simulation to False.

Longitudinal approximation

By default, Danton performs a comprehensive three-dimensional simulation of particle transport. However, at ultra-high energies, the deflection angles w.r.t. the initial direction of the neutrinos are typically minimal, and thus often disregarded. This approximation is activated by setting the longitudinal attribute of the simulation to True.

Target point

Danton’s approach is to assume that primary neutrinos always originate from the top of the atmosphere. In the forward Monte Carlo case, for the sake of convenience, the input position is considered to be a target point for the incoming neutrino, in lieu of its actual injection point. Prior to running the Monte Carlo, the neutrino is repositioned to the top of the atmosphere by extrapolating its track backwards from the target point.

Note

The relocation of primary neutrinos is purely geometrical. Thus, the simulated neutrinos may not intersect their target point, for example, if the longitudinal approximation is disabled.

Point estimate

A salient property of backward Monte Carlo is its ability to perform a point estimate for a final state of interest. To illustrate, Danton may compute the density of tau decays at a specific location on Earth, with a given direction and tau energy at decay. This is obtained through a Monte Carlo weighted average of the fluxes of backward sampled neutrinos. For instance,

simulation = danton.Simulation(mode="backward")
taus = danton.particles(N, # Create N identical particles.
    energy = 1E+10, # GeV
    altitude = 100, # m
    elevation = 1   # deg
)
primaries = simulation.run(taus).primaries # Backward propagate the taus.
density = sum(primaries["weight"] * flux(primaries["energy"])) / N # Weighted average.

Monte Carlo integration

In a typical scenario, the objective is to collect particles within a specified region of interest, e.g., one that matches the field of view of a detector. This requires performing a Monte Carlo integration over initial (final) states in the forward (backward) case. To elaborate, one would generate incoming neutrinos over the top of the atmosphere, in the forward case (or tau decay vertices over the region of interest, in the backward case).

In order to facilitate these generation procedures, Danton provides a ParticlesGenerator object, which is instantiated using the particles method of the Simulation class. The particles generator is configured using a builder pattern, and then triggered with the generate method. For instance, the following

particles = simulation.particles()       \
    .powerlaw(1E+06, 1E+12, exponent=-2) \
    .solid_angle(elevation=[-1, 1])      \
    .generate(N)

would result in the generation of N particles with energies distributed according to a \(1 / E^2\) power-law (between 1 PeV and 1 EeV), and across a solid angle encompassing elevation values between -1 and 1 deg.

Generation weight

The default setting is to weight particles by the inverse of their generation likelihood (\(\omega = 1 / \text{pdf}(\text{S})\), for a Monte Carlo state \(\text{S}\)). The weight argument of the particles method can be set to False to disable weighting. It should be noted that weighting can also be enabled or disabled individually at the level of each generator method.

Sample size

The generator may employ rejection sampling techniques (depending on its configuration), resulting in a Monte Carlo sample that is larger than the number of returned particles (N). In this case, along with the selected particles, the generate method also returns the sample size as a second argument.

Box region

The Box object allows one to specify a region of interest by means of a bounding box. This region may subsequently be targeted by a generator, e.g. when employing forward Monte Carlo. Conversely, in the backward case, the generator can be configured to sample decay vertices inside the box.

Monte Carlo history

Danton simulations are pseudo-random, relying on a Pseudo-Random Number Generator (PRNG) as the source of randomness (a Permuted Congruential Generator (PCG) namely). The generator is, by default, seeded from the operating system’s entropy. The seed value, in conjunction with the PRNG stream index, determines the subsequent fate of a Monte Carlo particle. This can be exploited to re-simulate a specific event, as Danton systematically logs the PRNG stream index in addition to the Monte Carlo particles data.

A typical use case is to re-simulate selected events by enabling a detailed log of Monte Carlo steps (by setting the record_steps attribute of the simulation to True).

Note

It is inadvisable to systematically record Monte Carlo steps for large samples due to the resulting memory overhead.