From e45344d49509e22b262f329618bbc1780d592b33 Mon Sep 17 00:00:00 2001 From: ralfulrich <ralf.ulrich@kit.edu> Date: Tue, 2 Apr 2019 11:27:16 +0200 Subject: [PATCH] Added SecondariesProcess, some wording -> style --- Framework/ProcessSequence/CMakeLists.txt | 1 + Framework/ProcessSequence/ContinuousProcess.h | 4 +- Framework/ProcessSequence/DecayProcess.h | 8 +- .../ProcessSequence/InteractionProcess.h | 4 +- Framework/ProcessSequence/ProcessSequence.h | 100 ++++++++++-------- Framework/ProcessSequence/SecondaryProcess.h | 48 +++++++++ .../ProcessSequence/testProcessSequence.cc | 59 +++++------ 7 files changed, 143 insertions(+), 81 deletions(-) create mode 100644 Framework/ProcessSequence/SecondaryProcess.h diff --git a/Framework/ProcessSequence/CMakeLists.txt b/Framework/ProcessSequence/CMakeLists.txt index 9ec9072da..ac49b2b38 100644 --- a/Framework/ProcessSequence/CMakeLists.txt +++ b/Framework/ProcessSequence/CMakeLists.txt @@ -13,6 +13,7 @@ set ( BaseProcess.h BoundaryCrossingProcess.h ContinuousProcess.h + SecondariesProcess.h InteractionProcess.h DecayProcess.h ProcessSequence.h diff --git a/Framework/ProcessSequence/ContinuousProcess.h b/Framework/ProcessSequence/ContinuousProcess.h index 9a3424a35..58d378aa3 100644 --- a/Framework/ProcessSequence/ContinuousProcess.h +++ b/Framework/ProcessSequence/ContinuousProcess.h @@ -33,8 +33,8 @@ namespace corsika::process { // here starts the interface part // -> enforce derived to implement DoContinuous... - template <typename Particle, typename Track, typename Stack> - EProcessReturn DoContinuous(Particle&, Track&, Stack&) const; + template <typename Particle, typename Track> + EProcessReturn DoContinuous(Particle&, Track&) const; // -> enforce derived to implement MaxStepLength... template <typename Particle, typename Track> diff --git a/Framework/ProcessSequence/DecayProcess.h b/Framework/ProcessSequence/DecayProcess.h index 82be572d6..16d63ae66 100644 --- a/Framework/ProcessSequence/DecayProcess.h +++ b/Framework/ProcessSequence/DecayProcess.h @@ -34,15 +34,15 @@ namespace corsika::process { /// here starts the interface-definition part // -> enforce derived to implement DoDecay... - template <typename Particle, typename Stack> - EProcessReturn DoDecay(Particle&, Stack&); + template <typename Particle> + EProcessReturn DoDecay(Particle&); template <typename Particle> corsika::units::si::TimeType GetLifetime(Particle& p); template <typename Particle> - corsika::units::si::InverseTimeType GetInverseLifetime(Particle& p) { - return 1. / GetRef().GetLifetime(p); + corsika::units::si::InverseTimeType GetInverseLifetime(Particle& vP) { + return 1. / GetRef().GetLifetime(vP); } }; diff --git a/Framework/ProcessSequence/InteractionProcess.h b/Framework/ProcessSequence/InteractionProcess.h index 723cbf02f..27ad9c222 100644 --- a/Framework/ProcessSequence/InteractionProcess.h +++ b/Framework/ProcessSequence/InteractionProcess.h @@ -35,8 +35,8 @@ namespace corsika::process { /// here starts the interface-definition part // -> enforce derived to implement DoInteraction... - template <typename P, typename S> - EProcessReturn DoInteraction(P&, S&); + template <typename Particle> + EProcessReturn DoInteraction(Particle&); template <typename Particle, typename Track> corsika::units::si::GrammageType GetInteractionLength(Particle& p, Track& t); diff --git a/Framework/ProcessSequence/ProcessSequence.h b/Framework/ProcessSequence/ProcessSequence.h index c697d09f4..bc1a5f779 100644 --- a/Framework/ProcessSequence/ProcessSequence.h +++ b/Framework/ProcessSequence/ProcessSequence.h @@ -18,6 +18,7 @@ #include <corsika/process/DecayProcess.h> #include <corsika/process/InteractionProcess.h> #include <corsika/process/ProcessReturn.h> +#include <corsika/process/SecondariesProcess.h> #include <corsika/units/PhysicalUnits.h> #include <cmath> @@ -89,78 +90,93 @@ namespace corsika::process { return ret; } - template <typename Particle, typename Track, typename Stack> - EProcessReturn DoContinuous(Particle& p, Track& t, Stack& s) { + template <typename TParticle, typename TTrack> + EProcessReturn DoContinuous(TParticle& vP, TTrack& vT) { EProcessReturn ret = EProcessReturn::eOk; if constexpr (std::is_base_of<ContinuousProcess<T1type>, T1type>::value || is_process_sequence<T1>::value) { - ret |= A.DoContinuous(p, t, s); + ret |= A.DoContinuous(vP, vT); } if constexpr (std::is_base_of<ContinuousProcess<T2type>, T2type>::value || is_process_sequence<T2>::value) { - ret |= B.DoContinuous(p, t, s); + ret |= B.DoContinuous(vP, vT); } return ret; } - template <typename Particle, typename Track> - corsika::units::si::LengthType MaxStepLength(Particle& p, Track& track) { + template <typename TSecondaries> + EProcessReturn DoSecondaries(TSecondaries& vS) { + EProcessReturn ret = EProcessReturn::eOk; + if constexpr (std::is_base_of<SecondariesProcess<T1type>, T1type>::value || + is_process_sequence<T1>::value) { + ret |= A.DoSecondaries(vS); + } + if constexpr (std::is_base_of<SecondariesProcess<T2type>, T2type>::value || + is_process_sequence<T2>::value) { + ret |= B.DoSecondaries(vS); + } + return ret; + } + + template <typename TParticle, typename TTrack> + corsika::units::si::LengthType MaxStepLength(TParticle& vP, TTrack& vTrack) { corsika::units::si::LengthType max_length = // if no other process in the sequence implements it std::numeric_limits<double>::infinity() * corsika::units::si::meter; if constexpr (std::is_base_of<ContinuousProcess<T1type>, T1type>::value || is_process_sequence<T1>::value) { - corsika::units::si::LengthType const len = A.MaxStepLength(p, track); + corsika::units::si::LengthType const len = A.MaxStepLength(vP, vTrack); max_length = std::min(max_length, len); } if constexpr (std::is_base_of<ContinuousProcess<T2type>, T2type>::value || is_process_sequence<T2>::value) { - corsika::units::si::LengthType const len = B.MaxStepLength(p, track); + corsika::units::si::LengthType const len = B.MaxStepLength(vP, vTrack); max_length = std::min(max_length, len); } return max_length; } - template <typename Particle, typename Track> - corsika::units::si::GrammageType GetTotalInteractionLength(Particle& p, Track& t) { - return 1. / GetInverseInteractionLength(p, t); + template <typename TParticle, typename TTrack> + corsika::units::si::GrammageType GetTotalInteractionLength(TParticle& vP, + TTrack& vT) { + return 1. / GetInverseInteractionLength(vP, vT); } - template <typename Particle, typename Track> - corsika::units::si::InverseGrammageType GetTotalInverseInteractionLength(Particle& p, - Track& t) { - return GetInverseInteractionLength(p, t); + template <typename TParticle, typename TTrack> + corsika::units::si::InverseGrammageType GetTotalInverseInteractionLength( + TParticle& vP, TTrack& vT) { + return GetInverseInteractionLength(vP, vT); } - template <typename Particle, typename Track> - corsika::units::si::InverseGrammageType GetInverseInteractionLength(Particle& p, - Track& t) { + template <typename TParticle, typename TTrack> + corsika::units::si::InverseGrammageType GetInverseInteractionLength(TParticle& vP, + TTrack& vT) { using namespace corsika::units::si; InverseGrammageType tot = 0 * meter * meter / gram; if constexpr (std::is_base_of<InteractionProcess<T1type>, T1type>::value || is_process_sequence<T1>::value) { - tot += A.GetInverseInteractionLength(p, t); + tot += A.GetInverseInteractionLength(vP, vT); } if constexpr (std::is_base_of<InteractionProcess<T2type>, T2type>::value || is_process_sequence<T2>::value) { - tot += B.GetInverseInteractionLength(p, t); + tot += B.GetInverseInteractionLength(vP, vT); } return tot; } - template <typename Particle, typename Track, typename Stack> + template <typename TParticle, typename TSecondaries, typename TTrack> EProcessReturn SelectInteraction( - Particle& vP, Track& vT, Stack& vS, + TParticle& vP, TSecondaries& vS, TTrack& vT, [[maybe_unused]] corsika::units::si::InverseGrammageType lambda_select, corsika::units::si::InverseGrammageType& lambda_inv_count) { if constexpr (is_process_sequence<T1type>::value) { // if A is a process sequence --> check inside const EProcessReturn ret = - A.SelectInteraction(vP, vT, vS, lambda_select, lambda_inv_count); + A.SelectInteraction(vP, vS, vT, lambda_select, lambda_inv_count); // if A did succeed, stop routine if (ret != EProcessReturn::eOk) { return ret; } } else if constexpr (std::is_base_of<InteractionProcess<T1type>, T1type>::value) { @@ -168,7 +184,7 @@ namespace corsika::process { lambda_inv_count += A.GetInverseInteractionLength(vP, vT); // check if we should execute THIS process and then EXIT if (lambda_select < lambda_inv_count) { - A.DoInteraction(vP, vS); + A.DoInteraction(vS); return EProcessReturn::eInteracted; } } // end branch A @@ -176,7 +192,7 @@ namespace corsika::process { if constexpr (is_process_sequence<T2>::value) { // if A is a process sequence --> check inside const EProcessReturn ret = - B.SelectInteraction(vP, vT, vS, lambda_select, lambda_inv_count); + B.SelectInteraction(vP, vS, vT, lambda_select, lambda_inv_count); // if A did succeed, stop routine if (ret != EProcessReturn::eOk) { return ret; } } else if constexpr (std::is_base_of<InteractionProcess<T2type>, T2type>::value) { @@ -184,25 +200,25 @@ namespace corsika::process { lambda_inv_count += B.GetInverseInteractionLength(vP, vT); // check if we should execute THIS process and then EXIT if (lambda_select < lambda_inv_count) { - B.DoInteraction(vP, vS); + B.DoInteraction(vS); return EProcessReturn::eInteracted; } } // end branch A return EProcessReturn::eOk; } - template <typename Particle> - corsika::units::si::TimeType GetTotalLifetime(Particle& p) { + template <typename TParticle> + corsika::units::si::TimeType GetTotalLifetime(TParticle& p) { return 1. / GetInverseLifetime(p); } - template <typename Particle> - corsika::units::si::InverseTimeType GetTotalInverseLifetime(Particle& p) { + template <typename TParticle> + corsika::units::si::InverseTimeType GetTotalInverseLifetime(TParticle& p) { return GetInverseLifetime(p); } - template <typename Particle> - corsika::units::si::InverseTimeType GetInverseLifetime(Particle& p) { + template <typename TParticle> + corsika::units::si::InverseTimeType GetInverseLifetime(TParticle& p) { using namespace corsika::units::si; corsika::units::si::InverseTimeType tot = 0 / second; @@ -219,38 +235,38 @@ namespace corsika::process { } // select decay process - template <typename Particle, typename Stack> + template <typename TParticle, typename TSecondaries> EProcessReturn SelectDecay( - Particle& p, Stack& s, + TParticle& vP, TSecondaries& vS, [[maybe_unused]] corsika::units::si::InverseTimeType decay_select, corsika::units::si::InverseTimeType& decay_inv_count) { if constexpr (is_process_sequence<T1>::value) { // if A is a process sequence --> check inside - const EProcessReturn ret = A.SelectDecay(p, s, decay_select, decay_inv_count); + const EProcessReturn ret = A.SelectDecay(vP, vS, decay_select, decay_inv_count); // if A did succeed, stop routine if (ret != EProcessReturn::eOk) { return ret; } } else if constexpr (std::is_base_of<DecayProcess<T1type>, T1type>::value) { // if this is not a ContinuousProcess --> evaluate probability - decay_inv_count += A.GetInverseLifetime(p); + decay_inv_count += A.GetInverseLifetime(vP); // check if we should execute THIS process and then EXIT if (decay_select < decay_inv_count) { // more pedagogical: rndm_select < // decay_inv_count / decay_inv_tot - A.DoDecay(p, s); + A.DoDecay(vS); return EProcessReturn::eDecayed; } } // end branch A if constexpr (is_process_sequence<T2>::value) { // if A is a process sequence --> check inside - const EProcessReturn ret = B.SelectDecay(p, s, decay_select, decay_inv_count); + const EProcessReturn ret = B.SelectDecay(vP, vS, decay_select, decay_inv_count); // if A did succeed, stop routine if (ret != EProcessReturn::eOk) { return ret; } } else if constexpr (std::is_base_of<DecayProcess<T2type>, T2type>::value) { // if this is not a ContinuousProcess --> evaluate probability - decay_inv_count += B.GetInverseLifetime(p); + decay_inv_count += B.GetInverseLifetime(vP); // check if we should execute THIS process and then EXIT if (decay_select < decay_inv_count) { - B.DoDecay(p, s); + B.DoDecay(vS); return EProcessReturn::eDecayed; } } // end branch B @@ -272,8 +288,8 @@ namespace corsika::process { typename P1, typename P2, typename std::enable_if<is_process<typename std::decay<P1>::type>::value && is_process<typename std::decay<P2>::type>::value>::type...> - inline auto operator<<(P1&& A, P2&& B) -> ProcessSequence<P1, P2> { - return ProcessSequence<P1, P2>(A.GetRef(), B.GetRef()); + inline auto operator<<(P1&& vA, P2&& vB) -> ProcessSequence<P1, P2> { + return ProcessSequence<P1, P2>(vA.GetRef(), vB.GetRef()); } /// marker to identify objectas ProcessSequence diff --git a/Framework/ProcessSequence/SecondaryProcess.h b/Framework/ProcessSequence/SecondaryProcess.h new file mode 100644 index 000000000..523ce4d88 --- /dev/null +++ b/Framework/ProcessSequence/SecondaryProcess.h @@ -0,0 +1,48 @@ + +/* + * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu + * + * See file AUTHORS for a list of contributors. + * + * This software is distributed under the terms of the GNU General Public + * Licence version 3 (GPL Version 3). See file LICENSE for a full version of + * the license. + */ + +#ifndef _include_corsika_secondariesprocess_h_ +#define _include_corsika_secondariesprocess_h_ + +#include <corsika/process/ProcessReturn.h> // for convenience +#include <corsika/setup/SetupTrajectory.h> +#include <corsika/units/PhysicalUnits.h> + +namespace corsika::process { + + /** + \class SecondariesProcess + + The structural base type of a process object in a + ProcessSequence. Both, the ProcessSequence and all its elements + are of type SecondariesProcess<T> + + */ + + template <typename derived> + struct SecondariesProcess { + + derived& GetRef() { return static_cast<derived&>(*this); } + const derived& GetRef() const { return static_cast<const derived&>(*this); } + + /// here starts the interface-definition part + // -> enforce derived to implement DoSecondaries... + template <typename Particle> + inline EProcessReturn DoSecondaries(Particle&); + }; + + // overwrite the default trait class, to mark BaseProcess<T> as useful process + template <class T> + std::true_type is_process_impl(const SecondariesProcess<T>* impl); + +} // namespace corsika::process + +#endif diff --git a/Framework/ProcessSequence/testProcessSequence.cc b/Framework/ProcessSequence/testProcessSequence.cc index 8fd46418a..583c6a1ad 100644 --- a/Framework/ProcessSequence/testProcessSequence.cc +++ b/Framework/ProcessSequence/testProcessSequence.cc @@ -39,8 +39,8 @@ public: assert(globalCount == fV); globalCount++; } - template <typename D, typename T, typename S> - inline EProcessReturn DoContinuous(D& d, T&, S&) const { + template <typename D, typename T> + inline EProcessReturn DoContinuous(D& d, T&) const { cout << "ContinuousProcess1::DoContinuous" << endl; for (int i = 0; i < nData; ++i) d.p[i] += 0.933; return EProcessReturn::eOk; @@ -58,8 +58,8 @@ public: assert(globalCount == fV); globalCount++; } - template <typename D, typename T, typename S> - inline EProcessReturn DoContinuous(D& d, T&, S&) const { + template <typename D, typename T> + inline EProcessReturn DoContinuous(D& d, T&) const { cout << "ContinuousProcess2::DoContinuous" << endl; for (int i = 0; i < nData; ++i) d.p[i] += 0.933; return EProcessReturn::eOk; @@ -95,8 +95,8 @@ public: assert(globalCount == fV); globalCount++; } - template <typename Particle, typename Stack> - inline EProcessReturn DoInteraction(Particle&, Stack&) const { + template <typename Particle> + inline EProcessReturn DoInteraction(Particle&) const { cout << "Process2::DoInteraction" << endl; return EProcessReturn::eOk; } @@ -118,8 +118,8 @@ public: assert(globalCount == fV); globalCount++; } - template <typename Particle, typename Stack> - inline EProcessReturn DoInteraction(Particle&, Stack&) const { + template <typename Particle> + inline EProcessReturn DoInteraction(Particle&) const { cout << "Process3::DoInteraction" << endl; return EProcessReturn::eOk; } @@ -141,14 +141,14 @@ public: assert(globalCount == fV); globalCount++; } - template <typename D, typename T, typename S> - inline EProcessReturn DoContinuous(D& d, T&, S&) const { + template <typename D, typename T> + inline EProcessReturn DoContinuous(D& d, T&) const { for (int i = 0; i < nData; ++i) { d.p[i] /= 1.2; } return EProcessReturn::eOk; } // inline double MinStepLength(D& d) { - template <typename Particle, typename Stack> - EProcessReturn DoInteraction(Particle&, Stack&) const { + template <typename Particle> + EProcessReturn DoInteraction(Particle&) const { return EProcessReturn::eOk; } }; @@ -168,8 +168,8 @@ public: TimeType GetLifetime(Particle&) const { return 1_s; } - template <typename Particle, typename Stack> - EProcessReturn DoDecay(Particle&, Stack&) const { + template <typename Particle> + EProcessReturn DoDecay(Particle&) const { return EProcessReturn::eOk; } }; @@ -177,7 +177,6 @@ public: struct DummyData { double p[nData] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; }; -struct DummyStack {}; struct DummyTrajectory {}; TEST_CASE("Process Sequence", "[Process Sequence]") { @@ -205,12 +204,13 @@ TEST_CASE("Process Sequence", "[Process Sequence]") { Process2 m2(1); Process3 m3(2); - DummyStack s; - DummyTrajectory t; + DummyData particle; + DummyTrajectory track; auto sequence2 = cp1 << m2 << m3; - GrammageType const tot = sequence2.GetTotalInteractionLength(s, t); - InverseGrammageType const tot_inv = sequence2.GetTotalInverseInteractionLength(s, t); + GrammageType const tot = sequence2.GetTotalInteractionLength(particle, track); + InverseGrammageType const tot_inv = + sequence2.GetTotalInverseInteractionLength(particle, track); cout << "lambda_tot=" << tot << "; lambda_tot_inv=" << tot_inv << endl; } @@ -220,11 +220,11 @@ TEST_CASE("Process Sequence", "[Process Sequence]") { Process3 m3(2); Decay1 d3(2); - DummyStack s; + DummyData particle; auto sequence2 = cp1 << m2 << m3 << d3; - TimeType const tot = sequence2.GetTotalLifetime(s); - InverseTimeType const tot_inv = sequence2.GetTotalInverseLifetime(s); + TimeType const tot = sequence2.GetTotalLifetime(particle); + InverseTimeType const tot_inv = sequence2.GetTotalInverseLifetime(particle); cout << "lambda_tot=" << tot << "; lambda_tot_inv=" << tot_inv << endl; } @@ -237,27 +237,24 @@ TEST_CASE("Process Sequence", "[Process Sequence]") { auto sequence2 = cp1 << m2 << m3 << cp2; - DummyData p; - DummyStack s; - DummyTrajectory t; + DummyData particle; + DummyTrajectory track; cout << "-->init sequence2" << endl; globalCount = 0; sequence2.Init(); cout << "-->docont" << endl; - sequence2.DoContinuous(p, t, s); + sequence2.DoContinuous(particle, track); cout << "-->dodisc" << endl; - // sequence2.DoInteraction(p, s); cout << "-->done" << endl; const int nLoop = 5; cout << "Running loop with n=" << nLoop << endl; - for (int i = 0; i < nLoop; ++i) { - sequence2.DoContinuous(p, t, s); - // sequence2.DoInteraction(p, s); + for (int i = 0; i < nLoop; ++i) { sequence2.DoContinuous(particle, track); } + for (int i = 0; i < nData; i++) { + cout << "data[" << i << "]=" << particle.p[i] << endl; } - for (int i = 0; i < nData; i++) { cout << "data[" << i << "]=" << p.p[i] << endl; } cout << "done" << endl; } } -- GitLab