em/xs/LivermorePEMicroXsCalculator.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/xs/LivermorePEMicroXsCalculator.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/PolyEvaluator.hh"
0013 #include "corecel/math/Quantity.hh"
0014 #include "celeritas/Quantities.hh"
0015 #include "celeritas/Types.hh"
0016 #include "celeritas/em/data/LivermorePEData.hh"
0017 #include "celeritas/grid/NonuniformGridCalculator.hh"
0018
0019 namespace celeritas
0020 {
0021 //---------------------------------------------------------------------------//
0022 /*!
0023  * Calculate photoelectric effect cross sections using the Livermore data.
0024  *
0025  * The Livermore photoelectric data is loaded from Geant4 low-energy EM data
0026  * files.
0027  */
0028 class LivermorePEMicroXsCalculator
0029 {
0030   public:
0031     //!@{
0032     //! \name Type aliases
0033     using ParamsRef = NativeCRef<LivermorePEData>;
0034     using Energy = RealQuantity<LivermoreSubshell::EnergyUnits>;
0035     using BarnXs = units::BarnXs;
0036     //!@}
0037
0038   public:
0039     // Construct with shared and state data
0040     inline CELER_FUNCTION
0041     LivermorePEMicroXsCalculator(ParamsRef const& shared, Energy energy);
0042
0043     // Compute cross section
0044     inline CELER_FUNCTION BarnXs operator()(ElementId el_id) const;
0045
0046   private:
0047     // Shared constant physics properties
0048     ParamsRef const& shared_;
0049     // Incident gamma energy
0050     Energy const inc_energy_;
0051 };
0052
0053 //---------------------------------------------------------------------------//
0054 // INLINE DEFINITIONS
0055 //---------------------------------------------------------------------------//
0056 /*!
0057  * Construct with shared and state data.
0058  */
0059 CELER_FUNCTION LivermorePEMicroXsCalculator::LivermorePEMicroXsCalculator(
0060     ParamsRef const& shared, Energy energy)
0061     : shared_(shared), inc_energy_(energy.value())
0062 {
0063 }
0064
0065 //---------------------------------------------------------------------------//
0066 /*!
0067  * Compute cross section.
0068  */
0069 CELER_FUNCTION
0070 auto LivermorePEMicroXsCalculator::operator()(ElementId el_id) const -> BarnXs
0071 {
0072     CELER_EXPECT(el_id);
0073     LivermoreElement const& el = shared_.xs.elements[el_id];
0074     auto const& shells = shared_.xs.shells[el.shells];
0075
0076     // In Geant4, if the incident gamma energy is below the lowest binding
0077     // energy, it is set to the binding energy so that the photoelectric cross
0078     // section is constant rather than zero for low energy gammas.
0079     Energy energy = max(inc_energy_, shells.back().binding_energy);
0080     real_type inv_energy = 1. / energy.value();
0081
0082     real_type result = 0.;
0083     if (energy >= el.thresh_lo)
0084     {
0085         // Fit parameters from the final shell are used to calculate the cross
0086         // section integrated over all subshells
0087         auto const& param = shells.back().param[energy < el.thresh_hi ? 0 : 1];
0088         PolyEvaluator<real_type, 5> eval_poly(param);
0089
0090         // Use the parameterization of the integrated subshell cross sections
0091         result = inv_energy * eval_poly(inv_energy);
0092     }
0093     else if (energy >= shells.front().binding_energy)
0094     {
0095         // Use tabulated cross sections above K-shell energy but below energy
0096         // limit for parameterization
0097         NonuniformGridCalculator calc_xs(el.xs_hi, shared_.xs.reals);
0098         result = ipow<3>(inv_energy) * calc_xs(energy.value());
0099     }
0100     else
0101     {
0102         CELER_ASSERT(el.xs_lo);
0103         // Use tabulated cross sections below K-shell energy
0104         NonuniformGridCalculator calc_xs(el.xs_lo, shared_.xs.reals);
0105         result = ipow<3>(inv_energy) * calc_xs(energy.value());
0106     }
0107     return BarnXs{result};
0108 }
0109
0110 //---------------------------------------------------------------------------//
0111 }  // namespace celeritas