From 0269572a7f633c4c974b14b2642659c933121588 Mon Sep 17 00:00:00 2001 From: haykh Date: Mon, 29 Jul 2024 09:56:03 -0400 Subject: [PATCH 01/21] patch for mpich send buffr --- src/framework/domain/metadomain.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/framework/domain/metadomain.cpp b/src/framework/domain/metadomain.cpp index b13475fd6..26a4f3168 100644 --- a/src/framework/domain/metadomain.cpp +++ b/src/framework/domain/metadomain.cpp @@ -390,7 +390,7 @@ namespace ntt { #if defined(MPI_ENABLED) auto dx_mins = std::vector(g_ndomains); dx_mins[g_mpi_rank] = dx_min; - MPI_Allgather(&dx_mins[g_mpi_rank], + MPI_Allgather(&dx_min, 1, mpi::get_type(), dx_mins.data(), From 35ebe5328b2aa7df5fe95f0f6e9927b0769411b7 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Sat, 3 Aug 2024 14:54:47 -0500 Subject: [PATCH 02/21] initial commit: modification to shock pgen to run magnetized shocks --- setups/srpic/shock/pgen.hpp | 57 ++++++++++++++++++++++++++++++++--- setups/srpic/shock/shock.toml | 4 +++ 2 files changed, 57 insertions(+), 4 deletions(-) diff --git a/setups/srpic/shock/pgen.hpp b/setups/srpic/shock/pgen.hpp index f07b99878..4a9cc3f09 100644 --- a/setups/srpic/shock/pgen.hpp +++ b/setups/srpic/shock/pgen.hpp @@ -14,6 +14,47 @@ namespace user { using namespace ntt; + template + struct InitFields + { + InitFields(real_t bmag, real_t btheta, real_t bphi, real_t bbeta) : + Bmag { bmag }, Btheta { btheta }, Bphi { bphi }, Bbeta { bbeta } {} + + // alternative: initialize magnetisation from simulation parameters as in Tristan? + // Bmag = math::sqrt(ppc0 * 0.5 * c * c * me * sigma); + + // magnetic field components + Inline auto bx1(const coord_t &x_Ph) const -> real_t + { + return Bmag * math::cos(Btheta / 180.0 * Kokkos::numbers::pi); + } + Inline auto bx2(const coord_t &x_Ph) const -> real_t + { + return Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); + } + Inline auto bx3(const coord_t &x_Ph) const -> real_t + { + return Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); + } + + // electric field components + Inline auto ex1(const coord_t &x_Ph) const -> real_t + { + return ZERO; + } + Inline auto ex2(const coord_t &x_Ph) const -> real_t + { + return -Bbeta * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); + } + Inline auto ex3(const coord_t &x_Ph) const -> real_t + { + return -Bbeta * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); + } + + private: + const real_t Btheta, Bphi, Bbeta, Bmag; + }; + template struct PGen : public arch::ProblemGenerator { // compatibility traits for the problem generator @@ -30,10 +71,18 @@ namespace user { const real_t drift_ux, temperature; - inline PGen(const SimulationParams& p, const Metadomain& m) - : arch::ProblemGenerator(p) - , drift_ux { p.template get("setup.drift_ux") } - , temperature { p.template get("setup.temperature") } {} + const real_t Btheta, Bphi, Bbeta, Bmag; + InitFields init_flds; + + inline PGen(const SimulationParams &p, const Metadomain &m) + : arch::ProblemGenerator { p } + , drift_ux { p.template get("setup.drift_ux") } + , temperature { p.template get("setup.temperature") } + , Bmag { p.template get("setup.Bmag", 0.0) } + , Btheta { p.template get("setup.Btheta", 0.0) } + , Bphi { p.template get("setup.Bphi", 0.0) } + , Bbeta { p.template get("setup.Bbeta", 0.0) } + , init_flds { Bmag, Btheta, Bphi, Bbeta } {} inline PGen() {} diff --git a/setups/srpic/shock/shock.toml b/setups/srpic/shock/shock.toml index f48edb2d6..90571631e 100644 --- a/setups/srpic/shock/shock.toml +++ b/setups/srpic/shock/shock.toml @@ -42,6 +42,10 @@ [setup] drift_ux = 0.1 temperature = 1e-3 + Bmag = 0.0 + Btheta = 0.0 + Bphi = 0.0 + Bbeta = 0.0 [output] interval_time = 0.1 From 51323608d56e7a2cb594d4d921c40f03309e63ae Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Thu, 5 Sep 2024 15:34:22 -0500 Subject: [PATCH 03/21] fix misunderstanding in setup --- setups/srpic/shock/pgen.hpp | 15 +++++++-------- 1 file changed, 7 insertions(+), 8 deletions(-) diff --git a/setups/srpic/shock/pgen.hpp b/setups/srpic/shock/pgen.hpp index 4a9cc3f09..c3771cde2 100644 --- a/setups/srpic/shock/pgen.hpp +++ b/setups/srpic/shock/pgen.hpp @@ -17,8 +17,8 @@ namespace user { template struct InitFields { - InitFields(real_t bmag, real_t btheta, real_t bphi, real_t bbeta) : - Bmag { bmag }, Btheta { btheta }, Bphi { bphi }, Bbeta { bbeta } {} + InitFields(real_t bmag, real_t btheta, real_t bphi, real_t drift_ux) : + Bmag { bmag }, Btheta { btheta }, Bphi { bphi }, Vx { drift_ux } {} // alternative: initialize magnetisation from simulation parameters as in Tristan? // Bmag = math::sqrt(ppc0 * 0.5 * c * c * me * sigma); @@ -44,15 +44,15 @@ namespace user { } Inline auto ex2(const coord_t &x_Ph) const -> real_t { - return -Bbeta * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); + return -Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); } Inline auto ex3(const coord_t &x_Ph) const -> real_t { - return -Bbeta * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); + return -Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); } private: - const real_t Btheta, Bphi, Bbeta, Bmag; + const real_t Btheta, Bphi, Vx, Bmag; }; template @@ -71,7 +71,7 @@ namespace user { const real_t drift_ux, temperature; - const real_t Btheta, Bphi, Bbeta, Bmag; + const real_t Btheta, Bphi, Bmag; InitFields init_flds; inline PGen(const SimulationParams &p, const Metadomain &m) @@ -81,8 +81,7 @@ namespace user { , Bmag { p.template get("setup.Bmag", 0.0) } , Btheta { p.template get("setup.Btheta", 0.0) } , Bphi { p.template get("setup.Bphi", 0.0) } - , Bbeta { p.template get("setup.Bbeta", 0.0) } - , init_flds { Bmag, Btheta, Bphi, Bbeta } {} + , init_flds { Bmag, Btheta, Bphi, drift_ux } {} inline PGen() {} From cdfd2859c1c5bf163478ac2f48a5ca90916611b1 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Wed, 11 Sep 2024 13:26:00 -0500 Subject: [PATCH 04/21] fix sign error --- setups/srpic/shock/pgen.hpp | 5 +---- 1 file changed, 1 insertion(+), 4 deletions(-) diff --git a/setups/srpic/shock/pgen.hpp b/setups/srpic/shock/pgen.hpp index c3771cde2..1194e7fed 100644 --- a/setups/srpic/shock/pgen.hpp +++ b/setups/srpic/shock/pgen.hpp @@ -20,9 +20,6 @@ namespace user { InitFields(real_t bmag, real_t btheta, real_t bphi, real_t drift_ux) : Bmag { bmag }, Btheta { btheta }, Bphi { bphi }, Vx { drift_ux } {} - // alternative: initialize magnetisation from simulation parameters as in Tristan? - // Bmag = math::sqrt(ppc0 * 0.5 * c * c * me * sigma); - // magnetic field components Inline auto bx1(const coord_t &x_Ph) const -> real_t { @@ -44,7 +41,7 @@ namespace user { } Inline auto ex2(const coord_t &x_Ph) const -> real_t { - return -Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); + return Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); } Inline auto ex3(const coord_t &x_Ph) const -> real_t { From e24078cadd4096d57e75c1dd14c9e92cc96f60aa Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Wed, 11 Sep 2024 14:07:43 -0500 Subject: [PATCH 05/21] fix sign error --- setups/srpic/shock/pgen.hpp | 5 +---- 1 file changed, 1 insertion(+), 4 deletions(-) diff --git a/setups/srpic/shock/pgen.hpp b/setups/srpic/shock/pgen.hpp index c3771cde2..1194e7fed 100644 --- a/setups/srpic/shock/pgen.hpp +++ b/setups/srpic/shock/pgen.hpp @@ -20,9 +20,6 @@ namespace user { InitFields(real_t bmag, real_t btheta, real_t bphi, real_t drift_ux) : Bmag { bmag }, Btheta { btheta }, Bphi { bphi }, Vx { drift_ux } {} - // alternative: initialize magnetisation from simulation parameters as in Tristan? - // Bmag = math::sqrt(ppc0 * 0.5 * c * c * me * sigma); - // magnetic field components Inline auto bx1(const coord_t &x_Ph) const -> real_t { @@ -44,7 +41,7 @@ namespace user { } Inline auto ex2(const coord_t &x_Ph) const -> real_t { - return -Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); + return Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); } Inline auto ex3(const coord_t &x_Ph) const -> real_t { From 0fc366c76db99099c6eb98ee903abc47e0cf711c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Wed, 11 Sep 2024 14:13:54 -0500 Subject: [PATCH 06/21] Added comment for `InitFields` --- setups/srpic/shock/pgen.hpp | 11 ++++++++++- 1 file changed, 10 insertions(+), 1 deletion(-) diff --git a/setups/srpic/shock/pgen.hpp b/setups/srpic/shock/pgen.hpp index 1194e7fed..715c222df 100644 --- a/setups/srpic/shock/pgen.hpp +++ b/setups/srpic/shock/pgen.hpp @@ -16,7 +16,16 @@ namespace user { template struct InitFields - { + { + /* + Sets up magnetic and electric field components for the simulation. + Must satisfy E = -v x B for Lorentz Force to be zero. + + @param bmag: magnetic field scaling + @param btheta: magnetic field polar angle + @param bphi: magnetic field azimuthal angle + @param drift_ux: drift velocity in the x direction + */ InitFields(real_t bmag, real_t btheta, real_t bphi, real_t drift_ux) : Bmag { bmag }, Btheta { btheta }, Bphi { bphi }, Vx { drift_ux } {} From c3b1018069bba7ec582b1d612d3d4c8b9b48ed5c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Mon, 23 Sep 2024 11:25:57 -0500 Subject: [PATCH 07/21] fix signs (again) --- setups/srpic/shock/pgen.hpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/setups/srpic/shock/pgen.hpp b/setups/srpic/shock/pgen.hpp index 715c222df..1fdd18faa 100644 --- a/setups/srpic/shock/pgen.hpp +++ b/setups/srpic/shock/pgen.hpp @@ -50,11 +50,11 @@ namespace user { } Inline auto ex2(const coord_t &x_Ph) const -> real_t { - return Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); + return -Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); } Inline auto ex3(const coord_t &x_Ph) const -> real_t { - return -Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); + return Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); } private: From d666d1ea0df93e38d048257a9fc48e7a24407504 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Mon, 23 Sep 2024 11:37:27 -0500 Subject: [PATCH 08/21] removed redundant parameter and added comments --- setups/srpic/shock/shock.toml | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/setups/srpic/shock/shock.toml b/setups/srpic/shock/shock.toml index 90571631e..e475ae097 100644 --- a/setups/srpic/shock/shock.toml +++ b/setups/srpic/shock/shock.toml @@ -42,10 +42,9 @@ [setup] drift_ux = 0.1 temperature = 1e-3 - Bmag = 0.0 - Btheta = 0.0 - Bphi = 0.0 - Bbeta = 0.0 + Bmag = 0.0 # set to 1.0 if magnetized shock is required + Btheta = 0.0 # magnetic field polar angle + Bphi = 0.0 # magnetic field azimuthal angle [output] interval_time = 0.1 From 0b65f557e56cbf2598c131f87cc8c8d9212179be Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Wed, 2 Oct 2024 10:10:44 -0500 Subject: [PATCH 09/21] added atmosphere bc to enforce initial magnetic field config at the boundaries --- setups/srpic/shock/pgen.hpp | 30 +++++++++++++++++++++++++++--- 1 file changed, 27 insertions(+), 3 deletions(-) diff --git a/setups/srpic/shock/pgen.hpp b/setups/srpic/shock/pgen.hpp index 1fdd18faa..999a7b608 100644 --- a/setups/srpic/shock/pgen.hpp +++ b/setups/srpic/shock/pgen.hpp @@ -57,9 +57,25 @@ namespace user { return Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); } - private: - const real_t Btheta, Bphi, Vx, Bmag; - }; + private: + const real_t Btheta, Bphi, Vx, Bmag; + }; + + template + struct DriveFields : public InitFields { + DriveFields(real_t bmag, real_t btheta, real_t bphi, real_t drift_ux) : + InitFields {bmag, btheta, bphi, drift_ux} {} + + /* Enforce resetting magnetic and electric field at the boundary + This avoids weird */ + using InitFields::bx1; + using InitFields::bx2; + using InitFields::bx3; + + using InitFields::ex1; + using InitFields::ex2; + using InitFields::ex3; + }; template struct PGen : public arch::ProblemGenerator { @@ -91,6 +107,14 @@ namespace user { inline PGen() {} + auto FieldDriver(real_t time) const -> DriveFields { + const real_t bmag = Bmag; + const real_t btheta = Btheta; + const real_t bphi = Bphi; + const real_t ux = drift_ux; + return DriveFields{bmag, btheta, bphi, ux}; + } + inline void InitPrtls(Domain& local_domain) { const auto energy_dist = arch::Maxwellian(local_domain.mesh.metric, local_domain.random_pool, From 96a295c62ed9ab8d1a5a5337b851cfc419869534 Mon Sep 17 00:00:00 2001 From: haykh Date: Mon, 29 Jul 2024 09:56:03 -0400 Subject: [PATCH 10/21] patch for mpich send buffr --- src/framework/domain/metadomain.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/framework/domain/metadomain.cpp b/src/framework/domain/metadomain.cpp index 5e66bc366..ca52a53cf 100644 --- a/src/framework/domain/metadomain.cpp +++ b/src/framework/domain/metadomain.cpp @@ -381,7 +381,7 @@ namespace ntt { #if defined(MPI_ENABLED) auto dx_mins = std::vector(g_ndomains); dx_mins[g_mpi_rank] = dx_min; - MPI_Allgather(&dx_mins[g_mpi_rank], + MPI_Allgather(&dx_min, 1, mpi::get_type(), dx_mins.data(), From 2344629a5210d54294cea1e35dfc0e1eae682b90 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Sat, 3 Aug 2024 14:54:47 -0500 Subject: [PATCH 11/21] initial commit: modification to shock pgen to run magnetized shocks --- setups/srpic/shock/pgen.hpp | 57 ++++++++++++++++++++++++++++++++--- setups/srpic/shock/shock.toml | 4 +++ 2 files changed, 57 insertions(+), 4 deletions(-) diff --git a/setups/srpic/shock/pgen.hpp b/setups/srpic/shock/pgen.hpp index f07b99878..4a9cc3f09 100644 --- a/setups/srpic/shock/pgen.hpp +++ b/setups/srpic/shock/pgen.hpp @@ -14,6 +14,47 @@ namespace user { using namespace ntt; + template + struct InitFields + { + InitFields(real_t bmag, real_t btheta, real_t bphi, real_t bbeta) : + Bmag { bmag }, Btheta { btheta }, Bphi { bphi }, Bbeta { bbeta } {} + + // alternative: initialize magnetisation from simulation parameters as in Tristan? + // Bmag = math::sqrt(ppc0 * 0.5 * c * c * me * sigma); + + // magnetic field components + Inline auto bx1(const coord_t &x_Ph) const -> real_t + { + return Bmag * math::cos(Btheta / 180.0 * Kokkos::numbers::pi); + } + Inline auto bx2(const coord_t &x_Ph) const -> real_t + { + return Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); + } + Inline auto bx3(const coord_t &x_Ph) const -> real_t + { + return Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); + } + + // electric field components + Inline auto ex1(const coord_t &x_Ph) const -> real_t + { + return ZERO; + } + Inline auto ex2(const coord_t &x_Ph) const -> real_t + { + return -Bbeta * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); + } + Inline auto ex3(const coord_t &x_Ph) const -> real_t + { + return -Bbeta * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); + } + + private: + const real_t Btheta, Bphi, Bbeta, Bmag; + }; + template struct PGen : public arch::ProblemGenerator { // compatibility traits for the problem generator @@ -30,10 +71,18 @@ namespace user { const real_t drift_ux, temperature; - inline PGen(const SimulationParams& p, const Metadomain& m) - : arch::ProblemGenerator(p) - , drift_ux { p.template get("setup.drift_ux") } - , temperature { p.template get("setup.temperature") } {} + const real_t Btheta, Bphi, Bbeta, Bmag; + InitFields init_flds; + + inline PGen(const SimulationParams &p, const Metadomain &m) + : arch::ProblemGenerator { p } + , drift_ux { p.template get("setup.drift_ux") } + , temperature { p.template get("setup.temperature") } + , Bmag { p.template get("setup.Bmag", 0.0) } + , Btheta { p.template get("setup.Btheta", 0.0) } + , Bphi { p.template get("setup.Bphi", 0.0) } + , Bbeta { p.template get("setup.Bbeta", 0.0) } + , init_flds { Bmag, Btheta, Bphi, Bbeta } {} inline PGen() {} diff --git a/setups/srpic/shock/shock.toml b/setups/srpic/shock/shock.toml index f48edb2d6..90571631e 100644 --- a/setups/srpic/shock/shock.toml +++ b/setups/srpic/shock/shock.toml @@ -42,6 +42,10 @@ [setup] drift_ux = 0.1 temperature = 1e-3 + Bmag = 0.0 + Btheta = 0.0 + Bphi = 0.0 + Bbeta = 0.0 [output] interval_time = 0.1 From a5f3485eb9e7f0cddeb45c3e6ebd032c72f1a7fa Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Thu, 5 Sep 2024 15:34:22 -0500 Subject: [PATCH 12/21] fix misunderstanding in setup --- setups/srpic/shock/pgen.hpp | 15 +++++++-------- 1 file changed, 7 insertions(+), 8 deletions(-) diff --git a/setups/srpic/shock/pgen.hpp b/setups/srpic/shock/pgen.hpp index 4a9cc3f09..c3771cde2 100644 --- a/setups/srpic/shock/pgen.hpp +++ b/setups/srpic/shock/pgen.hpp @@ -17,8 +17,8 @@ namespace user { template struct InitFields { - InitFields(real_t bmag, real_t btheta, real_t bphi, real_t bbeta) : - Bmag { bmag }, Btheta { btheta }, Bphi { bphi }, Bbeta { bbeta } {} + InitFields(real_t bmag, real_t btheta, real_t bphi, real_t drift_ux) : + Bmag { bmag }, Btheta { btheta }, Bphi { bphi }, Vx { drift_ux } {} // alternative: initialize magnetisation from simulation parameters as in Tristan? // Bmag = math::sqrt(ppc0 * 0.5 * c * c * me * sigma); @@ -44,15 +44,15 @@ namespace user { } Inline auto ex2(const coord_t &x_Ph) const -> real_t { - return -Bbeta * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); + return -Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); } Inline auto ex3(const coord_t &x_Ph) const -> real_t { - return -Bbeta * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); + return -Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); } private: - const real_t Btheta, Bphi, Bbeta, Bmag; + const real_t Btheta, Bphi, Vx, Bmag; }; template @@ -71,7 +71,7 @@ namespace user { const real_t drift_ux, temperature; - const real_t Btheta, Bphi, Bbeta, Bmag; + const real_t Btheta, Bphi, Bmag; InitFields init_flds; inline PGen(const SimulationParams &p, const Metadomain &m) @@ -81,8 +81,7 @@ namespace user { , Bmag { p.template get("setup.Bmag", 0.0) } , Btheta { p.template get("setup.Btheta", 0.0) } , Bphi { p.template get("setup.Bphi", 0.0) } - , Bbeta { p.template get("setup.Bbeta", 0.0) } - , init_flds { Bmag, Btheta, Bphi, Bbeta } {} + , init_flds { Bmag, Btheta, Bphi, drift_ux } {} inline PGen() {} From 6eb034f3e1637eb4f5e0d1c82747045c7a4ad8ad Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Wed, 11 Sep 2024 14:07:43 -0500 Subject: [PATCH 13/21] fix sign error --- setups/srpic/shock/pgen.hpp | 5 +---- 1 file changed, 1 insertion(+), 4 deletions(-) diff --git a/setups/srpic/shock/pgen.hpp b/setups/srpic/shock/pgen.hpp index c3771cde2..1194e7fed 100644 --- a/setups/srpic/shock/pgen.hpp +++ b/setups/srpic/shock/pgen.hpp @@ -20,9 +20,6 @@ namespace user { InitFields(real_t bmag, real_t btheta, real_t bphi, real_t drift_ux) : Bmag { bmag }, Btheta { btheta }, Bphi { bphi }, Vx { drift_ux } {} - // alternative: initialize magnetisation from simulation parameters as in Tristan? - // Bmag = math::sqrt(ppc0 * 0.5 * c * c * me * sigma); - // magnetic field components Inline auto bx1(const coord_t &x_Ph) const -> real_t { @@ -44,7 +41,7 @@ namespace user { } Inline auto ex2(const coord_t &x_Ph) const -> real_t { - return -Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); + return Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); } Inline auto ex3(const coord_t &x_Ph) const -> real_t { From 9d1b39b03a909fe3f7c234e85e1c43b7d6226ce3 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Wed, 11 Sep 2024 14:13:54 -0500 Subject: [PATCH 14/21] Added comment for `InitFields` --- setups/srpic/shock/pgen.hpp | 11 ++++++++++- 1 file changed, 10 insertions(+), 1 deletion(-) diff --git a/setups/srpic/shock/pgen.hpp b/setups/srpic/shock/pgen.hpp index 1194e7fed..715c222df 100644 --- a/setups/srpic/shock/pgen.hpp +++ b/setups/srpic/shock/pgen.hpp @@ -16,7 +16,16 @@ namespace user { template struct InitFields - { + { + /* + Sets up magnetic and electric field components for the simulation. + Must satisfy E = -v x B for Lorentz Force to be zero. + + @param bmag: magnetic field scaling + @param btheta: magnetic field polar angle + @param bphi: magnetic field azimuthal angle + @param drift_ux: drift velocity in the x direction + */ InitFields(real_t bmag, real_t btheta, real_t bphi, real_t drift_ux) : Bmag { bmag }, Btheta { btheta }, Bphi { bphi }, Vx { drift_ux } {} From 037c705a2c11649713e7b81e33778e37fdddc85c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Mon, 23 Sep 2024 11:25:57 -0500 Subject: [PATCH 15/21] fix signs (again) --- setups/srpic/shock/pgen.hpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/setups/srpic/shock/pgen.hpp b/setups/srpic/shock/pgen.hpp index 715c222df..1fdd18faa 100644 --- a/setups/srpic/shock/pgen.hpp +++ b/setups/srpic/shock/pgen.hpp @@ -50,11 +50,11 @@ namespace user { } Inline auto ex2(const coord_t &x_Ph) const -> real_t { - return Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); + return -Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); } Inline auto ex3(const coord_t &x_Ph) const -> real_t { - return -Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); + return Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); } private: From a613d0aa43bf932bf43e756d1a96c080ab9676b7 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Mon, 23 Sep 2024 11:37:27 -0500 Subject: [PATCH 16/21] removed redundant parameter and added comments --- setups/srpic/shock/shock.toml | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/setups/srpic/shock/shock.toml b/setups/srpic/shock/shock.toml index 90571631e..e475ae097 100644 --- a/setups/srpic/shock/shock.toml +++ b/setups/srpic/shock/shock.toml @@ -42,10 +42,9 @@ [setup] drift_ux = 0.1 temperature = 1e-3 - Bmag = 0.0 - Btheta = 0.0 - Bphi = 0.0 - Bbeta = 0.0 + Bmag = 0.0 # set to 1.0 if magnetized shock is required + Btheta = 0.0 # magnetic field polar angle + Bphi = 0.0 # magnetic field azimuthal angle [output] interval_time = 0.1 From 31a7d5b1105fdec1ebc48f4623b3a85fbf5513e0 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Wed, 2 Oct 2024 10:10:44 -0500 Subject: [PATCH 17/21] added atmosphere bc to enforce initial magnetic field config at the boundaries --- setups/srpic/shock/pgen.hpp | 30 +++++++++++++++++++++++++++--- 1 file changed, 27 insertions(+), 3 deletions(-) diff --git a/setups/srpic/shock/pgen.hpp b/setups/srpic/shock/pgen.hpp index 1fdd18faa..999a7b608 100644 --- a/setups/srpic/shock/pgen.hpp +++ b/setups/srpic/shock/pgen.hpp @@ -57,9 +57,25 @@ namespace user { return Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); } - private: - const real_t Btheta, Bphi, Vx, Bmag; - }; + private: + const real_t Btheta, Bphi, Vx, Bmag; + }; + + template + struct DriveFields : public InitFields { + DriveFields(real_t bmag, real_t btheta, real_t bphi, real_t drift_ux) : + InitFields {bmag, btheta, bphi, drift_ux} {} + + /* Enforce resetting magnetic and electric field at the boundary + This avoids weird */ + using InitFields::bx1; + using InitFields::bx2; + using InitFields::bx3; + + using InitFields::ex1; + using InitFields::ex2; + using InitFields::ex3; + }; template struct PGen : public arch::ProblemGenerator { @@ -91,6 +107,14 @@ namespace user { inline PGen() {} + auto FieldDriver(real_t time) const -> DriveFields { + const real_t bmag = Bmag; + const real_t btheta = Btheta; + const real_t bphi = Bphi; + const real_t ux = drift_ux; + return DriveFields{bmag, btheta, bphi, ux}; + } + inline void InitPrtls(Domain& local_domain) { const auto energy_dist = arch::Maxwellian(local_domain.mesh.metric, local_domain.random_pool, From 78412c85f7b20b592a6919679d747f5711aaf762 Mon Sep 17 00:00:00 2001 From: hayk Date: Mon, 4 Nov 2024 21:10:04 -0500 Subject: [PATCH 18/21] changed FieldBC to FIXED --- .gitignore | 1 + input.example.toml | 4 +- setups/srpic/magnetar/magnetar.toml | 2 +- setups/srpic/magnetosphere/magnetosphere.toml | 2 +- setups/srpic/monopole/monopole.toml | 2 +- setups/srpic/shock/pgen.hpp | 143 +++++++++--------- setups/srpic/shock/shock.py | 4 +- setups/srpic/shock/shock.toml | 8 +- src/engines/srpic.hpp | 40 +++-- src/global/enums.h | 32 ++-- src/global/utils/numeric.h | 4 + src/kernels/fields_bcs.hpp | 12 +- 12 files changed, 127 insertions(+), 127 deletions(-) diff --git a/.gitignore b/.gitignore index a1b05e751..9a167b9d5 100644 --- a/.gitignore +++ b/.gitignore @@ -51,6 +51,7 @@ venv/ # CMake testing files Testing/ +tags .clangd .schema.json *_old/ diff --git a/input.example.toml b/input.example.toml index 5ee34d65d..e541f6e2b 100644 --- a/input.example.toml +++ b/input.example.toml @@ -90,10 +90,10 @@ # Boundary conditions for fields: # @required # @type: 1/2/3-size array of string tuples, each of size 1 or 2 - # @valid: "PERIODIC", "ABSORB", "ATMOSPHERE", "CUSTOM", "HORIZON" + # @valid: "PERIODIC", "ABSORB", "FIXED", "CUSTOM", "HORIZON" # @example: [["CUSTOM", "ABSORB"]] (for 2D spherical [[rmin, rmax]]) # @note: When periodic in any of the directions, you should only set one value [..., ["PERIODIC"], ...] - # @note: In spherical, bondaries in theta/phi are set automatically (only specify bc @ [rmin, rmax]) [["ATMOSPHERE", "ABSORB"]] + # @note: In spherical, bondaries in theta/phi are set automatically (only specify bc @ [rmin, rmax]) [["FIXED", "ABSORB"]] # @note: In GR, the horizon boundary is set automatically (only specify bc @ rmax): [["ABSORB"]] fields = "" # Boundary conditions for fields: diff --git a/setups/srpic/magnetar/magnetar.toml b/setups/srpic/magnetar/magnetar.toml index 2a2260af5..cd9ff5695 100644 --- a/setups/srpic/magnetar/magnetar.toml +++ b/setups/srpic/magnetar/magnetar.toml @@ -11,7 +11,7 @@ metric = "qspherical" [grid.boundaries] - fields = [["ATMOSPHERE", "ABSORB"]] + fields = [["FIXED", "ABSORB"]] particles = [["ATMOSPHERE", "ABSORB"]] [grid.boundaries.absorb] diff --git a/setups/srpic/magnetosphere/magnetosphere.toml b/setups/srpic/magnetosphere/magnetosphere.toml index 34e04b02d..83ade6e48 100644 --- a/setups/srpic/magnetosphere/magnetosphere.toml +++ b/setups/srpic/magnetosphere/magnetosphere.toml @@ -11,7 +11,7 @@ metric = "qspherical" [grid.boundaries] - fields = [["ATMOSPHERE", "ABSORB"]] + fields = [["FIXED", "ABSORB"]] particles = [["ATMOSPHERE", "ABSORB"]] [grid.boundaries.absorb] diff --git a/setups/srpic/monopole/monopole.toml b/setups/srpic/monopole/monopole.toml index 169837489..322c15dd4 100644 --- a/setups/srpic/monopole/monopole.toml +++ b/setups/srpic/monopole/monopole.toml @@ -11,7 +11,7 @@ metric = "qspherical" [grid.boundaries] - fields = [["ATMOSPHERE", "ABSORB"]] + fields = [["FIXED", "ABSORB"]] particles = [["ATMOSPHERE", "ABSORB"]] [grid.boundaries.absorb] diff --git a/setups/srpic/shock/pgen.hpp b/setups/srpic/shock/pgen.hpp index 999a7b608..30929383f 100644 --- a/setups/srpic/shock/pgen.hpp +++ b/setups/srpic/shock/pgen.hpp @@ -5,6 +5,7 @@ #include "global.h" #include "arch/traits.h" +#include "utils/numeric.h" #include "archetypes/energy_dist.h" #include "archetypes/particle_injector.h" @@ -15,67 +16,66 @@ namespace user { using namespace ntt; template - struct InitFields - { - /* - Sets up magnetic and electric field components for the simulation. - Must satisfy E = -v x B for Lorentz Force to be zero. - - @param bmag: magnetic field scaling - @param btheta: magnetic field polar angle - @param bphi: magnetic field azimuthal angle - @param drift_ux: drift velocity in the x direction - */ - InitFields(real_t bmag, real_t btheta, real_t bphi, real_t drift_ux) : - Bmag { bmag }, Btheta { btheta }, Bphi { bphi }, Vx { drift_ux } {} - - // magnetic field components - Inline auto bx1(const coord_t &x_Ph) const -> real_t - { - return Bmag * math::cos(Btheta / 180.0 * Kokkos::numbers::pi); - } - Inline auto bx2(const coord_t &x_Ph) const -> real_t - { - return Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); - } - Inline auto bx3(const coord_t &x_Ph) const -> real_t - { - return Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); - } - - // electric field components - Inline auto ex1(const coord_t &x_Ph) const -> real_t - { - return ZERO; - } - Inline auto ex2(const coord_t &x_Ph) const -> real_t - { - return -Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::cos(Bphi / 180.0 * Kokkos::numbers::pi); - } - Inline auto ex3(const coord_t &x_Ph) const -> real_t - { - return Vx * Bmag * math::sin(Btheta / 180.0 * Kokkos::numbers::pi) * math::sin(Bphi / 180.0 * Kokkos::numbers::pi); - } - - private: - const real_t Btheta, Bphi, Vx, Bmag; - }; + struct InitFields { + /* + Sets up magnetic and electric field components for the simulation. + Must satisfy E = -v x B for Lorentz Force to be zero. + + @param bmag: magnetic field scaling + @param btheta: magnetic field polar angle + @param bphi: magnetic field azimuthal angle + @param drift_ux: drift velocity in the x direction + */ + InitFields(real_t bmag, real_t btheta, real_t bphi, real_t drift_ux) + : Bmag { bmag * static_cast(convert::deg2rad) } + , Btheta { btheta * static_cast(convert::deg2rad) } + , Bphi { bphi * static_cast(convert::deg2rad) } + , Vx { drift_ux } {} + + // magnetic field components + Inline auto bx1(const coord_t& x_Ph) const -> real_t { + return Bmag * math::cos(Btheta); + } + + Inline auto bx2(const coord_t& x_Ph) const -> real_t { + return Bmag * math::sin(Btheta) * math::sin(Bphi); + } + + Inline auto bx3(const coord_t& x_Ph) const -> real_t { + return Bmag * math::sin(Btheta) * math::cos(Bphi); + } + + // electric field components + Inline auto ex1(const coord_t& x_Ph) const -> real_t { + return ZERO; + } + + Inline auto ex2(const coord_t& x_Ph) const -> real_t { + return -Vx * Bmag * math::sin(Btheta) * math::cos(Bphi); + } + Inline auto ex3(const coord_t& x_Ph) const -> real_t { + return Vx * Bmag * math::sin(Btheta) * math::sin(Bphi); + } + + private: + const real_t Btheta, Bphi, Vx, Bmag; + }; + + /* Enforce resetting magnetic and electric field at the boundary */ template struct DriveFields : public InitFields { - DriveFields(real_t bmag, real_t btheta, real_t bphi, real_t drift_ux) : - InitFields {bmag, btheta, bphi, drift_ux} {} - - /* Enforce resetting magnetic and electric field at the boundary - This avoids weird */ - using InitFields::bx1; - using InitFields::bx2; - using InitFields::bx3; - - using InitFields::ex1; - using InitFields::ex2; - using InitFields::ex3; - }; + DriveFields(real_t bmag, real_t btheta, real_t bphi, real_t drift_ux) + : InitFields { bmag, btheta, bphi, drift_ux } {} + + using InitFields::bx1; + using InitFields::bx2; + using InitFields::bx3; + + using InitFields::ex1; + using InitFields::ex2; + using InitFields::ex3; + }; template struct PGen : public arch::ProblemGenerator { @@ -93,26 +93,22 @@ namespace user { const real_t drift_ux, temperature; - const real_t Btheta, Bphi, Bmag; + const real_t Btheta, Bphi, Bmag; InitFields init_flds; - inline PGen(const SimulationParams &p, const Metadomain &m) - : arch::ProblemGenerator { p } - , drift_ux { p.template get("setup.drift_ux") } - , temperature { p.template get("setup.temperature") } - , Bmag { p.template get("setup.Bmag", 0.0) } - , Btheta { p.template get("setup.Btheta", 0.0) } - , Bphi { p.template get("setup.Bphi", 0.0) } - , init_flds { Bmag, Btheta, Bphi, drift_ux } {} + inline PGen(const SimulationParams& p, const Metadomain& m) + : arch::ProblemGenerator { p } + , drift_ux { p.template get("setup.drift_ux") } + , temperature { p.template get("setup.temperature") } + , Bmag { p.template get("setup.Bmag", ZERO) } + , Btheta { p.template get("setup.Btheta", ZERO) } + , Bphi { p.template get("setup.Bphi", ZERO) } + , init_flds { Bmag, Btheta, Bphi, drift_ux } {} inline PGen() {} auto FieldDriver(real_t time) const -> DriveFields { - const real_t bmag = Bmag; - const real_t btheta = Btheta; - const real_t bphi = Bphi; - const real_t ux = drift_ux; - return DriveFields{bmag, btheta, bphi, ux}; + return DriveFields { Bmag, Btheta, Bphi, drift_ux }; } inline void InitPrtls(Domain& local_domain) { @@ -121,7 +117,8 @@ namespace user { temperature, -drift_ux, in::x1); - const auto injector = arch::UniformInjector( + + const auto injector = arch::UniformInjector( energy_dist, { 1, 2 }); arch::InjectUniform>( diff --git a/setups/srpic/shock/shock.py b/setups/srpic/shock/shock.py index 64224c728..dc1565572 100644 --- a/setups/srpic/shock/shock.py +++ b/setups/srpic/shock/shock.py @@ -2,7 +2,7 @@ import matplotlib.pyplot as plt import matplotlib as mpl -data = nt2r.Data("shock-03.h5") +data = nt2r.Data("shock.h5") def frame(ti, f): @@ -55,7 +55,7 @@ def frame(ti, f): axs = [fig.add_subplot(gs[i]) for i in range(len(quantities))] for ax, q in zip(axs, quantities): - q["compute"](f).coarsen(x=2, y=2).mean().plot( + q["compute"](f.isel(t=ti)).plot( ax=ax, cmap=q["cmap"], norm=q["norm"], diff --git a/setups/srpic/shock/shock.toml b/setups/srpic/shock/shock.toml index e475ae097..f8f5e81a7 100644 --- a/setups/srpic/shock/shock.toml +++ b/setups/srpic/shock/shock.toml @@ -11,7 +11,7 @@ metric = "minkowski" [grid.boundaries] - fields = [["CONDUCTOR", "ABSORB"], ["PERIODIC"]] + fields = [["FIXED", "ABSORB"], ["PERIODIC"]] particles = [["REFLECT", "ABSORB"], ["PERIODIC"]] [scales] @@ -42,9 +42,9 @@ [setup] drift_ux = 0.1 temperature = 1e-3 - Bmag = 0.0 # set to 1.0 if magnetized shock is required - Btheta = 0.0 # magnetic field polar angle - Bphi = 0.0 # magnetic field azimuthal angle + Bmag = 1.0 + Btheta = 0.0 + Bphi = 0.0 [output] interval_time = 0.1 diff --git a/src/engines/srpic.hpp b/src/engines/srpic.hpp index 78c8f371e..1747c5138 100644 --- a/src/engines/srpic.hpp +++ b/src/engines/srpic.hpp @@ -42,7 +42,6 @@ #include #include -#include #include namespace ntt { @@ -586,8 +585,8 @@ namespace ntt { if (domain.mesh.flds_bc_in(direction) == FldsBC::AXIS) { AxisFieldsIn(direction, domain, tags); } - } else if (m_metadomain.mesh().flds_bc_in(direction) == FldsBC::ATMOSPHERE) { - AtmosphereFieldsIn(direction, domain, tags); + } else if (m_metadomain.mesh().flds_bc_in(direction) == FldsBC::FIXED) { + FixedFieldsIn(direction, domain, tags); } else if (m_metadomain.mesh().flds_bc_in(direction) == FldsBC::CONDUCTOR) { if (domain.mesh.flds_bc_in(direction) == FldsBC::CONDUCTOR) { ConductorFieldsIn(direction, domain, tags); @@ -713,11 +712,11 @@ namespace ntt { } } - void AtmosphereFieldsIn(dir::direction_t direction, - domain_t& domain, - BCTags tags) { + void FixedFieldsIn(dir::direction_t direction, + domain_t& domain, + BCTags tags) { /** - * atmosphere boundaries + * fixed field boundaries */ if constexpr (traits::has_member::value) { const auto [sign, dim, xg_min, xg_max] = get_atm_extent(direction); @@ -759,9 +758,9 @@ namespace ntt { if (dim == in::x1) { if (sign > 0) { Kokkos::parallel_for( - "AtmosphereBCFields", + "FixedBCFields", range, - kernel::AtmosphereBoundaries_kernel( + kernel::FixedBoundaries_kernel( domain.fields.em, field_driver, domain.mesh.metric, @@ -769,9 +768,9 @@ namespace ntt { tags)); } else { Kokkos::parallel_for( - "AtmosphereBCFields", + "FixedBCFields", range, - kernel::AtmosphereBoundaries_kernel( + kernel::FixedBoundaries_kernel( domain.fields.em, field_driver, domain.mesh.metric, @@ -782,9 +781,9 @@ namespace ntt { if constexpr (M::Dim == Dim::_2D or M::Dim == Dim::_3D) { if (sign > 0) { Kokkos::parallel_for( - "AtmosphereBCFields", + "FixedBCFields", range, - kernel::AtmosphereBoundaries_kernel( + kernel::FixedBoundaries_kernel( domain.fields.em, field_driver, domain.mesh.metric, @@ -792,9 +791,9 @@ namespace ntt { tags)); } else { Kokkos::parallel_for( - "AtmosphereBCFields", + "FixedBCFields", range, - kernel::AtmosphereBoundaries_kernel( + kernel::FixedBoundaries_kernel( domain.fields.em, field_driver, domain.mesh.metric, @@ -808,9 +807,9 @@ namespace ntt { if constexpr (M::Dim == Dim::_3D) { if (sign > 0) { Kokkos::parallel_for( - "AtmosphereBCFields", + "FixedBCFields", range, - kernel::AtmosphereBoundaries_kernel( + kernel::FixedBoundaries_kernel( domain.fields.em, field_driver, domain.mesh.metric, @@ -818,9 +817,9 @@ namespace ntt { tags)); } else { Kokkos::parallel_for( - "AtmosphereBCFields", + "FixedBCFields", range, - kernel::AtmosphereBoundaries_kernel( + kernel::FixedBoundaries_kernel( domain.fields.em, field_driver, domain.mesh.metric, @@ -834,8 +833,7 @@ namespace ntt { raise::Error("Invalid dimension", HERE); } } else { - raise::Error("Field driver not implemented in PGEN for atmosphere BCs", - HERE); + raise::Error("Field driver not implemented in PGEN for fixed BCs", HERE); } } diff --git a/src/global/enums.h b/src/global/enums.h index 57822dec4..1946da4b8 100644 --- a/src/global/enums.h +++ b/src/global/enums.h @@ -8,7 +8,7 @@ * - enum ntt::SimEngine // SRPIC, GRPIC * - enum ntt::PrtlBC // periodic, absorb, atmosphere, custom, * reflect, horizon, axis, sync - * - enum ntt::FldsBC // periodic, absorb, atmosphere, custom, + * - enum ntt::FldsBC // periodic, absorb, fixed, custom, * conductor, horizon, axis, sync * - enum ntt::PrtlPusher // boris, vay, photon, none * - enum ntt::Cooling // synchrotron, none @@ -213,25 +213,25 @@ namespace ntt { static constexpr const char* label = "flds_bc"; enum type : uint8_t { - INVALID = 0, - PERIODIC = 1, - ABSORB = 2, - ATMOSPHERE = 3, - CUSTOM = 4, - CONDUCTOR = 5, - HORIZON = 6, - AXIS = 7, - SYNC = 8, // <- SYNC means synchronization with other domains + INVALID = 0, + PERIODIC = 1, + ABSORB = 2, + FIXED = 3, + CUSTOM = 4, + CONDUCTOR = 5, + HORIZON = 6, + AXIS = 7, + SYNC = 8, // <- SYNC means synchronization with other domains }; constexpr FldsBC(uint8_t c) : enums_hidden::BaseEnum { c } {} - static constexpr type variants[] = { PERIODIC, ABSORB, ATMOSPHERE, CUSTOM, - CONDUCTOR, HORIZON, AXIS, SYNC }; - static constexpr const char* lookup[] = { "periodic", "absorb", - "atmosphere", "custom", - "conductor", "horizon", - "axis", "sync" }; + static constexpr type variants[] = { PERIODIC, ABSORB, FIXED, CUSTOM, + CONDUCTOR, HORIZON, AXIS, SYNC }; + static constexpr const char* lookup[] = { "periodic", "absorb", + "fixed", "custom", + "conductor", "horizon", + "axis", "sync" }; static constexpr std::size_t total = sizeof(variants) / sizeof(variants[0]); }; diff --git a/src/global/utils/numeric.h b/src/global/utils/numeric.h index 0b09f6c11..719256d1d 100644 --- a/src/global/utils/numeric.h +++ b/src/global/utils/numeric.h @@ -91,4 +91,8 @@ namespace constant { inline constexpr double SQRT3 = 1.73205080756887729352; } // namespace constant +namespace convert { + inline constexpr double deg2rad = constant::PI / 180.0; +} // namespace convert + #endif // GLOBAL_UTILS_NUMERIC_H diff --git a/src/kernels/fields_bcs.hpp b/src/kernels/fields_bcs.hpp index e617010b4..1a9ffcff4 100644 --- a/src/kernels/fields_bcs.hpp +++ b/src/kernels/fields_bcs.hpp @@ -217,7 +217,7 @@ namespace kernel { }; template - struct AtmosphereBoundaries_kernel { + struct FixedBoundaries_kernel { static constexpr Dimension D = M::Dim; static constexpr bool defines_ex1 = traits::has_method::value; static constexpr bool defines_ex2 = traits::has_method::value; @@ -240,11 +240,11 @@ namespace kernel { const std::size_t i_edge; const bool setE, setB; - AtmosphereBoundaries_kernel(ndfield_t& Fld, - const I& finit, - const M& metric, - std::size_t i_edge, - BCTags tags) + FixedBoundaries_kernel(ndfield_t& Fld, + const I& finit, + const M& metric, + std::size_t i_edge, + BCTags tags) : Fld { Fld } , finit { finit } , metric { metric } From 98349950f7392364248c353906e6440369c322ca Mon Sep 17 00:00:00 2001 From: hayk Date: Tue, 5 Nov 2024 13:43:39 -0500 Subject: [PATCH 19/21] FIXED bc + ATM --- CMakeLists.txt | 6 +- input.example.toml | 6 +- setups/srpic/magnetar/magnetar.toml | 2 +- setups/srpic/magnetosphere/magnetosphere.toml | 2 +- setups/srpic/monopole/monopole.toml | 2 +- setups/srpic/shock/pgen.hpp | 42 ++-- src/engines/srpic.hpp | 206 ++++++++++-------- src/global/arch/traits.h | 3 + src/global/enums.h | 33 ++- src/global/tests/enums.cpp | 4 +- src/kernels/fields_bcs.hpp | 12 +- 11 files changed, 166 insertions(+), 152 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 4ee00d1b4..1a977c990 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -55,7 +55,9 @@ if(${DEBUG} STREQUAL "OFF") set(CMAKE_BUILD_TYPE Release CACHE STRING "CMake build type") - set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -DNDEBUG") + set(CMAKE_CXX_FLAGS + "${CMAKE_CXX_FLAGS} -DNDEBUG -Wno-unused-local-typedefs -Wno-unknown-cuda-version" + ) else() set(CMAKE_BUILD_TYPE Debug @@ -64,8 +66,6 @@ else() "${CMAKE_CXX_FLAGS} -DDEBUG -Wall -Wextra -Wno-unknown-pragmas") endif() -set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wno-unused-local-typedefs") - # options set(precisions "single" "double" diff --git a/input.example.toml b/input.example.toml index e541f6e2b..4bf005b0e 100644 --- a/input.example.toml +++ b/input.example.toml @@ -90,10 +90,10 @@ # Boundary conditions for fields: # @required # @type: 1/2/3-size array of string tuples, each of size 1 or 2 - # @valid: "PERIODIC", "ABSORB", "FIXED", "CUSTOM", "HORIZON" + # @valid: "PERIODIC", "ABSORB", "FIXED", "ATMOSPHERE", "CUSTOM", "HORIZON" # @example: [["CUSTOM", "ABSORB"]] (for 2D spherical [[rmin, rmax]]) - # @note: When periodic in any of the directions, you should only set one value [..., ["PERIODIC"], ...] - # @note: In spherical, bondaries in theta/phi are set automatically (only specify bc @ [rmin, rmax]) [["FIXED", "ABSORB"]] + # @note: When periodic in any of the directions, you should only set one value: [..., ["PERIODIC"], ...] + # @note: In spherical, bondaries in theta/phi are set automatically (only specify bc @ [rmin, rmax]): [["ATMOSPHERE", "ABSORB"]] # @note: In GR, the horizon boundary is set automatically (only specify bc @ rmax): [["ABSORB"]] fields = "" # Boundary conditions for fields: diff --git a/setups/srpic/magnetar/magnetar.toml b/setups/srpic/magnetar/magnetar.toml index cd9ff5695..2a2260af5 100644 --- a/setups/srpic/magnetar/magnetar.toml +++ b/setups/srpic/magnetar/magnetar.toml @@ -11,7 +11,7 @@ metric = "qspherical" [grid.boundaries] - fields = [["FIXED", "ABSORB"]] + fields = [["ATMOSPHERE", "ABSORB"]] particles = [["ATMOSPHERE", "ABSORB"]] [grid.boundaries.absorb] diff --git a/setups/srpic/magnetosphere/magnetosphere.toml b/setups/srpic/magnetosphere/magnetosphere.toml index 83ade6e48..34e04b02d 100644 --- a/setups/srpic/magnetosphere/magnetosphere.toml +++ b/setups/srpic/magnetosphere/magnetosphere.toml @@ -11,7 +11,7 @@ metric = "qspherical" [grid.boundaries] - fields = [["FIXED", "ABSORB"]] + fields = [["ATMOSPHERE", "ABSORB"]] particles = [["ATMOSPHERE", "ABSORB"]] [grid.boundaries.absorb] diff --git a/setups/srpic/monopole/monopole.toml b/setups/srpic/monopole/monopole.toml index 322c15dd4..169837489 100644 --- a/setups/srpic/monopole/monopole.toml +++ b/setups/srpic/monopole/monopole.toml @@ -11,7 +11,7 @@ metric = "qspherical" [grid.boundaries] - fields = [["FIXED", "ABSORB"]] + fields = [["ATMOSPHERE", "ABSORB"]] particles = [["ATMOSPHERE", "ABSORB"]] [grid.boundaries.absorb] diff --git a/setups/srpic/shock/pgen.hpp b/setups/srpic/shock/pgen.hpp index 30929383f..1eedb3a01 100644 --- a/setups/srpic/shock/pgen.hpp +++ b/setups/srpic/shock/pgen.hpp @@ -5,6 +5,7 @@ #include "global.h" #include "arch/traits.h" +#include "utils/error.h" #include "utils/numeric.h" #include "archetypes/energy_dist.h" @@ -33,28 +34,28 @@ namespace user { , Vx { drift_ux } {} // magnetic field components - Inline auto bx1(const coord_t& x_Ph) const -> real_t { + Inline auto bx1(const coord_t&) const -> real_t { return Bmag * math::cos(Btheta); } - Inline auto bx2(const coord_t& x_Ph) const -> real_t { + Inline auto bx2(const coord_t&) const -> real_t { return Bmag * math::sin(Btheta) * math::sin(Bphi); } - Inline auto bx3(const coord_t& x_Ph) const -> real_t { + Inline auto bx3(const coord_t&) const -> real_t { return Bmag * math::sin(Btheta) * math::cos(Bphi); } // electric field components - Inline auto ex1(const coord_t& x_Ph) const -> real_t { + Inline auto ex1(const coord_t&) const -> real_t { return ZERO; } - Inline auto ex2(const coord_t& x_Ph) const -> real_t { + Inline auto ex2(const coord_t&) const -> real_t { return -Vx * Bmag * math::sin(Btheta) * math::cos(Bphi); } - Inline auto ex3(const coord_t& x_Ph) const -> real_t { + Inline auto ex3(const coord_t&) const -> real_t { return Vx * Bmag * math::sin(Btheta) * math::sin(Bphi); } @@ -62,21 +63,6 @@ namespace user { const real_t Btheta, Bphi, Vx, Bmag; }; - /* Enforce resetting magnetic and electric field at the boundary */ - template - struct DriveFields : public InitFields { - DriveFields(real_t bmag, real_t btheta, real_t bphi, real_t drift_ux) - : InitFields { bmag, btheta, bphi, drift_ux } {} - - using InitFields::bx1; - using InitFields::bx2; - using InitFields::bx3; - - using InitFields::ex1; - using InitFields::ex2; - using InitFields::ex3; - }; - template struct PGen : public arch::ProblemGenerator { // compatibility traits for the problem generator @@ -107,8 +93,18 @@ namespace user { inline PGen() {} - auto FieldDriver(real_t time) const -> DriveFields { - return DriveFields { Bmag, Btheta, Bphi, drift_ux }; + auto FixField(const em& comp) const -> real_t { + if (comp == em::ex2) { + return init_flds.ex2({ ZERO }); + } else if (comp == em::ex3) { + return init_flds.ex3({ ZERO }); + } else if (comp == em::bx1) { + return init_flds.bx1({ ZERO }); + } else { + raise::Error("Other components should not be requested when BC is in X", + HERE); + return ZERO; + } } inline void InitPrtls(Domain& local_domain) { diff --git a/src/engines/srpic.hpp b/src/engines/srpic.hpp index 1747c5138..244a5f863 100644 --- a/src/engines/srpic.hpp +++ b/src/engines/srpic.hpp @@ -585,11 +585,11 @@ namespace ntt { if (domain.mesh.flds_bc_in(direction) == FldsBC::AXIS) { AxisFieldsIn(direction, domain, tags); } + } else if (m_metadomain.mesh().flds_bc_in(direction) == FldsBC::ATMOSPHERE) { + AtmosphereFieldsIn(direction, domain, tags); } else if (m_metadomain.mesh().flds_bc_in(direction) == FldsBC::FIXED) { - FixedFieldsIn(direction, domain, tags); - } else if (m_metadomain.mesh().flds_bc_in(direction) == FldsBC::CONDUCTOR) { - if (domain.mesh.flds_bc_in(direction) == FldsBC::CONDUCTOR) { - ConductorFieldsIn(direction, domain, tags); + if (domain.mesh.flds_bc_in(direction) == FldsBC::FIXED) { + FixedFieldsIn(direction, domain, tags); } } else if (m_metadomain.mesh().flds_bc_in(direction) == FldsBC::CUSTOM) { if (domain.mesh.flds_bc_in(direction) == FldsBC::CUSTOM) { @@ -718,6 +718,100 @@ namespace ntt { /** * fixed field boundaries */ + const auto sign = direction.get_sign(); + const auto dim = direction.get_dim(); + raise::ErrorIf(dim != in::x1 and M::CoordType != Coord::Cart, + "Fixed BCs only implemented for x1 in " + "non-cartesian coordinates", + HERE); + em normal_b_comp, tang_e_comp1, tang_e_comp2; + if (dim == in::x1) { + normal_b_comp = em::bx1; + tang_e_comp1 = em::ex2; + tang_e_comp2 = em::ex3; + } else if (dim == in::x2) { + normal_b_comp = em::bx2; + tang_e_comp1 = em::ex1; + tang_e_comp2 = em::ex3; + } else if (dim == in::x3) { + normal_b_comp = em::bx3; + tang_e_comp1 = em::ex1; + tang_e_comp2 = em::ex2; + } else { + raise::Error("Invalid dimension", HERE); + } + std::vector xi_min, xi_max; + const std::vector all_dirs { in::x1, in::x2, in::x3 }; + for (unsigned short d { 0 }; d < static_cast(M::Dim); ++d) { + const auto dd = all_dirs[d]; + if (dim == dd) { + if (sign > 0) { // + direction + xi_min.push_back(domain.mesh.n_all(dd) - N_GHOSTS); + xi_max.push_back(domain.mesh.n_all(dd)); + } else { // - direction + xi_min.push_back(0); + xi_max.push_back(N_GHOSTS); + } + } else { + xi_min.push_back(0); + xi_max.push_back(domain.mesh.n_all(dd)); + } + } + raise::ErrorIf(xi_min.size() != xi_max.size() or + xi_min.size() != static_cast(M::Dim), + "Invalid range size", + HERE); + std::vector comps; + if (tags & BC::E) { + comps.push_back(tang_e_comp1); + comps.push_back(tang_e_comp2); + } + if (tags & BC::B) { + comps.push_back(normal_b_comp); + } + if constexpr (traits::has_member::value) { + raise::Error("Field driver for fixed fields not implemented", HERE); + } else { + // if field driver not present, set fields to fixed values + for (const auto& comp : comps) { + auto value = ZERO; + if constexpr ( + traits::has_member::value) { + // if fix field function present, read from it + value = m_pgen.FixField((em)comp); + } + if constexpr (M::Dim == Dim::_1D) { + Kokkos::deep_copy(Kokkos::subview(domain.fields.em, + std::make_pair(xi_min[0], xi_max[0]), + comp), + value); + } else if constexpr (M::Dim == Dim::_2D) { + Kokkos::deep_copy(Kokkos::subview(domain.fields.em, + std::make_pair(xi_min[0], xi_max[0]), + std::make_pair(xi_min[1], xi_max[1]), + comp), + value); + } else if constexpr (M::Dim == Dim::_3D) { + Kokkos::deep_copy( + Kokkos::subview(domain.fields.em, + std::make_pair(xi_min[0], xi_max[0]), + std::make_pair(xi_min[1], xi_max[1]), + std::make_pair(xi_min[2], xi_max[2]), + comp), + value); + } else { + raise::Error("Invalid dimension", HERE); + } + } + } + } + + void AtmosphereFieldsIn(dir::direction_t direction, + domain_t& domain, + BCTags tags) { + /** + * atmosphere field boundaries + */ if constexpr (traits::has_member::value) { const auto [sign, dim, xg_min, xg_max] = get_atm_extent(direction); const auto dd = static_cast(dim); @@ -758,9 +852,9 @@ namespace ntt { if (dim == in::x1) { if (sign > 0) { Kokkos::parallel_for( - "FixedBCFields", + "AtmosphereBCFields", range, - kernel::FixedBoundaries_kernel( + kernel::EnforcedBoundaries_kernel( domain.fields.em, field_driver, domain.mesh.metric, @@ -768,9 +862,9 @@ namespace ntt { tags)); } else { Kokkos::parallel_for( - "FixedBCFields", + "AtmosphereBCFields", range, - kernel::FixedBoundaries_kernel( + kernel::EnforcedBoundaries_kernel( domain.fields.em, field_driver, domain.mesh.metric, @@ -781,9 +875,9 @@ namespace ntt { if constexpr (M::Dim == Dim::_2D or M::Dim == Dim::_3D) { if (sign > 0) { Kokkos::parallel_for( - "FixedBCFields", + "AtmosphereBCFields", range, - kernel::FixedBoundaries_kernel( + kernel::EnforcedBoundaries_kernel( domain.fields.em, field_driver, domain.mesh.metric, @@ -791,9 +885,9 @@ namespace ntt { tags)); } else { Kokkos::parallel_for( - "FixedBCFields", + "AtmosphereBCFields", range, - kernel::FixedBoundaries_kernel( + kernel::EnforcedBoundaries_kernel( domain.fields.em, field_driver, domain.mesh.metric, @@ -807,9 +901,9 @@ namespace ntt { if constexpr (M::Dim == Dim::_3D) { if (sign > 0) { Kokkos::parallel_for( - "FixedBCFields", + "AtmosphereBCFields", range, - kernel::FixedBoundaries_kernel( + kernel::EnforcedBoundaries_kernel( domain.fields.em, field_driver, domain.mesh.metric, @@ -817,9 +911,9 @@ namespace ntt { tags)); } else { Kokkos::parallel_for( - "FixedBCFields", + "AtmosphereBCFields", range, - kernel::FixedBoundaries_kernel( + kernel::EnforcedBoundaries_kernel( domain.fields.em, field_driver, domain.mesh.metric, @@ -833,86 +927,8 @@ namespace ntt { raise::Error("Invalid dimension", HERE); } } else { - raise::Error("Field driver not implemented in PGEN for fixed BCs", HERE); - } - } - - void ConductorFieldsIn(dir::direction_t direction, - domain_t& domain, - BCTags tags) { - const auto sign = direction.get_sign(); - const auto dim = direction.get_dim(); - raise::ErrorIf( - dim != in::x1 and M::CoordType != Coord::Cart, - "Conductor BCs only implemented for x1 in non-cartesian coordinates", - HERE); - em normal_b_comp, tang_e_comp1, tang_e_comp2; - if (dim == in::x1) { - normal_b_comp = em::bx1; - tang_e_comp1 = em::ex2; - tang_e_comp2 = em::ex3; - } else if (dim == in::x2) { - normal_b_comp = em::bx2; - tang_e_comp1 = em::ex1; - tang_e_comp2 = em::ex3; - } else if (dim == in::x3) { - normal_b_comp = em::bx3; - tang_e_comp1 = em::ex1; - tang_e_comp2 = em::ex2; - } else { - raise::Error("Invalid dimension", HERE); - } - std::vector xi_min, xi_max; - const std::vector all_dirs { in::x1, in::x2, in::x3 }; - for (unsigned short d { 0 }; d < static_cast(M::Dim); ++d) { - const auto dd = all_dirs[d]; - if (dim == dd) { - if (sign > 0) { // + direction - xi_min.push_back(domain.mesh.n_all(dd) - N_GHOSTS); - xi_max.push_back(domain.mesh.n_all(dd)); - } else { // - direction - xi_min.push_back(0); - xi_max.push_back(N_GHOSTS); - } - } else { - xi_min.push_back(0); - xi_max.push_back(domain.mesh.n_all(dd)); - } - } - raise::ErrorIf(xi_min.size() != xi_max.size() or - xi_min.size() != static_cast(M::Dim), - "Invalid range size", + raise::Error("Field driver not implemented in PGEN for atmosphere BCs", HERE); - std::vector comps; - if (tags & BC::E) { - comps.push_back(tang_e_comp1); - comps.push_back(tang_e_comp2); - } - if (tags & BC::B) { - comps.push_back(normal_b_comp); - } - for (const auto& comp : comps) { - if constexpr (M::Dim == Dim::_1D) { - Kokkos::deep_copy(Kokkos::subview(domain.fields.em, - std::make_pair(xi_min[0], xi_max[0]), - comp), - ZERO); - } else if constexpr (M::Dim == Dim::_2D) { - Kokkos::deep_copy(Kokkos::subview(domain.fields.em, - std::make_pair(xi_min[0], xi_max[0]), - std::make_pair(xi_min[1], xi_max[1]), - comp), - ZERO); - } else if constexpr (M::Dim == Dim::_3D) { - Kokkos::deep_copy(Kokkos::subview(domain.fields.em, - std::make_pair(xi_min[0], xi_max[0]), - std::make_pair(xi_min[1], xi_max[1]), - std::make_pair(xi_min[2], xi_max[2]), - comp), - ZERO); - } else { - raise::Error("Invalid dimension", HERE); - } } } diff --git a/src/global/arch/traits.h b/src/global/arch/traits.h index e915bdf1a..9fd40e201 100644 --- a/src/global/arch/traits.h +++ b/src/global/arch/traits.h @@ -96,6 +96,9 @@ namespace traits { template using field_driver_t = decltype(&T::FieldDriver); + template + using fix_field_t = decltype(&T::FixField); + template using custom_fields_t = decltype(&T::CustomFields); diff --git a/src/global/enums.h b/src/global/enums.h index 1946da4b8..283cb456d 100644 --- a/src/global/enums.h +++ b/src/global/enums.h @@ -8,8 +8,8 @@ * - enum ntt::SimEngine // SRPIC, GRPIC * - enum ntt::PrtlBC // periodic, absorb, atmosphere, custom, * reflect, horizon, axis, sync - * - enum ntt::FldsBC // periodic, absorb, fixed, custom, - * conductor, horizon, axis, sync + * - enum ntt::FldsBC // periodic, absorb, fixed, atmosphere, + * custom, horizon, axis, sync * - enum ntt::PrtlPusher // boris, vay, photon, none * - enum ntt::Cooling // synchrotron, none * - enum ntt::FldsID // e, dive, d, divd, b, h, j, @@ -213,25 +213,24 @@ namespace ntt { static constexpr const char* label = "flds_bc"; enum type : uint8_t { - INVALID = 0, - PERIODIC = 1, - ABSORB = 2, - FIXED = 3, - CUSTOM = 4, - CONDUCTOR = 5, - HORIZON = 6, - AXIS = 7, - SYNC = 8, // <- SYNC means synchronization with other domains + INVALID = 0, + PERIODIC = 1, + ABSORB = 2, + FIXED = 3, + ATMOSPHERE = 4, + CUSTOM = 5, + HORIZON = 6, + AXIS = 7, + SYNC = 8, // <- SYNC means synchronization with other domains }; constexpr FldsBC(uint8_t c) : enums_hidden::BaseEnum { c } {} - static constexpr type variants[] = { PERIODIC, ABSORB, FIXED, CUSTOM, - CONDUCTOR, HORIZON, AXIS, SYNC }; - static constexpr const char* lookup[] = { "periodic", "absorb", - "fixed", "custom", - "conductor", "horizon", - "axis", "sync" }; + static constexpr type variants[] = { PERIODIC, ABSORB, FIXED, ATMOSPHERE, + CUSTOM, HORIZON, AXIS, SYNC }; + static constexpr const char* lookup[] = { "periodic", "absorb", "fixed", + "atmosphere", "custom", "horizon", + "axis", "sync" }; static constexpr std::size_t total = sizeof(variants) / sizeof(variants[0]); }; diff --git a/src/global/tests/enums.cpp b/src/global/tests/enums.cpp index 1fc57398f..4d678e85e 100644 --- a/src/global/tests/enums.cpp +++ b/src/global/tests/enums.cpp @@ -61,8 +61,8 @@ auto main() -> int { enum_str_t all_simulation_engines = { "srpic", "grpic" }; enum_str_t all_particle_bcs = { "periodic", "absorb", "atmosphere", "custom", "reflect", "horizon", "axis", "sync" }; - enum_str_t all_fields_bcs = { "periodic", "absorb", "atmosphere", "custom", - "horizon", "conductor", "axis", "sync" }; + enum_str_t all_fields_bcs = { "periodic", "absorb", "fixed", "atmosphere", + "custom", "horizon", "axis", "sync" }; enum_str_t all_particle_pushers = { "boris", "vay", "photon", "none" }; enum_str_t all_coolings = { "synchrotron", "none" }; diff --git a/src/kernels/fields_bcs.hpp b/src/kernels/fields_bcs.hpp index 1a9ffcff4..2f2a458bb 100644 --- a/src/kernels/fields_bcs.hpp +++ b/src/kernels/fields_bcs.hpp @@ -217,7 +217,7 @@ namespace kernel { }; template - struct FixedBoundaries_kernel { + struct EnforcedBoundaries_kernel { static constexpr Dimension D = M::Dim; static constexpr bool defines_ex1 = traits::has_method::value; static constexpr bool defines_ex2 = traits::has_method::value; @@ -240,11 +240,11 @@ namespace kernel { const std::size_t i_edge; const bool setE, setB; - FixedBoundaries_kernel(ndfield_t& Fld, - const I& finit, - const M& metric, - std::size_t i_edge, - BCTags tags) + EnforcedBoundaries_kernel(ndfield_t& Fld, + const I& finit, + const M& metric, + std::size_t i_edge, + BCTags tags) : Fld { Fld } , finit { finit } , metric { metric } From 833f038e88e2b0b00bccdea2425127dc777d4e1a Mon Sep 17 00:00:00 2001 From: hayk Date: Wed, 6 Nov 2024 17:02:16 -0500 Subject: [PATCH 20/21] multifile output --- input.example.toml | 4 ++ src/framework/domain/output.cpp | 11 +++-- src/framework/parameters.cpp | 10 ++++- src/global/global.h | 11 +++++ src/output/CMakeLists.txt | 1 + src/output/writer.cpp | 73 ++++++++++++++++++++++++++------- src/output/writer.h | 10 +++-- 7 files changed, 98 insertions(+), 22 deletions(-) diff --git a/input.example.toml b/input.example.toml index 4bf005b0e..ce5a2079d 100644 --- a/input.example.toml +++ b/input.example.toml @@ -320,6 +320,10 @@ # @default: -1.0 (disabled) # @note: When `interval_time` < 0, the output is controlled by `interval`, otherwise by `interval_time` interval_time = "" + # Whether to output each timestep into separate files: + # @type: bool + # @default: true + separate_files = "" [output.fields] # Toggle for the field output: diff --git a/src/framework/domain/output.cpp b/src/framework/domain/output.cpp index c7cb6bb65..d88e593b5 100644 --- a/src/framework/domain/output.cpp +++ b/src/framework/domain/output.cpp @@ -65,7 +65,8 @@ namespace ntt { g_writer.init(ptr_adios, params.template get("output.format"), - params.template get("simulation.name")); + params.template get("simulation.name"), + params.template get("output.separate_files")); g_writer.defineMeshLayout(glob_shape_with_ghosts, off_ncells_with_ghosts, loc_shape_with_ghosts, @@ -215,8 +216,8 @@ namespace ntt { "local_domain is a placeholder", HERE); logger::Checkpoint("Writing output", HERE); - g_writer.beginWriting(current_step, current_time); if (write_fields) { + g_writer.beginWriting(WriteMode::Fields, current_step, current_time); const auto incl_ghosts = params.template get("output.debug.ghosts"); const auto dwn = params.template get>( "output.fields.downsampling"); @@ -467,9 +468,11 @@ namespace ntt { } g_writer.writeField(names, local_domain->fields.bckp, addresses); } + g_writer.endWriting(WriteMode::Fields); } // end shouldWrite("fields", step, time) if (write_particles) { + g_writer.beginWriting(WriteMode::Particles, current_step, current_time); const auto prtl_stride = params.template get( "output.particles.stride"); for (const auto& prtl : g_writer.speciesWriters()) { @@ -547,9 +550,11 @@ namespace ntt { g_writer.writeParticleQuantity(buff_x3, glob_tot, offset, prtl.name("X", 3)); } } + g_writer.endWriting(WriteMode::Particles); } // end shouldWrite("particles", step, time) if (write_spectra) { + g_writer.beginWriting(WriteMode::Spectra, current_step, current_time); const auto log_bins = params.template get( "output.spectra.log_bins"); const auto n_bins = params.template get( @@ -613,9 +618,9 @@ namespace ntt { g_writer.writeSpectrum(dn, spec.name()); } g_writer.writeSpectrumBins(energy, "sEbn"); + g_writer.endWriting(WriteMode::Spectra); } // end shouldWrite("spectra", step, time) - g_writer.endWriting(); return true; } diff --git a/src/framework/parameters.cpp b/src/framework/parameters.cpp index b667b5ac9..91a14ae09 100644 --- a/src/framework/parameters.cpp +++ b/src/framework/parameters.cpp @@ -463,6 +463,9 @@ namespace ntt { toml::find_or(toml_data, "output", "interval", defaults::output::interval)); set("output.interval_time", toml::find_or(toml_data, "output", "interval_time", -1.0)); + set("output.separate_files", + toml::find_or(toml_data, "output", "separate_files", true)); + promiseToDefine("output.fields.interval"); promiseToDefine("output.fields.interval_time"); promiseToDefine("output.fields.enable"); @@ -509,11 +512,16 @@ namespace ntt { set("output.fields.downsampling", field_dwn); // particles + auto all_specs = std::vector {}; + const auto nspec = get("particles.nspec"); + for (auto i = 0u; i < nspec; ++i) { + all_specs.push_back(static_cast(i + 1)); + } const auto prtl_out = toml::find_or(toml_data, "output", "particles", "species", - std::vector {}); + all_specs); set("output.particles.species", prtl_out); set("output.particles.stride", toml::find_or(toml_data, diff --git a/src/global/global.h b/src/global/global.h index ad524fb0e..d55d8e21e 100644 --- a/src/global/global.h +++ b/src/global/global.h @@ -249,6 +249,17 @@ namespace Comm { typedef int CommTags; +namespace WriteMode { + enum WriteModeTags_ { + None = 0, + Fields = 1 << 0, + Particles = 1 << 1, + Spectra = 1 << 2, + }; +} // namespace WriteMode + +typedef int WriteModeTags; + namespace BC { enum BCTags_ { None = 0, diff --git a/src/output/CMakeLists.txt b/src/output/CMakeLists.txt index e6dbcc03a..81333e9ff 100644 --- a/src/output/CMakeLists.txt +++ b/src/output/CMakeLists.txt @@ -30,6 +30,7 @@ add_library(ntt_output ${SOURCES}) set(libs ntt_global) add_dependencies(ntt_output ${libs}) target_link_libraries(ntt_output PUBLIC ${libs}) +target_link_libraries(ntt_output PRIVATE stdc++fs) target_include_directories( ntt_output diff --git a/src/output/writer.cpp b/src/output/writer.cpp index 4ba0ea14c..95965c864 100644 --- a/src/output/writer.cpp +++ b/src/output/writer.cpp @@ -18,6 +18,7 @@ #include #endif +#include #include #include @@ -25,9 +26,11 @@ namespace out { void Writer::init(adios2::ADIOS* ptr_adios, const std::string& engine, - const std::string& title) { - m_engine = engine; - p_adios = ptr_adios; + const std::string& title, + bool use_separate_files) { + m_separate_files = use_separate_files; + m_engine = engine; + p_adios = ptr_adios; raise::ErrorIf(p_adios == nullptr, "ADIOS pointer is null", HERE); @@ -36,7 +39,7 @@ namespace out { m_io.DefineVariable("Step"); m_io.DefineVariable("Time"); - m_fname = title + (m_engine == "hdf5" ? ".h5" : ".bp"); + m_fname = title; } void Writer::addTracker(const std::string& type, @@ -412,33 +415,75 @@ namespace out { m_writer.Put(vare, xe_h); } - void Writer::beginWriting(std::size_t tstep, long double time) { + void Writer::beginWriting(WriteModeTags write_mode, + std::size_t tstep, + long double time) { + raise::ErrorIf(write_mode == WriteMode::None, "None is not a valid mode", HERE); raise::ErrorIf(p_adios == nullptr, "ADIOS pointer is null", HERE); - p_adios->ExitComputationBlock(); - if (m_writing_mode) { + if (m_active_mode != WriteMode::None) { raise::Fatal("Already writing", HERE); } - m_writing_mode = true; + m_active_mode = write_mode; try { - m_writer = m_io.Open(m_fname, m_mode); + std::string filename; + const std::string ext = m_engine == "hdf5" ? "h5" : "bp"; + if (m_separate_files) { + std::string mode_str; + if (m_active_mode == WriteMode::Fields) { + mode_str = "fields"; + } else if (m_active_mode == WriteMode::Particles) { + mode_str = "particles"; + } else if (m_active_mode == WriteMode::Spectra) { + mode_str = "spectra"; + } else { + raise::Fatal("Unknown write mode", HERE); + } + CallOnce( + [](auto& main_path, auto& mode_path) { + const std::filesystem::path main { main_path }; + const std::filesystem::path mode { mode_path }; + if (!std::filesystem::exists(main_path)) { + std::filesystem::create_directory(main_path); + } + if (!std::filesystem::exists(main / mode)) { + std::filesystem::create_directory(main / mode); + } + }, + m_fname, + mode_str); + filename = fmt::format("%s/%s/%s.%08lu.%s", + m_fname.c_str(), + mode_str.c_str(), + mode_str.c_str(), + tstep, + ext.c_str()); + m_mode = adios2::Mode::Write; + } else { + filename = fmt::format("%s.%s", m_fname.c_str(), ext.c_str()); + m_mode = std::filesystem::exists(filename) ? adios2::Mode::Append + : adios2::Mode::Write; + } + m_writer = m_io.Open(filename, m_mode); } catch (std::exception& e) { raise::Fatal(e.what(), HERE); } - m_mode = adios2::Mode::Append; m_writer.BeginStep(); m_writer.Put(m_io.InquireVariable("Step"), &tstep); m_writer.Put(m_io.InquireVariable("Time"), &time); } - void Writer::endWriting() { + void Writer::endWriting(WriteModeTags write_mode) { + raise::ErrorIf(write_mode == WriteMode::None, "None is not a valid mode", HERE); raise::ErrorIf(p_adios == nullptr, "ADIOS pointer is null", HERE); - if (!m_writing_mode) { + if (m_active_mode == WriteMode::None) { raise::Fatal("Not writing", HERE); } - m_writing_mode = false; + if (m_active_mode != write_mode) { + raise::Fatal("Writing mode mismatch", HERE); + } + m_active_mode = WriteMode::None; m_writer.EndStep(); m_writer.Close(); - p_adios->EnterComputationBlock(); } template void Writer::writeField(const std::vector&, diff --git a/src/output/writer.h b/src/output/writer.h index 566da44b2..a8abf4b12 100644 --- a/src/output/writer.h +++ b/src/output/writer.h @@ -36,6 +36,8 @@ namespace out { adios2::Engine m_writer; adios2::Mode m_mode { adios2::Mode::Write }; + bool m_separate_files; + // global shape of the fields array to output std::vector m_flds_g_shape; // local corner of the fields array to output @@ -63,7 +65,7 @@ namespace out { std::vector m_prtl_writers; std::vector m_spectra_writers; - bool m_writing_mode { false }; + WriteModeTags m_active_mode { WriteMode::None }; public: Writer() {} @@ -72,7 +74,7 @@ namespace out { Writer(Writer&&) = default; - void init(adios2::ADIOS*, const std::string&, const std::string&); + void init(adios2::ADIOS*, const std::string&, const std::string&, bool); void setMode(adios2::Mode); @@ -106,8 +108,8 @@ namespace out { void writeSpectrum(const array_t&, const std::string&); void writeSpectrumBins(const array_t&, const std::string&); - void beginWriting(std::size_t, long double); - void endWriting(); + void beginWriting(WriteModeTags, std::size_t, long double); + void endWriting(WriteModeTags); /* getters -------------------------------------------------------------- */ auto fname() const -> const std::string& { From 7c940f7586b95173546106622297b85bb9c01cec Mon Sep 17 00:00:00 2001 From: haykh Date: Mon, 11 Nov 2024 19:37:57 -0500 Subject: [PATCH 21/21] 0 particle case in timers --- src/global/utils/timer.cpp | 11 ++++++----- 1 file changed, 6 insertions(+), 5 deletions(-) diff --git a/src/global/utils/timer.cpp b/src/global/utils/timer.cpp index b5f4408ca..249feb6f5 100644 --- a/src/global/utils/timer.cpp +++ b/src/global/utils/timer.cpp @@ -116,11 +116,12 @@ namespace timer { (timer.second / local_tot) * 100.0); timer_stats.insert( { name, - std::make_tuple(timer.second, - timer.second / static_cast(npart), - timer.second / static_cast(ncells), - pcent, - 0u) }); + std::make_tuple( + timer.second, + npart > 0 ? timer.second / static_cast(npart) : 0.0, + timer.second / static_cast(ncells), + pcent, + 0u) }); } timer_stats.insert({ "Total", std::make_tuple(local_tot, 0.0, 0.0, 100u, 0u) }); #endif