From ecd73694c9bf30c51c8cc136df766fc6f6f62bc6 Mon Sep 17 00:00:00 2001
From: Maximilian Reininghaus <maximilian.reininghaus@kit.edu>
Date: Thu, 17 Jan 2019 18:51:57 +0100
Subject: [PATCH] conversion function HEP to SI based on dim. analysis
---
Framework/Units/PhysicalConstants.h | 8 +++--
Framework/Units/PhysicalUnits.h | 45 +++++++++++++++++++++++++++++
Framework/Units/testUnits.cc | 31 ++++++++++++++++++++
3 files changed, 82 insertions(+), 2 deletions(-)
diff --git a/Framework/Units/PhysicalConstants.h b/Framework/Units/PhysicalConstants.h
index 95436959..ae47bd34 100644
--- a/Framework/Units/PhysicalConstants.h
+++ b/Framework/Units/PhysicalConstants.h
@@ -37,7 +37,7 @@ namespace corsika::units::constants {
phys::units::square(phys::units::second)};
// Avogadro constant
- constexpr quantity<dimensions<0, 0, 0, 0, 0, -1> > N_sub_A{Rep(6.02214199e+23L) / mole};
+ constexpr quantity<dimensions<0, 0, 0, 0, 0, -1>> N_sub_A{Rep(6.02214199e+23L) / mole};
// elementary charge
constexpr quantity<electric_charge_d> e{Rep(1.6021766208e-19L) * coulomb};
@@ -46,12 +46,16 @@ namespace corsika::units::constants {
// constexpr quantity<hepenergy_d> eV{e / coulomb * joule};
// Planck constant
- constexpr quantity<dimensions<2, 1, -1> > h{Rep(6.62606876e-34L) * joule * second};
+ constexpr quantity<dimensions<2, 1, -1>> h{Rep(6.62606876e-34L) * joule * second};
+ constexpr quantity<dimensions<2, 1, -1>> hBar{h / (2 * M_PI)};
// speed of light in a vacuum
constexpr quantity<speed_d> c{Rep(299792458L) * meter / second};
constexpr auto cSquared = c * c;
+ constexpr quantity<dimensions<1, 0, 0, 0, 0, 0, 0, 1>> hBarC{
+ Rep(1.973'269'78e-7L) * electronvolt * meter}; // from RPP 2018
+
// unified atomic mass unit
constexpr quantity<mass_d> u{Rep(1.6605402e-27L) * kilogram};
diff --git a/Framework/Units/PhysicalUnits.h b/Framework/Units/PhysicalUnits.h
index 262e55d1..ffdb2480 100644
--- a/Framework/Units/PhysicalUnits.h
+++ b/Framework/Units/PhysicalUnits.h
@@ -74,6 +74,51 @@ namespace corsika::units::si {
phys::units::quantity<phys::units::dimensions<0, 0, -1>, double>;
using InverseGrammageType =
phys::units::quantity<phys::units::dimensions<2, -1, 0>, double>;
+
+ namespace detail {
+ template <int N, typename T>
+ auto constexpr static_pow([[maybe_unused]] T x) {
+ if constexpr (N == 0) {
+ return 1;
+ } else if constexpr (N > 0) {
+ return x * static_pow<N - 1, T>(x);
+ } else {
+ return 1 / static_pow<-N, T>(x);
+ }
+ }
+ } // namespace detail
+
+ template <typename DimFrom, typename DimTo>
+ auto constexpr ConversionFactorHEPToSI() {
+ static_assert(DimFrom::dim1 == 0 && DimFrom::dim2 == 0 && DimFrom::dim3 == 0 &&
+ DimFrom::dim4 == 0 && DimFrom::dim5 == 0 && DimFrom::dim6 == 0 &&
+ DimFrom::dim7 == 0,
+ "must be a pure HEP type");
+
+ static_assert(
+ DimTo::dim4 == 0 && DimTo::dim5 == 0 && DimTo::dim6 == 0 && DimTo::dim7 == 0,
+ "conversion possible only into L, M, T dimensions");
+
+ int constexpr e = DimFrom::dim8; // HEP dim.
+
+ int constexpr l = DimTo::dim1; // SI length dim.
+ int constexpr m = DimTo::dim2; // SI mass dim.
+ int constexpr t = DimTo::dim3; // SI time dim.
+
+ int constexpr p = m;
+ int constexpr q = -m - t;
+ static_assert(q == l + e - 2 * m, "HEP/SI dimension mismatch!");
+
+ using namespace detail;
+ return static_pow<-e>(corsika::units::constants::hBarC) *
+ static_pow<p>(corsika::units::constants::hBar) *
+ static_pow<q>(corsika::units::constants::c);
+ }
+
+ template <typename DimFrom, typename DimTo>
+ auto constexpr ConvertHEPToSI(quantity<DimFrom> q) {
+ return ConversionFactorHEPToSI<DimFrom, DimTo>() * q;
+ }
} // end namespace corsika::units::si
/**
diff --git a/Framework/Units/testUnits.cc b/Framework/Units/testUnits.cc
index 38b3b299..ecf2ec51 100644
--- a/Framework/Units/testUnits.cc
+++ b/Framework/Units/testUnits.cc
@@ -115,4 +115,35 @@ TEST_CASE("PhysicalUnits", "[Units]") {
REQUIRE(farAway > 100000_m);
REQUIRE_FALSE(farAway < 1e19 * meter);
}
+
+ SECTION("static_pow") {
+ using namespace corsika::units::si::detail;
+ double x = 235.7913;
+ REQUIRE(1 == static_pow<0, double>(x));
+ REQUIRE(x == static_pow<1, double>(x));
+ REQUIRE(x * x == static_pow<2, double>(x));
+ REQUIRE(1 / x == static_pow<-1, double>(x));
+ REQUIRE(1 / x / x == static_pow<-2, double>(x));
+ }
+
+ SECTION("HEP/SI conversion") {
+ auto const invEnergy = 1 / 197.326978_MeV; // should be convertible to length or time
+
+ LengthType const length =
+ ConvertHEPToSI<decltype(invEnergy)::dimension_type, LengthType::dimension_type>(
+ invEnergy);
+ REQUIRE((length / 1_fm) == Approx(1));
+
+ TimeType const time =
+ ConvertHEPToSI<decltype(invEnergy)::dimension_type, TimeType::dimension_type>(
+ invEnergy);
+ REQUIRE((time / (1_fm / corsika::units::constants::c)) == Approx(1));
+
+ auto const protonMass = 938.272'081'3_MeV; // convertible to mass or SI energy
+ MassType protonMassSI =
+ ConvertHEPToSI<decltype(protonMass)::dimension_type, MassType::dimension_type>(
+ protonMass);
+ REQUIRE((protonMassSI / 1.672'621'898e-27_kg) == Approx(1));
+ REQUIRE((protonMassSI / (1.007'276 * corsika::units::constants::u)) == Approx(1));
+ }
}
--
GitLab