From 654a40e8276effea6b6eb5019cf9485737345c0b Mon Sep 17 00:00:00 2001
From: Maximilian Reininghaus <maximilian.reininghaus@kit.edu>
Date: Tue, 16 Apr 2019 18:14:51 -0300
Subject: [PATCH] very basic test with NN interactions
---
Processes/UrQMD/CMakeLists.txt | 5 +++
Processes/UrQMD/UrQMD.cc | 12 +++----
Processes/UrQMD/UrQMD.h | 11 +++---
Processes/UrQMD/testUrQMD.cc | 66 +++++++++++++++++++++++++++++++++-
4 files changed, 82 insertions(+), 12 deletions(-)
diff --git a/Processes/UrQMD/CMakeLists.txt b/Processes/UrQMD/CMakeLists.txt
index 4c8d8e2cb..897492b2e 100644
--- a/Processes/UrQMD/CMakeLists.txt
+++ b/Processes/UrQMD/CMakeLists.txt
@@ -63,6 +63,7 @@ target_link_libraries (
CORSIKAunits
CORSIKAgeometry
CORSIKArandom
+ CORSIKAsetup
)
target_include_directories (
@@ -90,6 +91,10 @@ add_executable (testUrQMD
target_link_libraries (
testUrQMD
ProcessUrQMD
+ CORSIKAsetup
+ CORSIKArandom
+ CORSIKAgeometry
+ CORSIKAunits
CORSIKAthirdparty # for catch2
)
CORSIKA_ADD_TEST(testUrQMD)
diff --git a/Processes/UrQMD/UrQMD.cc b/Processes/UrQMD/UrQMD.cc
index 72d773cd9..9acc9b4d1 100644
--- a/Processes/UrQMD/UrQMD.cc
+++ b/Processes/UrQMD/UrQMD.cc
@@ -17,7 +17,8 @@ using namespace corsika::units::si;
UrQMD::UrQMD() { iniurqmd_(); }
using SetupStack = corsika::setup::Stack;
-using SetupParticle = SetupStack::StackIterator;
+using SetupParticle = corsika::setup::Stack::StackIterator;
+using SetupProjectile = corsika::setup::StackView::StackIterator;
using SetupTrack = corsika::setup::Trajectory;
CrossSectionType UrQMD::GetCrossSection(
@@ -27,7 +28,6 @@ CrossSectionType UrQMD::GetCrossSection(
return 10_mb; // TODO: implement
}
-template <>
GrammageType UrQMD::GetInteractionLength(SetupParticle& particle, SetupTrack&) const {
auto const& mediumComposition =
particle.GetNode()->GetModelProperties().GetNuclearComposition();
@@ -41,9 +41,7 @@ GrammageType UrQMD::GetInteractionLength(SetupParticle& particle, SetupTrack&) c
weightedProdCrossSection;
}
-template <>
-corsika::process::EProcessReturn UrQMD::DoInteraction(SetupParticle& projectile,
- SetupStack&) {
+corsika::process::EProcessReturn UrQMD::DoInteraction(SetupProjectile& projectile) {
using namespace units::si;
auto const projectileCode = projectile.GetPID();
@@ -104,7 +102,7 @@ corsika::process::EProcessReturn UrQMD::DoInteraction(SetupParticle& projectile,
inputs_.spiso3[0] = iso3;
}
- rsys_.ebeam = (projectileEnergyLab - particles::GetMass(projectileCode)) * (1 / 1_GeV);
+ rsys_.ebeam = (projectileEnergyLab - projectile.GetMass()) * (1 / 1_GeV);
// initilazation regarding target
auto const& mediumComposition =
@@ -159,6 +157,8 @@ corsika::process::EProcessReturn UrQMD::DoInteraction(SetupParticle& projectile,
TimeType>{code, energy, momentum, projectilePosition, projectileTime});
}
+ std::cout << "UrQMD generated " << sys_.npart << " secondaries!" << std::endl;
+
if (sys_.npart > 0) // delete only in case of inelastic collision, otherwise keep
projectile.Delete();
}
diff --git a/Processes/UrQMD/UrQMD.h b/Processes/UrQMD/UrQMD.h
index 9693d6ecc..ae3c4f019 100644
--- a/Processes/UrQMD/UrQMD.h
+++ b/Processes/UrQMD/UrQMD.h
@@ -4,6 +4,8 @@
#include <corsika/particles/ParticleProperties.h>
#include <corsika/process/InteractionProcess.h>
#include <corsika/random/RNGManager.h>
+#include <corsika/setup/SetupStack.h>
+#include <corsika/setup/SetupTrajectory.h>
#include <corsika/units/PhysicalUnits.h>
#include <array>
@@ -13,16 +15,15 @@ namespace corsika::process::UrQMD {
class UrQMD : public corsika::process::InteractionProcess<UrQMD> {
public:
UrQMD();
-
- template <typename Particle, typename Track>
- corsika::units::si::GrammageType GetInteractionLength(Particle&, Track&) const;
+ corsika::units::si::GrammageType GetInteractionLength(
+ corsika::setup::Stack::StackIterator&, corsika::setup::Trajectory&) const;
corsika::units::si::CrossSectionType GetCrossSection(
corsika::particles::Code, corsika::particles::Code,
corsika::units::si::HEPEnergyType) const;
- template <typename Particle, typename Stack>
- corsika::process::EProcessReturn DoInteraction(Particle&, Stack&);
+ corsika::process::EProcessReturn DoInteraction(
+ corsika::setup::StackView::StackIterator&);
private:
corsika::random::RNG& fRNG =
diff --git a/Processes/UrQMD/testUrQMD.cc b/Processes/UrQMD/testUrQMD.cc
index 6d26eb036..2dad439de 100644
--- a/Processes/UrQMD/testUrQMD.cc
+++ b/Processes/UrQMD/testUrQMD.cc
@@ -9,11 +9,75 @@
*/
#include <corsika/process/urqmd/UrQMD.h>
+#include <corsika/random/RNGManager.h>
+
+#include <corsika/geometry/Point.h>
+#include <corsika/geometry/RootCoordinateSystem.h>
+#include <corsika/geometry/Vector.h>
+
+#include <corsika/units/PhysicalConstants.h>
+#include <corsika/units/PhysicalUnits.h>
+
+#include <corsika/particles/ParticleProperties.h>
+#include <corsika/setup/SetupStack.h>
+#include <corsika/setup/SetupTrajectory.h>
+
+#include <corsika/environment/Environment.h>
+#include <corsika/environment/HomogeneousMedium.h>
+#include <corsika/environment/NuclearComposition.h>
#define CATCH_CONFIG_MAIN // This tells Catch to provide a main() - only do this in one
// cpp file
#include <catch2/catch.hpp>
+using namespace corsika;
using namespace corsika::process::UrQMD;
+using namespace corsika::units::si;
+
+TEST_CASE("UrQMD") {
+ corsika::random::RNGManager::GetInstance().RegisterRandomStream("UrQMD");
+ UrQMD urqmd;
+
+ // setup environment, geometry
+ environment::Environment<environment::IMediumModel> env;
+ auto& universe = *(env.GetUniverse());
+ const geometry::CoordinateSystem& cs = env.GetCoordinateSystem();
+
+ auto theMedium =
+ environment::Environment<environment::IMediumModel>::CreateNode<geometry::Sphere>(
+ geometry::Point{cs, 0_m, 0_m, 0_m},
+ 1_km * std::numeric_limits<double>::infinity());
+
+ using MyHomogeneousModel = environment::HomogeneousMedium<environment::IMediumModel>;
+ theMedium->SetModelProperties<MyHomogeneousModel>(
+ 1_kg / (1_m * 1_m * 1_m),
+ environment::NuclearComposition(
+ std::vector<particles::Code>{particles::Code::Oxygen}, std::vector<float>{1.}));
+
+ auto const* nodePtr = theMedium.get();
+ universe.AddChild(std::move(theMedium));
+
+ geometry::Point const origin(cs, {0_m, 0_m, 0_m});
+ geometry::Vector<units::si::SpeedType::dimension_type> v(cs, 0_m / second, 0_m / second,
+ 1_m / second);
+ geometry::Line line(origin, v);
+ geometry::Trajectory<geometry::Line> track(line, 10_s);
+
+ auto constexpr mN = corsika::units::constants::nucleonMass;
+
+ setup::Stack stack;
+ const HEPEnergyType P0 = 50_GeV;
+ unsigned short constexpr A = 56, Z = 26;
+ HEPMomentumType E0 = sqrt(A * A * mN * mN + P0 * P0);
+ auto plab = corsika::stack::MomentumVector(cs, {P0, 0_GeV, 0_GeV});
+ auto particle =
+ stack.AddParticle(std::tuple<particles::Code, units::si::HEPEnergyType,
+ corsika::stack::MomentumVector, geometry::Point,
+ units::si::TimeType, unsigned short, unsigned short>{
+ particles::Code::Nucleus, E0, plab, origin, 0_ns, A, Z}); // iron
+ particle.SetNode(nodePtr);
+ corsika::stack::SecondaryView view(particle);
+ auto projectile = view.GetProjectile();
-TEST_CASE("UrQMD") { UrQMD proc; }
+ [[maybe_unused]] const process::EProcessReturn ret = urqmd.DoInteraction(projectile);
+}
--
GitLab