Add a powerlaw injection to the archetypes

README.md

@@ -12,6 +12,8 @@ Our [detailed documentation](https://entity-toolkit.github.io/) includes everyth
 
 ## Core developers (alphabetical)
 
+👀 __Yangyang Cai__ {[@StaticObserver](https://github.com/StaticObserver): GRPIC}
+
 💁‍♂️ __Alexander Chernoglazov__ {[@SChernoglazov](https://github.com/SChernoglazov): PIC}
 
 🍵 __Benjamin Crinquand__ {[@bcrinquand](https://github.com/bcrinquand): GRPIC, cubed-sphere}

dev/runners/Dockerfile.runner.cuda

@@ -59,14 +59,17 @@ ARG HOME=/home/$USER
 WORKDIR $HOME
 
 # gh runner
+ARG RUNNER_VERSION=2.317.0
 RUN mkdir actions-runner
 WORKDIR $HOME/actions-runner
 
-RUN --mount=type=secret,id=ghtoken \
-    curl -o actions-runner-linux-x64-2.317.0.tar.gz \
-         -L https://github.com/actions/runner/releases/download/v2.317.0/actions-runner-linux-x64-2.317.0.tar.gz && \
-    tar xzf ./actions-runner-linux-x64-2.317.0.tar.gz && \
-    sudo ./bin/installdependencies.sh && \
-    ./config.sh --url https://github.com/entity-toolkit/entity --token "$(sudo cat /run/secrets/ghtoken)" --labels nvidia-gpu
+RUN curl -o actions-runner-linux-x64-${RUNNER_VERSION}.tar.gz \
+         -L https://github.com/actions/runner/releases/download/v${RUNNER_VERSION}/actions-runner-linux-x64-${RUNNER_VERSION}.tar.gz && \
+    tar xzf ./actions-runner-linux-x64-${RUNNER_VERSION}.tar.gz && \
+    sudo ./bin/installdependencies.sh
 
-ENTRYPOINT ["./run.sh"]
+ADD start.sh start.sh
+RUN sudo chown $USER:$USER start.sh && \
+    sudo chmod +x start.sh
+
+ENTRYPOINT ["./start.sh"]

legacy/src/framework/utils/particle_injectors.hpp

@@ -165,7 +165,7 @@ namespace ntt {
    * @brief Volumetrically uniform particle injector parallelized over particles.
    * @tparam D dimension.
    * @tparam S simulation engine.
-   * @tparam EnDist energy distribution [default = ColdDist].
+   * @tparam EnDist energy distribution [default = Cold].
    *
    * @param params simulation parameters.
    * @param mblock meshblock.
@@ -174,7 +174,7 @@ namespace ntt {
    * @param region region to inject particles as a list of coordinates [optional].
    * @param time current time [optional].
    */
-  template <Dimension D, SimulationEngine S, template <Dimension, SimulationEngine> class EnDist = ColdDist>
+  template <Dimension D, SimulationEngine S, template <Dimension, SimulationEngine> class EnDist = Cold>
   inline void InjectUniform(const SimulationParams& params,
                             Meshblock<D, S>&        mblock,
                             const std::vector<int>& species,
@@ -613,8 +613,8 @@ namespace ntt {
    * @brief Particle injector parallelized by cells in a volume.
    * @tparam D dimension.
    * @tparam S simulation engine.
-   * @tparam EnDist energy distribution [default = ColdDist].
-   * @tparam SpDist spatial distribution [default = UniformDist].
+   * @tparam EnDist energy distribution [default = Cold].
+   * @tparam SpDist spatial distribution [default = Uniform].
    * @tparam InjCrit injection criterion [default = NoCriterion].
    *
    * @param params simulation parameters.
@@ -626,8 +626,8 @@ namespace ntt {
    */
   template <Dimension        D,
             SimulationEngine S,
-            template <Dimension, SimulationEngine> class EnDist  = ColdDist,
-            template <Dimension, SimulationEngine> class SpDist  = UniformDist,
+            template <Dimension, SimulationEngine> class EnDist  = Cold,
+            template <Dimension, SimulationEngine> class SpDist  = Uniform,
             template <Dimension, SimulationEngine> class InjCrit = NoCriterion>
   inline void InjectInVolume(const SimulationParams& params,
                              Meshblock<D, S>&        mblock,
@@ -928,7 +928,7 @@ namespace ntt {
    * @brief ... up to certain number density.
    * @tparam D dimension.
    * @tparam S simulation engine.
-   * @tparam EnDist energy distribution [default = ColdDist].
+   * @tparam EnDist energy distribution [default = Cold].
    * @tparam InjCrit injection criterion [default = NoCriterion].
    *
    * @param params simulation parameters.
@@ -940,7 +940,7 @@ namespace ntt {
    */
   template <Dimension        D,
             SimulationEngine S,
-            template <Dimension, SimulationEngine> class EnDist  = ColdDist,
+            template <Dimension, SimulationEngine> class EnDist  = Cold,
             template <Dimension, SimulationEngine> class InjCrit = NoCriterion>
   inline void InjectNonUniform(const SimulationParams&      params,
                                Meshblock<D, S>&             mblock,

setups/srpic/turbulence/pgen.hpp

@@ -185,30 +185,25 @@ namespace user {
         const auto injector = arch::UniformInjector<S, M, arch::Maxwellian>(
           energy_dist,
           { 1, 2 });
-        const real_t ndens = 1.0;
+        const real_t ndens = 0.9;
         arch::InjectUniform<S, M, decltype(injector)>(params,
                                                       local_domain,
                                                       injector,
                                                       ndens);
       }
 
       {
-        // const auto energy_dist = arch::Maxwellian<S, M>(local_domain.mesh.metric,
-        //                                                 local_domain.random_pool,
-        //                                                 temperature*10);        
-        // // const auto energy_dist = arch::Maxwellian<S, M>(local_domain.mesh.metric,
-        // //                                                 local_domain.random_pool,
-        // //                                                 temperature * 2,
-        // //                                                 10.0,
-        // //                                                 1);
-        // const auto injector = arch::UniformInjector<S, M, arch::Maxwellian>(
-        //   energy_dist,
-        //   { 1, 2 });
-        // const real_t ndens = 0.01;
-        // arch::InjectUniform<S, M, decltype(injector)>(params,
-        //                                               local_domain,
-        //                                               injector,
-        //                                               ndens);
+        const auto energy_dist = arch::PowerlawDist<S, M>(local_domain.mesh.metric,
+                                                        local_domain.random_pool,
+                                                        0.1, 100.0, -3.0);
+        const auto injector = arch::UniformInjector<S, M, arch::PowerlawDist>(
+          energy_dist,
+          { 1, 2 });
+        const real_t ndens = 0.1;
+        arch::InjectUniform<S, M, decltype(injector)>(params,
+                                                      local_domain,
+                                                      injector,
+                                                      ndens);
       }
     }
 

src/archetypes/energy_dist.h

@@ -3,7 +3,8 @@
  * @brief Defines an archetype for energy distributions
  * @implements
  *   - arch::EnergyDistribution<>
- *   - arch::ColdDist<> : arch::EnergyDistribution<>
+ *   - arch::Cold<> : arch::EnergyDistribution<>
+ *   - arch::Powerlaw<> : arch::EnergyDistribution<>
  *   - arch::Maxwellian<> : arch::EnergyDistribution<>
  * @namespaces:
  *   - arch::
@@ -55,8 +56,8 @@ namespace arch {
   };
 
   template <SimEngine::type S, class M>
-  struct ColdDist : public EnergyDistribution<S, M> {
-    ColdDist(const M& metric) : EnergyDistribution<S, M> { metric } {}
+  struct Cold : public EnergyDistribution<S, M> {
+    Cold(const M& metric) : EnergyDistribution<S, M> { metric } {}
 
     Inline void operator()(const coord_t<M::Dim>&,
                            vec_t<Dim::_3D>& v,
@@ -67,6 +68,63 @@ namespace arch {
     }
   };
 
+  template <SimEngine::type S, class M>
+  struct Powerlaw : public EnergyDistribution<S, M> {
+    using EnergyDistribution<S, M>::metric;
+
+    Powerlaw(const M&              metric,
+             random_number_pool_t& pool,
+             real_t                g_min,
+             real_t                g_max,
+             real_t                pl_ind)
+      : EnergyDistribution<S, M> { metric }
+      , pool { pool }
+      , g_min { g_min }
+      , g_max { g_max }
+      , pl_ind { pl_ind } {}
+
+    Inline void operator()(const coord_t<M::Dim>& x_Code,
+                           vec_t<Dim::_3D>&       v,
+                           unsigned short         sp = 0) const override {
+      auto rand_gen = pool.get_state();
+      auto rand_X1  = Random<real_t>(rand_gen);
+      auto rand_gam = ONE;
+
+      // Power-law distribution from uniform (see https://mathworld.wolfram.com/RandomNumber.html)
+      if (pl_ind != -ONE) {
+        rand_gam += math::pow(
+          math::pow(g_min, ONE + pl_ind) +
+            (-math::pow(g_min, ONE + pl_ind) + math::pow(g_max, ONE + pl_ind)) *
+              rand_X1,
+          ONE / (ONE + pl_ind));
+      } else {
+        rand_gam += math::pow(g_min, ONE - rand_X1) * math::pow(g_max, rand_X1);
+      }
+      auto rand_u  = math::sqrt(SQR(rand_gam) - ONE);
+      auto rand_X2 = Random<real_t>(rand_gen);
+      auto rand_X3 = Random<real_t>(rand_gen);
+      v[0]         = rand_u * (TWO * rand_X2 - ONE);
+      v[2]         = TWO * rand_u * math::sqrt(rand_X2 * (ONE - rand_X2));
+      v[1]         = v[2] * math::cos(constant::TWO_PI * rand_X3);
+      v[2]         = v[2] * math::sin(constant::TWO_PI * rand_X3);
+
+      if constexpr (S == SimEngine::GRPIC) {
+        // convert from the tetrad basis to covariant
+        vec_t<Dim::_3D> v_Hat;
+        v_Hat[0] = v[0];
+        v_Hat[1] = v[1];
+        v_Hat[2] = v[2];
+        metric.template transform<Idx::T, Idx::D>(x_Code, v_Hat, v);
+      }
+
+      pool.free_state(rand_gen);
+    }
+
+  private:
+    const real_t         g_min, g_max, pl_ind;
+    random_number_pool_t pool;
+  };
+
   template <SimEngine::type S, class M>
   struct Maxwellian : public EnergyDistribution<S, M> {
     using EnergyDistribution<S, M>::metric;

src/archetypes/particle_injector.h

@@ -142,7 +142,7 @@ namespace arch {
     };
 
     using energy_dist_t  = Maxwellian<S, M>;
-    using spatial_dist_t = ReplenishDist<S, M, TargetDensityProfile>;
+    using spatial_dist_t = Replenish<S, M, TargetDensityProfile>;
     static_assert(M::is_metric, "M must be a metric class");
     static constexpr bool      is_nonuniform_injector { true };
     static constexpr Dimension D { M::Dim };

src/archetypes/spatial_dist.h

@@ -3,8 +3,8 @@
  * @brief Spatial distribution class passed to injectors
  * @implements
  *   - arch::SpatialDistribution<>
- *   - arch::UniformDist<> : arch::SpatialDistribution<>
- *   - arch::ReplenishDist<> : arch::SpatialDistribution<>
+ *   - arch::Uniform<> : arch::SpatialDistribution<>
+ *   - arch::Replenish<> : arch::SpatialDistribution<>
  * @namespace
  *   - arch::
  * @note
@@ -41,28 +41,28 @@ namespace arch {
   };
 
   template <SimEngine::type S, class M>
-  struct UniformDist : public SpatialDistribution<S, M> {
-    UniformDist(const M& metric) : SpatialDistribution<S, M> { metric } {}
+  struct Uniform : public SpatialDistribution<S, M> {
+    Uniform(const M& metric) : SpatialDistribution<S, M> { metric } {}
 
     Inline auto operator()(const coord_t<M::Dim>&) const -> real_t override {
       return ONE;
     }
   };
 
   template <SimEngine::type S, class M, class T>
-  struct ReplenishDist : public SpatialDistribution<S, M> {
+  struct Replenish : public SpatialDistribution<S, M> {
     using SpatialDistribution<S, M>::metric;
     const ndfield_t<M::Dim, 6> density;
     const unsigned short       idx;
 
     const T      target_density;
     const real_t target_max_density;
 
-    ReplenishDist(const M&                    metric,
-                  const ndfield_t<M::Dim, 6>& density,
-                  unsigned short              idx,
-                  const T&                    target_density,
-                  real_t                      target_max_density)
+    Replenish(const M&                    metric,
+              const ndfield_t<M::Dim, 6>& density,
+              unsigned short              idx,
+              const T&                    target_density,
+              real_t                      target_max_density)
       : SpatialDistribution<S, M> { metric }
       , density { density }
       , idx { idx }

src/archetypes/tests/CMakeLists.txt

@@ -22,4 +22,5 @@ endfunction()
 
 gen_test(energy_dist)
 gen_test(spatial_dist)
-gen_test(field_setter)
+gen_test(field_setter)
+gen_test(powerlaw)