msc/detail/UrbanMscHelper.hh

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