msc/detail/UrbanMscHelper.hh

0001 //----------------------------------*-C++-*----------------------------------//
0002 // Copyright 2020-2024 UT-Battelle, LLC, and other Celeritas developers.
0003 // See the top-level COPYRIGHT file for details.
0004 // SPDX-License-Identifier: (Apache-2.0 OR MIT)
0005 //---------------------------------------------------------------------------//
0006 //! \file celeritas/em/msc/detail/UrbanMscHelper.hh
0007 //---------------------------------------------------------------------------//
0008 #pragma once
0009
0010 #include "corecel/Macros.hh"
0011 #include "corecel/Types.hh"
0012 #include "corecel/math/Algorithms.hh"
0013 #include "corecel/math/ArrayUtils.hh"
0014 #include "celeritas/Quantities.hh"
0015 #include "celeritas/Types.hh"
0016 #include "celeritas/em/data/UrbanMscData.hh"
0017 #include "celeritas/grid/EnergyLossCalculator.hh"
0018 #include "celeritas/grid/InverseRangeCalculator.hh"
0019 #include "celeritas/grid/RangeCalculator.hh"
0020 #include "celeritas/grid/ValueGridType.hh"
0021 #include "celeritas/grid/XsCalculator.hh"
0022 #include "celeritas/phys/ParticleTrackView.hh"
0023 #include "celeritas/phys/PhysicsTrackView.hh"
0024
0025 namespace celeritas
0026 {
0027 namespace detail
0028 {
0029 //---------------------------------------------------------------------------//
0030 /*!
0031  * This is a helper class for the UrbanMscStepLimit and UrbanMscScatter.
0032  *
0033  * NOTE: units are "native" units, listed here as CGS.
0034  *
0035  * \todo Refactor to UrbanMscTrackView .
0036  */
0037 class UrbanMscHelper
0038 {
0039   public:
0040     //!@{
0041     //! \name Type aliases
0042     using Energy = units::MevEnergy;
0043     using MaterialData = UrbanMscMaterialData;
0044     using UrbanMscRef = NativeCRef<UrbanMscData>;
0045     //!@}
0046
0047   public:
0048     // Construct with shared and state data
0049     inline CELER_FUNCTION UrbanMscHelper(UrbanMscRef const& shared,
0050                                          ParticleTrackView const& particle,
0051                                          PhysicsTrackView const& physics);
0052
0053     //// HELPER FUNCTIONS ////
0054
0055     //! The mean free path of the multiple scattering at the current energy
0056     //! [len]
0057     CELER_FUNCTION real_type msc_mfp() const { return lambda_; }
0058
0059     // The mean free path of the multiple scattering for a given energy [len]
0060     inline CELER_FUNCTION real_type calc_msc_mfp(Energy energy) const;
0061
0062     // TODO: the following methods are used only by MscStepLimit
0063
0064     // Calculate the energy corresponding to a given particle range
0065     inline CELER_FUNCTION Energy calc_inverse_range(real_type step) const;
0066
0067     //! Step limit scaling based on atomic number and particle type
0068     CELER_FUNCTION real_type scaled_zeff() const
0069     {
0070         return this->pmdata().scaled_zeff;
0071     }
0072
0073     // Maximum expected distance based on the track's range
0074     inline CELER_FUNCTION real_type max_step() const;
0075
0076     // The kinetic energy at the end of a given step length corrected by dedx
0077     inline CELER_FUNCTION Energy calc_end_energy(real_type step) const;
0078
0079   private:
0080     //// DATA ////
0081
0082     // References to external data
0083     UrbanMscRef const& shared_;
0084     ParticleTrackView const& particle_;
0085     PhysicsTrackView const& physics_;
0086
0087     // Precalculated mean free path (TODO: move to physics step view)
0088     real_type lambda_;  // [len]
0089
0090     // Data for this particle+material
0091     CELER_FUNCTION UrbanMscParMatData const& pmdata() const
0092     {
0093         return shared_.par_mat_data[shared_.at<UrbanMscParMatData>(
0094             physics_.material_id(), particle_.particle_id())];
0095     }
0096
0097     // Scaled cross section data for this particle+material
0098     CELER_FUNCTION XsGridData const& xs() const
0099     {
0100         return shared_.xs[shared_.at<XsGridData>(physics_.material_id(),
0101                                                  particle_.particle_id())];
0102     }
0103 };
0104
0105 //---------------------------------------------------------------------------//
0106 // INLINE DEFINITIONS
0107 //---------------------------------------------------------------------------//
0108 /*!
0109  * Construct with shared and state data.
0110  */
0111 CELER_FUNCTION
0112 UrbanMscHelper::UrbanMscHelper(UrbanMscRef const& shared,
0113                                ParticleTrackView const& particle,
0114                                PhysicsTrackView const& physics)
0115     : shared_(shared)
0116     , particle_(particle)
0117     , physics_(physics)
0118     , lambda_(this->calc_msc_mfp(particle_.energy()))
0119 {
0120     CELER_EXPECT(particle.particle_id() == shared_.ids.electron
0121                  || particle.particle_id() == shared_.ids.positron);
0122 }
0123
0124 //---------------------------------------------------------------------------//
0125 /*!
0126  * Calculate the mean free path of the msc for a given particle energy.
0127  */
0128 CELER_FUNCTION real_type UrbanMscHelper::calc_msc_mfp(Energy energy) const
0129 {
0130     CELER_EXPECT(energy > zero_quantity());
0131     XsCalculator calc_scaled_xs(this->xs(), shared_.reals);
0132
0133     real_type xsec = calc_scaled_xs(energy) / ipow<2>(energy.value());
0134     CELER_ENSURE(xsec >= 0 && 1 / xsec > 0);
0135     return 1 / xsec;
0136 }
0137
0138 //---------------------------------------------------------------------------//
0139 /*!
0140  * Calculate the energy corresponding to a given particle range.
0141  *
0142  * This is an exact value based on the range claculation. It can be used to
0143  * find the exact energy loss over a step.
0144  */
0145 CELER_FUNCTION auto
0146 UrbanMscHelper::calc_inverse_range(real_type step) const -> Energy
0147 {
0148     auto range_gid
0149         = physics_.value_grid(ValueGridType::range, physics_.eloss_ppid());
0150     auto range_to_energy
0151         = physics_.make_calculator<InverseRangeCalculator>(range_gid);
0152     return range_to_energy(step);
0153 }
0154
0155 //---------------------------------------------------------------------------//
0156 /*!
0157  * Maximum expected step length based on the track's range.
0158  */
0159 CELER_FUNCTION real_type UrbanMscHelper::max_step() const
0160 {
0161     return physics_.dedx_range() * this->pmdata().d_over_r;
0162 }
0163
0164 //---------------------------------------------------------------------------//
0165 /*!
0166  * Evaluate the kinetic energy at the end of a given msc step.
0167  */
0168 CELER_FUNCTION auto
0169 UrbanMscHelper::calc_end_energy(real_type step) const -> Energy
0170 {
0171     CELER_EXPECT(step <= physics_.dedx_range());
0172     real_type range = physics_.dedx_range();
0173     if (step <= range * shared_.params.dtrl())
0174     {
0175         auto eloss_gid = physics_.value_grid(ValueGridType::energy_loss,
0176                                              physics_.eloss_ppid());
0177         // Assume constant energy loss rate over the step
0178         real_type dedx = physics_.make_calculator<EnergyLossCalculator>(
0179             eloss_gid)(particle_.energy());
0180
0181         return particle_.energy() - Energy{step * dedx};
0182     }
0183     else
0184     {
0185         // Longer step is calculated exactly with inverse range
0186         return this->calc_inverse_range(range - step);
0187     }
0188 }
0189
0190 //---------------------------------------------------------------------------//
0191 }  // namespace detail
0192 }  // namespace celeritas