Requirement list for ProcessSequences

This is an as detailed as possible outline of the requirement list for physics processes and process-lists:

1. The current (first) version of the physics process machinery is entirely designed on static C++ design patterns. However, we should aim for a parallel drop-in replacement that is entirely dynamic. Both can be very useful for different applications. And to develop a dynamic version of the static one might be interesting to some collaborators.

2. In all scenarios, individual processes must be stateless (if possible, for fortran this won't work) in order to allow parallel execution.

3. All modifications of particles on the stack should be done by "processes"

4. The collection of all processes defines the air shower physics

5. All processes can (or even must) be collected in one single C++ container, this we call ProcessList, however, ProcessList is a templated type and we would typically replace it with auto, and two different ProcessLists will not be simply type-compatible.

6. The ProcessList is one argument to the air shower machine of CORSIKA, i.e. ProcessList list=...; Cascade EAS(list,...);

7. We want ProcessLists to be mergeable, combineable to create more complicated sequences starting from simpler ones, thus logically: auto list1 = /*ProcessList*/; auto list2= /*ProcessList*/; auto list3 = list1 + list2;

8. Internally, we must distinguish different types of contained Processes which is driven by the fact that particles need to be modified in different ways:

   1. continuous processes modify properties of one particle, they need access to the current trajectory-step of the particle

   2. stochastic processes produce new particles, AND remove the current beam particle. This is typical for collisions, or decays. They do not need trajectory data.

   3. (we also need continuous processes that can also produce extra secondaries, e.g. Bremsstrahlung, maybe Cherenkov photons, etc. - but this is maybe beyond what we need to define here)

   4. processes that act on new secondaries (SecondariesProcess), typical applications are:

      1. thinning, re-weighting
      2. low-energy cuts etc.

   5. survey processes that act on all currently existing particles (from time-to-time) to determine combined statistical properties, etc. (energy conservation, run-time, ETA, sanity checks)

   6. boundary crossing process: when a particle is moving from one volume into another there may be interface-physics there (refraction, absorption, reflection, etc.)

9. The detailed interface to all the above processes is slightly different and driven by the underlying physics.

10. In addition we need a machinery to switch between two "branches" in the sequence based on phase-space regions, the most simple use-case is switching for hadrons between low- and high- energy interaction model. But this can include small process-sequences by itself (switching between sub-process lists)

11. processes lists are ordered. They are executed in the order of appearance in the code.

12. in particular thinking about point 6.4.1/2 we have to either further break down the interface, or we define a way to let "SecondariesProcesses" act only in specific cases. Consider a thinning SecondariesProcess, which should act on the output of collisons but not on the output of decays etc. Thus, logically we would like to distinguish interaction << decay << thinning from (interaction << thinning) << decay

Further discussions (to be resolved and merged in the requirements list above):

Notes from HD:

• ProcessSequences should not be created with operator<<, but just like std::tuple, by passing the elements to ctor. ProcessSequence(process1, process2, ...). The operator<< is not good for initialization.

• There are two kinds of thinning, which are distinctly different. There is the thinning that drops low-energy secondaries. This should be handled by a process, not by the sequence. A ThinnedInteractionProcess would derive from the generic InteractionProcess and apply thinning to its output. For the ProcessSequence this is transparent. The second kind of thinning happens after the shower has been fully simulated. To safe disk space, the shower core is either removed completely or further thinned as a function of lateral distance to the shower axis. This should be done by the particle writer, it is also not a responsibility of the ProcessSequence.

  • RU: indeed, there a physics-driven thinning that must happen inside the cascade - and it is up to the user to decide what and how exactly he wants to do that. There is not one thinning algorithm, there are many options with different applicability. In addition there is output-thinning just to save disk space. This is not physics-driven.

  • RU: there is one aspect we should keep in mind: thinning will be applied to a longish list of processes in the cascade, e.g. high-energy hadronic model, low-energy hadronic model, nuclear hadronic model, eventually decays, e.m. pair production, e.m. bremsstrahlung, and maybe others that we don't think about right now. So changing the algorithm should be safe against making simple mistakes of inconsistency. I think ThinnedInteractionProcess<ThinAlgo, InteractionAlgo> might be good, but also ProcessSequence sub_list = (interactionAlog << thinAlgo);. In the former case they are applied to a particular Process in the latter case to a corresponding piece of ProcessSequence. I think this is completely analogous and just a matter of taste.

  • HD to RU: No, it is not a matter of taste. Passing thinAlgo to the sequence is bad, because then it is the responsibility of the sequence to pass the thinning algorithm to the processes which need it. The flow of information is less clear and you prevent users from using different thinning algorithms for different processes. This is breaking modularity. The Thinning algorithm has to be provided to the InteractionProcess.

  • DB: I am in favour of some kind of Wrapper to apply different thinning algorithms to a Process. This would reduce code duplication quite a bit and is still easy enough to understand. The difficulty here is that the interaction process needs to forward the information which particle was created to the thinning algorithm.

  • HD to DB: You can have such a wrapper, but then the underlying processes must have a common interface that the wrapper can use. Which brings me back to my other point, that InteractionProcess and DecayProcess and possibly a future CherenkovProcess should inherit from a common base process.

• Interactions and Decays are very similar and called by the ProcessSequence in essentially the same way. This common aspect of the two classes should be factored into a base class. InteractionProcess and DecayProcess remain distinct, but have a common ancestor, perhaps called InteractionOrDecayProcess.

  • RU: if we want a common anchestor, it probably should be DiscreteProcess. However, what exactly should it do? The decay has "GetLifetime(projectile)" "DoDecay(projectile)" while the interaction has "GetInteractionLength(projectile)", "DoInteraction(projectile, target)". Note: the "target" in the last call is still missing but should be added for physics reasons.

  • HD to RU: Ok, DiscreteProcess then. The common ancestor would have GetStepLength and MakeSecondaries or something similar. That GetStepLength means something different for an InteractionProcess and a DecayProcess can be explained in the comments/docs.

• CherenkovProcess also shares functionality with DecayProcess and InteractionProcess. This could also possibility factored into a common base class (but I don't know enough about Cherenkov production to see whether this makes sense). Regarding thinning, there should be a gain a ThinnedCherenkovProcess which inherits from the common CherenkovProcess. The thinning algorithm for Cherenkov photons is likely different from the one applied to interactions, so these two should not be forced together.

  • RU: Cherenkov is a very good candidate for ContinuousProcess. In fact, Cherenkov radiation IS part of the continuous energy loss calculations a la Bethe-Bloch. Radio emission, Cherenkov, dEdX are all ContinuousProcesses. At least currently Cherenkov photons (as well as radio photons/waves) are not put on the particle stack in C7, they are directly processed and written to file. But is is up to discussion to eventually change this in C8 and also pipe it through the stack (?). This might be highly inefficient since we would all the time create a large number of very short-lived object dynamically on the stack...

  • HD to RU: The Stack should hold on to allocated memory and not give it back, to minimize the very costly allocations and deallocations. Ideally, the Stack would be placed on the stack (I am talking about what people normally call stack in C++, the non-dynamically allocated memory used by a program), by implementing it based on std::array or similar. The particles placed on the Stack should be trivially constructible types, so that copying a particle is a matter of copying bits around. Under these conditions, using the Stack could be fast enough. If the photons have to be stored somewhere anyway, they should be put on the stack.

  • DB: Cherenkov photons are in physics terms a ContinuousProcess even if they there contribution is rather small. For most implementations the exact track must be known (similar to Radio). Cherenkov Photons do require some sort of storage/stack... (see Requirement DynamicStack)