Skip to content
Draft
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
580 changes: 580 additions & 0 deletions pgens/stripedwind_new/pgen.hpp

Large diffs are not rendered by default.

73 changes: 73 additions & 0 deletions pgens/stripedwind_new/stripedwind.toml
Original file line number Diff line number Diff line change
@@ -0,0 +1,73 @@
[simulation]
name = "stripedwind"
engine = "srpic"
runtime = 50.0

[simulation.domain]
decomposition = [1, -1]

[grid]
resolution = [4096, 128]
extent = [[0.0, 409.6], [0.0, 12.8]]

[grid.metric]
metric = "minkowski"

[grid.boundaries]
fields = [["CONDUCTOR", "MATCH"], ["PERIODIC"]]
particles = [["REFLECT", "ABSORB"], ["PERIODIC"]]


[scales]
larmor0 = 0.316227766017
skindepth0 = 1.0

[algorithms]
current_filters = 8

[algorithms.timestep]
CFL = 0.5

[particles]
ppc0 = 8.0

[[particles.species]]
label = "e-"
mass = 1.0
charge = -1.0
maxnpart = 8e7

[[particles.species]]
label = "e+"
mass = 1.0
charge = 1.0
maxnpart = 8e7

[setup]
drift_ux = 0.0 # speed towards the wall [c]
temperature = 0.0001 # temperature of maxwell distribution [kB T / (m_i c^2)]
temperature_ratio = 1.0 # temperature ratio of electrons to protons
Bmag = 1.0 # magnetic field strength as fraction of magnetisation
Btheta = 0.0 # magnetic field angle in the plane
Bphi = 0.0 # magnetic field angle out of plane
alpha_s = 0.1 # stripe averaged field (see Sironi & Spitkovsky 2011)
lambda_s = 20.0 # striped wind wavelength
delta_s = 1.0 # width of region where magnetic field flips
filling_fraction = 0.1 # fraction of the shock piston filled with plasma
injector_velocity = 0.99 # speed of injector [c]
injection_start = 0.0 # start time of moving injector
injector_interval = 1 # inject particles every 1 timestep

[output]
interval_time = 0.1
format = "BPFile"

[output.fields]
quantities = ["N_1", "N_2", "B", "E"]

[output.particles]
enable = true
stride = 10

[output.spectra]
enable = false
4 changes: 2 additions & 2 deletions src/archetypes/particle_injector.h
Original file line number Diff line number Diff line change
Expand Up @@ -254,7 +254,7 @@ namespace arch {
energy_dists.first,
energy_dists.second,
ONE / params.template get<real_t>("scales.V0"),
domain.random_pool()));
domain.random_pool));
domain.species[species.first - 1].set_npart(
domain.species[species.first - 1].npart() + nparticles);
domain.species[species.second - 1].set_npart(
Expand Down Expand Up @@ -372,7 +372,7 @@ namespace arch {
energy_dists.second,
spatial_dist,
ONE / params.template get<real_t>("scales.V0"),
domain.random_pool());
domain.random_pool);
Kokkos::parallel_for("InjectNonUniformNumberDensity",
cell_range,
injector_kernel);
Expand Down
4 changes: 2 additions & 2 deletions src/archetypes/utils.h
Original file line number Diff line number Diff line change
Expand Up @@ -53,11 +53,11 @@ namespace arch {
const auto temperature_2 = temperatures.second / mass_2;

const auto maxwellian_1 = arch::Maxwellian<S, M>(domain.mesh.metric,
domain.random_pool(),
domain.random_pool,
temperature_1,
drift_four_vels.first);
const auto maxwellian_2 = arch::Maxwellian<S, M>(domain.mesh.metric,
domain.random_pool(),
domain.random_pool,
temperature_2,
drift_four_vels.second);

Expand Down
21 changes: 15 additions & 6 deletions src/engines/engine_init.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -3,8 +3,14 @@

#include "arch/traits.h"

#include "metrics/kerr_schild.h"
#include "metrics/kerr_schild_0.h"
#include "metrics/minkowski.h"
#include "metrics/qkerr_schild.h"
#include "metrics/qspherical.h"
#include "metrics/spherical.h"

#include "archetypes/field_setter.h"
#include "framework/specialization_registry.h"

#include "engines/engine.hpp"

Expand Down Expand Up @@ -65,10 +71,13 @@ namespace ntt {
print_report();
}

#define ENGINE_INIT(S, M, D) template class Engine<S, M<D>>;

NTT_FOREACH_SPECIALIZATION(ENGINE_INIT)

#undef ENGINE_INIT
template class Engine<SimEngine::SRPIC, metric::Minkowski<Dim::_1D>>;
template class Engine<SimEngine::SRPIC, metric::Minkowski<Dim::_2D>>;
template class Engine<SimEngine::SRPIC, metric::Minkowski<Dim::_3D>>;
template class Engine<SimEngine::SRPIC, metric::Spherical<Dim::_2D>>;
template class Engine<SimEngine::SRPIC, metric::QSpherical<Dim::_2D>>;
template class Engine<SimEngine::GRPIC, metric::KerrSchild<Dim::_2D>>;
template class Engine<SimEngine::GRPIC, metric::KerrSchild0<Dim::_2D>>;
template class Engine<SimEngine::GRPIC, metric::QKerrSchild<Dim::_2D>>;

} // namespace ntt
27 changes: 14 additions & 13 deletions src/engines/engine_printer.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,12 @@
#include "utils/colors.h"
#include "utils/formatting.h"

#include "framework/specialization_registry.h"
#include "metrics/kerr_schild.h"
#include "metrics/kerr_schild_0.h"
#include "metrics/minkowski.h"
#include "metrics/qkerr_schild.h"
#include "metrics/qspherical.h"
#include "metrics/spherical.h"

#include "engines/engine.hpp"

Expand Down Expand Up @@ -301,12 +306,6 @@ namespace ntt {
add_param(report, 4, "Problem generator", "%s", pgen.c_str());
add_param(report, 4, "Engine", "%s", SimEngine(S).to_string());
add_param(report, 4, "Metric", "%s", Metric(M::MetricType).to_string());
#if SHAPE_ORDER == 0
add_param(report, 4, "Deposit", "%s", "zigzag");
#else
add_param(report, 4, "Deposit", "%s", "esirkepov");
add_param(report, 4, "Interpolation order", "%i", SHAPE_ORDER);
#endif
add_param(report, 4, "Timestep [dt]", "%.3e", dt);
add_param(report, 4, "Runtime", "%.3e [%d steps]", runtime, max_steps);
report += "\n";
Expand Down Expand Up @@ -490,11 +489,13 @@ namespace ntt {
}
}

#define ENGINE_PRINTER(S, M, D) \
template void Engine<S, M<D>>::print_report() const;

NTT_FOREACH_SPECIALIZATION(ENGINE_PRINTER)

#undef ENGINE_PRINTER
template void Engine<SimEngine::SRPIC, metric::Minkowski<Dim::_1D>>::print_report() const;
template void Engine<SimEngine::SRPIC, metric::Minkowski<Dim::_2D>>::print_report() const;
template void Engine<SimEngine::SRPIC, metric::Minkowski<Dim::_3D>>::print_report() const;
template void Engine<SimEngine::SRPIC, metric::Spherical<Dim::_2D>>::print_report() const;
template void Engine<SimEngine::SRPIC, metric::QSpherical<Dim::_2D>>::print_report() const;
template void Engine<SimEngine::GRPIC, metric::KerrSchild<Dim::_2D>>::print_report() const;
template void Engine<SimEngine::GRPIC, metric::KerrSchild0<Dim::_2D>>::print_report() const;
template void Engine<SimEngine::GRPIC, metric::QKerrSchild<Dim::_2D>>::print_report() const;

} // namespace ntt
21 changes: 15 additions & 6 deletions src/engines/engine_run.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -3,8 +3,14 @@
#include "arch/traits.h"
#include "utils/diag.h"

#include "metrics/kerr_schild.h"
#include "metrics/kerr_schild_0.h"
#include "metrics/minkowski.h"
#include "metrics/qkerr_schild.h"
#include "metrics/qspherical.h"
#include "metrics/spherical.h"

#include "framework/domain/domain.h"
#include "framework/specialization_registry.h"

#include "engines/engine.hpp"

Expand Down Expand Up @@ -137,10 +143,13 @@ namespace ntt {
}
}

#define ENGINE_RUN(S, M, D) template void Engine<S, M<D>>::run();

NTT_FOREACH_SPECIALIZATION(ENGINE_RUN)

#undef ENGINE_RUN
template void Engine<SimEngine::SRPIC, metric::Minkowski<Dim::_1D>>::run();
template void Engine<SimEngine::SRPIC, metric::Minkowski<Dim::_2D>>::run();
template void Engine<SimEngine::SRPIC, metric::Minkowski<Dim::_3D>>::run();
template void Engine<SimEngine::SRPIC, metric::Spherical<Dim::_2D>>::run();
template void Engine<SimEngine::SRPIC, metric::QSpherical<Dim::_2D>>::run();
template void Engine<SimEngine::GRPIC, metric::KerrSchild<Dim::_2D>>::run();
template void Engine<SimEngine::GRPIC, metric::KerrSchild0<Dim::_2D>>::run();
template void Engine<SimEngine::GRPIC, metric::QKerrSchild<Dim::_2D>>::run();

} // namespace ntt
30 changes: 0 additions & 30 deletions src/engines/engine_traits.h

This file was deleted.

4 changes: 2 additions & 2 deletions src/engines/grpic.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -94,9 +94,9 @@ namespace ntt {

void step_forward(timer::Timers& timers, domain_t& dom) override {
const auto fieldsolver_enabled = m_params.template get<bool>(
"algorithms.fieldsolver.enable");
"algorithms.toggles.fieldsolver");
const auto deposit_enabled = m_params.template get<bool>(
"algorithms.deposit.enable");
"algorithms.toggles.deposit");
const auto clear_interval = m_params.template get<std::size_t>(
"particles.clear_interval");

Expand Down
64 changes: 28 additions & 36 deletions src/engines/srpic.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -78,9 +78,9 @@ namespace ntt {

void step_forward(timer::Timers& timers, domain_t& dom) override {
const auto fieldsolver_enabled = m_params.template get<bool>(
"algorithms.fieldsolver.enable");
"algorithms.toggles.fieldsolver");
const auto deposit_enabled = m_params.template get<bool>(
"algorithms.deposit.enable");
"algorithms.toggles.deposit");
const auto clear_interval = m_params.template get<std::size_t>(
"particles.clear_interval");

Expand Down Expand Up @@ -348,7 +348,7 @@ namespace ntt {
: ZERO;
// cooling
const auto has_synchrotron = (cooling == Cooling::SYNCHROTRON);
const auto has_compton = (cooling == Cooling::COMPTON);
const auto has_compton = (cooling == Cooling::COMPTON);
const auto sync_grad = has_synchrotron
? m_params.template get<real_t>(
"algorithms.synchrotron.gamma_rad")
Expand All @@ -359,16 +359,17 @@ namespace ntt {
"scales.omegaB0") /
(SQR(sync_grad) * species.mass())
: ZERO;
const auto comp_grad = has_compton ? m_params.template get<real_t>(
"algorithms.compton.gamma_rad")
: ZERO;
const auto comp_coeff = has_compton ? (real_t)(0.1) * dt *
m_params.template get<real_t>(
"scales.omegaB0") /
(SQR(comp_grad) * species.mass())
: ZERO;
const auto comp_grad = has_compton
? m_params.template get<real_t>(
"algorithms.compton.gamma_rad")
: ZERO;
const auto comp_coeff = has_compton
? (real_t)(0.1) * dt *
m_params.template get<real_t>(
"scales.omegaB0") / (SQR(comp_grad) * species.mass())
: ZERO;
// toggle to indicate whether pgen defines the external force
bool has_extforce = false;
bool has_extforce = false;
if constexpr (traits::has_member<traits::pgen::ext_force_t, pgen_t>::value) {
has_extforce = true;
// toggle to indicate whether the ext force applies to current species
Expand Down Expand Up @@ -518,30 +519,9 @@ namespace ntt {
}
}

template <unsigned short O>
void deposit_with(const Particles<M::Dim, M::CoordType>& species,
const M& metric,
const scatter_ndfield_t<M::Dim, 3>& scatter_cur,
real_t dt) {
// clang-format off
Kokkos::parallel_for("CurrentsDeposit",
species.rangeActiveParticles(),
kernel::DepositCurrents_kernel<SimEngine::SRPIC, M, O>(
scatter_cur,
species.i1, species.i2, species.i3,
species.i1_prev, species.i2_prev, species.i3_prev,
species.dx1, species.dx2, species.dx3,
species.dx1_prev, species.dx2_prev, species.dx3_prev,
species.ux1, species.ux2, species.ux3,
species.phi, species.weight, species.tag,
metric, (real_t)(species.charge()), dt));
// clang-format on
}

void CurrentsDeposit(domain_t& domain) {
auto scatter_cur = Kokkos::Experimental::create_scatter_view(
domain.fields.cur);
auto shape_order = m_params.template get<int>("algorithms.deposit.order");
for (auto& species : domain.species) {
if ((species.pusher() == PrtlPusher::NONE) or (species.npart() == 0) or
cmp::AlmostZero_host(species.charge())) {
Expand All @@ -554,8 +534,20 @@ namespace ntt {
species.npart(),
(double)species.charge()),
HERE);

deposit_with<SHAPE_ORDER>(species, domain.mesh.metric, scatter_cur, dt);
// clang-format off
Kokkos::parallel_for("CurrentsDeposit",
species.rangeActiveParticles(),
kernel::DepositCurrents_kernel<SimEngine::SRPIC, M>(
scatter_cur,
species.i1, species.i2, species.i3,
species.i1_prev, species.i2_prev, species.i3_prev,
species.dx1, species.dx2, species.dx3,
species.dx1_prev, species.dx2_prev, species.dx3_prev,
species.ux1, species.ux2, species.ux3,
species.phi, species.weight, species.tag,
domain.mesh.metric,
(real_t)(species.charge()), dt));
// clang-format on
}
Kokkos::Experimental::contribute(domain.fields.cur, scatter_cur);
}
Expand Down Expand Up @@ -1291,7 +1283,7 @@ namespace ntt {
}

const auto maxwellian = arch::Maxwellian<S, M> { domain.mesh.metric,
domain.random_pool(),
domain.random_pool,
temp };

if (dim == in::x1) {
Expand Down
Loading