celeritas/mat/ElementSelector.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/mat/ElementSelector.hh
0007 //---------------------------------------------------------------------------//
0008 #pragma once
0009
0010 #include "corecel/Assert.hh"
0011 #include "corecel/Macros.hh"
0012 #include "corecel/Types.hh"
0013 #include "corecel/cont/Range.hh"
0014 #include "corecel/cont/Span.hh"
0015 #include "celeritas/Types.hh"
0016 #include "celeritas/random/distribution/GenerateCanonical.hh"
0017
0018 #include "MaterialView.hh"
0019
0020 namespace celeritas
0021 {
0022 //---------------------------------------------------------------------------//
0023 /*!
0024  * Make a weighted random selection of an element.
0025  *
0026  * The element chooser is for selecting an elemental component (atom) of a
0027  * material based on the microscopic cross section and the abundance fraction
0028  * of the element in the material.
0029  *
0030  * On construction, the element chooser uses the provided arguments to
0031  * precalculate all the microscopic cross sections in the given storage space.
0032  * The given function `calc_micro_xs` must accept a `ElementId` and return a
0033  * `real_type`, a non-negative microscopic cross section.
0034  *
0035  * The element chooser does \em not calculate macroscopic cross sections
0036  * because they're multiplied by fraction, not number density, and we only
0037  * care about the fractional abundances and cross section weighting.
0038  *
0039  * \code
0040     ElementSelector select_element(mat, calc_micro, storage);
0041     real_type total_macro_xs
0042         = select_element.material_micro_xs() * mat.number_density();
0043     ElementComponentId id = select_element(rng);
0044     real_type selected_micro_xs
0045         = select_element.elemental_micro_xs()[el.get()];
0046     ElementView el = mat.make_element_view(id);
0047     // use el.Z(), etc.
0048    \endcode
0049  *
0050  * Note that the units of the calculated microscopic cross section will be
0051  * identical to the units returned by `calc_micro_xs`. The macroscopic cross
0052  * section units (micro times \c mat.number_density() ) will be 1/len if and
0053  * only if calc_micro units are len^2.
0054  *
0055  * \todo Refactor to use Selector.
0056  */
0057 class ElementSelector
0058 {
0059   public:
0060     //!@{
0061     //! \name Type aliases
0062     using SpanReal = Span<real_type>;
0063     using SpanConstReal = LdgSpan<real_type const>;
0064     //!@}
0065
0066   public:
0067     // Construct with material, xs calculator, and storage.
0068     template<class MicroXsCalc>
0069     inline CELER_FUNCTION ElementSelector(MaterialView const& material,
0070                                           MicroXsCalc&& calc_micro_xs,
0071                                           SpanReal micro_xs_storage);
0072
0073     // Sample with the given RNG
0074     template<class Engine>
0075     inline CELER_FUNCTION ElementComponentId operator()(Engine& rng) const;
0076
0077     //! Weighted material microscopic cross section
0078     CELER_FUNCTION real_type material_micro_xs() const { return material_xs_; }
0079
0080     // Individual (unweighted) microscopic cross sections
0081     inline CELER_FUNCTION SpanConstReal elemental_micro_xs() const;
0082
0083   private:
0084     Span<MatElementComponent const> elements_;
0085     real_type material_xs_;
0086     real_type* elemental_xs_;
0087 };
0088
0089 //---------------------------------------------------------------------------//
0090 // INLINE DEFINITIONS
0091 //---------------------------------------------------------------------------//
0092 /*!
0093  * Construct with material, xs calculator, and storage.
0094  */
0095 template<class MicroXsCalc>
0096 CELER_FUNCTION ElementSelector::ElementSelector(MaterialView const& material,
0097                                                 MicroXsCalc&& calc_micro_xs,
0098                                                 SpanReal storage)
0099     : elements_(material.elements())
0100     , material_xs_(0)
0101     , elemental_xs_(storage.data())
0102 {
0103     CELER_EXPECT(!elements_.empty());
0104     CELER_EXPECT(storage.size() >= material.num_elements());
0105     for (auto i : range<size_type>(elements_.size()))
0106     {
0107         real_type const micro_xs
0108             = (calc_micro_xs(elements_[i].element)).value();
0109         CELER_ASSERT(micro_xs >= 0);
0110         real_type const frac = elements_[i].fraction;
0111
0112         elemental_xs_[i] = micro_xs;
0113         material_xs_ += micro_xs * frac;
0114     }
0115     CELER_ENSURE(material_xs_ >= 0);
0116 }
0117
0118 //---------------------------------------------------------------------------//
0119 /*!
0120  * Sample the element with the given RNG.
0121  *
0122  * To reduce register usage, this function starts with the cumulative sums and
0123  * counts backward.
0124  */
0125 template<class Engine>
0126 CELER_FUNCTION ElementComponentId ElementSelector::operator()(Engine& rng) const
0127 {
0128     real_type accum_xs = -material_xs_ * generate_canonical(rng);
0129     size_type i = 0;
0130     size_type imax = elements_.size() - 1;
0131     for (; i != imax; ++i)
0132     {
0133         accum_xs += elements_[i].fraction * elemental_xs_[i];
0134         if (accum_xs > 0)
0135             break;
0136     }
0137     return ElementComponentId{i};
0138 }
0139
0140 //---------------------------------------------------------------------------//
0141 /*!
0142  * Get individual microscopic cross sections.
0143  */
0144 CELER_FUNCTION auto ElementSelector::elemental_micro_xs() const -> SpanConstReal
0145 {
0146     return {elemental_xs_, elements_.size()};
0147 }
0148
0149 //---------------------------------------------------------------------------//
0150 }  // namespace celeritas