Moving window

CLAUDE.md

@@ -49,6 +49,7 @@ entity
 │   └── entity.cpp               #   main entry-point
 ├── tests                        # unit tests for all submodules
 ├── .clang-format                # code formatting guidelines for clang-format
+├── .clang-tidy                  # configurations for the clang-tidy
 ├── .gitattributes
 ├── .gitignore
 ├── .gitmodules
@@ -59,9 +60,7 @@ entity
 ├── CODE_OF_CONDUCT.md
 ├── LICENSE
 ├── README.md
-├── conda-entity-nompi.sh
 ├── dependencies.py              # deployment scripts on various machines
-├── docker-compose.yml
 └── input.example.toml           # most complete toml file with all possible input options
 ```
 

pgens/examples/README.md

@@ -6,11 +6,11 @@ Problem generators in this directory are just examples demonstrating how to use
 
   <img width="1024" alt="image" src="https://github.com/user-attachments/assets/d8499ab7-31dc-49e8-8381-34aa8a994be0" />
 
-- `custom_spatial_distribution`: example of using the non-uniform plasma injector with a custom spatial distirbution
+- `custom_spatial_distribution`: example of using the non-uniform plasma injector with a custom spatial distirbution.
 
   https://github.com/user-attachments/assets/bea0c290-e7e4-4ec7-b360-68ce76beab5b
 
-- `match_fix_field_boundaries`: example of setting matching and/or fixed (coordinate-independent) field boundaries for the electromagnetic fields
+- `match_fix_field_boundaries`: example of setting matching and/or fixed (coordinate-independent) field boundaries for the electromagnetic fields.
 
   https://github.com/user-attachments/assets/a1b9ea22-34ce-474f-a9b0-49789a2e52b3
 
@@ -22,13 +22,23 @@ Problem generators in this directory are just examples demonstrating how to use
 
   <img width="512" alt="image" src="https://github.com/user-attachments/assets/7126fcdb-7484-4695-86c1-cdd425cc4655" />
 
-- `atmosphere`: setting up a gravitationally bound "atmosphere" with a constant particle replenisher (the plot also highlights the importance of having a constant replenisher and a gravity force acting on the particles; both are enabled by default when using the `ATMOSPHERE` particle boundary conditions)
+- `atmosphere`: setting up a gravitationally bound "atmosphere" with a constant particle replenisher (the plot also highlights the importance of having a constant replenisher and a gravity force acting on the particles; both are enabled by default when using the `ATMOSPHERE` particle boundary conditions).
 
   <img width="1024" alt="atmosphere" src="https://github.com/user-attachments/assets/70ac2cff-6775-47a3-a394-bf6df2533344" />
 
-- `replenish_injector`: demonstration of how to use the replenish injector to periodically inject new plasma to a target density (both uniform and non-uniform) in `CustomPostStep`
+- `replenish_injector`: demonstration of how to use the replenish injector to periodically inject new plasma to a target density (both uniform and non-uniform) in `CustomPostStep`.
 
   https://github.com/user-attachments/assets/4955f05e-4795-439f-a4a5-8196e42e987b
 
   https://github.com/user-attachments/assets/ebc303ed-7b21-4f97-9822-d81124fdf962
 
+- `piston`: implementing a moving piston which continuously pushes the plasma.
+
+  https://github.com/user-attachments/assets/7d2dfc39-5de4-49cb-9d55-be8700cbcaee
+
+- `moving_window`: demonstrates the domain being updated to follow certain features in your simulation in one of the direction (jerkiness is due to updates taking place every `n`-th timestep, if you limit plotting to a specific physical interval, this goes away).
+
+  https://github.com/user-attachments/assets/00edb648-af80-4d37-b1c2-c49872227d47
+
+
+

pgens/examples/moving_window/moving_window.py

@@ -0,0 +1,40 @@
+import nt2
+import matplotlib.pyplot as plt
+
+
+def plot(t, data):
+    plt.rcParams["text.usetex"] = True
+    plt.rcParams["font.family"] = "serif"
+    plt.rcParams["figure.dpi"] = 250
+    ax = plt.gca()
+    data.fields.N_1_2.sel(t=t, method="nearest").plot(
+        vmin=0, vmax=1.5, ax=ax, cmap="inferno"
+    )
+    plt.gcf().axes[1].set_ylabel(r"$n_-+n_+$")
+
+    ax.annotate(
+        "",
+        xy=(0.28, 0.08),
+        xycoords="axes fraction",
+        xytext=(0.08, 0.08),
+        textcoords="axes fraction",
+        arrowprops=dict(arrowstyle="-|>", color="#e05c2e", lw=1.8),
+    )
+    ax.text(
+        0.095,
+        0.095,
+        r"$v_0$",
+        transform=ax.transAxes,
+        color="#e05c2e",
+        fontsize=12,
+        va="center",
+    )
+    ax.set(
+        xlabel=r"$x$",
+        ylabel=r"$y$",
+        title=rf"$t={{{t:.2f}}}$",
+    )
+
+
+moving_window = nt2.Data("moving_window")
+moving_window.makeMovie(plot, framerate=10)

pgens/examples/moving_window/moving_window.toml

@@ -0,0 +1,57 @@
+[simulation]
+  name    = "moving_window"
+  engine  = "srpic"
+  runtime = 2.0
+
+[grid]
+  resolution = [128, 128]
+  extent     = [[0.0, 1.0], [0.0, 1.0]]
+
+  [grid.metric]
+    metric = "minkowski"
+
+  [grid.boundaries]
+    fields    = [["PERIODIC"], ["PERIODIC"]]
+    particles = [["ABSORB", "ABSORB"], ["PERIODIC"]]
+
+[scales]
+  larmor0    = 0.01
+  skindepth0 = 0.1
+
+[algorithms]
+
+  [algorithms.timestep]
+    CFL = 0.5
+
+[particles]
+  ppc0 = 8.0
+
+  [[particles.species]]
+    label    = "e+"
+    mass     = 1.0
+    charge   = 1.0
+    maxnpart = 1e6
+
+  [[particles.species]]
+    label    = "e-"
+    mass     = 1.0
+    charge   = -1.0
+    maxnpart = 1e6
+
+[setup]
+  window_velocity = 0.2
+
+[output]
+  interval = 10
+
+  [output.fields]
+    quantities = ["N_1_2"]
+
+  [output.particles]
+    enable = false
+
+  [output.spectra]
+    enable = false
+
+[checkpoint]
+  keep = 0

pgens/examples/moving_window/pgen.hpp

@@ -0,0 +1,114 @@
+#ifndef PROBLEM_GENERATOR_H
+#define PROBLEM_GENERATOR_H
+
+#include "enums.h"
+#include "global.h"
+
+#include "arch/kokkos_aliases.h"
+#include "traits/pgen.h"
+
+#include "archetypes/energy_dist.h"
+#include "archetypes/moving_window.h"
+#include "archetypes/particle_injector.h"
+#include "framework/domain/domain.h"
+#include "framework/domain/metadomain.h"
+#include "framework/parameters/parameters.h"
+
+namespace user {
+  using namespace ntt;
+
+  template <Dimension D>
+  struct NonUniformDensity {
+    Inline auto operator()(const coord_t<D>& x_Ph) const -> real_t {
+      // example of a non-uniform target density that peaks at the center of the domain
+      real_t r2 { ZERO };
+      for (auto d = 0u; d < D; ++d) {
+        r2 += SQR(x_Ph[d] - 0.5);
+      }
+      return math::exp(-r2 / SQR(static_cast<real_t>(0.1)));
+      //                                                ^
+      //                        characteristic width of the density profile
+    }
+  };
+
+  template <SimEngine::type S, class M>
+  struct PGen {
+    static constexpr auto engines {
+      ::traits::pgen::compatible_with<SimEngine::SRPIC> {}
+    };
+    static constexpr auto metrics {
+      ::traits::pgen::compatible_with<Metric::Minkowski> {}
+    };
+    static constexpr auto dimensions {
+      ::traits::pgen::compatible_with<Dim::_2D, Dim::_3D> {}
+    };
+    const SimulationParams& params;
+    Metadomain<S, M>&       metadomain;
+
+    struct MovingWindow {
+      int    pos_i { -1 };
+      int    init_pos_i { -1 };
+      real_t pos_di { ZERO };
+      real_t vel_Cd { ZERO };
+
+      void init(const M& metric, real_t global_x, real_t velocity) {
+        const auto pos_Cd = metric.template convert<1, Crd::Ph, Crd::Cd>(global_x);
+        pos_i      = static_cast<int>(pos_Cd + 1) - 1;
+        init_pos_i = pos_i;
+        pos_di     = pos_Cd - static_cast<real_t>(pos_i);
+        vel_Cd = metric.template transform<1, Idx::XYZ, Idx::U>({}, velocity);
+      }
+
+      void update(real_t            dt,
+                  int               shift,
+                  Metadomain<S, M>& metadomain,
+                  Domain<S, M>&     domain) {
+        pos_di += vel_Cd * dt;
+        pos_i  += static_cast<int>(pos_di >= ONE);
+        pos_di -= static_cast<real_t>(pos_di >= ONE);
+        if ((pos_i - init_pos_i) >= shift) {
+          // move the window and all fields and particles in it
+          arch::MoveWindow<S, M, in::x1>(domain, metadomain, shift);
+          // update window index for next update
+          pos_i -= shift;
+        }
+      }
+    };
+
+    MovingWindow moving_window;
+    const real_t dt;
+
+    PGen(const SimulationParams& p, Metadomain<S, M>& m)
+      : params { p }
+      , metadomain { m }
+      , dt { params.template get<real_t>("algorithms.timestep.dt") } {}
+
+    void InitPrtls(Domain<S, M>& domain) {
+      const auto window_velocity = params.template get<real_t>(
+        "setup.window_velocity",
+        ZERO);
+      moving_window.init(domain.mesh.metric, ZERO, window_velocity);
+      // inject particles
+      const auto energy_dist = arch::energy_dist::Maxwellian<M::Dim, M::CoordType>(
+        domain.random_pool(),
+        0.001,
+        { window_velocity, ZERO, ZERO });
+
+      const auto density_profile = NonUniformDensity<M::Dim> {};
+      arch::InjectNonUniform<S, M, decltype(energy_dist), decltype(energy_dist), decltype(density_profile)>(
+        params,
+        domain,
+        { 1, 2 },
+        { energy_dist, energy_dist },
+        density_profile,
+        ONE);
+    }
+
+    void CustomPostStep(timestep_t /*step*/, simtime_t /*time*/, Domain<S, M>& domain) {
+      moving_window.update(dt, N_GHOSTS, metadomain, domain);
+    }
+  };
+
+} // namespace user
+
+#endif

pgens/examples/piston/pgen.hpp

@@ -0,0 +1,109 @@
+#ifndef PROBLEM_GENERATOR_H
+#define PROBLEM_GENERATOR_H
+
+#include "enums.h"
+#include "global.h"
+
+#include "arch/kokkos_aliases.h"
+#include "traits/pgen.h"
+
+#include "archetypes/energy_dist.h"
+#include "archetypes/particle_injector.h"
+#include "archetypes/piston.h"
+#include "framework/domain/domain.h"
+#include "framework/domain/metadomain.h"
+#include "framework/parameters/parameters.h"
+
+namespace user {
+  using namespace ntt;
+
+  template <Dimension D>
+  struct NonUniformDensity {
+    Inline auto operator()(const coord_t<D>& x_Ph) const -> real_t {
+      // example of a non-uniform target density that peaks at the center of the domain
+      real_t r2 { ZERO };
+      for (auto d = 0u; d < D; ++d) {
+        r2 += SQR(x_Ph[d] - 0.5);
+      }
+      return math::exp(-r2 / SQR(static_cast<real_t>(0.2)));
+      //                                                ^
+      //                        characteristic width of the density profile
+    }
+  };
+
+  template <SimEngine::type S, class M>
+  struct PGen {
+
+    static constexpr auto engines {
+      ::traits::pgen::compatible_with<SimEngine::SRPIC> {}
+    };
+    static constexpr auto metrics {
+      ::traits::pgen::compatible_with<Metric::Minkowski> {}
+    };
+    static constexpr auto dimensions {
+      ::traits::pgen::compatible_with<Dim::_2D, Dim::_3D> {}
+    };
+
+    const SimulationParams& params;
+    Metadomain<S, M>&       metadomain;
+    const real_t            piston_velocity;
+
+    PGen(const SimulationParams& p, Metadomain<S, M>& m)
+      : params { p }
+      , metadomain { m }
+      , piston_velocity { params.template get<real_t>("setup.piston_velocity",
+                                                      0.0) } {}
+
+    void InitPrtls(Domain<S, M>& domain) {
+      // inject particles
+      const auto energy_dist = arch::energy_dist::Maxwellian<M::Dim, M::CoordType>(
+        domain.random_pool(),
+        0.2); // <-- target temperature for injection
+
+      const auto density_profile = NonUniformDensity<M::Dim> {};
+      arch::InjectNonUniform<S, M, decltype(energy_dist), decltype(energy_dist), decltype(density_profile)>(
+        params,
+        domain,
+        { 1, 2 },
+        { energy_dist, energy_dist },
+        density_profile,
+        ONE);
+    }
+
+    struct CustomPrtlUpdate {
+      real_t x_piston; // current position of piston
+      real_t v_piston;
+      real_t global_xmax;
+      bool   is_left;
+      bool   massive;
+
+      Inline void operator()(index_t                          p,
+                             const kernel::sr::PusherContext& ctx,
+                             const kernel::sr::PusherBoundaries<M::Dim>&,
+                             const kernel::PusherArrays& particles,
+                             const M&                    metric) const {
+
+        real_t piston_pos;
+        if (x_piston > global_xmax) {
+          // piston has moved beyond the domain, set position to a large value so that no particles interact with it
+          piston_pos = global_xmax;
+        } else {
+          piston_pos = x_piston;
+        }
+
+        arch::Piston<M>(p, ctx.dt, particles, metric, piston_pos, v_piston, massive);
+      }
+    };
+
+    template <class DOM>
+    auto CustomParticleUpdate(simtime_t time, spidx_t /*sp*/, DOM& /*domain*/) const {
+      return CustomPrtlUpdate { piston_velocity * static_cast<real_t>(time),
+                                piston_velocity,
+                                metadomain.mesh().extent(in::x1).second,
+                                true,
+                                true };
+    };
+  };
+} // namespace user
+
+#endif