Line data    Source code

       1            0 : // Distributed under the MIT License.
       2              : // See LICENSE.txt for details.
       3              : 
       4              : #pragma once
       5              : 
       6              : #include <array>
       7              : #include <pup.h>
       8              : #include <string>
       9              : #include <vector>
      10              : 
      11              : #include "DataStructures/BoostMultiArray.hpp"
      12              : #include "DataStructures/Tensor/TypeAliases.hpp"
      13              : #include "NumericalAlgorithms/Interpolation/MultiLinearSpanInterpolation.hpp"
      14              : #include "Options/String.hpp"
      15              : #include "Utilities/Gsl.hpp"
      16              : #include "Utilities/Serialization/CharmPupable.hpp"
      17              : 
      18              : /// \cond
      19              : class DataVector;
      20              : /// \endcond
      21              : 
      22              : namespace Particles::MonteCarlo {
      23              : 
      24              : /// Class responsible for reading neutrino-matter interaction
      25              : /// tables.
      26              : template <size_t EnergyBins, size_t NeutrinoSpecies>
      27            1 : class NeutrinoInteractionTable : public PUP::able {
      28              :  public:
      29              :   /// Read table from disk and stores interaction rates.
      30            1 :   explicit NeutrinoInteractionTable(const std::string& filename);
      31              : 
      32            0 :   static constexpr Options::String help = {
      33              :       "Tabulated neutrino-matter interactions in NuLib format.\n"
      34              :       "Emissivity, absorption and scattering opacities are    \n"
      35              :       "tabulated as a function of density, eletron fraction   \n"
      36              :       "and temperature."};
      37              : 
      38            0 :   struct TableFilename {
      39            0 :     using type = std::string;
      40            0 :     static constexpr Options::String help{"File name of the NuLib table"};
      41              :   };
      42              : 
      43            0 :   using options = tmpl::list<TableFilename>;
      44              : 
      45              :   /// Explicit instantiation from table values, for tests
      46            1 :   NeutrinoInteractionTable(
      47              :       std::vector<double> table_data_,
      48              :       const std::array<double, EnergyBins>& table_neutrino_energies_,
      49              :       std::vector<double> table_log_density_,
      50              :       std::vector<double> table_log_temperature_,
      51              :       std::vector<double> table_electron_fraction_);
      52              : 
      53            0 :   explicit NeutrinoInteractionTable(CkMigrateMessage* msg) : PUP::able(msg) {}
      54              : 
      55              :   using PUP::able::register_constructor;
      56            0 :   void pup(PUP::er& p) override;
      57            0 :   WRAPPED_PUPable_decl_template(NeutrinoInteractionTable);
      58              : 
      59              :   /// Interpolate interaction rates to given values of density,
      60              :   /// temperature and electron fraction.
      61            1 :   void get_neutrino_matter_interactions(
      62              :       gsl::not_null<
      63              :           std::array<std::array<DataVector, EnergyBins>, NeutrinoSpecies>*>
      64              :           emissivity_in_cell,
      65              :       gsl::not_null<
      66              :           std::array<std::array<DataVector, EnergyBins>, NeutrinoSpecies>*>
      67              :           absorption_opacity,
      68              :       gsl::not_null<
      69              :           std::array<std::array<DataVector, EnergyBins>, NeutrinoSpecies>*>
      70              :           scattering_opacity,
      71              :       const Scalar<DataVector>& electron_fraction,
      72              :       const Scalar<DataVector>& rest_mass_density,
      73              :       const Scalar<DataVector>& temperature,
      74              :       const double& minimum_temperature) const;
      75              : 
      76            0 :   const std::array<double, EnergyBins>& get_neutrino_energies() const {
      77              :     return table_neutrino_energies;
      78              :   }
      79              : 
      80              :  private:
      81            0 :   void initialize_interpolator();
      82              : 
      83              :   // Stores emissivity, absorption_opacity, scattering_opacity
      84              :   // For each quantities, there are NeutrinoSpecies * EnergyBins
      85              :   // variables store as a function of log(density), log(temperature)
      86              :   // and electron fraction.
      87              :   // The indexing varies fastest in EnergyBins, then Species, then
      88              :   // log(density), then log(temperature), and finally Ye.
      89            0 :   std::vector<double> table_data{};
      90              :   // Central energy of each bin
      91            0 :   std::array<double, EnergyBins> table_neutrino_energies;
      92              :   // Table discretization
      93            0 :   std::vector<double> table_log_density{};
      94            0 :   std::vector<double> table_log_temperature{};
      95            0 :   std::vector<double> table_electron_fraction{};
      96              : 
      97              :   intrp::UniformMultiLinearSpanInterpolation<3,
      98              :                                              3 * EnergyBins * NeutrinoSpecies>
      99            0 :       interpolator_{};
     100              : 
     101            0 :   const double min_kappa = 1.e-70;
     102            0 :   const double max_kappa = 1.e70;
     103              : };
     104              : 
     105              : }  // namespace Particles::MonteCarlo