IAP GITLAB

Skip to content
Snippets Groups Projects
Commit c12607d1 authored by Maximilian Reininghaus's avatar Maximilian Reininghaus :vulcan: Committed by Maximilian Reininghaus
Browse files

changed some units

parent fb2adb40
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Environment/HomogeneousMedium.h
+ 4
− 4
View file @ c12607d1
......@@ -43,15 +43,15 @@ namespace corsika::environment {
corsika::units::si::GrammageType IntegratedGrammage(
corsika::geometry::Trajectory<corsika::geometry::Line> const& pTraj,
corsika::units::si::TimeType pTo) const override {
corsika::units::si::LengthType pTo) const override {
using namespace corsika::units::si;
return pTraj.ArcLength(0_s, pTo) * fDensity;
pTo * fDensity;
}
corsika::units::si::TimeType FromGrammage(
corsika::units::si::TimeType ArclengthFromGrammage(
corsika::geometry::Trajectory<corsika::geometry::Line> const& pTraj,
corsika::units::si::GrammageType pGrammage) const override {
return pTraj.TimeFromArclength(pGrammage / fDensity);
return pGrammage / fDensity;
}
};
......
Environment/IMediumModel.h
+ 2
− 2
View file @ c12607d1
......@@ -20,9 +20,9 @@ namespace corsika::environment {
// approach for now, only lines are supported
virtual corsika::units::si::GrammageType IntegratedGrammage(
corsika::geometry::Trajectory<corsika::geometry::Line> const&,
corsika::units::si::TimeType) const = 0;
corsika::units::si::LengthType) const = 0;
virtual corsika::units::si::TimeType FromGrammage(
virtual corsika::units::si::LengthType ArclengthFromGrammage(
corsika::geometry::Trajectory<corsika::geometry::Line> const&,
corsika::units::si::GrammageType) const = 0;
......
Framework/Cascade/Cascade.h
+ 14
− 13
View file @ c12607d1
......@@ -16,6 +16,7 @@
#include <corsika/random/RNGManager.h>
#include <corsika/setup/SetupTrajectory.h>
#include <corsika/units/PhysicalUnits.h>
#include <corsika/environment/Environment.h>
#include <type_traits>
......@@ -30,12 +31,12 @@ namespace corsika::cascade {
public:
Cascade(Tracking& tr, ProcessList& pl, Stack& stack)
: fTracking(tr)
, fProcesseList(pl)
, fProcessSequence(pl)
, fStack(stack) {}
void Init() {
fTracking.Init();
fProcesseList.Init();
fProcessSequence.Init();
fStack.Init();
}
......@@ -58,7 +59,7 @@ namespace corsika::cascade {
// determine combined total interaction length (inverse)
InverseLengthType const total_inv_lambda =
fProcesseList.GetTotalInverseInteractionLength(particle, step);
fProcessSequence.GetTotalInverseInteractionLength(particle, step);
// sample random exponential step length
std::exponential_distribution expDist(total_inv_lambda * 1_m);
......@@ -67,11 +68,11 @@ namespace corsika::cascade {
<< ", next_interact=" << next_interact << std::endl;
// determine the maximum geometric step length
LengthType const distance_max = fProcesseList.MaxStepLength(particle, step);
LengthType const distance_max = fProcessSequence.MaxStepLength(particle, step);
std::cout << "distance_max=" << distance_max << std::endl;
// determine combined total inverse decay time
InverseTimeType const total_inv_lifetime = fProcesseList.GetTotalInverseLifetime(particle);
InverseTimeType const total_inv_lifetime = fProcessSequence.GetTotalInverseLifetime(particle);
// sample random exponential decay time
std::exponential_distribution expDistDecay(total_inv_lifetime * 1_s);
......@@ -81,7 +82,7 @@ namespace corsika::cascade {
// convert next_step from grammage to length [m]
// Environment::GetDistance(step, next_step);
const GrammageType distance_interact = ;
const GrammageType distance_interact = fEnvironment.GetUniverse()->GetContainingNode(particle.GetPosition())->GetModelProperties().FromGrammage()
// ....
// convert next_decay from time to length [m]
......@@ -93,14 +94,14 @@ namespace corsika::cascade {
const double distance_decay_interact = std::min(next_decay, next_interact);
const double distance_next = std::min(distance_decay_interact, distance_max);
/// here the particle is actually moved along the trajectory to new position:
// here the particle is actually moved along the trajectory to new position:
// std::visit(corsika::setup::ParticleUpdate<Particle>{particle}, step);
particle.SetPosition(step.GetPosition(1));
// .... also update time, momentum, direction, ...
// apply all continuous processes on particle + track
corsika::process::EProcessReturn status =
fProcesseList.DoContinuous(particle, step, fStack);
fProcessSequence.DoContinuous(particle, step, fStack);
if (status == corsika::process::EProcessReturn::eParticleAbsorbed) {
// fStack.Delete(particle); // TODO: check if this is really needed
......@@ -114,18 +115,18 @@ namespace corsika::cascade {
if (distance_decay > distance_interact) {
std::cout << "collide" << std::endl;
const double actual_inv_length =
fProcesseList.GetTotalInverseInteractionLength(particle, step);
fProcessSequence.GetTotalInverseInteractionLength(particle, step);
const double sample_process = rmng() / (double)rmng.max();
double inv_lambda_count = 0;
fProcesseList.SelectInteraction(particle, fStack, actual_inv_length,
fProcessSequence.SelectInteraction(particle, fStack, actual_inv_length,
sample_process, inv_lambda_count);
} else {
std::cout << "decay" << std::endl;
const double actual_decay_time =
fProcesseList.GetTotalInverseLifetime(particle);
fProcessSequence.GetTotalInverseLifetime(particle);
const double sample_process = rmng() / (double)rmng.max();
double inv_decay_count = 0;
fProcesseList.SelectDecay(particle, fStack, actual_decay_time, sample_process,
fProcessSequence.SelectDecay(particle, fStack, actual_decay_time, sample_process,
inv_decay_count);
}
}
......@@ -134,7 +135,7 @@ namespace corsika::cascade {
private:
Tracking& fTracking;
ProcessList& fProcesseList;
ProcessList& fProcessSequence;
Stack& fStack;
corsika::environment::Environment const& fEnvironment;
corsika::random::RNG& fRNG =
......
Framework/ProcessSequence/ProcessSequence.h
+ 20
− 20
View file @ c12607d1
......@@ -160,7 +160,7 @@ namespace corsika::process {
// void Hello() const { detail::CallHello<T1,T2>::Call(A, B); }
template <typename Particle, typename Track, typename Stack>
inline EProcessReturn DoContinuous(Particle& p, Track& t, Stack& s) const {
EProcessReturn DoContinuous(Particle& p, Track& t, Stack& s) const {
EProcessReturn ret = EProcessReturn::eOk;
if constexpr (std::is_base_of<ContinuousProcess<T1>, T1>::value ||
is_process_sequence<T1>::value) {
......@@ -174,8 +174,8 @@ namespace corsika::process {
}
template <typename Particle, typename Track>
inline double MaxStepLength(Particle& p, Track& track) const {
double max_length = std::numeric_limits<double>::infinity();
LengthType MaxStepLength(Particle& p, Track& track) const {
LengthType max_length = std::numeric_limits<double>::infinity() * 1_m;
if constexpr (std::is_base_of<ContinuousProcess<T1>, T1>::value ||
is_process_sequence<T1>::value) {
max_length = std::min(max_length, A.MaxStepLength(p, track));
......@@ -188,18 +188,18 @@ namespace corsika::process {
}
template <typename Particle, typename Track>
inline double GetTotalInteractionLength(Particle& p, Track& t) const {
GrammageType GetTotalInteractionLength(Particle& p, Track& t) const {
return 1. / GetInverseInteractionLength(p, t);
}
template <typename Particle, typename Track>
inline double GetTotalInverseInteractionLength(Particle& p, Track& t) const {
InverseGrammageType GetTotalInverseInteractionLength(Particle& p, Track& t) const {
return GetInverseInteractionLength(p, t);
}
template <typename Particle, typename Track>
inline double GetInverseInteractionLength(Particle& p, Track& t) const {
double tot = 0;
InverseGrammageType GetInverseInteractionLength(Particle& p, Track& t) const {
InverseGrammageType tot = 0;
if constexpr (std::is_base_of<InteractionProcess<T1>, T1>::value ||
is_process_sequence<T1>::value) {
tot += A.GetInverseInteractionLength(p, t);
......@@ -213,9 +213,9 @@ namespace corsika::process {
template <typename Particle, typename Stack>
inline EProcessReturn SelectInteraction(Particle& p, Stack& s,
const double lambda_inv_tot,
const double rndm_select,
double& lambda_inv_count) const {
InverseGrammageType lambda_inv_tot,
InverseGrammageType rndm_select,
InverseGrammageType& lambda_inv_count) const {
if constexpr (is_process_sequence<T1>::value) {
// if A is a process sequence --> check inside
const EProcessReturn ret =
......@@ -226,7 +226,7 @@ namespace corsika::process {
// if this is not a ContinuousProcess --> evaluate probability
lambda_inv_count += A.GetInverseInteractionLength(p, s);
// check if we should execute THIS process and then EXIT
if (rndm_select < lambda_inv_count / lambda_inv_tot) {
if (rndm_select * lambda_inv_tot < lambda_inv_count) { // more pedagogical: rndm_select < lambda_inv_count / lambda_inv_tot
A.DoInteraction(p, s);
return EProcessReturn::eInteracted;
}
......@@ -251,18 +251,18 @@ namespace corsika::process {
}
template <typename Particle>
inline double GetTotalLifetime(Particle& p) const {
TimeType GetTotalLifetime(Particle& p) const {
return 1. / GetInverseLifetime(p);
}
template <typename Particle>
inline double GetTotalInverseLifetime(Particle& p) const {
InverseTimeType GetTotalInverseLifetime(Particle& p) const {
return GetInverseLifetime(p);
}
template <typename Particle>
inline double GetInverseLifetime(Particle& p) const {
double tot = 0;
InverseTimeType GetInverseLifetime(Particle& p) const {
InverseTimeType tot = 0;
if constexpr (std::is_base_of<DecayProcess<T1>, T1>::value ||
is_process_sequence<T1>::value) {
tot += A.GetInverseLifetime(p);
......@@ -276,9 +276,9 @@ namespace corsika::process {
// select decay process
template <typename Particle, typename Stack>
inline EProcessReturn SelectDecay(Particle& p, Stack& s, const double decay_inv_tot,
const double rndm_select,
double& decay_inv_count) const {
EProcessReturn SelectDecay(Particle& p, Stack& s, InverseTimeType decay_inv_tot,
InverseTimeType rndm_select,
InverseTimeType& decay_inv_count) const {
if constexpr (is_process_sequence<T1>::value) {
// if A is a process sequence --> check inside
const EProcessReturn ret =
......@@ -289,7 +289,7 @@ namespace corsika::process {
// if this is not a ContinuousProcess --> evaluate probability
decay_inv_count += A.GetInverseLifetime(p);
// check if we should execute THIS process and then EXIT
if (rndm_select < decay_inv_count / decay_inv_tot) {
if (rndm_select * decay_inv_tot < decay_inv_count) { // more pedagogical: rndm_select < decay_inv_count / decay_inv_tot
A.DoDecay(p, s);
return EProcessReturn::eDecayed;
}
......@@ -314,7 +314,7 @@ namespace corsika::process {
}
/// TODO the const_cast is not nice, think about the constness here
inline void Init() const {
void Init() const {
const_cast<T1*>(&A)->Init();
const_cast<T2*>(&B)->Init();
}
......
Framework/Units/PhysicalUnits.h
+ 1
− 0
View file @ c12607d1
......@@ -56,6 +56,7 @@ namespace corsika::units::si {
using CrossSectionType = phys::units::quantity<area_d, double>;
using InverseLengthType = phys::units::quantity<phys::units::dimensions<-1, 0, 0>, double>;
using InverseTimeType = phys::units::quantity<phys::units::dimensions<0, 0, -1>, double>;
using InverseGrammageType = phys::units::quantity<phys::units::dimensions<2, -1, 0>, double>;
} // end namespace corsika::units::si
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment

data privacy - impressum