celeritas/random/TabulatedElementSelector.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/random/TabulatedElementSelector.hh
0006 //---------------------------------------------------------------------------//
0007 #pragma once
0008
0009 #include "corecel/Assert.hh"
0010 #include "corecel/Macros.hh"
0011 #include "corecel/Types.hh"
0012 #include "corecel/random/distribution/GenerateCanonical.hh"
0013 #include "celeritas/Types.hh"
0014 #include "celeritas/grid/XsCalculator.hh"
0015 #include "celeritas/phys/PhysicsData.hh"
0016
0017 namespace celeritas
0018 {
0019 //---------------------------------------------------------------------------//
0020 /*!
0021  * Make a weighted random selection of an element.
0022  *
0023  * This selects an elemental component (atom) of a material based on the
0024  * precalculated cross section CDF tables of the elements in the material.
0025  * Unlike \c ElementSelector which calculates the microscopic cross sections on
0026  * the fly, this interpolates the values using tabulated CDF grids.
0027  */
0028 class TabulatedElementSelector
0029 {
0030   public:
0031     //!@{
0032     //! \name Type aliases
0033     using Energy = RealQuantity<XsGridRecord::EnergyUnits>;
0034     using GridValues
0035         = Collection<ValueGrid, Ownership::const_reference, MemSpace::native>;
0036     using GridIdValues
0037         = Collection<ValueGridId, Ownership::const_reference, MemSpace::native>;
0038     using Values
0039         = Collection<real_type, Ownership::const_reference, MemSpace::native>;
0040     //!@}
0041
0042   public:
0043     // Construct with xs CDF data for a particular model and material
0044     inline CELER_FUNCTION TabulatedElementSelector(ValueTable const& table,
0045                                                    GridValues const& grids,
0046                                                    GridIdValues const& ids,
0047                                                    Values const& reals,
0048                                                    Energy energy);
0049
0050     // Sample with the given RNG
0051     template<class Engine>
0052     inline CELER_FUNCTION ElementComponentId operator()(Engine& rng) const;
0053
0054   private:
0055     ValueTable const& table_;
0056     GridValues const& grids_;
0057     GridIdValues const& ids_;
0058     Values const& reals_;
0059     Energy const energy_;
0060 };
0061
0062 //---------------------------------------------------------------------------//
0063 // INLINE DEFINITIONS
0064 //---------------------------------------------------------------------------//
0065 /*!
0066  * Construct with xs CDF data for a particular model and material.
0067  */
0068 CELER_FUNCTION
0069 TabulatedElementSelector::TabulatedElementSelector(ValueTable const& table,
0070                                                    GridValues const& grids,
0071                                                    GridIdValues const& ids,
0072                                                    Values const& reals,
0073                                                    Energy energy)
0074     : table_(table), grids_(grids), ids_(ids), reals_(reals), energy_(energy)
0075 {
0076     CELER_EXPECT(table);
0077 }
0078
0079 //---------------------------------------------------------------------------//
0080 /*!
0081  * Sample the element with the given RNG.
0082  */
0083 template<class Engine>
0084 CELER_FUNCTION ElementComponentId
0085 TabulatedElementSelector::operator()(Engine& rng) const
0086 {
0087     size_type i = 0;
0088     real_type u = generate_canonical(rng);
0089     for (; i < table_.grids.size() - 1; ++i)
0090     {
0091         ValueGridId grid_id = ids_[table_.grids[i]];
0092         CELER_ASSERT(grid_id < grids_.size());
0093         XsCalculator calc_xs(grids_[grid_id], reals_);
0094         if (calc_xs(energy_) > u)
0095             break;
0096     }
0097     return ElementComponentId{i};
0098 }
0099
0100 //---------------------------------------------------------------------------//
0101 }  // namespace celeritas