update main

CMakeLists.txt

@@ -53,8 +53,8 @@ target_include_directories(singularity-opac::flags
                            INTERFACE
                            $<BUILD_INTERFACE:${PROJECT_SOURCE_DIR}>)
 
-include (SetupDeps)
 include (SetupOptions)
+include (SetupDeps)
 include (SetupCompilers)
 include (SetupFlags)
 

README.md

@@ -57,6 +57,21 @@ A number of options are avaialable for compiling:
 | --------------------------------- | ------- | ------------------------------------------------------------------------------------ |
 | SINGULARITY_BUILD_TESTS           | OFF     | Build test infrastructure.                                                           |
 | SINGULARITY_USE_HDF5              | ON      | Enables HDF5. Required for Spiner opacities.                                         |
+| SINGULARITY_KOKKOS_IN_TREE        | OFF     | Force cmake to use Kokkos source included in tree.                                   |
+
+### Loading ASCII Data
+
+Currently, the MeanOpacity class, defined in singularity-opac/photons/mean_opacity_photons.hpp, supports
+loading grey Rosseland and Planck opacity data in an ASCII format.  An example of this format is
+provided by singularity-opac/photons/example_ascii/kap_plaw.txt.  The 1st row of the header has the
+number of density points, NRho, then the number of temperature points, NT.  The 2nd (3rd) row of the header
+has min and max density (temperature) bounds.  These bounds are inclusive, so the opacity data in the file
+should have evaluations at these min and max values.  The rest of the ASCII file is a two-column table, where
+the 1st (2nd) column is Rosseland (Planck) opacity.  The number of rows in each column is NRhoxNT, where
+density is the slow index and temperature is the fast index (thus the row index = temperature index
++ NT x (density index), indexing from 0).  Each opacity is assumed to be evaluated on log-spaced
+density and temperature grids, where these grids are defined by NRho, NT, and the (again inclusive) min and
+max bounds the header.
 
 ## Copyright
 

cmake/SetupDeps.cmake

@@ -10,23 +10,32 @@ else()
   message(status "CUDA::toolkit provided by parent package")
 endif()
 
+#=======================================
+# Setup Kokkos
+# - provides Kokkos::kokkos
+#=======================================
+if (SINGULARITY_USE_KOKKOS)
+  if (NOT TARGET Kokkos::kokkos)
+    message(status "Kokkos::kokkos must be found")
+    if (SINGULARITY_KOKKOS_IN_TREE)
+      message(status "Using in-tree Kokkos")
+      add_subdirectory(utils/kokkos)
+    else()
+      message(status "Using system Kokkos if available")
+      find_package(Kokkos REQUIRED)
+    endif()
+  else()
+    message(status "Kokkos::kokkos provided by parent package")
+  endif()
+endif()
+
 #=======================================
 # Setup ports of call
 # - provides PortsofCall::PortsofCall
 #=======================================
 find_package(PortsofCall REQUIRED)
 target_link_libraries(singularity-opac::flags INTERFACE PortsofCall::PortsofCall)
 
-#=======================================
-# Setup Kokkos
-# - provides Kokkos::kokkos
-#=======================================
-if (NOT TARGET Kokkos::kokkos)
-  find_package(Kokkos QUIET)
-else()
-  message(status "Kokkos::kokkos provided by parent package")
-endif()
-
 #=======================================
 # Find HDF5
 # - [email protected]+ provides HDF5::HDF5, but

cmake/SetupOptions.cmake

@@ -17,11 +17,18 @@ option (SINGULARITY_USE_FMATH "Enable fast-math logarithms" ON)
 # Dependency options
 #=======================================
 # check for using in-tree third-party dependencies
-option (SINULARITY_KOKKOS_IN_TREE "Use in-tree dependencies" OFF)
-
 option (SINGULARITY_USE_KOKKOS "Use Kokkos for portability" OFF)
 option (SINGULARITY_USE_CUDA "Enable cuda support" OFF)
 option (SINGULARITY_USE_HDF5 "Pull in hdf5" OFF)
+cmake_dependent_option(SINULARITY_KOKKOS_IN_TREE
+  "Use in-tree dependencies" OFF
+  ${SINGULARITY_USE_KOKKOS} OFF)
+
+# Kill cmake's package registry because it can interfere
+if (SINGULARITY_KOKKOS_IN_TREE)
+  set(CMAKE_FIND_PACKAGE_NO_PACKAGE_REGISTRY ON CACHE BOOL "" FORCE)
+  set(CMAKE_FIND_PACKAGE_NO_SYSTEM_PACKAGE_REGISTRY ON CACHE BOOL "" FORCE)
+endif()
 
 # If the conditional is TRUE, it's the first default, else it's the
 # second.

singularity-opac/constants/constants.hpp

@@ -1,5 +1,5 @@
 // ======================================================================
-// © 2021. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract
 // 89233218CNA000001 for Los Alamos National Laboratory (LANL), which
 // is operated by Triad National Security, LLC for the U.S.
@@ -194,6 +194,47 @@ struct PhysicalConstants {
   static constexpr Real gA = axial_vector_coupling_constant;
 };
 
+struct RuntimePhysicalConstants {
+  template <typename T>
+  PORTABLE_INLINE_FUNCTION
+  RuntimePhysicalConstants(T pc)
+      : na(pc.na), alpha(pc.alpha), h(pc.h), hbar(pc.hbar), kb(pc.kb),
+        r_gas(pc.r_gas), qe(pc.qe), c(pc.c), g_newt(pc.g_newt),
+        g_accel(pc.g_accel), me(pc.me), mp(pc.mp), mn(pc.mn), amu(pc.amu),
+        sb(pc.sb), ar(pc.ar), faraday(pc.faraday), mu0(pc.mu0), eps0(pc.eps0),
+        eV(pc.eV), Fc(pc.Fc), nu_sigma0(pc.nu_sigma0), gA(pc.gA) {}
+
+  const Real na;
+  const Real alpha;
+  const Real h;
+  const Real hbar;
+  const Real kb;
+  const Real r_gas;
+  const Real qe;
+  const Real c;
+  const Real g_newt;
+  const Real g_accel;
+  const Real me;
+  const Real mp;
+  const Real mn;
+  const Real amu;
+  const Real sb;
+  const Real ar;
+  const Real faraday;
+  const Real mu0;
+  const Real eps0;
+  const Real eV;
+  const Real Fc;
+  const Real nu_sigma0;
+  const Real gA;
+};
+
+template <typename T>
+PORTABLE_INLINE_FUNCTION
+RuntimePhysicalConstants GetRuntimePhysicalConstants(T phys_constants) {
+  return RuntimePhysicalConstants(phys_constants);
+}
+
 using PhysicalConstantsUnity =
     PhysicalConstants<BaseUnity, UnitConversionDefault>;
 using PhysicalConstantsSI =

singularity-opac/neutrinos/brt_neutrinos.hpp

@@ -1,5 +1,5 @@
 // ======================================================================
-// © 2021. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract
 // 89233218CNA000001 for Los Alamos National Laboratory (LANL), which
 // is operated by Triad National Security, LLC for the U.S.
@@ -32,6 +32,8 @@ namespace neutrinos {
 template <typename ThermalDistribution, typename pc = PhysicalConstantsCGS>
 class BRTOpacity {
  public:
+  using PC = pc;
+
   BRTOpacity() = default;
   BRTOpacity(const ThermalDistribution &dist) : dist_(dist) {}
   BRTOpacity GetOnDevice() { return *this; }

singularity-opac/neutrinos/gray_opacity_neutrinos.hpp

@@ -1,5 +1,5 @@
 // ======================================================================
-// © 2021. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract
 // 89233218CNA000001 for Los Alamos National Laboratory (LANL), which
 // is operated by Triad National Security, LLC for the U.S.
@@ -30,6 +30,8 @@ namespace neutrinos {
 template <typename ThermalDistribution, typename pc = PhysicalConstantsCGS>
 class GrayOpacity {
  public:
+  using PC = pc;
+
   GrayOpacity() = default;
   GrayOpacity(const Real kappa) : kappa_(kappa) {}
   GrayOpacity(const ThermalDistribution &dist, const Real kappa)

singularity-opac/neutrinos/mean_neutrino_variant.hpp

@@ -1,5 +1,5 @@
 // ======================================================================
-// © 2021. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract
 // 89233218CNA000001 for Los Alamos National Laboratory (LANL), which
 // is operated by Triad National Security, LLC for the U.S.
@@ -69,6 +69,16 @@ class MeanVariant {
         opac_);
   }
 
+  PORTABLE_INLINE_FUNCTION RuntimePhysicalConstants
+  GetRuntimePhysicalConstants() const {
+    return mpark::visit(
+        [=](const auto &opac) {
+          using PC = typename std::decay_t<decltype(opac)>::PC;
+          return singularity::GetRuntimePhysicalConstants(PC());
+        },
+        opac_);
+  }
+
   PORTABLE_INLINE_FUNCTION Real PlanckMeanAbsorptionCoefficient(
       const Real rho, const Real temp, const Real Ye,
       const RadiationType type) const {
@@ -91,6 +101,12 @@ class MeanVariant {
   inline void Finalize() noexcept {
     return mpark::visit([](auto &opac) { return opac.Finalize(); }, opac_);
   }
+
+#ifdef SPINER_USE_HDF
+  void Save(const std::string &filename) const {
+    return mpark::visit([=](auto &opac) { return opac.Save(filename); }, opac_);
+  }
+#endif
 };
 
 } // namespace impl

singularity-opac/neutrinos/mean_opacity_neutrinos.hpp

@@ -1,5 +1,5 @@
 // ======================================================================
-// © 2022. Triad National Security, LLC. All rights reserved.  This
+// © 2022-2024. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract
 // 89233218CNA000001 for Los Alamos National Laboratory (LANL), which
 // is operated by Triad National Security, LLC for the U.S.
@@ -37,7 +37,6 @@ namespace impl {
 // TODO(BRR) Note: It is assumed that lambda is constant for all densities,
 // temperatures, and Ye
 
-template <typename pc = PhysicalConstantsCGS>
 class MeanOpacity {
 
  public:
@@ -134,6 +133,8 @@ class MeanOpacity {
                         const Real lTMax, const int NT, const Real YeMin,
                         const Real YeMax, const int NYe, Real lNuMin,
                         Real lNuMax, const int NNu, Real *lambda = nullptr) {
+    using PC = typename Opacity::PC;
+
     lkappaPlanck_.resize(NRho, NT, NYe, NEUTRINO_NTYPES);
     // index 0 is the species and is not interpolatable
     lkappaPlanck_.setRange(1, YeMin, YeMax, NYe);
@@ -159,8 +160,8 @@ class MeanOpacity {
             // Choose default temperature-specific frequency grid if frequency
             // grid not specified
             if (AUTOFREQ) {
-              lNuMin = toLog_(1.e-3 * pc::kb * fromLog_(lTMin) / pc::h);
-              lNuMax = toLog_(1.e3 * pc::kb * fromLog_(lTMax) / pc::h);
+              lNuMin = toLog_(1.e-3 * PC::kb * fromLog_(lTMin) / PC::h);
+              lNuMax = toLog_(1.e3 * PC::kb * fromLog_(lTMax) / PC::h);
             }
             const Real dlnu = (lNuMax - lNuMin) / (NNu - 1);
             // Integrate over frequency
@@ -219,10 +220,9 @@ class MeanOpacity {
 
 } // namespace impl
 
-using MeanOpacityScaleFree = impl::MeanOpacity<PhysicalConstantsUnity>;
-using MeanOpacityCGS = impl::MeanOpacity<PhysicalConstantsCGS>;
-using MeanOpacity = impl::MeanVariant<MeanOpacityScaleFree, MeanOpacityCGS,
-                                      MeanNonCGSUnits<MeanOpacityCGS>>;
+using MeanOpacityBase = impl::MeanOpacity;
+using MeanOpacity =
+    impl::MeanVariant<MeanOpacityBase, MeanNonCGSUnits<MeanOpacityBase>>;
 
 } // namespace neutrinos
 } // namespace singularity

singularity-opac/neutrinos/neutrino_variant.hpp

@@ -1,5 +1,5 @@
 // ======================================================================
-// © 2021. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract
 // 89233218CNA000001 for Los Alamos National Laboratory (LANL), which
 // is operated by Triad National Security, LLC for the U.S.
@@ -71,6 +71,16 @@ class Variant {
         opac_);
   }
 
+  PORTABLE_INLINE_FUNCTION RuntimePhysicalConstants
+  GetRuntimePhysicalConstants() const {
+    return mpark::visit(
+        [=](const auto &opac) {
+          using PC = typename std::decay_t<decltype(opac)>::PC;
+          return singularity::GetRuntimePhysicalConstants(PC());
+        },
+        opac_);
+  }
+
   // Directional absorption coefficient with units of 1/length
   // Signature should be at least
   // rho, temp, Ye, type, nu, lambda

singularity-opac/neutrinos/non_cgs_neutrinos.hpp

@@ -1,5 +1,5 @@
 // ======================================================================
-// © 2021. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract
 // 89233218CNA000001 for Los Alamos National Laboratory (LANL), which
 // is operated by Triad National Security, LLC for the U.S.
@@ -29,6 +29,8 @@ namespace neutrinos {
 template <typename Opac>
 class NonCGSUnits {
  public:
+  using PC = typename Opac::PC;
+
   NonCGSUnits() = default;
   NonCGSUnits(Opac &&opac, const Real time_unit, const Real mass_unit,
               const Real length_unit, const Real temp_unit)
@@ -47,6 +49,7 @@ class NonCGSUnits {
     return NonCGSUnits<Opac>(opac_.GetOnDevice(), time_unit_, mass_unit_,
                              length_unit_, temp_unit_);
   }
+
   inline void Finalize() noexcept { opac_.Finalize(); }
 
   PORTABLE_INLINE_FUNCTION
@@ -66,10 +69,11 @@ class NonCGSUnits {
   }
 
   template <typename FrequencyIndexer, typename DataIndexer>
-  PORTABLE_INLINE_FUNCTION void AbsorptionCoefficient(
-      const Real rho, const Real temp, const Real Ye, RadiationType type,
-      FrequencyIndexer &nu_bins, DataIndexer &coeffs, const int nbins,
-      Real *lambda = nullptr) const {
+  PORTABLE_INLINE_FUNCTION void
+  AbsorptionCoefficient(const Real rho, const Real temp, const Real Ye,
+                        RadiationType type, FrequencyIndexer &nu_bins,
+                        DataIndexer &coeffs, const int nbins,
+                        Real *lambda = nullptr) const {
     for (int i = 0; i < nbins; ++i) {
       nu_bins[i] *= freq_unit_;
     }
@@ -262,6 +266,12 @@ class MeanNonCGSUnits {
   PORTABLE_INLINE_FUNCTION
   int nlambda() const noexcept { return mean_opac_.nlambda(); }
 
+#ifdef SPINER_USE_HDF
+  void Save(const std::string &filename) const {
+    return mean_opac_.Save(filename);
+  }
+#endif
+
   PORTABLE_INLINE_FUNCTION
   Real PlanckMeanAbsorptionCoefficient(const Real rho, const Real temp,
                                        const Real Ye,

singularity-opac/neutrinos/spiner_opac_neutrinos.hpp

@@ -1,5 +1,5 @@
 // ======================================================================
-// © 2021. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract
 // 89233218CNA000001 for Los Alamos National Laboratory (LANL), which
 // is operated by Triad National Security, LLC for the U.S.
@@ -63,7 +63,9 @@ enum class DataStatus { Deallocated, OnDevice, OnHost };
 template <typename ThermalDistribution, typename pc = PhysicalConstantsCGS>
 class SpinerOpacity {
  public:
+  using PC = pc;
   using DataBox = Spiner::DataBox<Real>;
+
   static constexpr Real MEV = 1e6 * pc::eV;
   static constexpr Real EPS = 10.0 * std::numeric_limits<Real>::min();
   static constexpr Real Hz2MeV = pc::h / (1e6 * pc::eV);

singularity-opac/neutrinos/tophat_emissivity_neutrinos.hpp

@@ -1,5 +1,5 @@
 // ======================================================================
-// © 2021. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract
 // 89233218CNA000001 for Los Alamos National Laboratory (LANL), which
 // is operated by Triad National Security, LLC for the U.S.
@@ -32,6 +32,8 @@ namespace neutrinos {
 template <typename ThermalDistribution, typename pc = PhysicalConstantsCGS>
 class TophatEmissivity {
  public:
+  using PC = pc;
+
   TophatEmissivity(const Real C, const Real numin, const Real numax)
       : C_(C), numin_(numin), numax_(numax) {}
   TophatEmissivity(const ThermalDistribution &dist, const Real C,
@@ -46,6 +48,7 @@ class TophatEmissivity {
     printf("Tophat emissivity. C, numin, numax = %g, %g, %g\n", C_, numin_,
            numax_);
   }
+
   inline void Finalize() noexcept {}
 
   PORTABLE_INLINE_FUNCTION