From 1dc19e6452fca0aa0d5c2626c138df4a9fec7cc1 Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Mon, 29 Dec 2025 12:35:24 -0500 Subject: [PATCH 01/26] Add initial muon CDF data --- src/celeritas/inp/MucfPhysics.cc | 13 ++- src/celeritas/inp/MucfPhysics.hh | 37 ++++--- src/celeritas/inp/MucfPhysicsData.hh | 117 +++++++++++++++++++++++ test/celeritas/ext/GeantImporter.test.cc | 20 +++- 4 files changed, 162 insertions(+), 25 deletions(-) create mode 100644 src/celeritas/inp/MucfPhysicsData.hh diff --git a/src/celeritas/inp/MucfPhysics.cc b/src/celeritas/inp/MucfPhysics.cc index 9724a6e0b0..472fb54590 100644 --- a/src/celeritas/inp/MucfPhysics.cc +++ b/src/celeritas/inp/MucfPhysics.cc @@ -6,6 +6,8 @@ //---------------------------------------------------------------------------// #include "MucfPhysics.hh" +#include "MucfPhysicsData.hh" + namespace celeritas { namespace inp @@ -21,16 +23,13 @@ namespace inp MucfPhysics MucfPhysics::from_default() { MucfPhysics result; - - //! \todo Initialize hardcoded CDF data - //! \todo Initialize hardcoded cycle rate data - //! \todo Initialize hardcoded atom transfer data - //! \todo Initialize hardcoded spin flip data + result.muon_energy_cdf = mucf_muon_energy_cdf(); + result.cycle_rates = mucf_cycle_rates(); + result.atom_transfer = mucf_atom_transfer_rates(); + result.atom_spin_flip = mucf_atom_spin_flip_rates(); // Temporary test dummy data to verify correct import { - result.muon_energy_cdf = Grid::from_constant(1.0); - MucfCycleRate dt_cycle; dt_cycle.molecule = MucfMuonicMolecule::deuterium_tritium; diff --git a/src/celeritas/inp/MucfPhysics.hh b/src/celeritas/inp/MucfPhysics.hh index 88c67a90e4..7854b85e0a 100644 --- a/src/celeritas/inp/MucfPhysics.hh +++ b/src/celeritas/inp/MucfPhysics.hh @@ -68,6 +68,9 @@ struct MucfCycleRate struct MucfAtomTransferRate { //! \todo Implement + + //! True if data is assigned \todo Implement + explicit operator bool() const { return true; } }; //---------------------------------------------------------------------------// @@ -75,22 +78,27 @@ struct MucfAtomTransferRate * Muon-catalyzed fusion mean atom spin flip data. * * Spin flip rates are as a function of temperature, with each grid/table - * representing an atom pair combination and its spin (e.g. deuterium-tritium, - * spin 1). Ordering is important, thus same spin deuterium-tritium and - * tritium-deuterium have different tables, which leads to a different final - * spin flip rate for a given material definition. + * representing an atom pair combination and its spin (e.g. + * deuterium-tritium, spin 1). Ordering is important, thus same spin + * deuterium-tritium and tritium-deuterium have different tables, which + * leads to a different final spin flip rate for a given material + * definition. * * Each struct holds one of such combinations, with the full set of - * combinations being the \c vector in \c MucfPhysics . + * combinations being the \c vector in \c MucfPhysics + * . * * \note These grids are host-only, with only the final spin flip rate per - * state (which is just a \c real_type ) for each combination being needed in - * the stepping loop. This is because these rates are material dependent, and - * thus can be cached at model construction. + * state (which is just a \c real_type ) for each combination being needed + * in the stepping loop. This is because these rates are material + * dependent, and thus can be cached at model construction. */ struct MucfAtomSpinFlipRate { //! \todo Implement + + //! True if data is assigned \todo Implement + explicit operator bool() const { return true; } }; //---------------------------------------------------------------------------// @@ -98,21 +106,22 @@ struct MucfAtomSpinFlipRate * Muon-catalyzed fusion physics options and data import. * * Minimum requirements for muon-catalyzed fusion: - * - Muon energy CDF data, required for sampling the outgoing muCF muon, and + * - Muon energy CDF data, required for sampling the outgoing muCF muon, + * and * - Mean cycle rate data for dd, dt, and tt muonic molecules. * - * Muonic atom transfer and muonic atom spin flip are secondary effects and not - * required for muCF to function. + * Muonic atom transfer and muonic atom spin flip are secondary effects and + * not required for muCF to function. */ struct MucfPhysics { template using Vec = std::vector; - Grid muon_energy_cdf; //!< CDF for outgoing muCF muon + Grid muon_energy_cdf; //!< CDF for sampling the outgoing muCF muon Vec cycle_rates; //!< Mean cycle rates for muonic molecules - Vec atom_transfer; //!< Muon atom transfer rates - Vec atom_spin_flip; //!< Muon atom spin flip rates + Vec atom_transfer; //!< Muonic atom transfer rates + Vec atom_spin_flip; //!< Muonic atom spin flip rates //! Whether muon-catalyzed fusion physics is enabled explicit operator bool() const diff --git a/src/celeritas/inp/MucfPhysicsData.hh b/src/celeritas/inp/MucfPhysicsData.hh new file mode 100644 index 0000000000..20c023dc07 --- /dev/null +++ b/src/celeritas/inp/MucfPhysicsData.hh @@ -0,0 +1,117 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/inp/MucfPhysicsData.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "corecel/Macros.hh" +#include "corecel/inp/Grid.hh" + +#include "MucfPhysics.hh" + +namespace celeritas +{ +namespace inp +{ +//---------------------------------------------------------------------------// +/*! + * Static muon energy CDF data for muon-catalyzed fusion. + * + * \todo This data may be loaded from an external file in the future. + */ +static Grid mucf_muon_energy_cdf() +{ + Grid cdf; + cdf.interpolation.type = InterpolationType::cubic_spline; + + // Cumulative distribution data [unitless] + cdf.x = {0, + 0.04169381107491854, + 0.08664495114006499, + 0.14332247557003264, + 0.20456026058631915, + 0.2723127035830618, + 0.34136807817589576, + 0.41563517915309456, + 0.48990228013029324, + 0.5667752442996744, + 0.6306188925081434, + 0.6866449511400652, + 0.7309446254071662, + 0.7778501628664496, + 0.8104234527687297, + 0.8403908794788275, + 0.8618892508143323, + 0.8814332247557004, + 0.8970684039087949, + 0.903583061889251, + 1.0}; + + // Energy [keV] + cdf.y = {0, + 0.48850540675768084, + 0.8390389347819425, + 1.2521213482687141, + 1.7153033196164724, + 2.253638712180777, + 2.854653691809707, + 3.606073540073316, + 4.470346052913727, + 5.560291219507215, + 6.700556502915258, + 7.953772477101693, + 9.194596305637525, + 10.849180562221111, + 12.353474314071864, + 14.045888515617822, + 15.650634617544647, + 17.38079707555165, + 19.111008546659452, + 19.976130619913615, + 80.0}; + + CELER_ENSURE(cdf); + return cdf; +} + +//---------------------------------------------------------------------------// +/*! + * Static cycle rate data for muon-catalyzed fusion. + * + * \todo This data may be loaded from an external file in the future. + */ +static std::vector mucf_cycle_rates() +{ + //! \todo Implement + return {}; +}; + +//---------------------------------------------------------------------------// +/*! + * Static spin-flip data for muon-catalyzed fusion. + * + * \todo This data may be loaded from an external file in the future. + */ +static std::vector mucf_atom_spin_flip_rates() +{ + //! \todo Implement + return {}; +}; + +//---------------------------------------------------------------------------// +/*! + * Static atom transfer data for muon-catalyzed fusion. + * + * \todo This data may be loaded from an external file in the future. + */ +static std::vector mucf_atom_transfer_rates() +{ + //! \todo Implement + return {}; +}; + +//---------------------------------------------------------------------------// +} // namespace inp +} // namespace celeritas diff --git a/test/celeritas/ext/GeantImporter.test.cc b/test/celeritas/ext/GeantImporter.test.cc index 68022e0a32..a42578192a 100644 --- a/test/celeritas/ext/GeantImporter.test.cc +++ b/test/celeritas/ext/GeantImporter.test.cc @@ -2132,15 +2132,27 @@ TEST_F(OpticalSurfaces, surfaces) } //---------------------------------------------------------------------------// -TEST_F(MucfBox, run) +TEST_F(MucfBox, static_data) { auto const& mucf = this->imported_data().mucf_physics; EXPECT_TRUE(mucf); - static double const expected_muon_energy_cdf_y[] = {1, 1}; - EXPECT_EQ(2, mucf.muon_energy_cdf.x.size()); - EXPECT_VEC_EQ(expected_muon_energy_cdf_y, mucf.muon_energy_cdf.y); + auto average = [](inp::Grid::VecDbl const& data) -> double { + double sum{0}; + for (auto y : data) + { + sum += y; + } + return sum / data.size(); + }; + + static size_t const expected_muon_energy_cdf_size = 21; + EXPECT_EQ(expected_muon_energy_cdf_size, mucf.muon_energy_cdf.x.size()); + EXPECT_EQ(expected_muon_energy_cdf_size, mucf.muon_energy_cdf.y.size()); + EXPECT_SOFT_EQ(0.55157437567861023, average(mucf.muon_energy_cdf.x)); + EXPECT_SOFT_EQ(11.250286274435437, average(mucf.muon_energy_cdf.y)); + // Dummy data auto const& cycle_f0 = mucf.cycle_rates[0]; static double const expected_cycle_rate_f0_y[] = {2, 2}; EXPECT_TRUE(cycle_f0); From 23c85856634925d69f4313748c6110aed9f2892f Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Mon, 29 Dec 2025 12:45:40 -0500 Subject: [PATCH 02/26] Add process/model skeletons --- src/celeritas/mucf/model/DTMixMucfModel.hh | 75 ++++++++++++++++++++++ src/celeritas/mucf/process/MucfProcess.cc | 66 +++++++++++++++++++ src/celeritas/mucf/process/MucfProcess.hh | 59 +++++++++++++++++ 3 files changed, 200 insertions(+) create mode 100644 src/celeritas/mucf/model/DTMixMucfModel.hh create mode 100644 src/celeritas/mucf/process/MucfProcess.cc create mode 100644 src/celeritas/mucf/process/MucfProcess.hh diff --git a/src/celeritas/mucf/model/DTMixMucfModel.hh b/src/celeritas/mucf/model/DTMixMucfModel.hh new file mode 100644 index 0000000000..d6cd37e531 --- /dev/null +++ b/src/celeritas/mucf/model/DTMixMucfModel.hh @@ -0,0 +1,75 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/model/DTMixMucfModel.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "celeritas/mat/MaterialParams.hh" +#include "celeritas/mucf/data/DTMixMucfData.hh" +#include "celeritas/mucf/executor/DTMixMucfExecutor.hh" // IWYU pragma: associated +#include "celeritas/phys/InteractionApplier.hh" // IWYU pragma: associated +#include "celeritas/phys/Model.hh" +#include "celeritas/phys/ParticleParams.hh" + +namespace celeritas +{ +//---------------------------------------------------------------------------// +/*! + * Muon-catalyzed fusion model for dd, dt, and tt molecules. + * + * In this model the executor performs the full muon-catalyzed fusion workflow. + * It forms a muonic d or t atom, samples which muonic molecule will be + * produced, selects the channel, and calls the appropriate interactor. + * + * The full set of "actions" is as follows, and in this ordering: + * - Define muon decay time to compete with the rest of the execution + * - Form muonic atom and select its spin + * - May execute atom spin flip or atom transfer + * - Form muonic molecule and select its spin + * - Calculate mean cycle time (time it takes from atom formation to fusion) + * - Confirm if fusion happens or the if the muon should decay + * - Call appropriate Interactor: Muon decay, or one of the muCF interactors + * + * \note This is an at-rest model. + */ +class DTMixMucfModel final : public Model, public StaticConcreteAction +{ + public: + //!@{ + //! \name Type aliases + using HostRef = HostCRef; + using DeviceRef = DeviceCRef; + //!@} + + // Construct from model ID and other necessary data + DTMixMucfModel(ActionId id, + ParticleParams const& particles, + MaterialParams const& materials); + + // Particle types and energy ranges that this model applies to + SetApplicability applicability() const final; + + // Get the microscopic cross sections for the given particle and material + XsTable micro_xs(Applicability) const final; + + // Apply the interaction kernel on host + void step(CoreParams const&, CoreStateHost&) const final; + + // Apply the interaction kernel on device + void step(CoreParams const&, CoreStateDevice&) const final; + + //!@{ + //! Access model data + HostRef const& host_ref() const { return data_.host_ref(); } + DeviceRef const& device_ref() const { return data_.device_ref(); } + //!@} + + private: + // Host/device storage and reference + CollectionMirror data_; +}; + +//---------------------------------------------------------------------------// +} // namespace celeritas diff --git a/src/celeritas/mucf/process/MucfProcess.cc b/src/celeritas/mucf/process/MucfProcess.cc new file mode 100644 index 0000000000..27a9256e19 --- /dev/null +++ b/src/celeritas/mucf/process/MucfProcess.cc @@ -0,0 +1,66 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/process/MucfProcess.cc +//---------------------------------------------------------------------------// +#include "MucfProcess.hh" + +#include + +#include "celeritas/mucf/model/DTMixMucfModel.hh" +#include "celeritas/phys/Model.hh" + +namespace celeritas +{ +//---------------------------------------------------------------------------// +/*! + * Construct from host data. + */ +MucfProcess::MucfProcess(SPConstParticles particles, SPConstMaterials materials) + : particles_(particles), materials_(materials) +{ + //! \todo Fix ImportProcessClass + CELER_EXPECT(particles_); + CELER_EXPECT(materials_); + CELER_EXPECT(particles_->find(pdg::mu_minus())); +} + +//---------------------------------------------------------------------------// +/*! + * Construct the models associated with this process. + */ +auto MucfProcess::build_models(ActionIdIter start_id) const -> VecModel +{ + return {std::make_shared( + *start_id++, *particles_, *materials_)}; +} + +//---------------------------------------------------------------------------// +/*! + * Get the interaction cross sections for the given energy range. + */ +auto MucfProcess::macro_xs(Applicability applic) const -> XsGrid +{ + return {}; +} + +//---------------------------------------------------------------------------// +/*! + * Get the energy loss for the given energy range. + */ +auto MucfProcess::energy_loss(Applicability) const -> EnergyLossGrid +{ + return {}; +} +//---------------------------------------------------------------------------// +/*! + * Name of the process. + */ +std::string_view MucfProcess::label() const +{ + return "Muon-catalyzed fusion"; +} + +//---------------------------------------------------------------------------// +} // namespace celeritas diff --git a/src/celeritas/mucf/process/MucfProcess.hh b/src/celeritas/mucf/process/MucfProcess.hh new file mode 100644 index 0000000000..bb14033d79 --- /dev/null +++ b/src/celeritas/mucf/process/MucfProcess.hh @@ -0,0 +1,59 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/process/MucfProcess.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "celeritas/Types.hh" +#include "celeritas/mat/MaterialParams.hh" +#include "celeritas/phys/Process.hh" + +namespace celeritas +{ +//---------------------------------------------------------------------------// +/*! + * Muon-catalyzed fusion of muonic dd, dt, or tt molecules. + * + * \note This is an at-rest process. + */ +class MucfProcess final : public Process +{ + public: + //!@{ + //! \name Type aliases + using SPConstParticles = std::shared_ptr; + using SPConstMaterials = std::shared_ptr; + //!@} + + public: + // Construct from shared and imported data + explicit MucfProcess(SPConstParticles particles, + SPConstMaterials materials); + + // Construct the models associated with this process + VecModel build_models(ActionIdIter start_id) const final; + + // Get the interaction cross sections for the given energy range + XsGrid macro_xs(Applicability range) const final; + + // Get the energy loss for the given energy range + EnergyLossGrid energy_loss(Applicability range) const final; + + //! Whether the integral method can be used to sample interaction length + bool supports_integral_xs() const final { return true; } //! \todo Check + + //! Whether the process applies when the particle is stopped + bool applies_at_rest() const final { return true; } + + // Name of the process + std::string_view label() const final; + + private: + SPConstParticles particles_; + SPConstMaterials materials_; +}; + +//---------------------------------------------------------------------------// +} // namespace celeritas From 1b8df651f9924ccd470e53309e655dbe00b3d6d5 Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Mon, 29 Dec 2025 12:48:03 -0500 Subject: [PATCH 03/26] Draft host/device data --- src/celeritas/Types.hh | 13 +- src/celeritas/mucf/data/DTMixMucfData.hh | 156 +++++++++++++++++++++++ 2 files changed, 167 insertions(+), 2 deletions(-) create mode 100644 src/celeritas/mucf/data/DTMixMucfData.hh diff --git a/src/celeritas/Types.hh b/src/celeritas/Types.hh index c2bc81764a..a2d39abf08 100644 --- a/src/celeritas/Types.hh +++ b/src/celeritas/Types.hh @@ -98,6 +98,9 @@ using UniformGridId = OpaqueId; //! Opaque index of a cross section grid using XsGridId = OpaqueId; +//! Opaque index of a muCF material component +using MuCfMatCompId = OpaqueId; + //---------------------------------------------------------------------------// // ENUMERATIONS //---------------------------------------------------------------------------// @@ -228,7 +231,10 @@ enum class CylAxis }; //---------------------------------------------------------------------------// -//! Muon-catalyzed fusion atoms +/*! + * Muonic atom selection from material data. This is *not* intended to be used + * by the transport loop. + */ enum class MucfMuonicAtom { deuterium, @@ -237,7 +243,10 @@ enum class MucfMuonicAtom }; //---------------------------------------------------------------------------// -//! Muon-catalyzed fusion molecules +/*! + * Muonic molecule selection from material data. This is *not* intended to be + * used by the transport loop. + */ enum class MucfMuonicMolecule { deuterium_deuterium, diff --git a/src/celeritas/mucf/data/DTMixMucfData.hh b/src/celeritas/mucf/data/DTMixMucfData.hh new file mode 100644 index 0000000000..43a57c2508 --- /dev/null +++ b/src/celeritas/mucf/data/DTMixMucfData.hh @@ -0,0 +1,156 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/data/DTMixMucfData.hh +//---------------------------------------------------------------------------// +#pragma once + +#include +#include + +#include "corecel/Macros.hh" +#include "corecel/cont/EnumArray.hh" +#include "corecel/grid/NonuniformGridData.hh" +#include "corecel/io/Join.hh" +#include "celeritas/Types.hh" +#include "celeritas/phys/ParticleParams.hh" + +namespace celeritas +{ +//---------------------------------------------------------------------------// +/*! + * ParticleIds used by the \c DTMixMucfModel . + */ +struct MucfParticles +{ + //! Primary + ParticleId mu_minus; + + //!@{ + //! Elementary particles and nuclei + ParticleId neutron; + ParticleId proton; + ParticleId alpha; + ParticleId helium3; + //!@} + + //!@{ + //! Muonic atoms + ParticleId muonic_hydrogen; + ParticleId muonic_deuteron; + ParticleId muonic_triton; + ParticleId muonic_alpha; + ParticleId muonic_helium3; + //!@} + + //! Check whether all particles are assigned + CELER_FUNCTION explicit operator bool() const + { + return mu_minus && neutron && proton && alpha && helium3 + && muonic_hydrogen && muonic_alpha && muonic_triton + && muonic_helium3; + } + + //! Assign from ParticleParams (anonymous free function in the model.cc?) + static MucfParticles from_params(ParticleParams const& particles) + { +#define MP_ADD(MEMBER) \ + pids.MEMBER = particles.find(pdg::MEMBER()); \ + if (!pids.MEMBER) \ + { \ + missing.push_back({#MEMBER, pdg::MEMBER()}); \ + } + + using PairStrPdg = std::pair; + std::vector missing; + MucfParticles pids; + + MP_ADD(mu_minus); + MP_ADD(neutron); + MP_ADD(proton); + MP_ADD(alpha); + //! \todo Decide whether to implement these PDGs in PDGNumber.hh +#if 0 + MP_ADD(helium_3); + MP_ADD(muonic_hydrogen); + MP_ADD(muonic_alpha); + MP_ADD(muonic_deuteron); + MP_ADD(muonic_triton); + MP_ADD(muonic_helium3); +#endif + CELER_VALIDATE( + missing.empty(), + << "missing particles required for muon-catalyzed fusion: " + << join(missing.begin(), missing.end(), ", ", [](PairStrPdg const& p) { + return p.first + " (PDG " + + std::to_string(p.second.unchecked_get()) + ')'; + })); + + return pids; + +#undef MP_ADD + } +}; + +//---------------------------------------------------------------------------// +/*! + * Data for for the \c DTMixMucfModel . + */ +template +struct DTMixMucfData +{ + template + using Items = Collection; + template + using MatCycleItems = Collection; + template + using MaterialItems = Collection; + using Grid = NonuniformGridRecord; + using CycleTimesArray = EnumArray>; + + //! Particle IDs + MucfParticles particles; + + //! Muon CDF energy grid for sampling outgoing muCF muons + Grid muon_energy_cdf; //! \todo Verify energy unit + Items reals; + + //!@{ + //! Material-dependent data calculated at model construction + //! \c PhysMatId indexed by \c MuCfMatCompId + MaterialItems matcompid_to_matid; + //! Cycle times per material: [mat_comp_id][muonic_molecule][spin_index] + MatCycleItems cycle_times; //!< In [s] + //! \todo Add mean atom spin flip times + //! \todo Add mean atom transfer times + //!@} + + //! \todo Check whether the data are assigned + explicit CELER_FUNCTION operator bool() const + { + //! \todo Finalize implementation + return particles && muon_energy_cdf && !matcompid_to_matid.empty() + && !cycle_times.empty() + && (matcompid_to_matid.size() == cycle_times.size()); + } + + //! Assign from another set of data + template + DTMixMucfData& operator=(DTMixMucfData const& other) + { + CELER_EXPECT(other); + + //! \todo Finish implementation + this->particles = other.particles; + this->reals = other.reals; + this->muon_energy_cdf = other.muon_energy_cdf; + this->matcompid_to_matid = other.matcompid_to_matid; + this->cycle_times = other.cycle_times; + + return *this; + } +}; + +//---------------------------------------------------------------------------// +} // namespace celeritas From 1a4d2ffe2164e34ace004cbe2a4cef0afd68b6f6 Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Mon, 29 Dec 2025 12:57:06 -0500 Subject: [PATCH 04/26] Add model implementation; draft host/device data initialization --- src/celeritas/mucf/model/DTMixMucfModel.cc | 146 ++++++++++++ src/celeritas/mucf/model/DTMixMucfModel.cu | 33 +++ .../model/detail/DTMixMaterialCalculator.cc | 212 ++++++++++++++++++ .../model/detail/DTMixMaterialCalculator.hh | 142 ++++++++++++ 4 files changed, 533 insertions(+) create mode 100644 src/celeritas/mucf/model/DTMixMucfModel.cc create mode 100644 src/celeritas/mucf/model/DTMixMucfModel.cu create mode 100644 src/celeritas/mucf/model/detail/DTMixMaterialCalculator.cc create mode 100644 src/celeritas/mucf/model/detail/DTMixMaterialCalculator.hh diff --git a/src/celeritas/mucf/model/DTMixMucfModel.cc b/src/celeritas/mucf/model/DTMixMucfModel.cc new file mode 100644 index 0000000000..ac4c4e9e56 --- /dev/null +++ b/src/celeritas/mucf/model/DTMixMucfModel.cc @@ -0,0 +1,146 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/model/DTMixMucfModel.cc +//---------------------------------------------------------------------------// +#include "DTMixMucfModel.hh" + +#include + +#include "corecel/inp/Grid.hh" +#include "celeritas/global/ActionLauncher.hh" +#include "celeritas/global/TrackExecutor.hh" +#include "celeritas/grid/NonuniformGridBuilder.hh" +#include "celeritas/inp/MucfPhysicsData.hh" +#include "celeritas/mucf/executor/DTMixMucfExecutor.hh" // IWYU pragma: associated +#include "celeritas/phys/InteractionApplier.hh" // IWYU pragma: associated +#include "celeritas/phys/PDGNumber.hh" + +#include "detail/DTMixMaterialCalculator.hh" + +namespace celeritas +{ +//---------------------------------------------------------------------------// +/*! + * Construct from model ID and other necessary data. + * + * Most of the muon-catalyzed fusion data is static throughout the simulation, + * as it is only material-dependent (DT mixture and temperature). Therefore, + * most grids can be host-only and used to calculate final values, which are + * then cached and copied to device. The exception to this is the muon energy + * CDF grid, needed to sample the final state of the outgoing muon after a muCF + * interaction. + * + * \todo Correctly update \c ImportProcessClass and \c ImportModelClass . These + * operate under the assumption that there is a one-to-one equivalente between + * Geant4 and Celeritas. But for muCF, everything is done via one + * process/model/executor in Celeritas, whereas in Geant4 atom formation, spin + * flip, atom transfer, etc., are are all separate processes. + */ +DTMixMucfModel::DTMixMucfModel(ActionId id, + ParticleParams const& particles, + MaterialParams const& materials) + : StaticConcreteAction( + id, + "dt-mucf", + R"(interact by muon forming and fusing a dd, dt, or tt muonic molecule)") +{ + CELER_EXPECT(id); + + using VecCycleTimes + = std::vector; + using VecMatId = std::vector; + + // Initialize static muCF physics input data + //! \todo This may be replaced by user-provided data in the future + inp::MucfPhysics inp_data = inp::MucfPhysics::from_default(); + CELER_EXPECT(inp_data); + + // Load particle IDs + HostVal host_data; + host_data.particles = MucfParticles::from_params(particles); + + // Copy muon energy CDF data using NonuniformGridBuilder + NonuniformGridBuilder build_grid_record{&host_data.reals}; + host_data.muon_energy_cdf = build_grid_record(inp_data.muon_energy_cdf); + + // Calculate and cache quantities for all materials with dt mixtures + VecCycleTimes vec_cycle_times; + VecMatId vec_physmatid; + for (auto const& matid : range(materials.num_materials())) + { + auto const& mat_view = materials.get(PhysMatId{matid}); + + // Construct calculator for this material + detail::DTMixMaterialCalculator mat_calculator(mat_view); + if (!mat_calculator) + { + // Skip non-dt mixture materials + continue; + } + + // Store material ID and calculated data + //! \todo Store mean atom spin flip and transfer times + vec_physmatid.push_back(PhysMatId{matid}); + vec_cycle_times.push_back(mat_calculator.cycle_times()); + } + + make_builder(&host_data.matcompid_to_matid) + .insert_back(vec_physmatid.begin(), vec_physmatid.end()); + make_builder(&host_data.cycle_times) + .insert_back(vec_cycle_times.begin(), vec_cycle_times.end()); + + // Copy to device + data_ = CollectionMirror{std::move(host_data)}; + CELER_ENSURE(this->data_); +} + +//---------------------------------------------------------------------------// +/*! + * Particle types and energy ranges that this model applies to. + */ +auto DTMixMucfModel::applicability() const -> SetApplicability +{ + Applicability applic; + applic.particle = this->host_ref().particles.mu_minus; + // At-rest model + applic.lower = zero_quantity(); + applic.upper = zero_quantity(); + + return {applic}; +} + +//---------------------------------------------------------------------------// +/*! + * At-rest model does not require microscopic cross sections. + */ +auto DTMixMucfModel::micro_xs(Applicability) const -> XsTable +{ + return {}; +} + +//---------------------------------------------------------------------------// +/*! + * Interact with host data. + */ +void DTMixMucfModel::step(CoreParams const& params, CoreStateHost& state) const +{ + auto execute = make_action_track_executor( + params.ptr(), + state.ptr(), + this->action_id(), + InteractionApplier{DTMixMucfExecutor{this->host_ref()}}); + return launch_action(*this, params, state, execute); +} + +//---------------------------------------------------------------------------// +#if !CELER_USE_DEVICE +void DTMixMucfModel::step(CoreParams const&, CoreStateDevice&) const +{ + CELER_NOT_CONFIGURED("CUDA OR HIP"); +} +#endif + +//---------------------------------------------------------------------------// +} // namespace celeritas diff --git a/src/celeritas/mucf/model/DTMixMucfModel.cu b/src/celeritas/mucf/model/DTMixMucfModel.cu new file mode 100644 index 0000000000..4ac694e935 --- /dev/null +++ b/src/celeritas/mucf/model/DTMixMucfModel.cu @@ -0,0 +1,33 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/model/DTMixMucfModel.cu +//---------------------------------------------------------------------------// +#include "DTMixMucfModel.hh" + +#include "celeritas/global/ActionLauncher.device.hh" +#include "celeritas/global/CoreParams.hh" +#include "celeritas/global/CoreState.hh" +#include "celeritas/global/TrackExecutor.hh" +#include "celeritas/mucf/executor/DTMixMucfExecutor.hh" +#include "celeritas/phys/InteractionApplier.hh" + +namespace celeritas +{ +//---------------------------------------------------------------------------// +/*! + * Interact with device data. + */ +void DTMixMucfModel::step(CoreParams const& params, CoreStateDevice& state) const +{ + auto execute = make_action_track_executor( + params.ptr(), + state.ptr(), + this->action_id(), + InteractionApplier{DTMixMucfExecutor{this->device_ref()}}); + static ActionLauncher const launch_kernel(*this); + launch_kernel(*this, params, state, execute); +} +//---------------------------------------------------------------------------// +} // namespace celeritas diff --git a/src/celeritas/mucf/model/detail/DTMixMaterialCalculator.cc b/src/celeritas/mucf/model/detail/DTMixMaterialCalculator.cc new file mode 100644 index 0000000000..b4da247641 --- /dev/null +++ b/src/celeritas/mucf/model/detail/DTMixMaterialCalculator.cc @@ -0,0 +1,212 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/model/detail/DTMixMaterialCalculator.cc +//---------------------------------------------------------------------------// +#include "DTMixMaterialCalculator.hh" + +#include "corecel/Assert.hh" + +namespace celeritas +{ +namespace detail +{ +//---------------------------------------------------------------------------// +/*! + * Construct with material data. + * + * Calculates and caches material-dependent properties needed by the + * \c DTMixMucfModel . If the material does not contain deuterium and/or + * tritium the object's operator bool will return false. + */ +DTMixMaterialCalculator::DTMixMaterialCalculator(MaterialView const& material) + : material_(material) +{ + for (auto elcompid : range(material_.num_elements())) + { + auto const& element_view + = material_.element_record(ElementComponentId{elcompid}); + if (element_view.atomic_number() != AtomicNumber{1}) + { + // Skip non-hydrogen elements + continue; + } + + has_isotope_ = {false, false}; + for (auto el_comp : range(element_view.num_isotopes())) + { + auto iso_view + = element_view.isotope_record(IsotopeComponentId{el_comp}); + auto mass = iso_view.atomic_mass_number(); + if (mass == AtomicMassNumber{1}) + { + // Skip protium + continue; + } + + if (auto const atom = this->from_mass_number(mass); + atom < MucfMuonicAtom::size_) + { + // D and/or t isotopes found; calculate properties + has_isotope_[atom] = true; + lhd_densities_ = calc_lhd_densities(); + eq_densities_ = calc_equilibrium_densities(); + cycle_times_ = calc_cycle_times(element_view); + } + } + } +} + +//---------------------------------------------------------------------------// +/*! + * Convert dt mixture densities to units of liquid hydrogen density. + * + * Used during cycle time calculations. + */ +DTMixMaterialCalculator::LhdArray DTMixMaterialCalculator::calc_lhd_densities() +{ + LhdArray result; + + //! \todo Implement + + return result; +} + +//---------------------------------------------------------------------------// +/*! + * Calculate dt mixture densities after reaching thermodynamical + * equilibrium. + * + * Used during cycle time calculations. + */ +DTMixMaterialCalculator::EquilibriumArray +DTMixMaterialCalculator::calc_equilibrium_densities() +{ + EquilibriumArray result; + + //! \todo Implement + + return result; +} + +//---------------------------------------------------------------------------// +/*! + * Calculate fusion mean cycle times. + * + * This is designed to work with the user's material definition being either: + * - Single element, multiple isotopes (H element, with H, d, and t isotopes); + * or + * - Multiple elements, single isotope each (separate H, d, and t elements). + */ +DTMixMaterialCalculator::CycleTimesArray +DTMixMaterialCalculator::calc_cycle_times(ElementView const& element) +{ + CycleTimesArray result; + for (auto el_comp : range(element.num_isotopes())) + { + auto iso_view = element.isotope_record(IsotopeComponentId{el_comp}); + + // Select possible muonic atom based on the isotope/element mass number + auto atom = this->from_mass_number(iso_view.atomic_mass_number()); + switch (atom) + { + // Calculate cycle times for dd molecules + case MucfMuonicAtom::deuterium: { + result[MucfMuonicMolecule::deuterium_deuterium] + = this->calc_dd_cycle(); + if (has_isotope_[MucfMuonicAtom::tritium]) + { + // Calculate cycle times for dt molecules + result[MucfMuonicMolecule::deuterium_tritium] + = this->calc_dt_cycle(); + } + break; + } + // Calculate cycle times for tt molecules + case MucfMuonicAtom::tritium: { + result[MucfMuonicMolecule::tritium_tritium] + = this->calc_tt_cycle(); + break; + } + default: + CELER_ASSERT_UNREACHABLE(); + } + } + return result; +} + +//---------------------------------------------------------------------------// +/*! + * Calculate dd muonic molecules cycle times from material properties and grid + * data. + * + * Cycle times for dd molecules come from F = 0 and F = 1 spin states. + */ +Array DTMixMaterialCalculator::calc_dd_cycle() +{ + Array result; + + //! \todo Implement + + // Reactive states are F = 0 and F = 1 + CELER_ENSURE(result[0] >= 0 && result[1] >= 0); + return result; +} + +//---------------------------------------------------------------------------// +/*! + * Calculate dt muonic molecules cycle times from material properties and grid + * data. + * + * Cycle times for dt molecules come from F = 1/2 and F = 3/2 spin states. + */ +Array DTMixMaterialCalculator::calc_dt_cycle() +{ + Array result; + + //! \todo Implement + + // Reactive states are F = 1/2 and F = 3/2 + CELER_ENSURE(result[0] >= 0 && result[1] >= 0); + return result; +} + +//---------------------------------------------------------------------------// +/*! + * Calculate tt muonic molecules cycle times from material properties and grid + * data. + * + * Cycle times for tt molecules come only from the F = 1/2 spin state. + */ +Array DTMixMaterialCalculator::calc_tt_cycle() +{ + Array result; + + //! \todo Implement + + // Only F = 1/2 is reactive + CELER_ENSURE(result[0] >= 0 && result[1] == 0); + return result; +} + +//---------------------------------------------------------------------------// +/*! + * Return \c MucfMuonicAtom from a given atomic mass number. + */ +MucfMuonicAtom DTMixMaterialCalculator::from_mass_number(AtomicMassNumber mass) +{ + if (mass == AtomicMassNumber{2}) + { + return MucfMuonicAtom::deuterium; + } + if (mass == AtomicMassNumber{3}) + { + return MucfMuonicAtom::tritium; + } + return MucfMuonicAtom::size_; +} + +//---------------------------------------------------------------------------// +} // namespace detail +} // namespace celeritas diff --git a/src/celeritas/mucf/model/detail/DTMixMaterialCalculator.hh b/src/celeritas/mucf/model/detail/DTMixMaterialCalculator.hh new file mode 100644 index 0000000000..84a93430f0 --- /dev/null +++ b/src/celeritas/mucf/model/detail/DTMixMaterialCalculator.hh @@ -0,0 +1,142 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/model/detail/DTMixMaterialCalculator.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "corecel/cont/EnumArray.hh" +#include "celeritas/inp/MucfPhysics.hh" +#include "celeritas/mat/MaterialView.hh" +#include "celeritas/mucf/data/DTMixMucfData.hh" + +namespace celeritas +{ +namespace detail +{ +//---------------------------------------------------------------------------// +/*! + * Enum for safely accessing hydrogen isoprotologues. + * + * Hydrogen isoprotologue molecules are: + * - Homonuclear: \f$ ^2H \f$, \f$ ^2d \f$, and \f$ ^2t \f$ + * - Heteronuclear: hd, ht, and dt. + * + * \note Muon-catalyzed fusion data is only applicable to a material with + * concentrations in thermodynamic equilibrium. This equilibrium is calculated + * at model construction from the material temperature and its h, d, and t + * fractions. + */ +enum class MucfIsoprotologueMolecule +{ + protium_protium, + protium_deuterium, + protium_tritium, + deuterium_tritium, + tritium_tritium, + size_ +}; + +//---------------------------------------------------------------------------// +/*! + * Calculate material-dependent quantities for muon-catalyzed fusion. + * + * This class calculates all the muCF data that can be cached during model + * construction. + * Use its operator bool to store material data into \c DTMixMucfData : + * \code + for (auto matid : range(materials.num_materials())) + { + auto mat_view = materials.material_view(PhysMatId{matid}); + DTMixMaterialCalculator material_calculator(mat_view); + if (material_calculator) + { + // Valid d-t mixture material; Store data + } + } + * \endcode + */ +class DTMixMaterialCalculator +{ + public: + //!@{ + //! \name Type aliases + using CycleTimesArray = EnumArray>; + //!@} + + //! Construct with material data and calculate all quantities + DTMixMaterialCalculator(MaterialView const& material); + + //! Get mean cycle times + CycleTimesArray cycle_times() const { return cycle_times_; } + + //! Check if the material is valid for muon-catalyzed fusion + explicit operator bool() const { return !cycle_times_.empty(); } + + private: + using LhdArray = EnumArray; + using EquilibriumArray = EnumArray; + using AtomicMassNumber = AtomicNumber; + + //// DATA //// + + MaterialView material_; + LhdArray lhd_densities_; + EquilibriumArray eq_densities_; + EnumArray has_isotope_; + CycleTimesArray cycle_times_; + + //// LOCAL SCALARS //// + //! \todo Values are the same used by Acceleron and may need revisiting. + // { + // Atomic masses + static constexpr units::AmuMass protium() + { + return units::AmuMass{1.007825031898}; + } + + static constexpr units::AmuMass deuterium() + { + return units::AmuMass{2.014101777844}; + } + + static constexpr units::AmuMass tritium() + { + return units::AmuMass{3.016049281320}; + } + //} + + // Liquid hydrogen density (LHD) unit [1/cm^3] + static constexpr auto liquid_hydrogen_density() + { + return units::InvCcDensity{4.25e22}; + } + + //// HELPER FUNCTIONS //// + + // Calculate dt mixture densities in units of liquid hydrogen density + LhdArray calc_lhd_densities(); + + // Calculate thermal equilibrium densities + EquilibriumArray calc_equilibrium_densities(); + + // Return muonic atom from given atomic mass number + MucfMuonicAtom from_mass_number(AtomicMassNumber mass); + + // Calculate mean fusion cycle times for all reactive muonic molecules + CycleTimesArray calc_cycle_times(ElementView const& element); + + // Calculate mean fusion cycle times for dd muonic molecules + Array calc_dd_cycle(); + + // Calculate mean fusion cycle times for dt muonic molecules + Array calc_dt_cycle(); + + // Calculate mean fusion cycle times for tt muonic molecules + Array calc_tt_cycle(); +}; + +//---------------------------------------------------------------------------// +} // namespace detail +} // namespace celeritas From 7e9e52734b2e7db4503b951ecdafab80294f8412 Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Mon, 29 Dec 2025 12:59:26 -0500 Subject: [PATCH 05/26] Add mucf executor + interactors + helper classes skeletons --- .../mucf/executor/DTMixMucfExecutor.hh | 144 +++++++++++++++++ .../mucf/executor/detail/DDChannelSelector.hh | 73 +++++++++ .../mucf/executor/detail/DTChannelSelector.hh | 73 +++++++++ .../detail/DTMixMuonicAtomSelector.hh | 61 +++++++ .../detail/DTMixMuonicMoleculeSelector.hh | 66 ++++++++ .../executor/detail/MuonicAtomSpinSelector.hh | 76 +++++++++ .../detail/MuonicMoleculeSpinSelector.hh | 78 +++++++++ .../mucf/executor/detail/TTChannelSelector.hh | 74 +++++++++ .../mucf/interactor/DDMucfInteractor.hh | 152 ++++++++++++++++++ .../mucf/interactor/DTMucfInteractor.hh | 145 +++++++++++++++++ .../mucf/interactor/TTMucfInteractor.hh | 145 +++++++++++++++++ src/celeritas/phys/PDGNumber.hh | 4 + 12 files changed, 1091 insertions(+) create mode 100644 src/celeritas/mucf/executor/DTMixMucfExecutor.hh create mode 100644 src/celeritas/mucf/executor/detail/DDChannelSelector.hh create mode 100644 src/celeritas/mucf/executor/detail/DTChannelSelector.hh create mode 100644 src/celeritas/mucf/executor/detail/DTMixMuonicAtomSelector.hh create mode 100644 src/celeritas/mucf/executor/detail/DTMixMuonicMoleculeSelector.hh create mode 100644 src/celeritas/mucf/executor/detail/MuonicAtomSpinSelector.hh create mode 100644 src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh create mode 100644 src/celeritas/mucf/executor/detail/TTChannelSelector.hh create mode 100644 src/celeritas/mucf/interactor/DDMucfInteractor.hh create mode 100644 src/celeritas/mucf/interactor/DTMucfInteractor.hh create mode 100644 src/celeritas/mucf/interactor/TTMucfInteractor.hh diff --git a/src/celeritas/mucf/executor/DTMixMucfExecutor.hh b/src/celeritas/mucf/executor/DTMixMucfExecutor.hh new file mode 100644 index 0000000000..7abfc107c1 --- /dev/null +++ b/src/celeritas/mucf/executor/DTMixMucfExecutor.hh @@ -0,0 +1,144 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/executor/DTMixMucfExecutor.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "corecel/Assert.hh" +#include "celeritas/global/CoreTrackView.hh" +#include "celeritas/mat/ElementView.hh" +#include "celeritas/mucf/data/DTMixMucfData.hh" +#include "celeritas/mucf/interactor/DDMucfInteractor.hh" +#include "celeritas/mucf/interactor/DTMucfInteractor.hh" +#include "celeritas/mucf/interactor/TTMucfInteractor.hh" + +#include "detail/DDChannelSelector.hh" +#include "detail/DTChannelSelector.hh" +#include "detail/DTMixMuonicAtomSelector.hh" +#include "detail/DTMixMuonicMoleculeSelector.hh" +#include "detail/MuonicAtomSpinSelector.hh" +#include "detail/MuonicMoleculeSpinSelector.hh" +#include "detail/TTChannelSelector.hh" + +namespace celeritas +{ +//---------------------------------------------------------------------------// +struct DTMixMucfExecutor +{ + inline CELER_FUNCTION Interaction + operator()(celeritas::CoreTrackView const& track); + + NativeCRef data; +}; + +//---------------------------------------------------------------------------// +/*! + * Execute muon-catalyzed fusion for muonic dd, dt, or tt molecules. + */ +CELER_FUNCTION Interaction +DTMixMucfExecutor::operator()(celeritas::CoreTrackView const& track) +{ + auto phys_step_view = track.physics_step(); + auto elcomp_id = phys_step_view.element(); + CELER_ASSERT(elcomp_id); + + auto element = track.material().material_record().element_record(elcomp_id); + CELER_ASSERT(element.atomic_number() == AtomicNumber{1}); // Must be H + + auto rng = track.rng(); + + // Muon decay may compete against other "actions" in this executor + real_type const decay_len{}; //! \todo Set muon decay interaction length + + // Form d or t muonic atom + detail::DTMixMuonicAtomSelector form_atom; + auto muonic_atom = form_atom(rng); + + // Select atom spin via a helper class + detail::MuonicAtomSpinSelector select_atom_spin(muonic_atom); + + // { + // Competing at-rest processes which add to the total track time + //! \todo Muonic atom transfer + //! \todo Muonic atom spin flip + // } + + // Form dd, dt, or tt muonic molecule + detail::DTMixMuonicMoleculeSelector form_muonic_molecule; + auto muonic_molecule = form_muonic_molecule(rng); + + // Select molecule spin + detail::MuonicMoleculeSpinSelector select_molecule_spin(muonic_molecule); + auto const molecule_spin = select_molecule_spin(rng); + + // Load cycle time for the selected molecule + auto find_component = [&](PhysMatId matid) -> MuCfMatCompId { + for (auto i : range(data.matcompid_to_matid.size())) + { + auto const comp_id = MuCfMatCompId{i}; + if (data.matcompid_to_matid[comp_id] == matid) + { + return MuCfMatCompId{i}; + } + } + return MuCfMatCompId{}; + }; + + auto const mat_comp = find_component(track.material().material_id()); + CELER_ASSERT(mat_comp); + auto const cycle_time + = data.cycle_times[mat_comp][muonic_molecule][molecule_spin]; + CELER_ASSERT(cycle_time > 0); + + // Check if muon decays before fusion happens + real_type const mucf_len = cycle_time * track.sim().step_length(); + if (decay_len < mucf_len) + { + // Muon decays and halts the interaction + //! \todo Update track time and return muon decay interactor + } + + //! \todo Correct track time update? Or should be done in Interactors? + track.sim().add_time(cycle_time); + + // Fuse molecule and generate secondaries + //! \todo Maybe move the channel selectors into the interactors + auto allocate_secondaries = phys_step_view.make_secondary_allocator(); + Interaction result; + switch (muonic_molecule) + { + case MucfMuonicMolecule::deuterium_deuterium: { + // Return DD interaction + DDMucfInteractor interact( + data, detail::DDChannelSelector()(rng), allocate_secondaries); + result = interact(rng); + break; + } + case MucfMuonicMolecule::deuterium_tritium: { + // Return DT interaction + DTMucfInteractor interact( + data, detail::DTChannelSelector()(rng), allocate_secondaries); + result = interact(rng); + break; + } + case MucfMuonicMolecule::tritium_tritium: { + // Return TT interaction + TTMucfInteractor interact( + data, detail::TTChannelSelector()(rng), allocate_secondaries); + result = interact(rng); + break; + } + default: + CELER_ASSERT_UNREACHABLE(); + } + + //! \todo Muon stripping: strip muon from muonic atom secondaries + // May be added as a separate discrete process in the stepping loop + + return result; +} + +//---------------------------------------------------------------------------// +} // namespace celeritas diff --git a/src/celeritas/mucf/executor/detail/DDChannelSelector.hh b/src/celeritas/mucf/executor/detail/DDChannelSelector.hh new file mode 100644 index 0000000000..f2793b0979 --- /dev/null +++ b/src/celeritas/mucf/executor/detail/DDChannelSelector.hh @@ -0,0 +1,73 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/executor/detail/DDChannelSelector.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "celeritas/mucf/interactor/DDMucfInteractor.hh" + +namespace celeritas +{ +namespace detail +{ +//---------------------------------------------------------------------------// +/*! + * Select final channel for muonic dd molecules. + * + * This selection already accounts for sticking, as that is one of the possible + * outcomes. + */ +class DDChannelSelector +{ + public: + //!@{ + //! \name Type aliases + using Channel = DDMucfInteractor::Channel; + //!@} + + public: + //! Construct with args; \todo Update documentation + inline DDChannelSelector(/* args */); + + // Select fusion channel to be used by the interactor + template + inline CELER_FUNCTION Channel operator()(Engine&); +}; + +//---------------------------------------------------------------------------// +// INLINE DEFINITIONS +//---------------------------------------------------------------------------// +/*! + * Construct with args. + * + * \todo Update documentation + */ +DDChannelSelector::DDChannelSelector(/* args */) +{ + CELER_NOT_IMPLEMENTED("Mucf dd fusion channel selection"); +} + +//---------------------------------------------------------------------------// +/*! + * Return a sampled channel to be used as input in the dd muCF interactor. + * + * \sa celeritas::DDMucfInteractor + */ +template +inline CELER_FUNCTION DDChannelSelector::Channel +DDChannelSelector::operator()(Engine&) +{ + Channel result{Channel::size_}; + + //! \todo Implement + // Final channel selection already takes into account sticking. + + CELER_ENSURE(result < Channel::size_); + return result; +} + +//---------------------------------------------------------------------------// +} // namespace detail +} // namespace celeritas diff --git a/src/celeritas/mucf/executor/detail/DTChannelSelector.hh b/src/celeritas/mucf/executor/detail/DTChannelSelector.hh new file mode 100644 index 0000000000..28483dea53 --- /dev/null +++ b/src/celeritas/mucf/executor/detail/DTChannelSelector.hh @@ -0,0 +1,73 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/executor/detail/DTChannelSelector.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "celeritas/mucf/interactor/DTMucfInteractor.hh" + +namespace celeritas +{ +namespace detail +{ +//---------------------------------------------------------------------------// +/*! + * Select final channel for muonic dd molecules. + * + * This selection already accounts for sticking, as that is one of the possible + * outcomes. + */ +class DTChannelSelector +{ + public: + //!@{ + //! \name Type aliases + using Channel = DTMucfInteractor::Channel; + //!@} + + public: + //! Construct with args; \todo Update documentation + inline DTChannelSelector(/* args */); + + // Select fusion channel to be used by the interactor + template + inline CELER_FUNCTION Channel operator()(Engine&); +}; + +//---------------------------------------------------------------------------// +// INLINE DEFINITIONS +//---------------------------------------------------------------------------// +/*! + * Construct with args. + * + * \todo Update documentation + */ +DTChannelSelector::DTChannelSelector(/* args */) +{ + CELER_NOT_IMPLEMENTED("Mucf dt fusion channel selection"); +} + +//---------------------------------------------------------------------------// +/*! + * Return a sampled channel to be used as input in the dt muCF interactor. + * + * \sa celeritas::DTMucfInteractor + */ +template +inline CELER_FUNCTION DTChannelSelector::Channel +DTChannelSelector::operator()(Engine&) +{ + Channel result{Channel::size_}; + + //! \todo Implement + // Final channel selection already takes into account sticking. + + CELER_ENSURE(result < Channel::size_); + return result; +} + +//---------------------------------------------------------------------------// +} // namespace detail +} // namespace celeritas diff --git a/src/celeritas/mucf/executor/detail/DTMixMuonicAtomSelector.hh b/src/celeritas/mucf/executor/detail/DTMixMuonicAtomSelector.hh new file mode 100644 index 0000000000..a371141417 --- /dev/null +++ b/src/celeritas/mucf/executor/detail/DTMixMuonicAtomSelector.hh @@ -0,0 +1,61 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/executor/detail/DTMixMuonicAtomSelector.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "corecel/cont/EnumArray.hh" +#include "celeritas/Types.hh" + +namespace celeritas +{ +namespace detail +{ +//---------------------------------------------------------------------------// +/*! + * Select a muonic atom given the material information. + */ +class DTMixMuonicAtomSelector +{ + public: + //! Construct with args; \todo Update documentation + inline DTMixMuonicAtomSelector(/* args */); + + // Select muonic atom + template + inline CELER_FUNCTION MucfMuonicAtom operator()(Engine&); +}; + +//---------------------------------------------------------------------------// +// INLINE DEFINITIONS +//---------------------------------------------------------------------------// +/*! + * Construct with args. + + * \todo Update documentation + */ +DTMixMuonicAtomSelector::DTMixMuonicAtomSelector(/* args */) +{ + CELER_NOT_IMPLEMENTED("Mucf muonic atom selection"); +} + +//---------------------------------------------------------------------------// +/*! + * Return selected muonic atom. + */ +template +inline CELER_FUNCTION MucfMuonicAtom DTMixMuonicAtomSelector::operator()(Engine&) +{ + MucfMuonicAtom result{MucfMuonicAtom::size_}; + + //! \todo Implement + + CELER_ENSURE(result < MucfMuonicAtom::size_); + return result; +} + +//---------------------------------------------------------------------------// +} // namespace detail +} // namespace celeritas diff --git a/src/celeritas/mucf/executor/detail/DTMixMuonicMoleculeSelector.hh b/src/celeritas/mucf/executor/detail/DTMixMuonicMoleculeSelector.hh new file mode 100644 index 0000000000..6558191a00 --- /dev/null +++ b/src/celeritas/mucf/executor/detail/DTMixMuonicMoleculeSelector.hh @@ -0,0 +1,66 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/executor/detail/DTMixMuonicMoleculeSelector.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "celeritas/Types.hh" + +namespace celeritas +{ +namespace detail +{ +//---------------------------------------------------------------------------// +/*! + * Select a muonic molecule by calculating the interaction lengths of the + * possible molecule formations. + * + * This is the equivalent of Geant4's + * \c G4VRestProcess::AtRestGetPhysicalInteractionLength + */ +class DTMixMuonicMoleculeSelector +{ + public: + //! Construct with args; \todo Update documentation + inline DTMixMuonicMoleculeSelector(/* args */); + + // Select muonic molecule + template + inline CELER_FUNCTION MucfMuonicMolecule operator()(Engine&); +}; + +//---------------------------------------------------------------------------// +// INLINE DEFINITIONS +//---------------------------------------------------------------------------// +/*! + * Construct with args. + * + * \todo Update documentation + */ +DTMixMuonicMoleculeSelector::DTMixMuonicMoleculeSelector(/* args */) +{ + //! \todo Implement + CELER_NOT_IMPLEMENTED("Mucf muonic molecule selection"); +} + +//---------------------------------------------------------------------------// +/*! + * Return selected muonic molecule. + */ +template +inline CELER_FUNCTION MucfMuonicMolecule +DTMixMuonicMoleculeSelector::operator()(Engine&) +{ + MucfMuonicMolecule result{MucfMuonicMolecule::size_}; + + //! \todo Implement + + CELER_ENSURE(result < MucfMuonicMolecule::size_); + return result; +} + +//---------------------------------------------------------------------------// +} // namespace detail +} // namespace celeritas diff --git a/src/celeritas/mucf/executor/detail/MuonicAtomSpinSelector.hh b/src/celeritas/mucf/executor/detail/MuonicAtomSpinSelector.hh new file mode 100644 index 0000000000..853950c4bd --- /dev/null +++ b/src/celeritas/mucf/executor/detail/MuonicAtomSpinSelector.hh @@ -0,0 +1,76 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/executor/detail/MuonicAtomSpinSelector.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "corecel/Assert.hh" +#include "celeritas/Types.hh" + +namespace celeritas +{ +namespace detail +{ +//---------------------------------------------------------------------------// +/*! + * Select muonic atom spin, in units of \f$ \frac{\hbar}{2} \f$. + * + * Applicable to \f$ (p)_\mu \f$, \f$ (d)_\mu \f$, and \f$ (t)_\mu \f$. + */ +class MuonicAtomSpinSelector +{ + public: + // Construct with muonic atom type + inline CELER_FUNCTION MuonicAtomSpinSelector(MucfMuonicAtom atom); + + // Sample and return a spin value in units of hbar / 2 + template + inline CELER_FUNCTION size_type operator()(Engine&); + + private: + MucfMuonicAtom atom_; + //! \todo Add constant atom spin sampling rejection fractions +}; + +//---------------------------------------------------------------------------// +// INLINE DEFINITIONS +//---------------------------------------------------------------------------// +/*! + * Construct with muonic atom. + */ +MuonicAtomSpinSelector::MuonicAtomSpinSelector(MucfMuonicAtom atom) + : atom_(atom) +{ + CELER_EXPECT(atom_ < MucfMuonicAtom::size_); + CELER_NOT_IMPLEMENTED("Mucf muonic atom spin selection"); +} + +//---------------------------------------------------------------------------// +/*! + * Return selected spin, in units of \f$ \hbar / 2 \f$. + */ +template +CELER_FUNCTION size_type MuonicAtomSpinSelector::operator()(Engine&) +{ + size_type result{}; + + //! \todo Implement + + switch (atom_) + { + case MucfMuonicAtom::deuterium: + break; + case MucfMuonicAtom::tritium: + // Protium and tritium samplings are equivalent + break; + default: + CELER_ASSERT_UNREACHABLE(); + } + return result; +} + +//---------------------------------------------------------------------------// +} // namespace detail +} // namespace celeritas diff --git a/src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh b/src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh new file mode 100644 index 0000000000..51ddedb741 --- /dev/null +++ b/src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh @@ -0,0 +1,78 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "corecel/Assert.hh" +#include "celeritas/Types.hh" + +namespace celeritas +{ +namespace detail +{ +//---------------------------------------------------------------------------// +/*! + * Select muonic molecule spin, in units of \f$ \frac{\hbar}{2} \f$. + * + * Applicable to \f$ (dd)_\mu \f$, * \f$ (dt)_\mu \f$ , and \f$ (tt)_\mu \f$. + */ +class MuonicMoleculeSpinSelector +{ + public: + // Construct with muonic molecule type + inline CELER_FUNCTION + MuonicMoleculeSpinSelector(MucfMuonicMolecule molecule); + + // Sample and return a spin value in units of hbar / 2 + template + inline CELER_FUNCTION size_type operator()(Engine&); + + private: + MucfMuonicMolecule molecule_; +}; + +//---------------------------------------------------------------------------// +// INLINE DEFINITIONS +//---------------------------------------------------------------------------// +/*! + * Construct with muonic molecule. + */ +MuonicMoleculeSpinSelector::MuonicMoleculeSpinSelector( + MucfMuonicMolecule molecule) + : molecule_(molecule) +{ + CELER_EXPECT(molecule_ < MucfMuonicMolecule::size_); + CELER_NOT_IMPLEMENTED("Mucf muonic molecule spin selection"); +} + +//---------------------------------------------------------------------------// +/*! + * Return selected spin, in units of \f$ \hbar / 2 \f$. + */ +template +CELER_FUNCTION size_type MuonicMoleculeSpinSelector::operator()(Engine&) +{ + size_type result{}; + + //! \todo Implement + + switch (molecule_) + { + case MucfMuonicMolecule::deuterium_deuterium: + break; + case MucfMuonicMolecule::deuterium_tritium: + break; + case MucfMuonicMolecule::tritium_tritium: + break; + default: + CELER_ASSERT_UNREACHABLE(); + } + return result; +} + +//---------------------------------------------------------------------------// +} // namespace detail +} // namespace celeritas diff --git a/src/celeritas/mucf/executor/detail/TTChannelSelector.hh b/src/celeritas/mucf/executor/detail/TTChannelSelector.hh new file mode 100644 index 0000000000..1ded52adb7 --- /dev/null +++ b/src/celeritas/mucf/executor/detail/TTChannelSelector.hh @@ -0,0 +1,74 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/executor/detail/TTChannelSelectionHelper.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "celeritas/mucf/interactor/TTMucfInteractor.hh" + +namespace celeritas +{ +namespace detail +{ +//---------------------------------------------------------------------------// +/*! + * Select final channel for muonic dd molecules. + * + * This selection already accounts for sticking, as that is one of the possible + * outcomes. + */ +class TTChannelSelector +{ + public: + //!@{ + //! \name Type aliases + using Channel = TTMucfInteractor::Channel; + //!@} + + public: + //! Construct with args; \todo Update documentation + inline TTChannelSelector(/* args */); + + // Select fusion channel to be used by the interactor + template + inline CELER_FUNCTION Channel operator()(Engine&); +}; + +//---------------------------------------------------------------------------// +// INLINE DEFINITIONS +//---------------------------------------------------------------------------// +/*! + * Construct with args. + * + * \todo Update documentation + */ +TTChannelSelector::TTChannelSelector(/* args */) +{ + //! \todo Implement + CELER_NOT_IMPLEMENTED("Mucf tt fusion channel selection"); +} + +//---------------------------------------------------------------------------// +/*! + * Return a sampled channel to be used as input in the tt muCF interactor. + * + * \sa celeritas::TTMucfInteractor + */ +template +inline CELER_FUNCTION TTChannelSelector::Channel +TTChannelSelector::operator()(Engine&) +{ + Channel result{Channel::size_}; + + //! \todo Implement + // Final channel selection already takes into account sticking. + + CELER_ENSURE(result < Channel::size_); + return result; +} + +//---------------------------------------------------------------------------// +} // namespace detail +} // namespace celeritas diff --git a/src/celeritas/mucf/interactor/DDMucfInteractor.hh b/src/celeritas/mucf/interactor/DDMucfInteractor.hh new file mode 100644 index 0000000000..48aae54812 --- /dev/null +++ b/src/celeritas/mucf/interactor/DDMucfInteractor.hh @@ -0,0 +1,152 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/interactor/DDMucfInteractor.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "corecel/Macros.hh" +#include "corecel/data/StackAllocator.hh" +#include "celeritas/mat/ElementView.hh" +#include "celeritas/mat/MaterialView.hh" +#include "celeritas/mucf/data/DTMixMucfData.hh" +#include "celeritas/phys/Interaction.hh" +#include "celeritas/phys/ParticleTrackView.hh" +#include "celeritas/phys/Secondary.hh" + +namespace celeritas +{ +//---------------------------------------------------------------------------// +/*! + * Muon-catalyzed fusion of \f$ (dd)_\mu \f$ molecules. + * + * This is an \em at-rest interaction. + */ +class DDMucfInteractor +{ + public: + //! \todo Implement muonichydrogen3_proton (\f$ (^3\text{H})_\mu + p \f$) + enum class Channel + { + helium3_muon_neutron, //!< \f$ ^3\text{He} + \mu + n \f$ + muonichelium3_neutron, //!< \f$ (^3\text{He})_\mu + n \f$ + hydrogen3_muon_proton, //!< \f$ ^3\text{H} + \mu + p \f$ + size_ + }; + + // Construct from shared and state data + inline CELER_FUNCTION + DDMucfInteractor(NativeCRef const& data, + Channel const channel, + StackAllocator& allocate); + + // Sample an interaction with the given RNG + template + inline CELER_FUNCTION Interaction operator()(Engine& rng); + + private: + // Shared constant physics properties + NativeCRef const& data_; + // Selected fusion channel + Channel channel_{Channel::size_}; + // Allocate space for secondary particles + StackAllocator& allocate_; + + // Get number of secondaries for each channel + inline CELER_FUNCTION size_type num_secondaries() const; + + // Sample Interaction secondaries + template + inline CELER_FUNCTION Span + sample_secondaries(Secondary* secondaries /*, other args */, Engine&); +}; + +//---------------------------------------------------------------------------// +// INLINE DEFINITIONS +//---------------------------------------------------------------------------// +/*! + * Construct with shared and state data. + */ +DDMucfInteractor::DDMucfInteractor(NativeCRef const& data, + Channel const channel, + StackAllocator& allocate) + : data_(data), channel_(channel), allocate_(allocate) +{ + CELER_EXPECT(data_); + CELER_EXPECT(channel < Channel::size_); +} + +//---------------------------------------------------------------------------// +/*! + * Sample a dd muonic molecule fusion. + */ +template +CELER_FUNCTION Interaction DDMucfInteractor::operator()(Engine& rng) +{ + // Allocate space for the final fusion channel + Secondary* secondaries = allocate_(this->num_secondaries()); + if (secondaries == nullptr) + { + // Failed to allocate space for secondaries + return Interaction::from_failure(); + } + + // Kill primary and generate secondaries + Interaction result = Interaction::from_absorption(); + result.secondaries = this->sample_secondaries(secondaries, rng); + + return result; +} + +//---------------------------------------------------------------------------// +/*! + * Return number of secondaries from each fusion channel. + */ +CELER_FUNCTION size_type DDMucfInteractor::num_secondaries() const +{ + switch (channel_) + { + case Channel::helium3_muon_neutron: + return 3; + case Channel::muonichelium3_neutron: + return 2; + case Channel::hydrogen3_muon_proton: + return 3; + default: + CELER_ASSERT_UNREACHABLE(); + } +} + +//---------------------------------------------------------------------------// +/*! + * Sample the secondaries of the selected channel. + * + * Since secondaries come from an at rest interaction, their final state is + * a simple combination of random direction + momentum conservation + */ +template +CELER_FUNCTION Span +DDMucfInteractor::sample_secondaries(Secondary* secondaries /*, other args */, + Engine&) +{ + switch (channel_) + { + case Channel::helium3_muon_neutron: + //! \todo Assign secondaries + break; + case Channel::muonichelium3_neutron: + //! \todo Assign secondaries + break; + case Channel::hydrogen3_muon_proton: + //! \todo Assign secondaries + break; + default: + CELER_ASSERT_UNREACHABLE(); + } + + return Span{secondaries, this->num_secondaries()}; +} + +//---------------------------------------------------------------------------// +} // namespace celeritas diff --git a/src/celeritas/mucf/interactor/DTMucfInteractor.hh b/src/celeritas/mucf/interactor/DTMucfInteractor.hh new file mode 100644 index 0000000000..5af9a5f3f8 --- /dev/null +++ b/src/celeritas/mucf/interactor/DTMucfInteractor.hh @@ -0,0 +1,145 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/interactor/DTMucfInteractor.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "corecel/Macros.hh" +#include "corecel/data/StackAllocator.hh" +#include "celeritas/mat/ElementView.hh" +#include "celeritas/mat/MaterialView.hh" +#include "celeritas/mucf/data/DTMixMucfData.hh" +#include "celeritas/phys/Interaction.hh" +#include "celeritas/phys/ParticleTrackView.hh" +#include "celeritas/phys/Secondary.hh" + +namespace celeritas +{ +//---------------------------------------------------------------------------// +/*! + * Muon-catalyzed fusion of \f$ (dt)_\mu \f$ molecules. + * + * This is an \em at-rest interaction. + */ +class DTMucfInteractor +{ + public: + enum class Channel + { + alpha_muon_neutron, //!< \f$ \alpha + \mu + n \f$ + muonicalpha_neutron, //!< \f$ (\alpha)_\mu + n \f$ + size_ + }; + + // Construct from shared and state data + inline CELER_FUNCTION + DTMucfInteractor(NativeCRef const& data, + Channel const channel, + StackAllocator& allocate); + + // Sample an interaction with the given RNG + template + inline CELER_FUNCTION Interaction operator()(Engine& rng); + + private: + // Shared constant physics properties + NativeCRef const& data_; + // Selected fusion channel + Channel channel_{Channel::size_}; + // Allocate space for secondary particles + StackAllocator& allocate_; + + // Get number of secondaries for each channel + inline CELER_FUNCTION size_type num_secondaries() const; + + // Sample Interaction secondaries + template + inline CELER_FUNCTION Span + sample_secondaries(Secondary* secondaries /*, other args */, Engine&); +}; + +//---------------------------------------------------------------------------// +// INLINE DEFINITIONS +//---------------------------------------------------------------------------// +/*! + * Construct with shared and state data. + */ +DTMucfInteractor::DTMucfInteractor(NativeCRef const& data, + Channel const channel, + StackAllocator& allocate) + : data_(data), channel_(channel), allocate_(allocate) +{ + CELER_EXPECT(data_); + CELER_EXPECT(channel_ < Channel::size_); +} + +//---------------------------------------------------------------------------// +/*! + * Sample a dt muonic molecule fusion. + */ +template +CELER_FUNCTION Interaction DTMucfInteractor::operator()(Engine& rng) +{ + // Allocate space for the final fusion channel + Secondary* secondaries = allocate_(this->num_secondaries()); + if (secondaries == nullptr) + { + // Failed to allocate space for secondaries + return Interaction::from_failure(); + } + + // Kill primary and generate secondaries + Interaction result = Interaction::from_absorption(); + result.secondaries = this->sample_secondaries(secondaries, rng); + + return result; +} + +//---------------------------------------------------------------------------// +/*! + * Return number of secondaries from each fusion channel. + */ +CELER_FUNCTION size_type DTMucfInteractor::num_secondaries() const +{ + switch (channel_) + { + case Channel::alpha_muon_neutron: + return 3; + case Channel::muonicalpha_neutron: + return 2; + default: + CELER_ASSERT_UNREACHABLE(); + } +} + +//---------------------------------------------------------------------------// +/*! + * Sample the secondaries of the selected channel. + * + * Since secondaries come from an at rest interaction, their final state is + * a simple combination of random direction + momentum conservation + */ +template +CELER_FUNCTION Span +DTMucfInteractor::sample_secondaries(Secondary* secondaries /*, other args */, + Engine&) +{ + switch (channel_) + { + case Channel::alpha_muon_neutron: + //! \todo Assign secondaries + break; + case Channel::muonicalpha_neutron: + //! \todo Assign secondaries + break; + default: + CELER_ASSERT_UNREACHABLE(); + } + + return Span{secondaries, this->num_secondaries()}; +} + +//---------------------------------------------------------------------------// +} // namespace celeritas diff --git a/src/celeritas/mucf/interactor/TTMucfInteractor.hh b/src/celeritas/mucf/interactor/TTMucfInteractor.hh new file mode 100644 index 0000000000..97ade2512f --- /dev/null +++ b/src/celeritas/mucf/interactor/TTMucfInteractor.hh @@ -0,0 +1,145 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/interactor/TTMucfInteractor.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "corecel/Macros.hh" +#include "corecel/data/StackAllocator.hh" +#include "celeritas/mat/ElementView.hh" +#include "celeritas/mat/MaterialView.hh" +#include "celeritas/mucf/data/DTMixMucfData.hh" +#include "celeritas/phys/Interaction.hh" +#include "celeritas/phys/ParticleTrackView.hh" +#include "celeritas/phys/Secondary.hh" + +namespace celeritas +{ +//---------------------------------------------------------------------------// +/*! + * Muon-catalyzed fusion of \f$ (tt)_\mu \f$ molecules. + * + * This is an \em at-rest interaction. + */ +class TTMucfInteractor +{ + public: + enum class Channel + { + alpha_muon_neutron_neutron, //!< \f$ \alpha + \mu + n + n \f$ + muonicalpha_neutron_neutron, //!< \f$ (\alpha)_\mu + n + n \f$ + size_ + }; + + // Construct from shared and state data + inline CELER_FUNCTION + TTMucfInteractor(NativeCRef const& data, + Channel const channel, + StackAllocator& allocate); + + // Sample an interaction with the given RNG + template + inline CELER_FUNCTION Interaction operator()(Engine& rng); + + private: + // Shared constant physics properties + NativeCRef const& data_; + // Selected fusion channel + Channel channel_{Channel::size_}; + // Allocate space for secondary particles + StackAllocator& allocate_; + + // Get number of secondaries for each channel + inline CELER_FUNCTION size_type num_secondaries() const; + + // Sample Interaction secondaries + template + inline CELER_FUNCTION Span + sample_secondaries(Secondary* secondaries /*, other args */, Engine&); +}; + +//---------------------------------------------------------------------------// +// INLINE DEFINITIONS +//---------------------------------------------------------------------------// +/*! + * Construct with shared and state data. + */ +TTMucfInteractor::TTMucfInteractor(NativeCRef const& data, + Channel const channel, + StackAllocator& allocate) + : data_(data), channel_(channel), allocate_(allocate) +{ + CELER_EXPECT(data_); + CELER_EXPECT(channel_ < Channel::size_); +} + +//---------------------------------------------------------------------------// +/*! + * Sample a dt muonic molecule fusion. + */ +template +CELER_FUNCTION Interaction TTMucfInteractor::operator()(Engine& rng) +{ + // Allocate space for the final fusion channel + Secondary* secondaries = allocate_(this->num_secondaries()); + if (secondaries == nullptr) + { + // Failed to allocate space for secondaries + return Interaction::from_failure(); + } + + // Kill primary and generate secondaries + Interaction result = Interaction::from_absorption(); + result.secondaries = this->sample_secondaries(secondaries, rng); + + return result; +} + +//---------------------------------------------------------------------------// +/*! + * Return number of secondaries from each fusion channel. + */ +CELER_FUNCTION size_type TTMucfInteractor::num_secondaries() const +{ + switch (channel_) + { + case Channel::alpha_muon_neutron_neutron: + return 4; + case Channel::muonicalpha_neutron_neutron: + return 3; + default: + CELER_ASSERT_UNREACHABLE(); + } +} + +//---------------------------------------------------------------------------// +/*! + * Sample the secondaries of the selected channel. + * + * Since secondaries come from an at rest interaction, their final state is + * a simple combination of random direction + momentum conservation + */ +template +CELER_FUNCTION Span +TTMucfInteractor::sample_secondaries(Secondary* secondaries /*, other args */, + Engine&) +{ + switch (channel_) + { + case Channel::alpha_muon_neutron_neutron: + //! \todo Assign secondaries + break; + case Channel::muonicalpha_neutron_neutron: + //! \todo Assign secondaries + break; + default: + CELER_ASSERT_UNREACHABLE(); + } + + return Span{secondaries, this->num_secondaries()}; +} + +//---------------------------------------------------------------------------// +} // namespace celeritas diff --git a/src/celeritas/phys/PDGNumber.hh b/src/celeritas/phys/PDGNumber.hh index 44755468c1..15baf0edbf 100644 --- a/src/celeritas/phys/PDGNumber.hh +++ b/src/celeritas/phys/PDGNumber.hh @@ -146,6 +146,10 @@ CELER_DEFINE_PDGNUMBER(anti_he3, -1000020030) CELER_DEFINE_PDGNUMBER(alpha, 1000020040) CELER_DEFINE_PDGNUMBER(anti_alpha, -1000020040) +// Muonic nuclei +CELER_DEFINE_PDGNUMBER(muonic_deuteron, 2000010020) +CELER_DEFINE_PDGNUMBER(muonic_triton, 2000010030) + //!@} #undef CELER_DEFINE_PDGNUMBER From 2f874685b36d0c69ec2d03ed62f8bc37355e2e78 Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Mon, 29 Dec 2025 13:01:16 -0500 Subject: [PATCH 06/26] Make sure it builds --- src/celeritas/CMakeLists.txt | 3 +++ src/celeritas/mucf/process/MucfProcess.hh | 1 + 2 files changed, 4 insertions(+) diff --git a/src/celeritas/CMakeLists.txt b/src/celeritas/CMakeLists.txt index 8d7cc5269d..c4a2365e3f 100644 --- a/src/celeritas/CMakeLists.txt +++ b/src/celeritas/CMakeLists.txt @@ -100,6 +100,8 @@ list(APPEND SOURCES mat/MaterialParams.cc mat/MaterialParamsOutput.cc mat/detail/Utils.cc + mucf/model/detail/DTMixMaterialCalculator.cc + mucf/process/MucfProcess.cc neutron/model/CascadeOptions.cc neutron/model/CascadeOptionsIO.json.cc neutron/process/NeutronElasticProcess.cc @@ -375,6 +377,7 @@ celeritas_polysource(em/model/CoulombScatteringModel) celeritas_polysource(geo/detail/BoundaryAction) celeritas_polysource(global/detail/KillActive) celeritas_polysource(global/detail/TrackSlotUtils) +celeritas_polysource(mucf/model/DTMixMucfModel) celeritas_polysource(neutron/model/ChipsNeutronElasticModel) celeritas_polysource(neutron/model/NeutronInelasticModel) celeritas_polysource(optical/model/AbsorptionModel) diff --git a/src/celeritas/mucf/process/MucfProcess.hh b/src/celeritas/mucf/process/MucfProcess.hh index bb14033d79..146ec3d2c7 100644 --- a/src/celeritas/mucf/process/MucfProcess.hh +++ b/src/celeritas/mucf/process/MucfProcess.hh @@ -8,6 +8,7 @@ #include "celeritas/Types.hh" #include "celeritas/mat/MaterialParams.hh" +#include "celeritas/phys/ParticleParams.hh" #include "celeritas/phys/Process.hh" namespace celeritas From f517bbb440f300912853c831e4479b743bd91ca5 Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Mon, 29 Dec 2025 13:01:23 -0500 Subject: [PATCH 07/26] Expand documentation --- doc/implementation/mucf-physics.rst | 27 +++++++++++++++++++++++++++ 1 file changed, 27 insertions(+) diff --git a/doc/implementation/mucf-physics.rst b/doc/implementation/mucf-physics.rst index ec6d81f84e..c2b1823c7e 100644 --- a/doc/implementation/mucf-physics.rst +++ b/doc/implementation/mucf-physics.rst @@ -114,3 +114,30 @@ fusion data is constructed when the ``G4MuonMinusAtomicCapture`` process is registered. .. todo:: Add process/model/executor details + +Code implementation +=================== + +The ``MucfProcess`` process only has the model ``DTMixMucfModel`` attached to +it, responsible for deuterium-tritium mixtures. It can simulate materials from +near absolute zero to 1500 kelvin, and it is an *at rest* model that encompasses +the full cycle---atom formation, molecule formation, and fusion. + +.. note:: Only reactive channels are implemented. + +.. doxygenclass:: celeritas::MucfProcess + +The main cycle is managed by the ``DTMixMucfExecutor``, with the +Interactors reserved for sampling final states of the outgoing secondaries. + +.. doxygenclass:: celeritas::DTMixMucfModel +.. doxygenclass:: celeritas::DTMixMucfExecutor + +Most of the data is material-dependent, and thus can be calculated and cached +during model construction. This is done by the ``DTMixMaterialCalculator``. + +.. doxygenclass:: celeritas::detail::DTMixMaterialCalculator +.. doxygenclass:: celeritas::DTMixMucfExecutor +.. doxygenclass:: celeritas::DDMucfInteractor +.. doxygenclass:: celeritas::DTMucfInteractor +.. doxygenclass:: celeritas::TTMucfInteractor From dc206341daaade4e05cd24bae3ddcf7630a6d71b Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Mon, 29 Dec 2025 13:17:07 -0500 Subject: [PATCH 08/26] Fixup --- src/celeritas/mucf/process/MucfProcess.cc | 2 +- src/celeritas/mucf/process/MucfProcess.hh | 4 ++-- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/src/celeritas/mucf/process/MucfProcess.cc b/src/celeritas/mucf/process/MucfProcess.cc index 27a9256e19..cde972edfb 100644 --- a/src/celeritas/mucf/process/MucfProcess.cc +++ b/src/celeritas/mucf/process/MucfProcess.cc @@ -40,7 +40,7 @@ auto MucfProcess::build_models(ActionIdIter start_id) const -> VecModel /*! * Get the interaction cross sections for the given energy range. */ -auto MucfProcess::macro_xs(Applicability applic) const -> XsGrid +auto MucfProcess::macro_xs(Applicability) const -> XsGrid { return {}; } diff --git a/src/celeritas/mucf/process/MucfProcess.hh b/src/celeritas/mucf/process/MucfProcess.hh index 146ec3d2c7..565b4c590c 100644 --- a/src/celeritas/mucf/process/MucfProcess.hh +++ b/src/celeritas/mucf/process/MucfProcess.hh @@ -37,10 +37,10 @@ class MucfProcess final : public Process VecModel build_models(ActionIdIter start_id) const final; // Get the interaction cross sections for the given energy range - XsGrid macro_xs(Applicability range) const final; + XsGrid macro_xs(Applicability) const final; // Get the energy loss for the given energy range - EnergyLossGrid energy_loss(Applicability range) const final; + EnergyLossGrid energy_loss(Applicability) const final; //! Whether the integral method can be used to sample interaction length bool supports_integral_xs() const final { return true; } //! \todo Check From 91fdbbf3762cfc92da36ade39bf957f18725a4ea Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Fri, 2 Jan 2026 15:28:11 -0500 Subject: [PATCH 09/26] Address most comments: - Move types to new mucf/Types.hh - Move MucfParticles::from_params to model.cc as a free function - Move all static data to MucfPhysics.cc - Refactor/Remove types in DTMixMucfData - Remove CELER_EXPECT(mu_minus) during MucfProcess construction --- src/celeritas/Types.hh | 28 ----- src/celeritas/inp/MucfPhysics.cc | 99 ++++++++++++++- src/celeritas/inp/MucfPhysics.hh | 2 +- src/celeritas/inp/MucfPhysicsData.hh | 117 ------------------ src/celeritas/mucf/Types.hh | 48 +++++++ src/celeritas/mucf/data/DTMixMucfData.hh | 53 +------- .../mucf/executor/detail/DTChannelSelector.hh | 2 +- .../mucf/executor/detail/TTChannelSelector.hh | 2 +- .../mucf/interactor/DDMucfInteractor.hh | 2 +- .../mucf/interactor/DTMucfInteractor.hh | 2 +- .../mucf/interactor/TTMucfInteractor.hh | 2 +- src/celeritas/mucf/model/DTMixMucfModel.cc | 55 +++++++- src/celeritas/mucf/process/MucfProcess.cc | 1 - 13 files changed, 209 insertions(+), 204 deletions(-) delete mode 100644 src/celeritas/inp/MucfPhysicsData.hh create mode 100644 src/celeritas/mucf/Types.hh diff --git a/src/celeritas/Types.hh b/src/celeritas/Types.hh index a2d39abf08..c2cd0e9e63 100644 --- a/src/celeritas/Types.hh +++ b/src/celeritas/Types.hh @@ -98,9 +98,6 @@ using UniformGridId = OpaqueId; //! Opaque index of a cross section grid using XsGridId = OpaqueId; -//! Opaque index of a muCF material component -using MuCfMatCompId = OpaqueId; - //---------------------------------------------------------------------------// // ENUMERATIONS //---------------------------------------------------------------------------// @@ -230,31 +227,6 @@ enum class CylAxis size_ }; -//---------------------------------------------------------------------------// -/*! - * Muonic atom selection from material data. This is *not* intended to be used - * by the transport loop. - */ -enum class MucfMuonicAtom -{ - deuterium, - tritium, - size_ -}; - -//---------------------------------------------------------------------------// -/*! - * Muonic molecule selection from material data. This is *not* intended to be - * used by the transport loop. - */ -enum class MucfMuonicMolecule -{ - deuterium_deuterium, - deuterium_tritium, - tritium_tritium, - size_ -}; - //---------------------------------------------------------------------------// // HELPER STRUCTS //---------------------------------------------------------------------------// diff --git a/src/celeritas/inp/MucfPhysics.cc b/src/celeritas/inp/MucfPhysics.cc index 472fb54590..32a7e37ff1 100644 --- a/src/celeritas/inp/MucfPhysics.cc +++ b/src/celeritas/inp/MucfPhysics.cc @@ -6,12 +6,107 @@ //---------------------------------------------------------------------------// #include "MucfPhysics.hh" -#include "MucfPhysicsData.hh" - namespace celeritas { namespace inp { +//---------------------------------------------------------------------------// +/*! + * Static muon energy CDF data for muon-catalyzed fusion. + * + * \todo This data may be loaded from an external file in the future. + */ +static Grid mucf_muon_energy_cdf() +{ + Grid cdf; + cdf.interpolation.type = InterpolationType::cubic_spline; + + // Cumulative distribution data [unitless] + cdf.x = {0, + 0.04169381107491854, + 0.08664495114006499, + 0.14332247557003264, + 0.20456026058631915, + 0.2723127035830618, + 0.34136807817589576, + 0.41563517915309456, + 0.48990228013029324, + 0.5667752442996744, + 0.6306188925081434, + 0.6866449511400652, + 0.7309446254071662, + 0.7778501628664496, + 0.8104234527687297, + 0.8403908794788275, + 0.8618892508143323, + 0.8814332247557004, + 0.8970684039087949, + 0.903583061889251, + 1.0}; + + // Energy [keV] + cdf.y = {0, + 0.48850540675768084, + 0.8390389347819425, + 1.2521213482687141, + 1.7153033196164724, + 2.253638712180777, + 2.854653691809707, + 3.606073540073316, + 4.470346052913727, + 5.560291219507215, + 6.700556502915258, + 7.953772477101693, + 9.194596305637525, + 10.849180562221111, + 12.353474314071864, + 14.045888515617822, + 15.650634617544647, + 17.38079707555165, + 19.111008546659452, + 19.976130619913615, + 80.0}; + + CELER_ENSURE(cdf); + return cdf; +} + +//---------------------------------------------------------------------------// +/*! + * Static cycle rate data for muon-catalyzed fusion. + * + * \todo This data may be loaded from an external file in the future. + */ +static std::vector mucf_cycle_rates() +{ + //! \todo Implement + return {}; +}; + +//---------------------------------------------------------------------------// +/*! + * Static spin-flip data for muon-catalyzed fusion. + * + * \todo This data may be loaded from an external file in the future. + */ +static std::vector mucf_atom_spin_flip_rates() +{ + //! \todo Implement + return {}; +}; + +//---------------------------------------------------------------------------// +/*! + * Static atom transfer data for muon-catalyzed fusion. + * + * \todo This data may be loaded from an external file in the future. + */ +static std::vector mucf_atom_transfer_rates() +{ + //! \todo Implement + return {}; +}; + //---------------------------------------------------------------------------// /*! * Construct hardcoded muon-catalyzed fusion physics data. diff --git a/src/celeritas/inp/MucfPhysics.hh b/src/celeritas/inp/MucfPhysics.hh index 7854b85e0a..c55f13fd95 100644 --- a/src/celeritas/inp/MucfPhysics.hh +++ b/src/celeritas/inp/MucfPhysics.hh @@ -10,7 +10,7 @@ #include #include "corecel/inp/Grid.hh" -#include "celeritas/Types.hh" +#include "celeritas/mucf/Types.hh" namespace celeritas { diff --git a/src/celeritas/inp/MucfPhysicsData.hh b/src/celeritas/inp/MucfPhysicsData.hh deleted file mode 100644 index 20c023dc07..0000000000 --- a/src/celeritas/inp/MucfPhysicsData.hh +++ /dev/null @@ -1,117 +0,0 @@ -//------------------------------- -*- C++ -*- -------------------------------// -// Copyright Celeritas contributors: see top-level COPYRIGHT file for details -// SPDX-License-Identifier: (Apache-2.0 OR MIT) -//---------------------------------------------------------------------------// -//! \file celeritas/inp/MucfPhysicsData.hh -//---------------------------------------------------------------------------// -#pragma once - -#include "corecel/Macros.hh" -#include "corecel/inp/Grid.hh" - -#include "MucfPhysics.hh" - -namespace celeritas -{ -namespace inp -{ -//---------------------------------------------------------------------------// -/*! - * Static muon energy CDF data for muon-catalyzed fusion. - * - * \todo This data may be loaded from an external file in the future. - */ -static Grid mucf_muon_energy_cdf() -{ - Grid cdf; - cdf.interpolation.type = InterpolationType::cubic_spline; - - // Cumulative distribution data [unitless] - cdf.x = {0, - 0.04169381107491854, - 0.08664495114006499, - 0.14332247557003264, - 0.20456026058631915, - 0.2723127035830618, - 0.34136807817589576, - 0.41563517915309456, - 0.48990228013029324, - 0.5667752442996744, - 0.6306188925081434, - 0.6866449511400652, - 0.7309446254071662, - 0.7778501628664496, - 0.8104234527687297, - 0.8403908794788275, - 0.8618892508143323, - 0.8814332247557004, - 0.8970684039087949, - 0.903583061889251, - 1.0}; - - // Energy [keV] - cdf.y = {0, - 0.48850540675768084, - 0.8390389347819425, - 1.2521213482687141, - 1.7153033196164724, - 2.253638712180777, - 2.854653691809707, - 3.606073540073316, - 4.470346052913727, - 5.560291219507215, - 6.700556502915258, - 7.953772477101693, - 9.194596305637525, - 10.849180562221111, - 12.353474314071864, - 14.045888515617822, - 15.650634617544647, - 17.38079707555165, - 19.111008546659452, - 19.976130619913615, - 80.0}; - - CELER_ENSURE(cdf); - return cdf; -} - -//---------------------------------------------------------------------------// -/*! - * Static cycle rate data for muon-catalyzed fusion. - * - * \todo This data may be loaded from an external file in the future. - */ -static std::vector mucf_cycle_rates() -{ - //! \todo Implement - return {}; -}; - -//---------------------------------------------------------------------------// -/*! - * Static spin-flip data for muon-catalyzed fusion. - * - * \todo This data may be loaded from an external file in the future. - */ -static std::vector mucf_atom_spin_flip_rates() -{ - //! \todo Implement - return {}; -}; - -//---------------------------------------------------------------------------// -/*! - * Static atom transfer data for muon-catalyzed fusion. - * - * \todo This data may be loaded from an external file in the future. - */ -static std::vector mucf_atom_transfer_rates() -{ - //! \todo Implement - return {}; -}; - -//---------------------------------------------------------------------------// -} // namespace inp -} // namespace celeritas diff --git a/src/celeritas/mucf/Types.hh b/src/celeritas/mucf/Types.hh new file mode 100644 index 0000000000..ed64826410 --- /dev/null +++ b/src/celeritas/mucf/Types.hh @@ -0,0 +1,48 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/Types.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "corecel/OpaqueId.hh" + +namespace celeritas +{ +//---------------------------------------------------------------------------// +// ENUMERATIONS +//---------------------------------------------------------------------------// +/*! + * Muonic atom selection from material data. This is *not* intended to be used + * by the transport loop. + */ +enum class MucfMuonicAtom +{ + deuterium, + tritium, + size_ +}; + +//---------------------------------------------------------------------------// +/*! + * Muonic molecule selection from material data. This is *not* intended to be + * used by the transport loop. + */ +enum class MucfMuonicMolecule +{ + deuterium_deuterium, + deuterium_tritium, + tritium_tritium, + size_ +}; + +//---------------------------------------------------------------------------// +// TYPE ALIASES +//---------------------------------------------------------------------------// + +//! Opaque index of a muCF material component +using MuCfMatCompId = OpaqueId; + +//---------------------------------------------------------------------------// +} // namespace celeritas diff --git a/src/celeritas/mucf/data/DTMixMucfData.hh b/src/celeritas/mucf/data/DTMixMucfData.hh index 43a57c2508..4074765d79 100644 --- a/src/celeritas/mucf/data/DTMixMucfData.hh +++ b/src/celeritas/mucf/data/DTMixMucfData.hh @@ -14,6 +14,7 @@ #include "corecel/grid/NonuniformGridData.hh" #include "corecel/io/Join.hh" #include "celeritas/Types.hh" +#include "celeritas/mucf/Types.hh" #include "celeritas/phys/ParticleParams.hh" namespace celeritas @@ -22,7 +23,7 @@ namespace celeritas /*! * ParticleIds used by the \c DTMixMucfModel . */ -struct MucfParticles +struct MucfParticleIds { //! Primary ParticleId mu_minus; @@ -51,46 +52,6 @@ struct MucfParticles && muonic_hydrogen && muonic_alpha && muonic_triton && muonic_helium3; } - - //! Assign from ParticleParams (anonymous free function in the model.cc?) - static MucfParticles from_params(ParticleParams const& particles) - { -#define MP_ADD(MEMBER) \ - pids.MEMBER = particles.find(pdg::MEMBER()); \ - if (!pids.MEMBER) \ - { \ - missing.push_back({#MEMBER, pdg::MEMBER()}); \ - } - - using PairStrPdg = std::pair; - std::vector missing; - MucfParticles pids; - - MP_ADD(mu_minus); - MP_ADD(neutron); - MP_ADD(proton); - MP_ADD(alpha); - //! \todo Decide whether to implement these PDGs in PDGNumber.hh -#if 0 - MP_ADD(helium_3); - MP_ADD(muonic_hydrogen); - MP_ADD(muonic_alpha); - MP_ADD(muonic_deuteron); - MP_ADD(muonic_triton); - MP_ADD(muonic_helium3); -#endif - CELER_VALIDATE( - missing.empty(), - << "missing particles required for muon-catalyzed fusion: " - << join(missing.begin(), missing.end(), ", ", [](PairStrPdg const& p) { - return p.first + " (PDG " - + std::to_string(p.second.unchecked_get()) + ')'; - })); - - return pids; - -#undef MP_ADD - } }; //---------------------------------------------------------------------------// @@ -103,17 +64,15 @@ struct DTMixMucfData template using Items = Collection; template - using MatCycleItems = Collection; - template using MaterialItems = Collection; - using Grid = NonuniformGridRecord; + using GridRecord = NonuniformGridRecord; using CycleTimesArray = EnumArray>; //! Particle IDs - MucfParticles particles; + MucfParticleIds particles; //! Muon CDF energy grid for sampling outgoing muCF muons - Grid muon_energy_cdf; //! \todo Verify energy unit + GridRecord muon_energy_cdf; //! \todo Verify energy unit Items reals; //!@{ @@ -121,7 +80,7 @@ struct DTMixMucfData //! \c PhysMatId indexed by \c MuCfMatCompId MaterialItems matcompid_to_matid; //! Cycle times per material: [mat_comp_id][muonic_molecule][spin_index] - MatCycleItems cycle_times; //!< In [s] + MaterialItems cycle_times; //!< In [s] //! \todo Add mean atom spin flip times //! \todo Add mean atom transfer times //!@} diff --git a/src/celeritas/mucf/executor/detail/DTChannelSelector.hh b/src/celeritas/mucf/executor/detail/DTChannelSelector.hh index 28483dea53..418571083a 100644 --- a/src/celeritas/mucf/executor/detail/DTChannelSelector.hh +++ b/src/celeritas/mucf/executor/detail/DTChannelSelector.hh @@ -14,7 +14,7 @@ namespace detail { //---------------------------------------------------------------------------// /*! - * Select final channel for muonic dd molecules. + * Select final channel for muonic dt molecules. * * This selection already accounts for sticking, as that is one of the possible * outcomes. diff --git a/src/celeritas/mucf/executor/detail/TTChannelSelector.hh b/src/celeritas/mucf/executor/detail/TTChannelSelector.hh index 1ded52adb7..d69bab7bfd 100644 --- a/src/celeritas/mucf/executor/detail/TTChannelSelector.hh +++ b/src/celeritas/mucf/executor/detail/TTChannelSelector.hh @@ -14,7 +14,7 @@ namespace detail { //---------------------------------------------------------------------------// /*! - * Select final channel for muonic dd molecules. + * Select final channel for muonic tt molecules. * * This selection already accounts for sticking, as that is one of the possible * outcomes. diff --git a/src/celeritas/mucf/interactor/DDMucfInteractor.hh b/src/celeritas/mucf/interactor/DDMucfInteractor.hh index 48aae54812..a659758f10 100644 --- a/src/celeritas/mucf/interactor/DDMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/DDMucfInteractor.hh @@ -38,7 +38,7 @@ class DDMucfInteractor // Construct from shared and state data inline CELER_FUNCTION DDMucfInteractor(NativeCRef const& data, - Channel const channel, + Channel channel, StackAllocator& allocate); // Sample an interaction with the given RNG diff --git a/src/celeritas/mucf/interactor/DTMucfInteractor.hh b/src/celeritas/mucf/interactor/DTMucfInteractor.hh index 5af9a5f3f8..94498de7f0 100644 --- a/src/celeritas/mucf/interactor/DTMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/DTMucfInteractor.hh @@ -36,7 +36,7 @@ class DTMucfInteractor // Construct from shared and state data inline CELER_FUNCTION DTMucfInteractor(NativeCRef const& data, - Channel const channel, + Channel channel, StackAllocator& allocate); // Sample an interaction with the given RNG diff --git a/src/celeritas/mucf/interactor/TTMucfInteractor.hh b/src/celeritas/mucf/interactor/TTMucfInteractor.hh index 97ade2512f..c37cda8733 100644 --- a/src/celeritas/mucf/interactor/TTMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/TTMucfInteractor.hh @@ -36,7 +36,7 @@ class TTMucfInteractor // Construct from shared and state data inline CELER_FUNCTION TTMucfInteractor(NativeCRef const& data, - Channel const channel, + Channel channel, StackAllocator& allocate); // Sample an interaction with the given RNG diff --git a/src/celeritas/mucf/model/DTMixMucfModel.cc b/src/celeritas/mucf/model/DTMixMucfModel.cc index ac4c4e9e56..b21e045710 100644 --- a/src/celeritas/mucf/model/DTMixMucfModel.cc +++ b/src/celeritas/mucf/model/DTMixMucfModel.cc @@ -8,11 +8,11 @@ #include +#include "corecel/OpaqueIdIO.hh" #include "corecel/inp/Grid.hh" #include "celeritas/global/ActionLauncher.hh" #include "celeritas/global/TrackExecutor.hh" #include "celeritas/grid/NonuniformGridBuilder.hh" -#include "celeritas/inp/MucfPhysicsData.hh" #include "celeritas/mucf/executor/DTMixMucfExecutor.hh" // IWYU pragma: associated #include "celeritas/phys/InteractionApplier.hh" // IWYU pragma: associated #include "celeritas/phys/PDGNumber.hh" @@ -21,6 +21,56 @@ namespace celeritas { +namespace +{ +//---------------------------------------------------------------------------// +/*! + * Assign particle IDs from \c ParticleParams . + */ +static MucfParticleIds from_params(ParticleParams const& particles) +{ + using PairStrPdg = std::pair; + std::vector missing; + MucfParticleIds result; + +#define MP_ADD(MEMBER) \ + result.MEMBER = particles.find(pdg::MEMBER()); \ + if (!result.MEMBER) \ + { \ + missing.push_back({#MEMBER, pdg::MEMBER()}); \ + } + + MP_ADD(mu_minus); + MP_ADD(neutron); + MP_ADD(proton); + MP_ADD(alpha); + + //! \todo Decide whether to implement these PDGs in PDGNumber.hh +#if 0 + MP_ADD(helium_3); + MP_ADD(muonic_hydrogen); + MP_ADD(muonic_deuteron); + MP_ADD(muonic_triton); + MP_ADD(muonic_alpha); + MP_ADD(muonic_helium3); +#endif + + CELER_VALIDATE(missing.empty(), + << "missing particles required for muon-catalyzed fusion: " + << join_stream(missing.begin(), + missing.end(), + ", ", + [](std::ostream& os, PairStrPdg const& p) { + os << p.first << " (PDG " + << p.second.unchecked_get() << ')'; + })); + return result; + +#undef MP_ADD +} +//---------------------------------------------------------------------------// +} // namespace + //---------------------------------------------------------------------------// /*! * Construct from model ID and other necessary data. @@ -57,9 +107,8 @@ DTMixMucfModel::DTMixMucfModel(ActionId id, inp::MucfPhysics inp_data = inp::MucfPhysics::from_default(); CELER_EXPECT(inp_data); - // Load particle IDs HostVal host_data; - host_data.particles = MucfParticles::from_params(particles); + host_data.particles = from_params(particles); // Copy muon energy CDF data using NonuniformGridBuilder NonuniformGridBuilder build_grid_record{&host_data.reals}; diff --git a/src/celeritas/mucf/process/MucfProcess.cc b/src/celeritas/mucf/process/MucfProcess.cc index cde972edfb..b3b3dc5d8c 100644 --- a/src/celeritas/mucf/process/MucfProcess.cc +++ b/src/celeritas/mucf/process/MucfProcess.cc @@ -23,7 +23,6 @@ MucfProcess::MucfProcess(SPConstParticles particles, SPConstMaterials materials) //! \todo Fix ImportProcessClass CELER_EXPECT(particles_); CELER_EXPECT(materials_); - CELER_EXPECT(particles_->find(pdg::mu_minus())); } //---------------------------------------------------------------------------// From 67be430ec278043cce795d55e47eb80f57cdb395 Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Mon, 5 Jan 2026 14:50:43 -0500 Subject: [PATCH 10/26] Move component id getter function into DTMixMucfData --- src/celeritas/mucf/data/DTMixMucfData.hh | 19 +++++++++++++++++++ .../mucf/executor/DTMixMucfExecutor.hh | 19 ++++--------------- 2 files changed, 23 insertions(+), 15 deletions(-) diff --git a/src/celeritas/mucf/data/DTMixMucfData.hh b/src/celeritas/mucf/data/DTMixMucfData.hh index 4074765d79..b53c0b1ab0 100644 --- a/src/celeritas/mucf/data/DTMixMucfData.hh +++ b/src/celeritas/mucf/data/DTMixMucfData.hh @@ -94,6 +94,25 @@ struct DTMixMucfData && (matcompid_to_matid.size() == cycle_times.size()); } + //! Get the material componend ID from a given \c PhysMatId . + CELER_FUNCTION MuCfMatCompId material_component_id(PhysMatId matid) const + { + CELER_EXPECT(matid); + CELER_EXPECT(this); + CELER_EXPECT(!this->matcompid_to_matid.empty()); + + for (auto i : range(this->matcompid_to_matid.size())) + { + if (auto const comp_id = MuCfMatCompId{i}; + this->matcompid_to_matid[comp_id] == matid) + { + return comp_id; + } + } + // Component ID not found + return MuCfMatCompId{}; + } + //! Assign from another set of data template DTMixMucfData& operator=(DTMixMucfData const& other) diff --git a/src/celeritas/mucf/executor/DTMixMucfExecutor.hh b/src/celeritas/mucf/executor/DTMixMucfExecutor.hh index 7abfc107c1..ce594cbbd5 100644 --- a/src/celeritas/mucf/executor/DTMixMucfExecutor.hh +++ b/src/celeritas/mucf/executor/DTMixMucfExecutor.hh @@ -74,22 +74,11 @@ DTMixMucfExecutor::operator()(celeritas::CoreTrackView const& track) auto const molecule_spin = select_molecule_spin(rng); // Load cycle time for the selected molecule - auto find_component = [&](PhysMatId matid) -> MuCfMatCompId { - for (auto i : range(data.matcompid_to_matid.size())) - { - auto const comp_id = MuCfMatCompId{i}; - if (data.matcompid_to_matid[comp_id] == matid) - { - return MuCfMatCompId{i}; - } - } - return MuCfMatCompId{}; - }; - - auto const mat_comp = find_component(track.material().material_id()); - CELER_ASSERT(mat_comp); + auto const mat_comp_id + = data.material_component_id(track.material().material_id()); + CELER_ASSERT(mat_comp_id); auto const cycle_time - = data.cycle_times[mat_comp][muonic_molecule][molecule_spin]; + = data.cycle_times[mat_comp_id][muonic_molecule][molecule_spin]; CELER_ASSERT(cycle_time > 0); // Check if muon decays before fusion happens From 0842322a8e1a34e7262db1038eb524526545b21b Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Mon, 5 Jan 2026 14:54:55 -0500 Subject: [PATCH 11/26] Assign all current PDG numbers --- src/celeritas/mucf/data/DTMixMucfData.hh | 9 ++++----- src/celeritas/mucf/model/DTMixMucfModel.cc | 8 ++++---- 2 files changed, 8 insertions(+), 9 deletions(-) diff --git a/src/celeritas/mucf/data/DTMixMucfData.hh b/src/celeritas/mucf/data/DTMixMucfData.hh index b53c0b1ab0..d6104ddea6 100644 --- a/src/celeritas/mucf/data/DTMixMucfData.hh +++ b/src/celeritas/mucf/data/DTMixMucfData.hh @@ -33,7 +33,7 @@ struct MucfParticleIds ParticleId neutron; ParticleId proton; ParticleId alpha; - ParticleId helium3; + ParticleId he3; //!@} //!@{ @@ -42,15 +42,14 @@ struct MucfParticleIds ParticleId muonic_deuteron; ParticleId muonic_triton; ParticleId muonic_alpha; - ParticleId muonic_helium3; + ParticleId muonic_he3; //!@} //! Check whether all particles are assigned CELER_FUNCTION explicit operator bool() const { - return mu_minus && neutron && proton && alpha && helium3 - && muonic_hydrogen && muonic_alpha && muonic_triton - && muonic_helium3; + return mu_minus && neutron && proton && alpha && he3 && muonic_hydrogen + && muonic_alpha && muonic_triton && muonic_he3; } }; diff --git a/src/celeritas/mucf/model/DTMixMucfModel.cc b/src/celeritas/mucf/model/DTMixMucfModel.cc index b21e045710..7328920b3a 100644 --- a/src/celeritas/mucf/model/DTMixMucfModel.cc +++ b/src/celeritas/mucf/model/DTMixMucfModel.cc @@ -44,15 +44,15 @@ static MucfParticleIds from_params(ParticleParams const& particles) MP_ADD(neutron); MP_ADD(proton); MP_ADD(alpha); + MP_ADD(he3); + MP_ADD(muonic_deuteron); + MP_ADD(muonic_triton); //! \todo Decide whether to implement these PDGs in PDGNumber.hh #if 0 - MP_ADD(helium_3); MP_ADD(muonic_hydrogen); - MP_ADD(muonic_deuteron); - MP_ADD(muonic_triton); MP_ADD(muonic_alpha); - MP_ADD(muonic_helium3); + MP_ADD(muonic_he3); #endif CELER_VALIDATE(missing.empty(), From 7123f3ed44531d1c6d29b72fe365cfb5af5bd1c2 Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Fri, 9 Jan 2026 14:25:46 -0500 Subject: [PATCH 12/26] Address comments: - Fix import data - Refactor the material calculator to behave as a material inserter - Interactors: define number of secondaries as EnumArrays - Rename internal material id to MucfMatId - Remove find MucfMatId function from DTMixMucfData and make it a lambda in the executor --- src/celeritas/CMakeLists.txt | 2 +- src/celeritas/inp/MucfPhysics.cc | 123 +++++++-------- src/celeritas/inp/MucfPhysics.hh | 50 ++++-- src/celeritas/mucf/Types.hh | 25 ++- src/celeritas/mucf/data/DTMixMucfData.hh | 34 +---- .../mucf/executor/DTMixMucfExecutor.hh | 22 ++- .../mucf/interactor/DDMucfInteractor.hh | 32 +--- .../mucf/interactor/DTMucfInteractor.hh | 29 +--- .../mucf/interactor/TTMucfInteractor.hh | 29 +--- src/celeritas/mucf/model/DTMixMucfModel.cc | 32 +--- .../model/detail/DTMixMaterialCalculator.hh | 142 ------------------ ...lCalculator.cc => MucfMaterialInserter.cc} | 117 +++++++++------ .../mucf/model/detail/MucfMaterialInserter.hh | 73 +++++++++ 13 files changed, 324 insertions(+), 386 deletions(-) delete mode 100644 src/celeritas/mucf/model/detail/DTMixMaterialCalculator.hh rename src/celeritas/mucf/model/detail/{DTMixMaterialCalculator.cc => MucfMaterialInserter.cc} (66%) create mode 100644 src/celeritas/mucf/model/detail/MucfMaterialInserter.hh diff --git a/src/celeritas/CMakeLists.txt b/src/celeritas/CMakeLists.txt index c4a2365e3f..04dec8e4e3 100644 --- a/src/celeritas/CMakeLists.txt +++ b/src/celeritas/CMakeLists.txt @@ -100,7 +100,7 @@ list(APPEND SOURCES mat/MaterialParams.cc mat/MaterialParamsOutput.cc mat/detail/Utils.cc - mucf/model/detail/DTMixMaterialCalculator.cc + mucf/model/detail/MucfMaterialInserter.cc mucf/process/MucfProcess.cc neutron/model/CascadeOptions.cc neutron/model/CascadeOptionsIO.json.cc diff --git a/src/celeritas/inp/MucfPhysics.cc b/src/celeritas/inp/MucfPhysics.cc index 32a7e37ff1..de87bbcba8 100644 --- a/src/celeritas/inp/MucfPhysics.cc +++ b/src/celeritas/inp/MucfPhysics.cc @@ -10,62 +10,66 @@ namespace celeritas { namespace inp { +namespace +{ //---------------------------------------------------------------------------// /*! - * Static muon energy CDF data for muon-catalyzed fusion. - * - * \todo This data may be loaded from an external file in the future. + * Muon energy CDF data for muon-catalyzed fusion. */ -static Grid mucf_muon_energy_cdf() +Grid mucf_muon_energy_cdf() { Grid cdf; cdf.interpolation.type = InterpolationType::cubic_spline; // Cumulative distribution data [unitless] - cdf.x = {0, - 0.04169381107491854, - 0.08664495114006499, - 0.14332247557003264, - 0.20456026058631915, - 0.2723127035830618, - 0.34136807817589576, - 0.41563517915309456, - 0.48990228013029324, - 0.5667752442996744, - 0.6306188925081434, - 0.6866449511400652, - 0.7309446254071662, - 0.7778501628664496, - 0.8104234527687297, - 0.8403908794788275, - 0.8618892508143323, - 0.8814332247557004, - 0.8970684039087949, - 0.903583061889251, - 1.0}; + cdf.x = { + 0, + 0.04169381107491854, + 0.08664495114006499, + 0.14332247557003264, + 0.20456026058631915, + 0.2723127035830618, + 0.34136807817589576, + 0.41563517915309456, + 0.48990228013029324, + 0.5667752442996744, + 0.6306188925081434, + 0.6866449511400652, + 0.7309446254071662, + 0.7778501628664496, + 0.8104234527687297, + 0.8403908794788275, + 0.8618892508143323, + 0.8814332247557004, + 0.8970684039087949, + 0.903583061889251, + 1.0, + }; // Energy [keV] - cdf.y = {0, - 0.48850540675768084, - 0.8390389347819425, - 1.2521213482687141, - 1.7153033196164724, - 2.253638712180777, - 2.854653691809707, - 3.606073540073316, - 4.470346052913727, - 5.560291219507215, - 6.700556502915258, - 7.953772477101693, - 9.194596305637525, - 10.849180562221111, - 12.353474314071864, - 14.045888515617822, - 15.650634617544647, - 17.38079707555165, - 19.111008546659452, - 19.976130619913615, - 80.0}; + cdf.y = { + 0, + 0.48850540675768084, + 0.8390389347819425, + 1.2521213482687141, + 1.7153033196164724, + 2.253638712180777, + 2.854653691809707, + 3.606073540073316, + 4.470346052913727, + 5.560291219507215, + 6.700556502915258, + 7.953772477101693, + 9.194596305637525, + 10.849180562221111, + 12.353474314071864, + 14.045888515617822, + 15.650634617544647, + 17.38079707555165, + 19.111008546659452, + 19.976130619913615, + 80.0, + }; CELER_ENSURE(cdf); return cdf; @@ -73,39 +77,35 @@ static Grid mucf_muon_energy_cdf() //---------------------------------------------------------------------------// /*! - * Static cycle rate data for muon-catalyzed fusion. - * - * \todo This data may be loaded from an external file in the future. + * Cycle rate data for muon-catalyzed fusion. */ -static std::vector mucf_cycle_rates() +std::vector mucf_cycle_rates() { //! \todo Implement return {}; -}; +} //---------------------------------------------------------------------------// /*! - * Static spin-flip data for muon-catalyzed fusion. - * - * \todo This data may be loaded from an external file in the future. + * Spin-flip data for muon-catalyzed fusion. */ -static std::vector mucf_atom_spin_flip_rates() +std::vector mucf_atom_spin_flip_rates() { //! \todo Implement return {}; -}; +} //---------------------------------------------------------------------------// /*! - * Static atom transfer data for muon-catalyzed fusion. - * - * \todo This data may be loaded from an external file in the future. + * Atom transfer data for muon-catalyzed fusion. */ -static std::vector mucf_atom_transfer_rates() +std::vector mucf_atom_transfer_rates() { //! \todo Implement return {}; -}; +} +//---------------------------------------------------------------------------// +} // namespace //---------------------------------------------------------------------------// /*! @@ -118,6 +118,7 @@ static std::vector mucf_atom_transfer_rates() MucfPhysics MucfPhysics::from_default() { MucfPhysics result; + result.scalars = MucfScalars::from_default(); result.muon_energy_cdf = mucf_muon_energy_cdf(); result.cycle_rates = mucf_cycle_rates(); result.atom_transfer = mucf_atom_transfer_rates(); diff --git a/src/celeritas/inp/MucfPhysics.hh b/src/celeritas/inp/MucfPhysics.hh index c55f13fd95..820f4fa358 100644 --- a/src/celeritas/inp/MucfPhysics.hh +++ b/src/celeritas/inp/MucfPhysics.hh @@ -10,12 +10,39 @@ #include #include "corecel/inp/Grid.hh" +#include "celeritas/Quantities.hh" #include "celeritas/mucf/Types.hh" namespace celeritas { namespace inp { +//---------------------------------------------------------------------------// +/*! + * Muon-catalyzed fusion scalars. + * + * Default values are the same used by Acceleron. + */ +struct MucfScalars +{ + // Atomic masses + units::AmuMass protium; //!< Protium atomic mass [AMU] + units::AmuMass deuterium; //!< Deuterium atomic mass [AMU] + units::AmuMass tritium; //!< Tritium atomic mass [AMU] + units::InvCcDensity liquid_hydrogen_density; //!< LHD unit [1/cm^3] + + //! Initialize with hardcoded values + static MucfScalars from_default() + { + MucfScalars result; + result.protium = units::AmuMass{1.007825031898}; + result.deuterium = units::AmuMass{2.014101777844}; + result.tritium = units::AmuMass{3.016049281320}; + result.liquid_hydrogen_density = units::InvCcDensity{4.25e22}; + return result; + } +}; + //---------------------------------------------------------------------------// /*! * Muon-catalyzed fusion mean cycle rate data. @@ -78,20 +105,18 @@ struct MucfAtomTransferRate * Muon-catalyzed fusion mean atom spin flip data. * * Spin flip rates are as a function of temperature, with each grid/table - * representing an atom pair combination and its spin (e.g. - * deuterium-tritium, spin 1). Ordering is important, thus same spin - * deuterium-tritium and tritium-deuterium have different tables, which - * leads to a different final spin flip rate for a given material - * definition. + * representing an atom pair combination and its spin (e.g. deuterium-tritium, + * spin 1). Ordering is important, thus same spin deuterium-tritium and + * tritium-deuterium have different tables, which leads to a different final + * spin flip rate for a given material definition. * * Each struct holds one of such combinations, with the full set of - * combinations being the \c vector in \c MucfPhysics - * . + * combinations being the \c vector in \c MucfPhysics . * * \note These grids are host-only, with only the final spin flip rate per - * state (which is just a \c real_type ) for each combination being needed - * in the stepping loop. This is because these rates are material - * dependent, and thus can be cached at model construction. + * state (which is just a \c real_type ) for each combination being needed in + * the stepping loop. This is because these rates are material dependent, and + * thus can be cached at model construction. */ struct MucfAtomSpinFlipRate { @@ -110,14 +135,15 @@ struct MucfAtomSpinFlipRate * and * - Mean cycle rate data for dd, dt, and tt muonic molecules. * - * Muonic atom transfer and muonic atom spin flip are secondary effects and - * not required for muCF to function. + * Muonic atom transfer and muonic atom spin flip are secondary effects and not + * required for muCF to function. */ struct MucfPhysics { template using Vec = std::vector; + MucfScalars scalars; Grid muon_energy_cdf; //!< CDF for sampling the outgoing muCF muon Vec cycle_rates; //!< Mean cycle rates for muonic molecules Vec atom_transfer; //!< Muonic atom transfer rates diff --git a/src/celeritas/mucf/Types.hh b/src/celeritas/mucf/Types.hh index ed64826410..0a7c410dda 100644 --- a/src/celeritas/mucf/Types.hh +++ b/src/celeritas/mucf/Types.hh @@ -37,12 +37,35 @@ enum class MucfMuonicMolecule size_ }; +//---------------------------------------------------------------------------// +/*! + * Enum for safely accessing hydrogen isoprotologues. + * + * Hydrogen isoprotologue molecules are: + * - Homonuclear: \f$ ^2H \f$, \f$ ^2d \f$, and \f$ ^2t \f$ + * - Heteronuclear: hd, ht, and dt. + * + * \note Muon-catalyzed fusion data is only applicable to a material with + * concentrations in thermodynamic equilibrium. This equilibrium is calculated + * at model construction from the material temperature and its h, d, and t + * fractions. + */ +enum class MucfIsoprotologueMolecule +{ + protium_protium, + protium_deuterium, + protium_tritium, + deuterium_tritium, + tritium_tritium, + size_ +}; + //---------------------------------------------------------------------------// // TYPE ALIASES //---------------------------------------------------------------------------// //! Opaque index of a muCF material component -using MuCfMatCompId = OpaqueId; +using MuCfMatId = OpaqueId; //---------------------------------------------------------------------------// } // namespace celeritas diff --git a/src/celeritas/mucf/data/DTMixMucfData.hh b/src/celeritas/mucf/data/DTMixMucfData.hh index d6104ddea6..63affe3867 100644 --- a/src/celeritas/mucf/data/DTMixMucfData.hh +++ b/src/celeritas/mucf/data/DTMixMucfData.hh @@ -63,7 +63,7 @@ struct DTMixMucfData template using Items = Collection; template - using MaterialItems = Collection; + using MaterialItems = Collection; using GridRecord = NonuniformGridRecord; using CycleTimesArray = EnumArray>; @@ -76,40 +76,20 @@ struct DTMixMucfData //!@{ //! Material-dependent data calculated at model construction - //! \c PhysMatId indexed by \c MuCfMatCompId - MaterialItems matcompid_to_matid; + //! \c PhysMatId indexed by \c MuCfMatId + MaterialItems mucfmatid_to_matid; //! Cycle times per material: [mat_comp_id][muonic_molecule][spin_index] MaterialItems cycle_times; //!< In [s] //! \todo Add mean atom spin flip times //! \todo Add mean atom transfer times //!@} - //! \todo Check whether the data are assigned + //! Check whether the data are assigned explicit CELER_FUNCTION operator bool() const { - //! \todo Finalize implementation - return particles && muon_energy_cdf && !matcompid_to_matid.empty() + return particles && muon_energy_cdf && !mucfmatid_to_matid.empty() && !cycle_times.empty() - && (matcompid_to_matid.size() == cycle_times.size()); - } - - //! Get the material componend ID from a given \c PhysMatId . - CELER_FUNCTION MuCfMatCompId material_component_id(PhysMatId matid) const - { - CELER_EXPECT(matid); - CELER_EXPECT(this); - CELER_EXPECT(!this->matcompid_to_matid.empty()); - - for (auto i : range(this->matcompid_to_matid.size())) - { - if (auto const comp_id = MuCfMatCompId{i}; - this->matcompid_to_matid[comp_id] == matid) - { - return comp_id; - } - } - // Component ID not found - return MuCfMatCompId{}; + && (mucfmatid_to_matid.size() == cycle_times.size()); } //! Assign from another set of data @@ -122,7 +102,7 @@ struct DTMixMucfData this->particles = other.particles; this->reals = other.reals; this->muon_energy_cdf = other.muon_energy_cdf; - this->matcompid_to_matid = other.matcompid_to_matid; + this->mucfmatid_to_matid = other.mucfmatid_to_matid; this->cycle_times = other.cycle_times; return *this; diff --git a/src/celeritas/mucf/executor/DTMixMucfExecutor.hh b/src/celeritas/mucf/executor/DTMixMucfExecutor.hh index ce594cbbd5..6988fb3cdf 100644 --- a/src/celeritas/mucf/executor/DTMixMucfExecutor.hh +++ b/src/celeritas/mucf/executor/DTMixMucfExecutor.hh @@ -73,12 +73,26 @@ DTMixMucfExecutor::operator()(celeritas::CoreTrackView const& track) detail::MuonicMoleculeSpinSelector select_molecule_spin(muonic_molecule); auto const molecule_spin = select_molecule_spin(rng); + // Find muCF material ID from PhysMatId + auto find = [&](PhysMatId matid) -> MuCfMatId { + CELER_EXPECT(matid); + for (auto i : range(data.mucfmatid_to_matid.size())) + { + if (auto const comp_id = MuCfMatId{i}; + data.mucfmatid_to_matid[comp_id] == matid) + { + return comp_id; + } + } + // MuCF material ID not found + return MuCfMatId{}; + }; + // Load cycle time for the selected molecule - auto const mat_comp_id - = data.material_component_id(track.material().material_id()); - CELER_ASSERT(mat_comp_id); + auto const mucf_matid = find(track.material().material_id()); + CELER_ASSERT(mucf_matid); auto const cycle_time - = data.cycle_times[mat_comp_id][muonic_molecule][molecule_spin]; + = data.cycle_times[mucf_matid][muonic_molecule][molecule_spin]; CELER_ASSERT(cycle_time > 0); // Check if muon decays before fusion happens diff --git a/src/celeritas/mucf/interactor/DDMucfInteractor.hh b/src/celeritas/mucf/interactor/DDMucfInteractor.hh index a659758f10..d5686b3008 100644 --- a/src/celeritas/mucf/interactor/DDMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/DDMucfInteractor.hh @@ -52,9 +52,12 @@ class DDMucfInteractor Channel channel_{Channel::size_}; // Allocate space for secondary particles StackAllocator& allocate_; - - // Get number of secondaries for each channel - inline CELER_FUNCTION size_type num_secondaries() const; + // Number of secondaries per channel + static constexpr EnumArray num_secondaries_{ + 3, // helium3_muon_neutron + 2, // muonichelium3_neutron + 3 // hydrogen3_muon_proton + }; // Sample Interaction secondaries template @@ -85,7 +88,7 @@ template CELER_FUNCTION Interaction DDMucfInteractor::operator()(Engine& rng) { // Allocate space for the final fusion channel - Secondary* secondaries = allocate_(this->num_secondaries()); + Secondary* secondaries = allocate_(num_secondaries_[channel_]); if (secondaries == nullptr) { // Failed to allocate space for secondaries @@ -99,25 +102,6 @@ CELER_FUNCTION Interaction DDMucfInteractor::operator()(Engine& rng) return result; } -//---------------------------------------------------------------------------// -/*! - * Return number of secondaries from each fusion channel. - */ -CELER_FUNCTION size_type DDMucfInteractor::num_secondaries() const -{ - switch (channel_) - { - case Channel::helium3_muon_neutron: - return 3; - case Channel::muonichelium3_neutron: - return 2; - case Channel::hydrogen3_muon_proton: - return 3; - default: - CELER_ASSERT_UNREACHABLE(); - } -} - //---------------------------------------------------------------------------// /*! * Sample the secondaries of the selected channel. @@ -145,7 +129,7 @@ DDMucfInteractor::sample_secondaries(Secondary* secondaries /*, other args */, CELER_ASSERT_UNREACHABLE(); } - return Span{secondaries, this->num_secondaries()}; + return Span{secondaries, num_secondaries_[channel_]}; } //---------------------------------------------------------------------------// diff --git a/src/celeritas/mucf/interactor/DTMucfInteractor.hh b/src/celeritas/mucf/interactor/DTMucfInteractor.hh index 94498de7f0..93a981811f 100644 --- a/src/celeritas/mucf/interactor/DTMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/DTMucfInteractor.hh @@ -50,9 +50,11 @@ class DTMucfInteractor Channel channel_{Channel::size_}; // Allocate space for secondary particles StackAllocator& allocate_; - - // Get number of secondaries for each channel - inline CELER_FUNCTION size_type num_secondaries() const; + // Number of secondaries per channel + static constexpr EnumArray num_secondaries_{ + 3, // alpha_muon_neutron + 2 // muonicalpha_neutron + }; // Sample Interaction secondaries template @@ -83,7 +85,7 @@ template CELER_FUNCTION Interaction DTMucfInteractor::operator()(Engine& rng) { // Allocate space for the final fusion channel - Secondary* secondaries = allocate_(this->num_secondaries()); + Secondary* secondaries = allocate_(num_secondaries_[channel_]); if (secondaries == nullptr) { // Failed to allocate space for secondaries @@ -97,23 +99,6 @@ CELER_FUNCTION Interaction DTMucfInteractor::operator()(Engine& rng) return result; } -//---------------------------------------------------------------------------// -/*! - * Return number of secondaries from each fusion channel. - */ -CELER_FUNCTION size_type DTMucfInteractor::num_secondaries() const -{ - switch (channel_) - { - case Channel::alpha_muon_neutron: - return 3; - case Channel::muonicalpha_neutron: - return 2; - default: - CELER_ASSERT_UNREACHABLE(); - } -} - //---------------------------------------------------------------------------// /*! * Sample the secondaries of the selected channel. @@ -138,7 +123,7 @@ DTMucfInteractor::sample_secondaries(Secondary* secondaries /*, other args */, CELER_ASSERT_UNREACHABLE(); } - return Span{secondaries, this->num_secondaries()}; + return Span{secondaries, num_secondaries_[channel_]}; } //---------------------------------------------------------------------------// diff --git a/src/celeritas/mucf/interactor/TTMucfInteractor.hh b/src/celeritas/mucf/interactor/TTMucfInteractor.hh index c37cda8733..30588024c2 100644 --- a/src/celeritas/mucf/interactor/TTMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/TTMucfInteractor.hh @@ -50,9 +50,11 @@ class TTMucfInteractor Channel channel_{Channel::size_}; // Allocate space for secondary particles StackAllocator& allocate_; - - // Get number of secondaries for each channel - inline CELER_FUNCTION size_type num_secondaries() const; + // Number of secondaries per channel + static constexpr EnumArray num_secondaries_{ + 4, // alpha_muon_neutron_neutron + 3 // muonicalpha_neutron_neutron + }; // Sample Interaction secondaries template @@ -83,7 +85,7 @@ template CELER_FUNCTION Interaction TTMucfInteractor::operator()(Engine& rng) { // Allocate space for the final fusion channel - Secondary* secondaries = allocate_(this->num_secondaries()); + Secondary* secondaries = allocate_(num_secondaries_[channel_]); if (secondaries == nullptr) { // Failed to allocate space for secondaries @@ -97,23 +99,6 @@ CELER_FUNCTION Interaction TTMucfInteractor::operator()(Engine& rng) return result; } -//---------------------------------------------------------------------------// -/*! - * Return number of secondaries from each fusion channel. - */ -CELER_FUNCTION size_type TTMucfInteractor::num_secondaries() const -{ - switch (channel_) - { - case Channel::alpha_muon_neutron_neutron: - return 4; - case Channel::muonicalpha_neutron_neutron: - return 3; - default: - CELER_ASSERT_UNREACHABLE(); - } -} - //---------------------------------------------------------------------------// /*! * Sample the secondaries of the selected channel. @@ -138,7 +123,7 @@ TTMucfInteractor::sample_secondaries(Secondary* secondaries /*, other args */, CELER_ASSERT_UNREACHABLE(); } - return Span{secondaries, this->num_secondaries()}; + return Span{secondaries, num_secondaries_[channel_]}; } //---------------------------------------------------------------------------// diff --git a/src/celeritas/mucf/model/DTMixMucfModel.cc b/src/celeritas/mucf/model/DTMixMucfModel.cc index 7328920b3a..0ad0e84372 100644 --- a/src/celeritas/mucf/model/DTMixMucfModel.cc +++ b/src/celeritas/mucf/model/DTMixMucfModel.cc @@ -13,11 +13,12 @@ #include "celeritas/global/ActionLauncher.hh" #include "celeritas/global/TrackExecutor.hh" #include "celeritas/grid/NonuniformGridBuilder.hh" +#include "celeritas/inp/MucfPhysics.hh" #include "celeritas/mucf/executor/DTMixMucfExecutor.hh" // IWYU pragma: associated #include "celeritas/phys/InteractionApplier.hh" // IWYU pragma: associated #include "celeritas/phys/PDGNumber.hh" -#include "detail/DTMixMaterialCalculator.hh" +#include "detail/MucfMaterialInserter.hh" namespace celeritas { @@ -98,12 +99,7 @@ DTMixMucfModel::DTMixMucfModel(ActionId id, { CELER_EXPECT(id); - using VecCycleTimes - = std::vector; - using VecMatId = std::vector; - - // Initialize static muCF physics input data - //! \todo This may be replaced by user-provided data in the future + // Initialize muCF physics input data inp::MucfPhysics inp_data = inp::MucfPhysics::from_default(); CELER_EXPECT(inp_data); @@ -115,31 +111,17 @@ DTMixMucfModel::DTMixMucfModel(ActionId id, host_data.muon_energy_cdf = build_grid_record(inp_data.muon_energy_cdf); // Calculate and cache quantities for all materials with dt mixtures - VecCycleTimes vec_cycle_times; - VecMatId vec_physmatid; + detail::MucfMaterialInserter insert(&host_data); for (auto const& matid : range(materials.num_materials())) { auto const& mat_view = materials.get(PhysMatId{matid}); - - // Construct calculator for this material - detail::DTMixMaterialCalculator mat_calculator(mat_view); - if (!mat_calculator) + if (insert(mat_view)) { - // Skip non-dt mixture materials - continue; + CELER_LOG(debug) << "Added material ID " << mat_view.material_id() + << " as a muCF d-t mixture"; } - - // Store material ID and calculated data - //! \todo Store mean atom spin flip and transfer times - vec_physmatid.push_back(PhysMatId{matid}); - vec_cycle_times.push_back(mat_calculator.cycle_times()); } - make_builder(&host_data.matcompid_to_matid) - .insert_back(vec_physmatid.begin(), vec_physmatid.end()); - make_builder(&host_data.cycle_times) - .insert_back(vec_cycle_times.begin(), vec_cycle_times.end()); - // Copy to device data_ = CollectionMirror{std::move(host_data)}; CELER_ENSURE(this->data_); diff --git a/src/celeritas/mucf/model/detail/DTMixMaterialCalculator.hh b/src/celeritas/mucf/model/detail/DTMixMaterialCalculator.hh deleted file mode 100644 index 84a93430f0..0000000000 --- a/src/celeritas/mucf/model/detail/DTMixMaterialCalculator.hh +++ /dev/null @@ -1,142 +0,0 @@ -//------------------------------- -*- C++ -*- -------------------------------// -// Copyright Celeritas contributors: see top-level COPYRIGHT file for details -// SPDX-License-Identifier: (Apache-2.0 OR MIT) -//---------------------------------------------------------------------------// -//! \file celeritas/mucf/model/detail/DTMixMaterialCalculator.hh -//---------------------------------------------------------------------------// -#pragma once - -#include "corecel/cont/EnumArray.hh" -#include "celeritas/inp/MucfPhysics.hh" -#include "celeritas/mat/MaterialView.hh" -#include "celeritas/mucf/data/DTMixMucfData.hh" - -namespace celeritas -{ -namespace detail -{ -//---------------------------------------------------------------------------// -/*! - * Enum for safely accessing hydrogen isoprotologues. - * - * Hydrogen isoprotologue molecules are: - * - Homonuclear: \f$ ^2H \f$, \f$ ^2d \f$, and \f$ ^2t \f$ - * - Heteronuclear: hd, ht, and dt. - * - * \note Muon-catalyzed fusion data is only applicable to a material with - * concentrations in thermodynamic equilibrium. This equilibrium is calculated - * at model construction from the material temperature and its h, d, and t - * fractions. - */ -enum class MucfIsoprotologueMolecule -{ - protium_protium, - protium_deuterium, - protium_tritium, - deuterium_tritium, - tritium_tritium, - size_ -}; - -//---------------------------------------------------------------------------// -/*! - * Calculate material-dependent quantities for muon-catalyzed fusion. - * - * This class calculates all the muCF data that can be cached during model - * construction. - * Use its operator bool to store material data into \c DTMixMucfData : - * \code - for (auto matid : range(materials.num_materials())) - { - auto mat_view = materials.material_view(PhysMatId{matid}); - DTMixMaterialCalculator material_calculator(mat_view); - if (material_calculator) - { - // Valid d-t mixture material; Store data - } - } - * \endcode - */ -class DTMixMaterialCalculator -{ - public: - //!@{ - //! \name Type aliases - using CycleTimesArray = EnumArray>; - //!@} - - //! Construct with material data and calculate all quantities - DTMixMaterialCalculator(MaterialView const& material); - - //! Get mean cycle times - CycleTimesArray cycle_times() const { return cycle_times_; } - - //! Check if the material is valid for muon-catalyzed fusion - explicit operator bool() const { return !cycle_times_.empty(); } - - private: - using LhdArray = EnumArray; - using EquilibriumArray = EnumArray; - using AtomicMassNumber = AtomicNumber; - - //// DATA //// - - MaterialView material_; - LhdArray lhd_densities_; - EquilibriumArray eq_densities_; - EnumArray has_isotope_; - CycleTimesArray cycle_times_; - - //// LOCAL SCALARS //// - //! \todo Values are the same used by Acceleron and may need revisiting. - // { - // Atomic masses - static constexpr units::AmuMass protium() - { - return units::AmuMass{1.007825031898}; - } - - static constexpr units::AmuMass deuterium() - { - return units::AmuMass{2.014101777844}; - } - - static constexpr units::AmuMass tritium() - { - return units::AmuMass{3.016049281320}; - } - //} - - // Liquid hydrogen density (LHD) unit [1/cm^3] - static constexpr auto liquid_hydrogen_density() - { - return units::InvCcDensity{4.25e22}; - } - - //// HELPER FUNCTIONS //// - - // Calculate dt mixture densities in units of liquid hydrogen density - LhdArray calc_lhd_densities(); - - // Calculate thermal equilibrium densities - EquilibriumArray calc_equilibrium_densities(); - - // Return muonic atom from given atomic mass number - MucfMuonicAtom from_mass_number(AtomicMassNumber mass); - - // Calculate mean fusion cycle times for all reactive muonic molecules - CycleTimesArray calc_cycle_times(ElementView const& element); - - // Calculate mean fusion cycle times for dd muonic molecules - Array calc_dd_cycle(); - - // Calculate mean fusion cycle times for dt muonic molecules - Array calc_dt_cycle(); - - // Calculate mean fusion cycle times for tt muonic molecules - Array calc_tt_cycle(); -}; - -//---------------------------------------------------------------------------// -} // namespace detail -} // namespace celeritas diff --git a/src/celeritas/mucf/model/detail/DTMixMaterialCalculator.cc b/src/celeritas/mucf/model/detail/MucfMaterialInserter.cc similarity index 66% rename from src/celeritas/mucf/model/detail/DTMixMaterialCalculator.cc rename to src/celeritas/mucf/model/detail/MucfMaterialInserter.cc index b4da247641..416008db44 100644 --- a/src/celeritas/mucf/model/detail/DTMixMaterialCalculator.cc +++ b/src/celeritas/mucf/model/detail/MucfMaterialInserter.cc @@ -2,9 +2,9 @@ // Copyright Celeritas contributors: see top-level COPYRIGHT file for details // SPDX-License-Identifier: (Apache-2.0 OR MIT) //---------------------------------------------------------------------------// -//! \file celeritas/mucf/model/detail/DTMixMaterialCalculator.cc +//! \file celeritas/mucf/model/detail/MucfMaterialInserter.cc //---------------------------------------------------------------------------// -#include "DTMixMaterialCalculator.hh" +#include "MucfMaterialInserter.hh" #include "corecel/Assert.hh" @@ -14,26 +14,41 @@ namespace detail { //---------------------------------------------------------------------------// /*! - * Construct with material data. + * Construct with \c DTMixMucfModel model data. + */ +MucfMaterialInserter::MucfMaterialInserter(HostVal* host_data) + : mucfmatid_to_matid_(&host_data->mucfmatid_to_matid) + , cycle_times_(&host_data->cycle_times) +{ + CELER_EXPECT(host_data); +} + +//---------------------------------------------------------------------------// +/*! + * Insert material information if applicable. * * Calculates and caches material-dependent properties needed by the * \c DTMixMucfModel . If the material does not contain deuterium and/or - * tritium the object's operator bool will return false. + * tritium the operator will return false. */ -DTMixMaterialCalculator::DTMixMaterialCalculator(MaterialView const& material) - : material_(material) +bool MucfMaterialInserter::operator()(MaterialView const& material) { - for (auto elcompid : range(material_.num_elements())) + LhdArray lhd_densities; + EquilibriumArray eq_densities; + CycleTimesArray cycle_times; + + for (auto elcompid : range(material.num_elements())) { auto const& element_view - = material_.element_record(ElementComponentId{elcompid}); + = material.element_record(ElementComponentId{elcompid}); if (element_view.atomic_number() != AtomicNumber{1}) { // Skip non-hydrogen elements continue; } - has_isotope_ = {false, false}; + // Found hydrogen; Check for d and t isotopes + IsotopeChecker has_isotope{false, false}; for (auto el_comp : range(element_view.num_isotopes())) { auto iso_view @@ -45,17 +60,44 @@ DTMixMaterialCalculator::DTMixMaterialCalculator(MaterialView const& material) continue; } - if (auto const atom = this->from_mass_number(mass); - atom < MucfMuonicAtom::size_) + auto const atom = this->from_mass_number(mass); + if (atom < MucfMuonicAtom::size_) { - // D and/or t isotopes found; calculate properties - has_isotope_[atom] = true; - lhd_densities_ = calc_lhd_densities(); - eq_densities_ = calc_equilibrium_densities(); - cycle_times_ = calc_cycle_times(element_view); + // Mark d or t isotope as found + has_isotope[atom] = true; } } + + if (!has_isotope[MucfMuonicAtom::deuterium] + && !has_isotope[MucfMuonicAtom::tritium]) + { + // No deuterium or tritium found; skip material + return false; + } + + // Calculate and insert material-dependent muCF properties + mucfmatid_to_matid_.push_back(material.material_id()); + cycle_times_.push_back(calc_cycle_times(element_view, has_isotope)); + //! \todo Store mean atom spin flip and transfer times } + return true; +} + +//---------------------------------------------------------------------------// +/*! + * Return \c MucfMuonicAtom from a given atomic mass number. + */ +MucfMuonicAtom MucfMaterialInserter::from_mass_number(AtomicMassNumber mass) +{ + if (mass == AtomicMassNumber{2}) + { + return MucfMuonicAtom::deuterium; + } + if (mass == AtomicMassNumber{3}) + { + return MucfMuonicAtom::tritium; + } + return MucfMuonicAtom::size_; } //---------------------------------------------------------------------------// @@ -64,7 +106,8 @@ DTMixMaterialCalculator::DTMixMaterialCalculator(MaterialView const& material) * * Used during cycle time calculations. */ -DTMixMaterialCalculator::LhdArray DTMixMaterialCalculator::calc_lhd_densities() +MucfMaterialInserter::LhdArray +MucfMaterialInserter::calc_lhd_densities(ElementView const&) { LhdArray result; @@ -80,8 +123,8 @@ DTMixMaterialCalculator::LhdArray DTMixMaterialCalculator::calc_lhd_densities() * * Used during cycle time calculations. */ -DTMixMaterialCalculator::EquilibriumArray -DTMixMaterialCalculator::calc_equilibrium_densities() +MucfMaterialInserter::EquilibriumArray +MucfMaterialInserter::calc_equilibrium_densities(ElementView const&) { EquilibriumArray result; @@ -99,8 +142,9 @@ DTMixMaterialCalculator::calc_equilibrium_densities() * or * - Multiple elements, single isotope each (separate H, d, and t elements). */ -DTMixMaterialCalculator::CycleTimesArray -DTMixMaterialCalculator::calc_cycle_times(ElementView const& element) +MucfMaterialInserter::CycleTimesArray +MucfMaterialInserter::calc_cycle_times(ElementView const& element, + IsotopeChecker const& has_isotope) { CycleTimesArray result; for (auto el_comp : range(element.num_isotopes())) @@ -114,19 +158,19 @@ DTMixMaterialCalculator::calc_cycle_times(ElementView const& element) // Calculate cycle times for dd molecules case MucfMuonicAtom::deuterium: { result[MucfMuonicMolecule::deuterium_deuterium] - = this->calc_dd_cycle(); - if (has_isotope_[MucfMuonicAtom::tritium]) + = this->calc_dd_cycle(element); + if (has_isotope[MucfMuonicAtom::tritium]) { // Calculate cycle times for dt molecules result[MucfMuonicMolecule::deuterium_tritium] - = this->calc_dt_cycle(); + = this->calc_dt_cycle(element); } break; } // Calculate cycle times for tt molecules case MucfMuonicAtom::tritium: { result[MucfMuonicMolecule::tritium_tritium] - = this->calc_tt_cycle(); + = this->calc_tt_cycle(element); break; } default: @@ -143,7 +187,7 @@ DTMixMaterialCalculator::calc_cycle_times(ElementView const& element) * * Cycle times for dd molecules come from F = 0 and F = 1 spin states. */ -Array DTMixMaterialCalculator::calc_dd_cycle() +Array MucfMaterialInserter::calc_dd_cycle(ElementView const&) { Array result; @@ -161,7 +205,7 @@ Array DTMixMaterialCalculator::calc_dd_cycle() * * Cycle times for dt molecules come from F = 1/2 and F = 3/2 spin states. */ -Array DTMixMaterialCalculator::calc_dt_cycle() +Array MucfMaterialInserter::calc_dt_cycle(ElementView const&) { Array result; @@ -179,7 +223,7 @@ Array DTMixMaterialCalculator::calc_dt_cycle() * * Cycle times for tt molecules come only from the F = 1/2 spin state. */ -Array DTMixMaterialCalculator::calc_tt_cycle() +Array MucfMaterialInserter::calc_tt_cycle(ElementView const&) { Array result; @@ -190,23 +234,6 @@ Array DTMixMaterialCalculator::calc_tt_cycle() return result; } -//---------------------------------------------------------------------------// -/*! - * Return \c MucfMuonicAtom from a given atomic mass number. - */ -MucfMuonicAtom DTMixMaterialCalculator::from_mass_number(AtomicMassNumber mass) -{ - if (mass == AtomicMassNumber{2}) - { - return MucfMuonicAtom::deuterium; - } - if (mass == AtomicMassNumber{3}) - { - return MucfMuonicAtom::tritium; - } - return MucfMuonicAtom::size_; -} - //---------------------------------------------------------------------------// } // namespace detail } // namespace celeritas diff --git a/src/celeritas/mucf/model/detail/MucfMaterialInserter.hh b/src/celeritas/mucf/model/detail/MucfMaterialInserter.hh new file mode 100644 index 0000000000..79d6109541 --- /dev/null +++ b/src/celeritas/mucf/model/detail/MucfMaterialInserter.hh @@ -0,0 +1,73 @@ +//------------------------------- -*- C++ -*- -------------------------------// +// Copyright Celeritas contributors: see top-level COPYRIGHT file for details +// SPDX-License-Identifier: (Apache-2.0 OR MIT) +//---------------------------------------------------------------------------// +//! \file celeritas/mucf/model/detail/MucfMaterialInserter.hh +//---------------------------------------------------------------------------// +#pragma once + +#include "corecel/data/CollectionBuilder.hh" +#include "celeritas/mat/MaterialView.hh" +#include "celeritas/mucf/Types.hh" +#include "celeritas/mucf/data/DTMixMucfData.hh" + +namespace celeritas +{ +namespace detail +{ +//---------------------------------------------------------------------------// +/*! + * Helper class to calculate and insert muCF material-dependent data into + * \c DTMixMucfData . + */ +class MucfMaterialInserter +{ + public: + // Construct with muCF host data + explicit MucfMaterialInserter(HostVal* host_data); + + // Insert material if it is a valid d-t mixture + bool operator()(MaterialView const& material); + + private: + //// DATA //// + + using CycleTimesArray = EnumArray>; + using LhdArray = EnumArray; + using EquilibriumArray = EnumArray; + using AtomicMassNumber = AtomicNumber; + using MucfIsotope = MucfMuonicAtom; + using IsotopeChecker = EnumArray; + + // DTMixMucfModel host data references populated by operator() + CollectionBuilder mucfmatid_to_matid_; + CollectionBuilder cycle_times_; + + //// HELPER FUNCTIONS //// + + // Return muonic atom from given atomic mass number + MucfMuonicAtom from_mass_number(AtomicMassNumber mass); + + // Calculate dt mixture densities in units of liquid hydrogen density + LhdArray calc_lhd_densities(ElementView const&); + + // Calculate thermal equilibrium densities + EquilibriumArray calc_equilibrium_densities(ElementView const&); + + // Calculate mean fusion cycle times for all reactive muonic molecules + CycleTimesArray calc_cycle_times(ElementView const& element, + IsotopeChecker const& has_isotope); + + // Calculate mean fusion cycle times for dd muonic molecules + Array calc_dd_cycle(ElementView const&); + + // Calculate mean fusion cycle times for dt muonic molecules + Array calc_dt_cycle(ElementView const&); + + // Calculate mean fusion cycle times for tt muonic molecules + Array calc_tt_cycle(ElementView const&); +}; + +//---------------------------------------------------------------------------// +} // namespace detail +} // namespace celeritas From 31a72c5f0a106248b500020ff3c06e40e59f964c Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Tue, 13 Jan 2026 11:04:39 -0500 Subject: [PATCH 13/26] Minor improvements --- .../mucf/executor/DTMixMucfExecutor.hh | 1 + .../mucf/model/detail/MucfMaterialInserter.cc | 28 +++++++++++-------- .../mucf/model/detail/MucfMaterialInserter.hh | 6 +++- 3 files changed, 22 insertions(+), 13 deletions(-) diff --git a/src/celeritas/mucf/executor/DTMixMucfExecutor.hh b/src/celeritas/mucf/executor/DTMixMucfExecutor.hh index 6988fb3cdf..030d2c2bf9 100644 --- a/src/celeritas/mucf/executor/DTMixMucfExecutor.hh +++ b/src/celeritas/mucf/executor/DTMixMucfExecutor.hh @@ -74,6 +74,7 @@ DTMixMucfExecutor::operator()(celeritas::CoreTrackView const& track) auto const molecule_spin = select_molecule_spin(rng); // Find muCF material ID from PhysMatId + // Make this a View if ever used beyond this executor auto find = [&](PhysMatId matid) -> MuCfMatId { CELER_EXPECT(matid); for (auto i : range(data.mucfmatid_to_matid.size())) diff --git a/src/celeritas/mucf/model/detail/MucfMaterialInserter.cc b/src/celeritas/mucf/model/detail/MucfMaterialInserter.cc index 416008db44..f859297ab9 100644 --- a/src/celeritas/mucf/model/detail/MucfMaterialInserter.cc +++ b/src/celeritas/mucf/model/detail/MucfMaterialInserter.cc @@ -33,10 +33,6 @@ MucfMaterialInserter::MucfMaterialInserter(HostVal* host_data) */ bool MucfMaterialInserter::operator()(MaterialView const& material) { - LhdArray lhd_densities; - EquilibriumArray eq_densities; - CycleTimesArray cycle_times; - for (auto elcompid : range(material.num_elements())) { auto const& element_view @@ -75,9 +71,14 @@ bool MucfMaterialInserter::operator()(MaterialView const& material) return false; } - // Calculate and insert material-dependent muCF properties + // Temporary data needed to calculate model data, such as cycle times + lhd_densities_ = this->calc_lhd_densities(element_view); + equilibrium_densities_ = this->calc_equilibrium_densities(element_view); + + // Calculate and insert muCF material data into model data mucfmatid_to_matid_.push_back(material.material_id()); - cycle_times_.push_back(calc_cycle_times(element_view, has_isotope)); + cycle_times_.push_back( + this->calc_cycle_times(element_view, has_isotope)); //! \todo Store mean atom spin flip and transfer times } return true; @@ -187,9 +188,10 @@ MucfMaterialInserter::calc_cycle_times(ElementView const& element, * * Cycle times for dd molecules come from F = 0 and F = 1 spin states. */ -Array MucfMaterialInserter::calc_dd_cycle(ElementView const&) +MucfMaterialInserter::MoleculeCycles +MucfMaterialInserter::calc_dd_cycle(ElementView const&) { - Array result; + MoleculeCycles result; //! \todo Implement @@ -205,9 +207,10 @@ Array MucfMaterialInserter::calc_dd_cycle(ElementView const&) * * Cycle times for dt molecules come from F = 1/2 and F = 3/2 spin states. */ -Array MucfMaterialInserter::calc_dt_cycle(ElementView const&) +MucfMaterialInserter::MoleculeCycles +MucfMaterialInserter::calc_dt_cycle(ElementView const&) { - Array result; + MoleculeCycles result; //! \todo Implement @@ -223,9 +226,10 @@ Array MucfMaterialInserter::calc_dt_cycle(ElementView const&) * * Cycle times for tt molecules come only from the F = 1/2 spin state. */ -Array MucfMaterialInserter::calc_tt_cycle(ElementView const&) +MucfMaterialInserter::MoleculeCycles +MucfMaterialInserter::calc_tt_cycle(ElementView const&) { - Array result; + MoleculeCycles result; //! \todo Implement diff --git a/src/celeritas/mucf/model/detail/MucfMaterialInserter.hh b/src/celeritas/mucf/model/detail/MucfMaterialInserter.hh index 79d6109541..ad7c307b9b 100644 --- a/src/celeritas/mucf/model/detail/MucfMaterialInserter.hh +++ b/src/celeritas/mucf/model/detail/MucfMaterialInserter.hh @@ -32,7 +32,8 @@ class MucfMaterialInserter private: //// DATA //// - using CycleTimesArray = EnumArray>; + using MoleculeCycles = Array; + using CycleTimesArray = EnumArray; using LhdArray = EnumArray; using EquilibriumArray = EnumArray; using AtomicMassNumber = AtomicNumber; @@ -42,6 +43,9 @@ class MucfMaterialInserter // DTMixMucfModel host data references populated by operator() CollectionBuilder mucfmatid_to_matid_; CollectionBuilder cycle_times_; + // Temporary quantities needed for calculating the model data + LhdArray lhd_densities_; + EquilibriumArray equilibrium_densities_; //// HELPER FUNCTIONS //// From 38667c506e6be97a57a6c46d2f75802938e6bc15 Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Tue, 13 Jan 2026 11:54:28 -0500 Subject: [PATCH 14/26] Update documentation --- doc/implementation/mucf-physics.rst | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/doc/implementation/mucf-physics.rst b/doc/implementation/mucf-physics.rst index c2b1823c7e..3671aa0428 100644 --- a/doc/implementation/mucf-physics.rst +++ b/doc/implementation/mucf-physics.rst @@ -134,9 +134,9 @@ Interactors reserved for sampling final states of the outgoing secondaries. .. doxygenclass:: celeritas::DTMixMucfExecutor Most of the data is material-dependent, and thus can be calculated and cached -during model construction. This is done by the ``DTMixMaterialCalculator``. +during model construction. This is done by the ``MucfMaterialInserter``. -.. doxygenclass:: celeritas::detail::DTMixMaterialCalculator +.. doxygenclass:: celeritas::detail::MucfMaterialInserter .. doxygenclass:: celeritas::DTMixMucfExecutor .. doxygenclass:: celeritas::DDMucfInteractor .. doxygenclass:: celeritas::DTMucfInteractor From 88acf69fc941747e70376ed6b90b6dbe3732354f Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Wed, 14 Jan 2026 10:06:11 -0500 Subject: [PATCH 15/26] Fixup --- doc/CMakeLists.txt | 5 +++++ doc/implementation/mucf-physics.rst | 1 - 2 files changed, 5 insertions(+), 1 deletion(-) diff --git a/doc/CMakeLists.txt b/doc/CMakeLists.txt index eb2d012e88..d8eb1d6382 100644 --- a/doc/CMakeLists.txt +++ b/doc/CMakeLists.txt @@ -263,6 +263,11 @@ file(GLOB _DOXYGEN_SOURCE "${PROJECT_SOURCE_DIR}/src/celeritas/em/msc/detail/*.hh" "${PROJECT_SOURCE_DIR}/src/celeritas/em/process/*.hh" "${PROJECT_SOURCE_DIR}/src/celeritas/em/xs/*.hh" + "${PROJECT_SOURCE_DIR}/src/celeritas/mucf/data/*.hh" + "${PROJECT_SOURCE_DIR}/src/celeritas/mucf/executor/*.hh" + "${PROJECT_SOURCE_DIR}/src/celeritas/mucf/interactor/*.hh" + "${PROJECT_SOURCE_DIR}/src/celeritas/mucf/model/*.hh" + "${PROJECT_SOURCE_DIR}/src/celeritas/mucf/process/*.hh" "${PROJECT_SOURCE_DIR}/src/celeritas/neutron/interactor/*.hh" "${PROJECT_SOURCE_DIR}/src/celeritas/neutron/interactor/detail/*.hh" "${PROJECT_SOURCE_DIR}/src/celeritas/neutron/model/*.hh" diff --git a/doc/implementation/mucf-physics.rst b/doc/implementation/mucf-physics.rst index 3671aa0428..511cba14a4 100644 --- a/doc/implementation/mucf-physics.rst +++ b/doc/implementation/mucf-physics.rst @@ -136,7 +136,6 @@ Interactors reserved for sampling final states of the outgoing secondaries. Most of the data is material-dependent, and thus can be calculated and cached during model construction. This is done by the ``MucfMaterialInserter``. -.. doxygenclass:: celeritas::detail::MucfMaterialInserter .. doxygenclass:: celeritas::DTMixMucfExecutor .. doxygenclass:: celeritas::DDMucfInteractor .. doxygenclass:: celeritas::DTMucfInteractor From 39b329eaf344aaa0a1dfb10cf1b2353a87aa4619 Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Wed, 14 Jan 2026 15:48:35 -0500 Subject: [PATCH 16/26] Improve documentation --- doc/CMakeLists.txt | 1 + doc/implementation/mucf-physics.rst | 43 ++++++++++--------- .../mucf/interactor/DDMucfInteractor.hh | 5 ++- .../mucf/interactor/DTMucfInteractor.hh | 4 +- .../mucf/interactor/TTMucfInteractor.hh | 4 +- 5 files changed, 33 insertions(+), 24 deletions(-) diff --git a/doc/CMakeLists.txt b/doc/CMakeLists.txt index d8eb1d6382..ad297309b8 100644 --- a/doc/CMakeLists.txt +++ b/doc/CMakeLists.txt @@ -267,6 +267,7 @@ file(GLOB _DOXYGEN_SOURCE "${PROJECT_SOURCE_DIR}/src/celeritas/mucf/executor/*.hh" "${PROJECT_SOURCE_DIR}/src/celeritas/mucf/interactor/*.hh" "${PROJECT_SOURCE_DIR}/src/celeritas/mucf/model/*.hh" + "${PROJECT_SOURCE_DIR}/src/celeritas/mucf/model/detail/*.hh" "${PROJECT_SOURCE_DIR}/src/celeritas/mucf/process/*.hh" "${PROJECT_SOURCE_DIR}/src/celeritas/neutron/interactor/*.hh" "${PROJECT_SOURCE_DIR}/src/celeritas/neutron/interactor/detail/*.hh" diff --git a/doc/implementation/mucf-physics.rst b/doc/implementation/mucf-physics.rst index 511cba14a4..75e13de3bf 100644 --- a/doc/implementation/mucf-physics.rst +++ b/doc/implementation/mucf-physics.rst @@ -72,8 +72,8 @@ sticking factor and the fusion cycle time are the main conditions that define how many fusion cycles a muon can undergo. The fusion cycle time depends on the d-t mixture, its temperature, and on the final spin of the molecule. Only muonic molecules where the total spin :math:`F = I_N \pm 1/2` is on, or has a -projection onto the total angular momentum J = 1 are reactive. The spin states -of the three possible muonic molecules are summarized in table +projection onto the total angular momentum :math:`J = 1` are reactive. The spin +states of the three possible muonic molecules are summarized in table :numref:`muon_spin_states`. @@ -108,35 +108,36 @@ enabling the ``mucf_physics`` option in Geant4 integration ------------------ -For integration interfaces, if ``mucf_physics`` option in -:cpp:class:`celeritas::ext::GeantPhysicsOptions` is enabled, the muon-catalyzed -fusion data is constructed when the ``G4MuonMinusAtomicCapture`` process is -registered. - -.. todo:: Add process/model/executor details +For integration interfaces, enabling the ``mucf_physics`` option in +:cpp:class:`celeritas::ext::GeantPhysicsOptions` will check if the +``G4MuonMinusAtomicCapture`` process is registered in the Geant4's Physics List. +If the process is present, the :cpp:class:`celeritas::inp::MucfPhysics` will be +populated, and the :cpp:class:`celeritas::MucfProcess` will be initialized. Code implementation =================== -The ``MucfProcess`` process only has the model ``DTMixMucfModel`` attached to -it, responsible for deuterium-tritium mixtures. It can simulate materials from -near absolute zero to 1500 kelvin, and it is an *at rest* model that encompasses -the full cycle---atom formation, molecule formation, and fusion. - -.. note:: Only reactive channels are implemented. +The :cpp:class:`celeritas::MucfProcess` process has only the +:cpp:class:`celeritas::DTMixMucfModel` attached to it, responsible for +deuterium-tritium mixtures. It can simulate materials from near absolute zero to +1500 kelvin. It is an *at rest* model that encompasses the full cycle---atom +formation, molecule formation, and fusion. .. doxygenclass:: celeritas::MucfProcess +.. doxygenclass:: celeritas::DTMixMucfModel -The main cycle is managed by the ``DTMixMucfExecutor``, with the -Interactors reserved for sampling final states of the outgoing secondaries. +Most of the data is material-dependent, being calculated and cached during model +construction. All of the cached quantities are calculated and added to +host/device data via :cpp:class:`celeritas::detail::MucfMaterialInserter`. -.. doxygenclass:: celeritas::DTMixMucfModel -.. doxygenclass:: celeritas::DTMixMucfExecutor +.. doxygenclass:: celeritas::detail::MucfMaterialInserter + +The main cycle is managed by the model's +:cpp:class:`celeritas::DTMixMucfExecutor`, with the Interactors reserved for +sampling final states of the outgoing secondaries. -Most of the data is material-dependent, and thus can be calculated and cached -during model construction. This is done by the ``MucfMaterialInserter``. +.. note:: Only reactive channels are implemented. -.. doxygenclass:: celeritas::DTMixMucfExecutor .. doxygenclass:: celeritas::DDMucfInteractor .. doxygenclass:: celeritas::DTMucfInteractor .. doxygenclass:: celeritas::TTMucfInteractor diff --git a/src/celeritas/mucf/interactor/DDMucfInteractor.hh b/src/celeritas/mucf/interactor/DDMucfInteractor.hh index d5686b3008..1f074097b6 100644 --- a/src/celeritas/mucf/interactor/DDMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/DDMucfInteractor.hh @@ -21,7 +21,10 @@ namespace celeritas /*! * Muon-catalyzed fusion of \f$ (dd)_\mu \f$ molecules. * - * This is an \em at-rest interaction. + * Fusion channels: + * - \f$ ^3\text{He} + \mu + n \f$ + * - \f$ (^3\text{He})_\mu + n \f$ + * - \f$ ^3\text{H} + \mu + p \f$ */ class DDMucfInteractor { diff --git a/src/celeritas/mucf/interactor/DTMucfInteractor.hh b/src/celeritas/mucf/interactor/DTMucfInteractor.hh index 93a981811f..32b413b54d 100644 --- a/src/celeritas/mucf/interactor/DTMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/DTMucfInteractor.hh @@ -21,7 +21,9 @@ namespace celeritas /*! * Muon-catalyzed fusion of \f$ (dt)_\mu \f$ molecules. * - * This is an \em at-rest interaction. + * Fusion channels: + * - \f$ \alpha + \mu + n \f$ + * - \f$ (\alpha)_\mu + n \f$ */ class DTMucfInteractor { diff --git a/src/celeritas/mucf/interactor/TTMucfInteractor.hh b/src/celeritas/mucf/interactor/TTMucfInteractor.hh index 30588024c2..0532c1b5bf 100644 --- a/src/celeritas/mucf/interactor/TTMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/TTMucfInteractor.hh @@ -21,7 +21,9 @@ namespace celeritas /*! * Muon-catalyzed fusion of \f$ (tt)_\mu \f$ molecules. * - * This is an \em at-rest interaction. + * Fusion channels: + * - \f$ \alpha + \mu + n + n \f$ + * - \f$ (\alpha)_\mu + n + n \f$ */ class TTMucfInteractor { From f7347639ff3dd4a3b7a1ca0b10ed08befd6810fd Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Wed, 14 Jan 2026 16:20:23 -0500 Subject: [PATCH 17/26] Fix `RootInterfaceLinkDef` --- src/celeritas/ext/RootInterfaceLinkDef.h | 8 ++++++-- .../data/four-steel-slabs.root-dump.json | 19 +++++++++++++++++++ 2 files changed, 25 insertions(+), 2 deletions(-) diff --git a/src/celeritas/ext/RootInterfaceLinkDef.h b/src/celeritas/ext/RootInterfaceLinkDef.h index 3e1ac8f498..0425bc3110 100644 --- a/src/celeritas/ext/RootInterfaceLinkDef.h +++ b/src/celeritas/ext/RootInterfaceLinkDef.h @@ -51,10 +51,11 @@ #pragma link C++ class celeritas::inp::Grid+; #pragma link C++ class celeritas::inp::GridReflection+; #pragma link C++ class celeritas::inp::Interpolation+; -#pragma link C++ class celeritas::inp::MucfPhysics+; -#pragma link C++ class celeritas::inp::MucfCycleRate+; #pragma link C++ class celeritas::inp::MucfAtomTransferRate+; #pragma link C++ class celeritas::inp::MucfAtomSpinFlipRate+; +#pragma link C++ class celeritas::inp::MucfCycleRate+; +#pragma link C++ class celeritas::inp::MucfPhysics+; +#pragma link C++ class celeritas::inp::MucfScalars+; #pragma link C++ class celeritas::inp::MuPairProductionEnergyTransferTable+; #pragma link C++ class celeritas::inp::NoRoughness+; #pragma link C++ class celeritas::inp::OpticalPhysics+; @@ -82,6 +83,9 @@ // Quantities #pragma link C++ class celeritas::Quantity+; #pragma link C++ class celeritas::Quantity+; +#pragma link C++ class celeritas::Quantity+; +#pragma link C++ class celeritas::Quantity+; + // Event data used by Geant4/Celeritas offloading applications #pragma link C++ class celeritas::EventData+; diff --git a/test/celeritas/data/four-steel-slabs.root-dump.json b/test/celeritas/data/four-steel-slabs.root-dump.json index 22948401fa..cba069f76c 100644 --- a/test/celeritas/data/four-steel-slabs.root-dump.json +++ b/test/celeritas/data/four-steel-slabs.root-dump.json @@ -569,6 +569,25 @@ }, "mucf_physics" : { "_typename" : "celeritas::inp::MucfPhysics", + "scalars" : { + "_typename" : "celeritas::inp::MucfScalars", + "protium" : { + "_typename" : "celeritas::Quantity", + "value_" : 0 + }, + "deuterium" : { + "_typename" : "celeritas::Quantity", + "value_" : 0 + }, + "tritium" : { + "_typename" : "celeritas::Quantity", + "value_" : 0 + }, + "liquid_hydrogen_density" : { + "_typename" : "celeritas::Quantity", + "value_" : 0 + } + }, "muon_energy_cdf" : { "_typename" : "celeritas::inp::Grid", "x" : [], From cea2b35d34e05ec475ddc6309d84b946155ece6c Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Thu, 15 Jan 2026 15:52:51 -0500 Subject: [PATCH 18/26] Address comments: - Fix integral_xs bool - Fwd declare material and particle params - Remove trailing includes --- src/celeritas/mucf/model/DTMixMucfModel.cc | 2 ++ src/celeritas/mucf/model/DTMixMucfModel.hh | 7 +++---- src/celeritas/mucf/process/MucfProcess.cc | 2 ++ src/celeritas/mucf/process/MucfProcess.hh | 7 ++++--- 4 files changed, 11 insertions(+), 7 deletions(-) diff --git a/src/celeritas/mucf/model/DTMixMucfModel.cc b/src/celeritas/mucf/model/DTMixMucfModel.cc index 0ad0e84372..d4cee2333a 100644 --- a/src/celeritas/mucf/model/DTMixMucfModel.cc +++ b/src/celeritas/mucf/model/DTMixMucfModel.cc @@ -14,9 +14,11 @@ #include "celeritas/global/TrackExecutor.hh" #include "celeritas/grid/NonuniformGridBuilder.hh" #include "celeritas/inp/MucfPhysics.hh" +#include "celeritas/mat/MaterialParams.hh" #include "celeritas/mucf/executor/DTMixMucfExecutor.hh" // IWYU pragma: associated #include "celeritas/phys/InteractionApplier.hh" // IWYU pragma: associated #include "celeritas/phys/PDGNumber.hh" +#include "celeritas/phys/ParticleParams.hh" #include "detail/MucfMaterialInserter.hh" diff --git a/src/celeritas/mucf/model/DTMixMucfModel.hh b/src/celeritas/mucf/model/DTMixMucfModel.hh index d6cd37e531..893d6c2e5e 100644 --- a/src/celeritas/mucf/model/DTMixMucfModel.hh +++ b/src/celeritas/mucf/model/DTMixMucfModel.hh @@ -6,15 +6,14 @@ //---------------------------------------------------------------------------// #pragma once -#include "celeritas/mat/MaterialParams.hh" #include "celeritas/mucf/data/DTMixMucfData.hh" -#include "celeritas/mucf/executor/DTMixMucfExecutor.hh" // IWYU pragma: associated -#include "celeritas/phys/InteractionApplier.hh" // IWYU pragma: associated #include "celeritas/phys/Model.hh" -#include "celeritas/phys/ParticleParams.hh" namespace celeritas { +class MaterialParams; +class ParticleParams; + //---------------------------------------------------------------------------// /*! * Muon-catalyzed fusion model for dd, dt, and tt molecules. diff --git a/src/celeritas/mucf/process/MucfProcess.cc b/src/celeritas/mucf/process/MucfProcess.cc index b3b3dc5d8c..cb2947e0a6 100644 --- a/src/celeritas/mucf/process/MucfProcess.cc +++ b/src/celeritas/mucf/process/MucfProcess.cc @@ -8,8 +8,10 @@ #include +#include "celeritas/mat/MaterialParams.hh" #include "celeritas/mucf/model/DTMixMucfModel.hh" #include "celeritas/phys/Model.hh" +#include "celeritas/phys/ParticleParams.hh" namespace celeritas { diff --git a/src/celeritas/mucf/process/MucfProcess.hh b/src/celeritas/mucf/process/MucfProcess.hh index 565b4c590c..23949c2acc 100644 --- a/src/celeritas/mucf/process/MucfProcess.hh +++ b/src/celeritas/mucf/process/MucfProcess.hh @@ -7,12 +7,13 @@ #pragma once #include "celeritas/Types.hh" -#include "celeritas/mat/MaterialParams.hh" -#include "celeritas/phys/ParticleParams.hh" #include "celeritas/phys/Process.hh" namespace celeritas { +class MaterialParams; +class ParticleParams; + //---------------------------------------------------------------------------// /*! * Muon-catalyzed fusion of muonic dd, dt, or tt molecules. @@ -43,7 +44,7 @@ class MucfProcess final : public Process EnergyLossGrid energy_loss(Applicability) const final; //! Whether the integral method can be used to sample interaction length - bool supports_integral_xs() const final { return true; } //! \todo Check + bool supports_integral_xs() const final { return false; } //! Whether the process applies when the particle is stopped bool applies_at_rest() const final { return true; } From 97cf7297a97eb75124cfaaf2b269e0399c43baa9 Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Fri, 16 Jan 2026 09:32:58 -0500 Subject: [PATCH 19/26] Fix missing CELER_FUNCTION --- src/celeritas/mucf/executor/detail/DDChannelSelector.hh | 7 +++---- src/celeritas/mucf/executor/detail/DTChannelSelector.hh | 7 +++---- .../mucf/executor/detail/DTMixMuonicAtomSelector.hh | 6 +++--- .../mucf/executor/detail/DTMixMuonicMoleculeSelector.hh | 5 +++-- .../mucf/executor/detail/MuonicAtomSpinSelector.hh | 1 + .../mucf/executor/detail/MuonicMoleculeSpinSelector.hh | 2 +- src/celeritas/mucf/executor/detail/TTChannelSelector.hh | 7 +++---- src/celeritas/mucf/interactor/DDMucfInteractor.hh | 1 + src/celeritas/mucf/interactor/DTMucfInteractor.hh | 1 + src/celeritas/mucf/interactor/TTMucfInteractor.hh | 3 ++- 10 files changed, 21 insertions(+), 19 deletions(-) diff --git a/src/celeritas/mucf/executor/detail/DDChannelSelector.hh b/src/celeritas/mucf/executor/detail/DDChannelSelector.hh index f2793b0979..05394c8019 100644 --- a/src/celeritas/mucf/executor/detail/DDChannelSelector.hh +++ b/src/celeritas/mucf/executor/detail/DDChannelSelector.hh @@ -29,7 +29,7 @@ class DDChannelSelector public: //! Construct with args; \todo Update documentation - inline DDChannelSelector(/* args */); + inline CELER_FUNCTION DDChannelSelector(/* args */); // Select fusion channel to be used by the interactor template @@ -44,7 +44,7 @@ class DDChannelSelector * * \todo Update documentation */ -DDChannelSelector::DDChannelSelector(/* args */) +CELER_FUNCTION DDChannelSelector::DDChannelSelector(/* args */) { CELER_NOT_IMPLEMENTED("Mucf dd fusion channel selection"); } @@ -56,8 +56,7 @@ DDChannelSelector::DDChannelSelector(/* args */) * \sa celeritas::DDMucfInteractor */ template -inline CELER_FUNCTION DDChannelSelector::Channel -DDChannelSelector::operator()(Engine&) +CELER_FUNCTION DDChannelSelector::Channel DDChannelSelector::operator()(Engine&) { Channel result{Channel::size_}; diff --git a/src/celeritas/mucf/executor/detail/DTChannelSelector.hh b/src/celeritas/mucf/executor/detail/DTChannelSelector.hh index 418571083a..47e097cdfd 100644 --- a/src/celeritas/mucf/executor/detail/DTChannelSelector.hh +++ b/src/celeritas/mucf/executor/detail/DTChannelSelector.hh @@ -29,7 +29,7 @@ class DTChannelSelector public: //! Construct with args; \todo Update documentation - inline DTChannelSelector(/* args */); + inline CELER_FUNCTION DTChannelSelector(/* args */); // Select fusion channel to be used by the interactor template @@ -44,7 +44,7 @@ class DTChannelSelector * * \todo Update documentation */ -DTChannelSelector::DTChannelSelector(/* args */) +CELER_FUNCTION DTChannelSelector::DTChannelSelector(/* args */) { CELER_NOT_IMPLEMENTED("Mucf dt fusion channel selection"); } @@ -56,8 +56,7 @@ DTChannelSelector::DTChannelSelector(/* args */) * \sa celeritas::DTMucfInteractor */ template -inline CELER_FUNCTION DTChannelSelector::Channel -DTChannelSelector::operator()(Engine&) +CELER_FUNCTION DTChannelSelector::Channel DTChannelSelector::operator()(Engine&) { Channel result{Channel::size_}; diff --git a/src/celeritas/mucf/executor/detail/DTMixMuonicAtomSelector.hh b/src/celeritas/mucf/executor/detail/DTMixMuonicAtomSelector.hh index a371141417..f193109e0c 100644 --- a/src/celeritas/mucf/executor/detail/DTMixMuonicAtomSelector.hh +++ b/src/celeritas/mucf/executor/detail/DTMixMuonicAtomSelector.hh @@ -21,7 +21,7 @@ class DTMixMuonicAtomSelector { public: //! Construct with args; \todo Update documentation - inline DTMixMuonicAtomSelector(/* args */); + inline CELER_FUNCTION DTMixMuonicAtomSelector(/* args */); // Select muonic atom template @@ -36,7 +36,7 @@ class DTMixMuonicAtomSelector * \todo Update documentation */ -DTMixMuonicAtomSelector::DTMixMuonicAtomSelector(/* args */) +CELER_FUNCTION DTMixMuonicAtomSelector::DTMixMuonicAtomSelector(/* args */) { CELER_NOT_IMPLEMENTED("Mucf muonic atom selection"); } @@ -46,7 +46,7 @@ DTMixMuonicAtomSelector::DTMixMuonicAtomSelector(/* args */) * Return selected muonic atom. */ template -inline CELER_FUNCTION MucfMuonicAtom DTMixMuonicAtomSelector::operator()(Engine&) +CELER_FUNCTION MucfMuonicAtom DTMixMuonicAtomSelector::operator()(Engine&) { MucfMuonicAtom result{MucfMuonicAtom::size_}; diff --git a/src/celeritas/mucf/executor/detail/DTMixMuonicMoleculeSelector.hh b/src/celeritas/mucf/executor/detail/DTMixMuonicMoleculeSelector.hh index 6558191a00..cc4c7bd30b 100644 --- a/src/celeritas/mucf/executor/detail/DTMixMuonicMoleculeSelector.hh +++ b/src/celeritas/mucf/executor/detail/DTMixMuonicMoleculeSelector.hh @@ -24,7 +24,7 @@ class DTMixMuonicMoleculeSelector { public: //! Construct with args; \todo Update documentation - inline DTMixMuonicMoleculeSelector(/* args */); + inline CELER_FUNCTION DTMixMuonicMoleculeSelector(/* args */); // Select muonic molecule template @@ -39,6 +39,7 @@ class DTMixMuonicMoleculeSelector * * \todo Update documentation */ +CELER_FUNCTION DTMixMuonicMoleculeSelector::DTMixMuonicMoleculeSelector(/* args */) { //! \todo Implement @@ -50,7 +51,7 @@ DTMixMuonicMoleculeSelector::DTMixMuonicMoleculeSelector(/* args */) * Return selected muonic molecule. */ template -inline CELER_FUNCTION MucfMuonicMolecule +CELER_FUNCTION MucfMuonicMolecule DTMixMuonicMoleculeSelector::operator()(Engine&) { MucfMuonicMolecule result{MucfMuonicMolecule::size_}; diff --git a/src/celeritas/mucf/executor/detail/MuonicAtomSpinSelector.hh b/src/celeritas/mucf/executor/detail/MuonicAtomSpinSelector.hh index 853950c4bd..15ff2c7fe4 100644 --- a/src/celeritas/mucf/executor/detail/MuonicAtomSpinSelector.hh +++ b/src/celeritas/mucf/executor/detail/MuonicAtomSpinSelector.hh @@ -40,6 +40,7 @@ class MuonicAtomSpinSelector /*! * Construct with muonic atom. */ +CELER_FUNCTION MuonicAtomSpinSelector::MuonicAtomSpinSelector(MucfMuonicAtom atom) : atom_(atom) { diff --git a/src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh b/src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh index 51ddedb741..ee4b5777bd 100644 --- a/src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh +++ b/src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh @@ -40,7 +40,7 @@ class MuonicMoleculeSpinSelector /*! * Construct with muonic molecule. */ -MuonicMoleculeSpinSelector::MuonicMoleculeSpinSelector( +CELER_FUNCTION MuonicMoleculeSpinSelector::MuonicMoleculeSpinSelector( MucfMuonicMolecule molecule) : molecule_(molecule) { diff --git a/src/celeritas/mucf/executor/detail/TTChannelSelector.hh b/src/celeritas/mucf/executor/detail/TTChannelSelector.hh index d69bab7bfd..1f7c0b40ea 100644 --- a/src/celeritas/mucf/executor/detail/TTChannelSelector.hh +++ b/src/celeritas/mucf/executor/detail/TTChannelSelector.hh @@ -29,7 +29,7 @@ class TTChannelSelector public: //! Construct with args; \todo Update documentation - inline TTChannelSelector(/* args */); + inline CELER_FUNCTION TTChannelSelector(/* args */); // Select fusion channel to be used by the interactor template @@ -44,7 +44,7 @@ class TTChannelSelector * * \todo Update documentation */ -TTChannelSelector::TTChannelSelector(/* args */) +CELER_FUNCTION TTChannelSelector::TTChannelSelector(/* args */) { //! \todo Implement CELER_NOT_IMPLEMENTED("Mucf tt fusion channel selection"); @@ -57,8 +57,7 @@ TTChannelSelector::TTChannelSelector(/* args */) * \sa celeritas::TTMucfInteractor */ template -inline CELER_FUNCTION TTChannelSelector::Channel -TTChannelSelector::operator()(Engine&) +CELER_FUNCTION TTChannelSelector::Channel TTChannelSelector::operator()(Engine&) { Channel result{Channel::size_}; diff --git a/src/celeritas/mucf/interactor/DDMucfInteractor.hh b/src/celeritas/mucf/interactor/DDMucfInteractor.hh index 1f074097b6..578e6f8b66 100644 --- a/src/celeritas/mucf/interactor/DDMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/DDMucfInteractor.hh @@ -74,6 +74,7 @@ class DDMucfInteractor /*! * Construct with shared and state data. */ +CELER_FUNCTION DDMucfInteractor::DDMucfInteractor(NativeCRef const& data, Channel const channel, StackAllocator& allocate) diff --git a/src/celeritas/mucf/interactor/DTMucfInteractor.hh b/src/celeritas/mucf/interactor/DTMucfInteractor.hh index 32b413b54d..87de422df8 100644 --- a/src/celeritas/mucf/interactor/DTMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/DTMucfInteractor.hh @@ -70,6 +70,7 @@ class DTMucfInteractor /*! * Construct with shared and state data. */ +CELER_FUNCTION DTMucfInteractor::DTMucfInteractor(NativeCRef const& data, Channel const channel, StackAllocator& allocate) diff --git a/src/celeritas/mucf/interactor/TTMucfInteractor.hh b/src/celeritas/mucf/interactor/TTMucfInteractor.hh index 0532c1b5bf..a9d0a03121 100644 --- a/src/celeritas/mucf/interactor/TTMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/TTMucfInteractor.hh @@ -70,6 +70,7 @@ class TTMucfInteractor /*! * Construct with shared and state data. */ +CELER_FUNCTION TTMucfInteractor::TTMucfInteractor(NativeCRef const& data, Channel const channel, StackAllocator& allocate) @@ -81,7 +82,7 @@ TTMucfInteractor::TTMucfInteractor(NativeCRef const& data, //---------------------------------------------------------------------------// /*! - * Sample a dt muonic molecule fusion. + * Sample a tt muonic molecule fusion. */ template CELER_FUNCTION Interaction TTMucfInteractor::operator()(Engine& rng) From d3b635478f54b78675ff4ec33519f8fb08d10eac Mon Sep 17 00:00:00 2001 From: Stefano Tognini Date: Mon, 19 Jan 2026 14:54:59 -0500 Subject: [PATCH 20/26] Fix device data --- src/celeritas/mucf/interactor/DDMucfInteractor.hh | 2 +- src/celeritas/mucf/interactor/DTMucfInteractor.hh | 2 +- src/celeritas/mucf/interactor/TTMucfInteractor.hh | 2 +- 3 files changed, 3 insertions(+), 3 deletions(-) diff --git a/src/celeritas/mucf/interactor/DDMucfInteractor.hh b/src/celeritas/mucf/interactor/DDMucfInteractor.hh index 578e6f8b66..3bafdd7163 100644 --- a/src/celeritas/mucf/interactor/DDMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/DDMucfInteractor.hh @@ -56,7 +56,7 @@ class DDMucfInteractor // Allocate space for secondary particles StackAllocator& allocate_; // Number of secondaries per channel - static constexpr EnumArray num_secondaries_{ + EnumArray num_secondaries_{ 3, // helium3_muon_neutron 2, // muonichelium3_neutron 3 // hydrogen3_muon_proton diff --git a/src/celeritas/mucf/interactor/DTMucfInteractor.hh b/src/celeritas/mucf/interactor/DTMucfInteractor.hh index 87de422df8..06d02c85c7 100644 --- a/src/celeritas/mucf/interactor/DTMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/DTMucfInteractor.hh @@ -53,7 +53,7 @@ class DTMucfInteractor // Allocate space for secondary particles StackAllocator& allocate_; // Number of secondaries per channel - static constexpr EnumArray num_secondaries_{ + EnumArray num_secondaries_{ 3, // alpha_muon_neutron 2 // muonicalpha_neutron }; diff --git a/src/celeritas/mucf/interactor/TTMucfInteractor.hh b/src/celeritas/mucf/interactor/TTMucfInteractor.hh index a9d0a03121..a44a078875 100644 --- a/src/celeritas/mucf/interactor/TTMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/TTMucfInteractor.hh @@ -53,7 +53,7 @@ class TTMucfInteractor // Allocate space for secondary particles StackAllocator& allocate_; // Number of secondaries per channel - static constexpr EnumArray num_secondaries_{ + EnumArray num_secondaries_{ 4, // alpha_muon_neutron_neutron 3 // muonicalpha_neutron_neutron }; From 64ece57c59852edaa39cf54bb2d1e40d77d5892b Mon Sep 17 00:00:00 2001 From: Seth R Johnson Date: Tue, 20 Jan 2026 07:17:17 -0500 Subject: [PATCH 21/26] Fix tidy warnings --- .../detail/MuonicMoleculeSpinSelector.hh | 16 ++-------------- .../mucf/interactor/DDMucfInteractor.hh | 16 ++-------------- .../mucf/interactor/DTMucfInteractor.hh | 14 ++------------ .../mucf/interactor/TTMucfInteractor.hh | 13 ++----------- 4 files changed, 8 insertions(+), 51 deletions(-) diff --git a/src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh b/src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh index ee4b5777bd..9e737a293a 100644 --- a/src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh +++ b/src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh @@ -57,20 +57,8 @@ CELER_FUNCTION size_type MuonicMoleculeSpinSelector::operator()(Engine&) { size_type result{}; - //! \todo Implement - - switch (molecule_) - { - case MucfMuonicMolecule::deuterium_deuterium: - break; - case MucfMuonicMolecule::deuterium_tritium: - break; - case MucfMuonicMolecule::tritium_tritium: - break; - default: - CELER_ASSERT_UNREACHABLE(); - } - return result; + //! \todo switch on molecule_ + CELER_ASSERT_UNREACHABLE(); } //---------------------------------------------------------------------------// diff --git a/src/celeritas/mucf/interactor/DDMucfInteractor.hh b/src/celeritas/mucf/interactor/DDMucfInteractor.hh index 3bafdd7163..b75fa71c02 100644 --- a/src/celeritas/mucf/interactor/DDMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/DDMucfInteractor.hh @@ -118,20 +118,8 @@ CELER_FUNCTION Span DDMucfInteractor::sample_secondaries(Secondary* secondaries /*, other args */, Engine&) { - switch (channel_) - { - case Channel::helium3_muon_neutron: - //! \todo Assign secondaries - break; - case Channel::muonichelium3_neutron: - //! \todo Assign secondaries - break; - case Channel::hydrogen3_muon_proton: - //! \todo Assign secondaries - break; - default: - CELER_ASSERT_UNREACHABLE(); - } + // TODO: switch on channel_ + CELER_ASSERT_UNREACHABLE(); return Span{secondaries, num_secondaries_[channel_]}; } diff --git a/src/celeritas/mucf/interactor/DTMucfInteractor.hh b/src/celeritas/mucf/interactor/DTMucfInteractor.hh index 06d02c85c7..568bf3ba0e 100644 --- a/src/celeritas/mucf/interactor/DTMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/DTMucfInteractor.hh @@ -114,18 +114,8 @@ CELER_FUNCTION Span DTMucfInteractor::sample_secondaries(Secondary* secondaries /*, other args */, Engine&) { - switch (channel_) - { - case Channel::alpha_muon_neutron: - //! \todo Assign secondaries - break; - case Channel::muonicalpha_neutron: - //! \todo Assign secondaries - break; - default: - CELER_ASSERT_UNREACHABLE(); - } - + // TODO: switch on channel_ + CELER_ASSERT_UNREACHABLE(); return Span{secondaries, num_secondaries_[channel_]}; } diff --git a/src/celeritas/mucf/interactor/TTMucfInteractor.hh b/src/celeritas/mucf/interactor/TTMucfInteractor.hh index a44a078875..d0239fc561 100644 --- a/src/celeritas/mucf/interactor/TTMucfInteractor.hh +++ b/src/celeritas/mucf/interactor/TTMucfInteractor.hh @@ -114,17 +114,8 @@ CELER_FUNCTION Span TTMucfInteractor::sample_secondaries(Secondary* secondaries /*, other args */, Engine&) { - switch (channel_) - { - case Channel::alpha_muon_neutron_neutron: - //! \todo Assign secondaries - break; - case Channel::muonicalpha_neutron_neutron: - //! \todo Assign secondaries - break; - default: - CELER_ASSERT_UNREACHABLE(); - } + // TODO: switch on channel_ + CELER_ASSERT_UNREACHABLE(); return Span{secondaries, num_secondaries_[channel_]}; } From 1fc65c96c01a872fb6650c4e8800b11adeea830a Mon Sep 17 00:00:00 2001 From: Seth R Johnson Date: Tue, 20 Jan 2026 07:23:36 -0500 Subject: [PATCH 22/26] Use double-precision for input to fix root errors, and un-inline from_default --- src/celeritas/inp/MucfPhysics.cc | 14 ++++++++++++++ src/celeritas/inp/MucfPhysics.hh | 27 ++++++++++++--------------- 2 files changed, 26 insertions(+), 15 deletions(-) diff --git a/src/celeritas/inp/MucfPhysics.cc b/src/celeritas/inp/MucfPhysics.cc index de87bbcba8..9a298eaa7b 100644 --- a/src/celeritas/inp/MucfPhysics.cc +++ b/src/celeritas/inp/MucfPhysics.cc @@ -107,6 +107,20 @@ std::vector mucf_atom_transfer_rates() //---------------------------------------------------------------------------// } // namespace +//---------------------------------------------------------------------------// +/*! + * Initialize with hardcoded values. + */ +MucfScalars MucfScalars::from_default() +{ + MucfScalars result; + result.protium = AmuMass{1.007825031898}; + result.deuterium = AmuMass{2.014101777844}; + result.tritium = AmuMass{3.016049281320}; + result.liquid_hydrogen_density = InvCcDensity{4.25e22}; + return result; +} + //---------------------------------------------------------------------------// /*! * Construct hardcoded muon-catalyzed fusion physics data. diff --git a/src/celeritas/inp/MucfPhysics.hh b/src/celeritas/inp/MucfPhysics.hh index 820f4fa358..b07430e1a3 100644 --- a/src/celeritas/inp/MucfPhysics.hh +++ b/src/celeritas/inp/MucfPhysics.hh @@ -10,9 +10,11 @@ #include #include "corecel/inp/Grid.hh" -#include "celeritas/Quantities.hh" +#include "corecel/math/Quantity.hh" #include "celeritas/mucf/Types.hh" +#include "UnitTypes.hh" + namespace celeritas { namespace inp @@ -25,22 +27,17 @@ namespace inp */ struct MucfScalars { + using AmuMass = Quantity; + using InvCcDensity = Quantity; + // Atomic masses - units::AmuMass protium; //!< Protium atomic mass [AMU] - units::AmuMass deuterium; //!< Deuterium atomic mass [AMU] - units::AmuMass tritium; //!< Tritium atomic mass [AMU] - units::InvCcDensity liquid_hydrogen_density; //!< LHD unit [1/cm^3] + AmuMass protium; //!< Protium atomic mass [AMU] + AmuMass deuterium; //!< Deuterium atomic mass [AMU] + AmuMass tritium; //!< Tritium atomic mass [AMU] + InvCcDensity liquid_hydrogen_density; //!< LHD unit [1/cm^3] - //! Initialize with hardcoded values - static MucfScalars from_default() - { - MucfScalars result; - result.protium = units::AmuMass{1.007825031898}; - result.deuterium = units::AmuMass{2.014101777844}; - result.tritium = units::AmuMass{3.016049281320}; - result.liquid_hydrogen_density = units::InvCcDensity{4.25e22}; - return result; - } + // Initialize with hardcoded values + static MucfScalars from_default(); }; //---------------------------------------------------------------------------// From e9e141dd8f854296d632431c7f42b9ee807f5d20 Mon Sep 17 00:00:00 2001 From: Seth R Johnson Date: Tue, 20 Jan 2026 07:25:33 -0500 Subject: [PATCH 23/26] Add operator bool --- src/celeritas/inp/MucfPhysics.hh | 10 +++++++++- 1 file changed, 9 insertions(+), 1 deletion(-) diff --git a/src/celeritas/inp/MucfPhysics.hh b/src/celeritas/inp/MucfPhysics.hh index b07430e1a3..66dde63063 100644 --- a/src/celeritas/inp/MucfPhysics.hh +++ b/src/celeritas/inp/MucfPhysics.hh @@ -36,6 +36,14 @@ struct MucfScalars AmuMass tritium; //!< Tritium atomic mass [AMU] InvCcDensity liquid_hydrogen_density; //!< LHD unit [1/cm^3] + //! Whether scalars have been defined + explicit operator bool() const + { + return protium > zero_quantity() && deuterium > zero_quantity() + && tritium > zero_quantity() + && liquid_hydrogen_density > zero_quantity(); + } + // Initialize with hardcoded values static MucfScalars from_default(); }; @@ -149,7 +157,7 @@ struct MucfPhysics //! Whether muon-catalyzed fusion physics is enabled explicit operator bool() const { - return muon_energy_cdf && !cycle_rates.empty(); + return scalars && muon_energy_cdf && !cycle_rates.empty(); } //! Construct hardcoded muon-catalyzed fusion physics data From c1bc9e5b07749ea62f383f37ebb11b484829f83f Mon Sep 17 00:00:00 2001 From: Seth R Johnson Date: Tue, 20 Jan 2026 07:28:27 -0500 Subject: [PATCH 24/26] fixup! Use double-precision for input to fix root errors, and un-inline from_default --- src/celeritas/inp/MucfPhysics.hh | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/src/celeritas/inp/MucfPhysics.hh b/src/celeritas/inp/MucfPhysics.hh index 66dde63063..401e58b525 100644 --- a/src/celeritas/inp/MucfPhysics.hh +++ b/src/celeritas/inp/MucfPhysics.hh @@ -11,10 +11,9 @@ #include "corecel/inp/Grid.hh" #include "corecel/math/Quantity.hh" +#include "celeritas/UnitTypes.hh" #include "celeritas/mucf/Types.hh" -#include "UnitTypes.hh" - namespace celeritas { namespace inp From f67069a26a9c9afc79f99da98d4f5abdbaf43b61 Mon Sep 17 00:00:00 2001 From: Seth R Johnson Date: Tue, 20 Jan 2026 07:29:49 -0500 Subject: [PATCH 25/26] fixup! Fix tidy warnings --- .../mucf/executor/detail/MuonicMoleculeSpinSelector.hh | 2 -- 1 file changed, 2 deletions(-) diff --git a/src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh b/src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh index 9e737a293a..67ae8d4e29 100644 --- a/src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh +++ b/src/celeritas/mucf/executor/detail/MuonicMoleculeSpinSelector.hh @@ -55,8 +55,6 @@ CELER_FUNCTION MuonicMoleculeSpinSelector::MuonicMoleculeSpinSelector( template CELER_FUNCTION size_type MuonicMoleculeSpinSelector::operator()(Engine&) { - size_type result{}; - //! \todo switch on molecule_ CELER_ASSERT_UNREACHABLE(); } From ef30f0cd88dfa5d88a2937f246722b3606d81da4 Mon Sep 17 00:00:00 2001 From: Seth R Johnson Date: Tue, 20 Jan 2026 08:27:00 -0500 Subject: [PATCH 26/26] fixup! Fix tidy warnings --- .../executor/detail/MuonicAtomSpinSelector.hh | 17 ++--------------- 1 file changed, 2 insertions(+), 15 deletions(-) diff --git a/src/celeritas/mucf/executor/detail/MuonicAtomSpinSelector.hh b/src/celeritas/mucf/executor/detail/MuonicAtomSpinSelector.hh index 15ff2c7fe4..6d22817345 100644 --- a/src/celeritas/mucf/executor/detail/MuonicAtomSpinSelector.hh +++ b/src/celeritas/mucf/executor/detail/MuonicAtomSpinSelector.hh @@ -55,21 +55,8 @@ MuonicAtomSpinSelector::MuonicAtomSpinSelector(MucfMuonicAtom atom) template CELER_FUNCTION size_type MuonicAtomSpinSelector::operator()(Engine&) { - size_type result{}; - - //! \todo Implement - - switch (atom_) - { - case MucfMuonicAtom::deuterium: - break; - case MucfMuonicAtom::tritium: - // Protium and tritium samplings are equivalent - break; - default: - CELER_ASSERT_UNREACHABLE(); - } - return result; + //! \todo switch on atom_ + CELER_ASSERT_UNREACHABLE(); } //---------------------------------------------------------------------------//