From e8698189eebb651ca664508aa5270f1b2bbe4aad Mon Sep 17 00:00:00 2001 From: ralfulrich <ralf.ulrich@kit.edu> Date: Wed, 13 Oct 2021 16:27:26 +0200 Subject: [PATCH] cascadeequations process --- corsika/detail/framework/core/Cascade.inl | 5 +++- .../framework/process/ProcessSequence.inl | 24 +++++++++++++++++++ corsika/detail/modules/conex/CONEXhybrid.inl | 4 ++-- corsika/framework/process/ProcessSequence.hpp | 7 +++++- corsika/modules/conex/CONEXhybrid.hpp | 9 ++++--- examples/hybrid_MC.cpp | 2 -- tests/modules/testCONEX.cpp | 2 +- 7 files changed, 43 insertions(+), 10 deletions(-) diff --git a/corsika/detail/framework/core/Cascade.inl b/corsika/detail/framework/core/Cascade.inl index 2bb692d06..64ca0a3aa 100644 --- a/corsika/detail/framework/core/Cascade.inl +++ b/corsika/detail/framework/core/Cascade.inl @@ -31,10 +31,12 @@ namespace corsika { setNodes(); // put each particle on stack in correct environment volume while (!stack_.isEmpty()) { + + sequence_.initCascadeEquations(); + while (!stack_.isEmpty()) { CORSIKA_LOG_TRACE("Stack: {}", stack_.asString()); count_++; - auto pNext = stack_.getNextParticle(); CORSIKA_LOG_TRACE( @@ -46,6 +48,7 @@ namespace corsika { step(pNext); sequence_.doStack(stack_); } + // do cascade equations, which can put new particles on Stack, // thus, the double loop sequence_.doCascadeEquations(stack_); diff --git a/corsika/detail/framework/process/ProcessSequence.inl b/corsika/detail/framework/process/ProcessSequence.inl index 78a28da60..4d381abb1 100644 --- a/corsika/detail/framework/process/ProcessSequence.inl +++ b/corsika/detail/framework/process/ProcessSequence.inl @@ -381,6 +381,30 @@ namespace corsika { } } + template <typename TProcess1, typename TProcess2, int IndexStart, int IndexProcess1, + int IndexProcess2> + void ProcessSequence<TProcess1, TProcess2, IndexStart, IndexProcess1, + IndexProcess2>::initCascadeEquations() { + + if constexpr (is_process_v<process1_type>) { // to protect from further compiler + // errors if process1_type is invalid + if constexpr ((process1_type::is_process_sequence && + !process1_type::is_switch_process_sequence) || + is_cascade_equations_process_v<process1_type>) { + A_.initCascadeEquations(); + } + } + + if constexpr (is_process_v<process2_type>) { // to protect from further compiler + // errors if process2_type is invalid + if constexpr ((process2_type::is_process_sequence && + !process2_type::is_switch_process_sequence) || + is_cascade_equations_process_v<process2_type>) { + B_.initCascadeEquations(); + } + } + } // namespace corsika + template <typename TProcess1, typename TProcess2, int IndexStart, int IndexProcess1, int IndexProcess2> template <typename TSecondaryView> diff --git a/corsika/detail/modules/conex/CONEXhybrid.inl b/corsika/detail/modules/conex/CONEXhybrid.inl index 47a072f8a..be172091f 100644 --- a/corsika/detail/modules/conex/CONEXhybrid.inl +++ b/corsika/detail/modules/conex/CONEXhybrid.inl @@ -99,7 +99,7 @@ namespace corsika { configPath.c_str(), configPath.size()); } - inline void CONEXhybrid::init() { + inline void CONEXhybrid::initCascadeEquations() { double eprima = primaryEnergy_ / 1_GeV; @@ -234,7 +234,7 @@ namespace corsika { } template <typename TStack> - inline void CONEXhybrid::doCascadeEquations(TStack& stack) { + inline void CONEXhybrid::doCascadeEquations(TStack&) { ::conex::conexcascade_(); diff --git a/corsika/framework/process/ProcessSequence.hpp b/corsika/framework/process/ProcessSequence.hpp index 432b78a28..407f40349 100644 --- a/corsika/framework/process/ProcessSequence.hpp +++ b/corsika/framework/process/ProcessSequence.hpp @@ -211,11 +211,16 @@ namespace corsika { void doStack(TStack& stack); /** - Execute the CascadeEquationsProcess-es in the ProcessSequence + * Execute the CascadeEquationsProcess-es in the ProcessSequence. */ template <typename TStack> void doCascadeEquations(TStack& stack); + /** + * Init the CascadeEquationsProcess-es in the ProcessSequence. + */ + void initCascadeEquations(); + /** * Calculate the maximum allowed length of the next tracking step, based on all * ContinuousProcess-es diff --git a/corsika/modules/conex/CONEXhybrid.hpp b/corsika/modules/conex/CONEXhybrid.hpp index cbc73b7af..0d3878e47 100644 --- a/corsika/modules/conex/CONEXhybrid.hpp +++ b/corsika/modules/conex/CONEXhybrid.hpp @@ -8,11 +8,13 @@ #pragma once +#include <corsika/framework/process/SecondariesProcess.hpp> +#include <corsika/framework/process/CascadeEquationsProcess.hpp> + #include <corsika/framework/core/PhysicalUnits.hpp> #include <corsika/framework/core/ParticleProperties.hpp> #include <corsika/framework/geometry/Point.hpp> #include <corsika/framework/geometry/Vector.hpp> -#include <corsika/framework/process/SecondariesProcess.hpp> #include <corsika/media/ShowerAxis.hpp> #include <corsika/modules/conex/CONEX_f.hpp> @@ -23,7 +25,8 @@ namespace corsika { LengthType constexpr earthRadius{6371315 * meter}; } // namespace conex - class CONEXhybrid : public SecondariesProcess<CONEXhybrid> { + class CONEXhybrid : public CascadeEquationsProcess<CONEXhybrid>, + public SecondariesProcess<CONEXhybrid> { public: CONEXhybrid(Point center, ShowerAxis const& showerAxis, LengthType groundDist, @@ -41,7 +44,7 @@ namespace corsika { * init currently needs to be called to initializa a new * event. All tables are cleared, etc. */ - void init(); + void initCascadeEquations(); /** * Cascade equations are solved basoned on the data in the tables diff --git a/examples/hybrid_MC.cpp b/examples/hybrid_MC.cpp index 04a126366..dabed46e9 100644 --- a/examples/hybrid_MC.cpp +++ b/examples/hybrid_MC.cpp @@ -291,8 +291,6 @@ int main(int argc, char** argv) { cut.reset(); eLoss.reset(); - conex_model.solveCE(); - auto const hists = sibyllCounted.getHistogram() + sibyllNucCounted.getHistogram() + urqmdCounted.getHistogram(); diff --git a/tests/modules/testCONEX.cpp b/tests/modules/testCONEX.cpp index baf808aed..985043d4d 100644 --- a/tests/modules/testCONEX.cpp +++ b/tests/modules/testCONEX.cpp @@ -103,7 +103,7 @@ TEST_CASE("CONEXSourceCut") { [[maybe_unused]] corsika::sibyll::NuclearInteraction sibyllNuc(sibyll, env); CONEXhybrid conex(center, showerAxis, t, injectionHeight, E0, get_PDG(Code::Proton)); - conex.init(); + conex.initCascadeEquations(); HEPEnergyType const Eem{1_PeV}; auto const momentum = showerAxis.getDirection() * Eem; -- GitLab