Skip to content

Commit 2162d55

Browse files
haykhhaykh
committed
replenish injectors
1 parent 2e45b5c commit 2162d55

8 files changed

Lines changed: 343 additions & 3 deletions

File tree

‎pgens/examples/README.md‎

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -25,3 +25,5 @@ Problem generators in this directory are just examples demonstrating how to use
2525
- `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)
2626

2727
<img width="1024" alt="atmosphere" src="https://github.com/user-attachments/assets/70ac2cff-6775-47a3-a394-bf6df2533344" />
28+
29+
- `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`

‎pgens/examples/atmosphere/atmosphere.toml‎

Lines changed: 0 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -62,6 +62,3 @@
6262

6363
[checkpoint]
6464
keep = 0
65-
66-
[diagnostics]
67-
colored_stdout = true

‎pgens/examples/replenish_injector/nonuniform_replenish.mp4‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,3 @@
1+
version https://git-lfs.github.com/spec/v1
2+
oid sha256:ec180a8026198809970dd325a07e332a431019ceae357495c702e200674b4d6a
3+
size 15809969

‎pgens/examples/replenish_injector/nonuniform_replenish.toml‎

Lines changed: 57 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,57 @@
1+
[simulation]
2+
name = "nonuniform_replenish"
3+
engine = "srpic"
4+
runtime = 2.0
5+
6+
[grid]
7+
resolution = [128, 128]
8+
extent = [[0.0, 1.0], [0.0, 1.0]]
9+
10+
[grid.metric]
11+
metric = "minkowski"
12+
13+
[grid.boundaries]
14+
fields = [["PERIODIC"], ["PERIODIC"]]
15+
particles = [["ABSORB", "ABSORB"], ["PERIODIC"]]
16+
17+
[scales]
18+
larmor0 = 0.01
19+
skindepth0 = 0.1
20+
21+
[algorithms]
22+
23+
[algorithms.timestep]
24+
CFL = 0.5
25+
26+
[particles]
27+
ppc0 = 8.0
28+
29+
[[particles.species]]
30+
label = "e+"
31+
mass = 1.0
32+
charge = 1.0
33+
maxnpart = 1e6
34+
35+
[[particles.species]]
36+
label = "e-"
37+
mass = 1.0
38+
charge = -1.0
39+
maxnpart = 1e6
40+
41+
[setup]
42+
target_density = "nonuniform"
43+
44+
[output]
45+
interval = 10
46+
47+
[output.fields]
48+
quantities = ["N_1_2"]
49+
50+
[output.particles]
51+
enable = false
52+
53+
[output.spectra]
54+
enable = false
55+
56+
[checkpoint]
57+
keep = 0

‎pgens/examples/replenish_injector/pgen.hpp‎

Lines changed: 155 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,155 @@
1+
#ifndef PROBLEM_GENERATOR_H
2+
#define PROBLEM_GENERATOR_H
3+
4+
#include "enums.h"
5+
#include "global.h"
6+
7+
#include "archetypes/energy_dist.h"
8+
#include "archetypes/particle_injector.h"
9+
#include "archetypes/problem_generator.h"
10+
#include "archetypes/spatial_dist.h"
11+
#include "archetypes/traits.h"
12+
#include "framework/domain/metadomain.h"
13+
#include "kernels/particle_moments.hpp"
14+
15+
namespace user {
16+
using namespace ntt;
17+
18+
template <Dimension D>
19+
struct EMFields {
20+
Inline auto ex2(const coord_t<D>&) const -> real_t {
21+
return 0.8;
22+
}
23+
24+
Inline auto bx3(const coord_t<D>&) const -> real_t {
25+
return ONE;
26+
}
27+
};
28+
29+
template <Dimension D>
30+
struct NonUniformTargetDensity {
31+
Inline auto operator()(const coord_t<D>& x_Ph) const -> real_t {
32+
// example of a non-uniform target density that peaks at the center of the domain
33+
real_t r2 { ZERO };
34+
for (auto d = 0u; d < D; ++d) {
35+
r2 += SQR(x_Ph[d] - 0.5);
36+
}
37+
return std::exp(
38+
-r2 / SQR(0.2)); // <-- characteristic width of the density profile
39+
}
40+
};
41+
42+
template <SimEngine::type S, class M>
43+
struct PGen : public arch::ProblemGenerator<S, M> {
44+
45+
static constexpr auto engines {
46+
arch::traits::pgen::compatible_with<SimEngine::SRPIC>::value
47+
};
48+
static constexpr auto metrics {
49+
arch::traits::pgen::compatible_with<Metric::Minkowski>::value
50+
};
51+
static constexpr auto dimensions {
52+
arch::traits::pgen::compatible_with<Dim::_2D, Dim::_3D>::value
53+
};
54+
55+
using arch::ProblemGenerator<S, M>::D;
56+
using arch::ProblemGenerator<S, M>::C;
57+
using arch::ProblemGenerator<S, M>::params;
58+
59+
EMFields<D> init_flds;
60+
61+
const std::string target_density;
62+
63+
inline PGen(const SimulationParams& p, const Metadomain<S, M>& metadomain)
64+
: arch::ProblemGenerator<S, M> { p }
65+
, target_density { params.template get<std::string>(
66+
"setup.target_density") } {}
67+
68+
void CustomPostStep(timestep_t step, simtime_t time, Domain<S, M>& domain) {
69+
if (step % 100u != 0u) {
70+
return;
71+
}
72+
// perform replenishment and injection every 100 timesteps
73+
74+
{ // compute density of species #1 and #2
75+
76+
// saves the density to domain.fields.buff(:,:,:,0)
77+
const auto ni2 = domain.mesh.n_active(in::x2);
78+
const auto inv_n0 = ONE / params.template get<real_t>("scales.n0");
79+
80+
auto scatter_buff = Kokkos::Experimental::create_scatter_view(
81+
domain.fields.buff);
82+
Kokkos::deep_copy(domain.fields.buff, ZERO);
83+
for (const auto sp : std::vector<spidx_t> { 1, 2 }) {
84+
const auto& prtl_spec = domain.species[sp - 1];
85+
Kokkos::parallel_for("ComputeDensity",
86+
prtl_spec.rangeActiveParticles(),
87+
kernel::ParticleMoments_kernel<S, M, FldsID::N, 3>(
88+
{},
89+
scatter_buff,
90+
0u,
91+
prtl_spec.i1,
92+
prtl_spec.i2,
93+
prtl_spec.i3,
94+
prtl_spec.dx1,
95+
prtl_spec.dx2,
96+
prtl_spec.dx3,
97+
prtl_spec.ux1,
98+
prtl_spec.ux2,
99+
prtl_spec.ux3,
100+
prtl_spec.phi,
101+
prtl_spec.weight,
102+
prtl_spec.tag,
103+
prtl_spec.mass(),
104+
prtl_spec.charge(),
105+
false,
106+
domain.mesh.metric,
107+
domain.mesh.flds_bc(),
108+
ni2,
109+
inv_n0,
110+
0u));
111+
}
112+
Kokkos::Experimental::contribute(domain.fields.buff, scatter_buff);
113+
}
114+
115+
const auto energy_dist = arch::Maxwellian<S, M>(
116+
domain.mesh.metric,
117+
domain.random_pool(),
118+
0.2); // <-- target temperature for injection
119+
if (target_density == "uniform") {
120+
// pass the computed density to the replenisher
121+
const auto replenish_sdist = arch::ReplenishUniform<S, M, 3>(
122+
domain.mesh.metric,
123+
domain.fields.buff,
124+
0u, // <-- index in buff where the density is stored
125+
ONE); // <-- target density for replenishment
126+
arch::InjectNonUniform<S, M, decltype(energy_dist), decltype(energy_dist), decltype(replenish_sdist)>(
127+
params,
128+
domain,
129+
{ 1, 2 },
130+
{ energy_dist, energy_dist },
131+
replenish_sdist,
132+
ONE);
133+
} else {
134+
const auto target_density_profile = NonUniformTargetDensity<D> {};
135+
const auto replenish_sdist =
136+
arch::Replenish<S, M, 3, decltype(target_density_profile)>(
137+
domain.mesh.metric,
138+
domain.fields.buff,
139+
0u, // <-- index in buff where the density is stored
140+
target_density_profile,
141+
ONE); // <-- target density for replenishment
142+
arch::InjectNonUniform<S, M, decltype(energy_dist), decltype(energy_dist), decltype(replenish_sdist)>(
143+
params,
144+
domain,
145+
{ 1, 2 },
146+
{ energy_dist, energy_dist },
147+
replenish_sdist,
148+
ONE);
149+
}
150+
}
151+
};
152+
153+
} // namespace user
154+
155+
#endif

‎pgens/examples/replenish_injector/replenish_injector.py‎

Lines changed: 66 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,66 @@
1+
import nt2
2+
import matplotlib.pyplot as plt
3+
4+
5+
def plot(t, data):
6+
plt.rcParams["text.usetex"] = True
7+
plt.rcParams["font.family"] = "serif"
8+
plt.rcParams["figure.dpi"] = 250
9+
ax = plt.gca()
10+
data.fields.N_1_2.sel(t=t, method="nearest").plot(
11+
vmin=0, vmax=1.5, ax=ax, cmap="inferno"
12+
)
13+
plt.gcf().axes[1].set_ylabel(r"$n_-+n_+$")
14+
15+
ax.annotate(
16+
"",
17+
xy=(0.08, 0.75),
18+
xycoords="axes fraction",
19+
xytext=(0.08, 0.55),
20+
textcoords="axes fraction",
21+
arrowprops=dict(arrowstyle="-|>", color="#e05c2e", lw=1.8),
22+
)
23+
ax.text(
24+
0.095,
25+
0.65,
26+
r"$E_y$",
27+
transform=ax.transAxes,
28+
color="#e05c2e",
29+
fontsize=12,
30+
va="center",
31+
)
32+
33+
circle = plt.Circle(
34+
(0.08, 0.35),
35+
0.035,
36+
transform=ax.transAxes,
37+
fill=False,
38+
color="#2e80e0",
39+
lw=1.8,
40+
clip_on=False,
41+
)
42+
ax.add_patch(circle)
43+
ax.plot(0.08, 0.35, ".", transform=ax.transAxes, color="#2e80e0", markersize=6)
44+
ax.text(
45+
0.125,
46+
0.35,
47+
r"$B_z$",
48+
transform=ax.transAxes,
49+
color="#2e80e0",
50+
fontsize=12,
51+
va="center",
52+
)
53+
ax.set(
54+
xlabel=r"$x$",
55+
ylabel=r"$y$",
56+
title=rf"$t={{{t:.2f}}}$",
57+
xlim=(0, 1),
58+
ylim=(0, 1),
59+
)
60+
61+
62+
data_uniform_replenish = nt2.Data("uniform_replenish")
63+
data_uniform_replenish.makeMovie(plot, framerate=10)
64+
65+
data_nonuniform_replenish = nt2.Data("nonuniform_replenish")
66+
data_nonuniform_replenish.makeMovie(plot, framerate=10)

‎pgens/examples/replenish_injector/uniform_replenish.mp4‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,3 @@
1+
version https://git-lfs.github.com/spec/v1
2+
oid sha256:72668cf43d686e4ee9c9085243c41dbcc9e5845f1c097d70ca7a8263d26548b3
3+
size 33877523

‎pgens/examples/replenish_injector/uniform_replenish.toml‎

Lines changed: 57 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,57 @@
1+
[simulation]
2+
name = "uniform_replenish"
3+
engine = "srpic"
4+
runtime = 2.0
5+
6+
[grid]
7+
resolution = [128, 128]
8+
extent = [[0.0, 1.0], [0.0, 1.0]]
9+
10+
[grid.metric]
11+
metric = "minkowski"
12+
13+
[grid.boundaries]
14+
fields = [["PERIODIC"], ["PERIODIC"]]
15+
particles = [["ABSORB", "ABSORB"], ["PERIODIC"]]
16+
17+
[scales]
18+
larmor0 = 0.01
19+
skindepth0 = 0.1
20+
21+
[algorithms]
22+
23+
[algorithms.timestep]
24+
CFL = 0.5
25+
26+
[particles]
27+
ppc0 = 8.0
28+
29+
[[particles.species]]
30+
label = "e+"
31+
mass = 1.0
32+
charge = 1.0
33+
maxnpart = 1e6
34+
35+
[[particles.species]]
36+
label = "e-"
37+
mass = 1.0
38+
charge = -1.0
39+
maxnpart = 1e6
40+
41+
[setup]
42+
target_density = "uniform"
43+
44+
[output]
45+
interval = 10
46+
47+
[output.fields]
48+
quantities = ["N_1_2"]
49+
50+
[output.particles]
51+
enable = false
52+
53+
[output.spectra]
54+
enable = false
55+
56+
[checkpoint]
57+
keep = 0

0 commit comments

Comments
 (0)