-
Notifications
You must be signed in to change notification settings - Fork 35
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
Add Wentzel macro xs calculator and fix
a_sq_factor (#1274)
* Add WentzelVI cross section calculator * Fix calculation of a_sq_factor * Add MSC test base and WentzelVI test * Rename calculator * Default Coulomb parameters to zero
- Loading branch information
1 parent
679baf7
commit c87c167
Showing
8 changed files
with
503 additions
and
122 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,109 @@ | ||
| //----------------------------------*-C++-*----------------------------------// | ||
| // Copyright 2024 UT-Battelle, LLC, and other Celeritas developers. | ||
| // See the top-level COPYRIGHT file for details. | ||
| // SPDX-License-Identifier: (Apache-2.0 OR MIT) | ||
| //---------------------------------------------------------------------------// | ||
| //! \file celeritas/em/xs/WentzelMacroXsCalculator.hh | ||
| //---------------------------------------------------------------------------// | ||
| #pragma once | ||
|
|
||
| #include "corecel/Assert.hh" | ||
| #include "corecel/Macros.hh" | ||
| #include "corecel/Types.hh" | ||
| #include "celeritas/em/data/WentzelOKVIData.hh" | ||
| #include "celeritas/em/data/WentzelVIMscData.hh" | ||
| #include "celeritas/mat/MaterialView.hh" | ||
| #include "celeritas/phys/ParticleTrackView.hh" | ||
|
|
||
| #include "WentzelHelper.hh" | ||
|
|
||
| namespace celeritas | ||
| { | ||
| //---------------------------------------------------------------------------// | ||
| /*! | ||
| * Calculate the total cross section for the Wentzel VI MSC model. | ||
| * | ||
| * \note This performs the same calculation of the total cross section (\c | ||
| * xtsec) as the Geant4 method | ||
| * G4WentzelVIModel::ComputeTransportXSectionPerVolume. | ||
| */ | ||
| class WentzelMacroXsCalculator | ||
| { | ||
| public: | ||
| //!@{ | ||
| //! \name Type aliases | ||
| using Energy = units::MevEnergy; | ||
| using XsUnits = units::Native; // [1/len] | ||
| //!@} | ||
|
|
||
| public: | ||
| // Construct with particle, material, and precalculatad Wentzel data | ||
| inline CELER_FUNCTION | ||
| WentzelMacroXsCalculator(ParticleTrackView const& particle, | ||
| MaterialView const& material, | ||
| NativeCRef<WentzelVIMscData> const& data, | ||
| NativeCRef<WentzelOKVIData> const& wentzel, | ||
| Energy cutoff); | ||
|
|
||
| // Compute the total cross section for the given angle | ||
| inline CELER_FUNCTION real_type operator()(real_type cos_theta) const; | ||
|
|
||
| private: | ||
| ParticleTrackView const& particle_; | ||
| MaterialView const& material_; | ||
| NativeCRef<WentzelOKVIData> const& wentzel_; | ||
| CoulombIds const& ids_; | ||
| Energy cutoff_; | ||
| }; | ||
|
|
||
| //---------------------------------------------------------------------------// | ||
| // INLINE DEFINITIONS | ||
| //---------------------------------------------------------------------------// | ||
| /*! | ||
| * Construct with shared model and material data. | ||
| */ | ||
| CELER_FUNCTION | ||
| WentzelMacroXsCalculator::WentzelMacroXsCalculator( | ||
| ParticleTrackView const& particle, | ||
| MaterialView const& material, | ||
| NativeCRef<WentzelVIMscData> const& data, | ||
| NativeCRef<WentzelOKVIData> const& wentzel, | ||
| Energy cutoff) | ||
| : particle_(particle) | ||
| , material_(material) | ||
| , wentzel_(wentzel) | ||
| , ids_(data.ids) | ||
| , cutoff_(cutoff) | ||
| { | ||
| } | ||
|
|
||
| //---------------------------------------------------------------------------// | ||
| /*! | ||
| * Compute the total cross section for the given angle. | ||
| */ | ||
| CELER_FUNCTION real_type | ||
| WentzelMacroXsCalculator::operator()(real_type cos_theta) const | ||
| { | ||
| real_type result = 0; | ||
|
|
||
| for (auto elcomp_id : range(ElementComponentId(material_.num_elements()))) | ||
| { | ||
| AtomicNumber z = material_.make_element_view(elcomp_id).atomic_number(); | ||
| WentzelHelper helper(particle_, material_, z, wentzel_, ids_, cutoff_); | ||
|
|
||
| real_type cos_thetamax = helper.cos_thetamax_nuclear(); | ||
| if (cos_thetamax < cos_theta) | ||
| { | ||
| result += material_.elements()[elcomp_id.get()].fraction | ||
| * (helper.calc_xs_nuclear(cos_theta, cos_thetamax) | ||
| + helper.calc_xs_electron(cos_theta, cos_thetamax)); | ||
| } | ||
| } | ||
| result *= material_.number_density(); | ||
|
|
||
| CELER_ENSURE(result >= 0); | ||
| return result; | ||
| } | ||
|
|
||
| //---------------------------------------------------------------------------// | ||
| } // namespace celeritas |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,98 @@ | ||
| //----------------------------------*-C++-*----------------------------------// | ||
| // Copyright 2024 UT-Battelle, LLC, and other Celeritas developers. | ||
| // See the top-level COPYRIGHT file for details. | ||
| // SPDX-License-Identifier: (Apache-2.0 OR MIT) | ||
| //---------------------------------------------------------------------------// | ||
| //! \file celeritas/em/MscTestBase.cc | ||
| //---------------------------------------------------------------------------// | ||
| #include "MscTestBase.hh" | ||
|
|
||
| namespace celeritas | ||
| { | ||
| namespace test | ||
| { | ||
| //---------------------------------------------------------------------------// | ||
| /*! | ||
| * Initialize states on construction. | ||
| */ | ||
| MscTestBase::MscTestBase() | ||
| { | ||
| size_type state_size = 1; | ||
| physics_state_ = StateStore<PhysicsStateData>(this->physics()->host_ref(), | ||
| state_size); | ||
| particle_state_ = StateStore<ParticleStateData>( | ||
| this->particle()->host_ref(), state_size); | ||
| geo_state_ | ||
| = StateStore<GeoStateData>(this->geometry()->host_ref(), state_size); | ||
| sim_state_ = StateStore<SimStateData>(state_size); | ||
| } | ||
|
|
||
| //---------------------------------------------------------------------------// | ||
| /*! | ||
| * Default destructor. | ||
| */ | ||
| MscTestBase::~MscTestBase() = default; | ||
|
|
||
| //---------------------------------------------------------------------------// | ||
| /*! | ||
| * Access particle state data. | ||
| */ | ||
| ParticleTrackView | ||
| MscTestBase::make_par_view(PDGNumber pdg, MevEnergy energy) const | ||
| { | ||
| CELER_EXPECT(pdg); | ||
| CELER_EXPECT(energy > zero_quantity()); | ||
| auto pid = this->particle()->find(pdg); | ||
| CELER_ASSERT(pid); | ||
|
|
||
| ParticleTrackView par{ | ||
| this->particle()->host_ref(), particle_state_.ref(), TrackSlotId{0}}; | ||
| ParticleTrackView::Initializer_t init; | ||
| init.particle_id = pid; | ||
| init.energy = energy; | ||
| par = init; | ||
| return par; | ||
| } | ||
|
|
||
| //---------------------------------------------------------------------------// | ||
| /*! | ||
| * Access particle state data. | ||
| */ | ||
| PhysicsTrackView | ||
| MscTestBase::make_phys_view(ParticleTrackView const& par, | ||
| std::string const& matname, | ||
| HostCRef<PhysicsParamsData> const& host_ref) const | ||
| { | ||
| auto mid = this->material()->find_material(matname); | ||
| CELER_ASSERT(mid); | ||
|
|
||
| // Initialize physics | ||
| PhysicsTrackView phys_view( | ||
| host_ref, physics_state_.ref(), par.particle_id(), mid, TrackSlotId{0}); | ||
| phys_view = PhysicsTrackInitializer{}; | ||
|
|
||
| // Calculate and store the energy loss (dedx) range limit | ||
| auto ppid = phys_view.eloss_ppid(); | ||
| auto grid_id = phys_view.value_grid(ValueGridType::range, ppid); | ||
| auto calc_range = phys_view.make_calculator<RangeCalculator>(grid_id); | ||
| real_type range = calc_range(par.energy()); | ||
| phys_view.dedx_range(range); | ||
|
|
||
| return phys_view; | ||
| } | ||
|
|
||
| //---------------------------------------------------------------------------// | ||
| /*! | ||
| * Access geometry state data. | ||
| */ | ||
| GeoTrackView MscTestBase::make_geo_view(real_type r) const | ||
| { | ||
| GeoTrackView geo_view( | ||
| this->geometry()->host_ref(), geo_state_.ref(), TrackSlotId{0}); | ||
| geo_view = {{r, r, r}, Real3{0, 0, 1}}; | ||
| return geo_view; | ||
| } | ||
|
|
||
| //---------------------------------------------------------------------------// | ||
| } // namespace test | ||
| } // namespace celeritas |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,96 @@ | ||
| //----------------------------------*-C++-*----------------------------------// | ||
| // Copyright 2024 UT-Battelle, LLC, and other Celeritas developers. | ||
| // See the top-level COPYRIGHT file for details. | ||
| // SPDX-License-Identifier: (Apache-2.0 OR MIT) | ||
| //---------------------------------------------------------------------------// | ||
| //! \file celeritas/em/MscTestBase.hh | ||
| //---------------------------------------------------------------------------// | ||
| #pragma once | ||
|
|
||
| #include <random> | ||
|
|
||
| #include "corecel/data/CollectionStateStore.hh" | ||
| #include "celeritas/RootTestBase.hh" | ||
| #include "celeritas/geo/GeoData.hh" | ||
| #include "celeritas/geo/GeoParams.hh" | ||
| #include "celeritas/geo/GeoTrackView.hh" | ||
| #include "celeritas/grid/RangeCalculator.hh" | ||
| #include "celeritas/mat/MaterialParams.hh" | ||
| #include "celeritas/phys/PDGNumber.hh" | ||
| #include "celeritas/phys/ParticleData.hh" | ||
| #include "celeritas/phys/ParticleParams.hh" | ||
| #include "celeritas/phys/ParticleTrackView.hh" | ||
| #include "celeritas/phys/PhysicsParams.hh" | ||
| #include "celeritas/phys/PhysicsTrackView.hh" | ||
| #include "celeritas/track/SimData.hh" | ||
| #include "celeritas/track/SimParams.hh" | ||
| #include "celeritas/track/SimTrackView.hh" | ||
|
|
||
| #include "DiagnosticRngEngine.hh" | ||
| #include "Test.hh" | ||
|
|
||
| namespace celeritas | ||
| { | ||
| namespace test | ||
| { | ||
| //---------------------------------------------------------------------------// | ||
| /*! | ||
| * Test harness base class for multiple scattering models. | ||
| */ | ||
| class MscTestBase : public RootTestBase | ||
| { | ||
| public: | ||
| //!@{ | ||
| //! \name Type aliases | ||
| using RandomEngine = DiagnosticRngEngine<std::mt19937>; | ||
| using MevEnergy = units::MevEnergy; | ||
| //!@} | ||
|
|
||
| public: | ||
| //!@{ | ||
| //! Initialize and destroy | ||
| MscTestBase(); | ||
| ~MscTestBase(); | ||
| //!@} | ||
|
|
||
| std::string_view geometry_basename() const final | ||
| { | ||
| return "four-steel-slabs"; | ||
| } | ||
|
|
||
| SPConstTrackInit build_init() override { CELER_ASSERT_UNREACHABLE(); } | ||
| SPConstAction build_along_step() override { CELER_ASSERT_UNREACHABLE(); } | ||
|
|
||
| //! Get random number generator with clean counter | ||
| RandomEngine& rng() | ||
| { | ||
| rng_.reset_count(); | ||
| return rng_; | ||
| } | ||
|
|
||
| // Access particle track data | ||
| ParticleTrackView make_par_view(PDGNumber pdg, MevEnergy energy) const; | ||
|
|
||
| // Access physics track data | ||
| PhysicsTrackView | ||
| make_phys_view(ParticleTrackView const& par, | ||
| std::string const& matname, | ||
| HostCRef<PhysicsParamsData> const& host_ref) const; | ||
|
|
||
| // Access geometry track data | ||
| GeoTrackView make_geo_view(real_type r) const; | ||
|
|
||
| private: | ||
| template<template<Ownership, MemSpace> class S> | ||
| using StateStore = CollectionStateStore<S, MemSpace::host>; | ||
|
|
||
| StateStore<PhysicsStateData> physics_state_; | ||
| StateStore<ParticleStateData> particle_state_; | ||
| StateStore<GeoStateData> geo_state_; | ||
| StateStore<SimStateData> sim_state_; | ||
| RandomEngine rng_; | ||
| }; | ||
|
|
||
| //---------------------------------------------------------------------------// | ||
| } // namespace test | ||
| } // namespace celeritas |
Oops, something went wrong.